#### /mj/HalfOctant4.pl

https://bitbucket.org/maphew/cahill-keyes
Perl | 1202 lines | 858 code | 91 blank | 253 comment | 179 complexity | 2d06a06c705e6d799fdf81a217d8321a MD5 | raw file
```   1#! /usr/bin/perl -w
2#	Program HalfOctant4. (made from HalfOctant3.odt on 2010-04-03)
3#
4#	In this version, all parallels are equidistant from equator to pole, along meridians.
5# There are no polar circular arcs. In other HalfOctant#.odt versions, dP was a hash, %dP,
6# with parallel spacings for polar, supple, temperate and torrid zones. In this version, it is
7# only an array, @dP, with a single parallel spacing value for each meridian.
8
9use Math::Trig;
10
11# Some constants for use by Subroutine MJtoG, which does coordinate axis transformation.
12# Angle of rotation is 60�.	Point G is (10000,0) in MJ, and (5000,0) in G
15\$yTranslate = 10000 * \$sin60;
16
17# Multidimensional arrays are done with hashes, with the keys used for i,j indeces.
18
19# Given input
20\$xM = 0;		\$yM = 0;	# Point M is the origin of the axes
21\$MA = 940;				# Indent of octant's apex, with respect to triangle's
22\$xG = 10000;	\$yG = 0;	# Point G, at center of base of octant
23\$DeltaFrigid = 104;			# Distance for each degree of latitude close to the pole
24
25# Other points and lengths of special interest
26\$xA = \$MA;		\$yA = 0;						# Point A, at the indent
27\$xN = \$xG;		\$yN = \$xG * tan (deg2rad(30));			# Point N
28(\$flag, \$xB, \$yB) = LineIntersection (\$xM,\$yM,30,\$xA,\$yA,45);	# Point B
29if (\$flag ne "Intersect") { die "\$flag while calculating point B.\n"; }
30\$AG = \$xG - \$xA;
31\$AB = Length (\$xA, \$yA, \$xB, \$yB);
32\$EF = \$AB;						# Sometimes I prefer AB, others EF
33\$MB = Length (\$xM, \$yM, \$xB, \$yB);
34\$MN = Length (\$xM, \$yM, \$xN, \$yN);
35\$xD = Interpolate (\$MB, \$MN, \$xN, \$yM);	# D is away from N as B is away from M
36\$yD = Interpolate (\$MB, \$MN, \$yN, \$yM);
37\$xF = \$xG;
38\$yF = \$yN - \$MB;
39\$xE = \$xN - \$MA * sin (deg2rad(30));
40\$yE = \$yN - \$MA * cos (deg2rad(30));
41\$FG = \$yF;
42\$BD = Length (\$xB, \$yB, \$xD, \$yD);
43\$EFG = \$EF + \$FG;
44
45# Calculate polar start, two joints, and equatorial end, for each meridian; hashes %xJ, %yJ.
46# Also calculate a few other values for other hashes, on the way.
47
48\$DeltaMEq = \$EFG / 45;	# 45 meridian spacings along equator for half an octant
49
50foreach \$m (0 .. 45) {	# Loop for each meridian.
51	\$M = \$m . ',' ;	# part of hash key (array row number)
52
53	# Start point of minor meridians will be modified to be at parallel 85� after that
54	# parallel has been calculated. For now, all meridians are set to start at point A.
55	\$xJ {\$M . 0} = \$xA;		\$yJ {\$M . 0} = \$yA;
56
57	# End point of meridians, that is, equidistant intersections with equator;
58	# These are also the meridian crossings for the zero-th parallel
59	\$distance = \$m * \$DeltaMEq;
60	if (\$distance <= \$yF) {
61		\$xJ {\$M . 3} = \$xP {\$M . 0} = \$xG;
62		\$yJ {\$M . 3} = \$yP {\$M . 0} = \$distance;
63	} else {	# Past joint at point F
64		\$extra = \$distance - \$yF;
65		\$xJ {\$M . 3} = \$xP {\$M . 0} = Interpolate (\$extra, \$EF, \$xF, \$xE);
66		\$yJ {\$M . 3} = \$yP {\$M . 0} = Interpolate (\$extra, \$EF, \$yF, \$yE);
67	}
68
69	# Calculate both joints of meridians
70
71	if (\$m == 0) {	# Straight, centre of octant meridian
72		# This is a good place to calculate a bunch of things for centre line, AG.
73
74		# No joints; give coordinates of point G for the joints
75		@xJ {\$M . 1, \$M . 2} = (\$xG, \$xG);
76		@yJ {\$M . 1, \$M . 2} = (\$yG, \$yG);
77		# Lenght of meridian segments; this whole meridian is a single segment
78		@L {\$M . 0, \$M . 1, \$M . 2, \$M . 3, \$M . 4} = (\$AG, 0, 0, \$AG, \$AG);
79		# May as well calculate parallel intersections along this straight line
80		\$dP [\$m] = \$L{\$M . 4} / 90;
81		foreach \$p (1 .. 89) {
82			\$xP {\$M . \$p} = \$xG - \$p * \$dP[\$m];
83			\$yP {\$M . \$p} = \$yG;
84		}
85	} else {	# Meridians other than 0�.
86		# Calculate joints by line intersection from A and M
87		(\$flag, \$xJ {\$M . 1}, \$yJ {\$M . 1}) =
88			LineIntersection (\$xM, \$yM, 2*\$m/3, \$xA, \$yA, \$m);
89		if (\$flag ne "Intersect") { die "\$flag while calculating 1st joint at \$m.\n"; }
90		(\$flag, \$xJ {\$M . 2}, \$yJ {\$M . 2}) =
91			LineIntersection (\$xM, \$yM, 2*\$m/3, \$xJ {\$M . 3}, \$yJ {\$M . 3}, \$m/3);
92		if (\$flag ne "Intersect") { die "\$flag while calculating 2nd joint at \$m.\n"; }
93		# Calculate length of meridian segments: These are calculated from point A
94		# to equator, even if this meridian starts at parallel 85�.
95		\$L {\$M . 0} = Length (\$xA, \$yA, \$xJ {\$M . 1}, \$yJ {\$M . 1});
96		\$L {\$M . 1} = Length (\$xJ {\$M . 1}, \$yJ {\$M . 1}, \$xJ {\$M . 2}, \$yJ {\$M . 2});
97		\$L {\$M . 2} = Length (\$xJ {\$M . 2}, \$yJ {\$M . 2}, \$xJ {\$M . 3}, \$yJ {\$M . 3});
98		\$L {\$M . 3} = \$L {\$M . 1} + \$L {\$M . 2};	# Torrid + Temperate
99		\$L {\$M . 4} = \$L {\$M . 3} + \$L {\$M . 0};	# Total meridian length from A
100		# Calculate parallel intersections
101		\$dP[\$m] = \$L {\$M . 4} / 90;	# Distance between parallels
102		foreach \$p (1 .. 89) {
103			\$distance = \$p * \$dP [\$m];
104			if (\$distance <= \$L {\$M . 2} ) {		# Torrid segment
105				\$xP{\$M.\$p} =
106					Interpolate (\$distance, \$L{\$M.2}, \$xJ{\$M.3}, \$xJ{\$M.2});
107				\$yP{\$M.\$p} =
108					Interpolate (\$distance, \$L{\$M.2}, \$yJ{\$M.3}, \$yJ{\$M.2});
109			} elsif (\$distance <= \$L {\$M . 3} ) {	# Temperate seg., past one joint
110				\$extra = \$distance - \$L {\$M . 2};
111				\$xP{\$M.\$p} =
112					Interpolate (\$extra, \$L{\$M.1}, \$xJ{\$M.2}, \$xJ{\$M.1});
113				\$yP{\$M.\$p} =
114					Interpolate (\$extra, \$L{\$M.1}, \$yJ{\$M.2}, \$yJ{\$M.1});
115			} else {					# Frigid segment, past two joints
116				\$extra = \$distance - \$L {\$M . 3};
117				\$xP{\$M.\$p} = Interpolate (\$extra, \$L{\$M.0}, \$xJ{\$M.1}, \$xA);
118				\$yP{\$M.\$p} = Interpolate (\$extra, \$L{\$M.0}, \$yJ{\$M.1}, \$yA);
119			}
120		}	# End of foreach \$p loop, that is, end of parallel calculation
121	}	# End of calculating meridian joints
122}	# End of loop for each meridian.
123
124#	Start point of every meridian was previously set to be point A. However, the drawn
125#	line for the minor meridians will start at parallel 85�; this could not be set previously
126#	because that parallel still hadn't been calculated.
127foreach \$m (0 .. 45) {
128	\$M = \$m . ',' ;
129	unless (\$m % 5 == 0) {
130		\$xJ {\$M . 0} = \$xP {\$M . 85};
131		\$yJ {\$M . 0} = \$yP {\$M . 85};
132	}
133}
134
135# Output all the values.
136\$PointFormat = "%3s \t%11.4f \t%11.4f\n";
137\$LengthFormat = "%s \t%11.4f\n";
138\$HashKeyValue = "%4s \t%11.4f\n";
139
140\$Skip = "Yes";
141unless (\$Skip) {	# Option to skip printing hashes in spreadsheet readable csv format
142	open (CSV, ">Macros4\/Hashes4.csv");
143	print CSV "\ndP\n";		# @dP is an array, not a hash; can't use sub PrintHash
144	foreach \$i (0 .. 45) { printf CSV "\t%10d",\$i; } print CSV "\n";
145	print CSV " ";
146	foreach \$m (0 .. 45) {
147		if (defined (\$dP[\$m]) ) {printf CSV "\t%11.4f",\$dP[\$m];
148		} else { printf CSV "\t undef"; }
149	}
150	print CSV "\n";
151	print CSV "\nL\n"; PrintHash (4,%L);
152	print CSV "\nxJ\n"; PrintHash (3,%xJ);
153	print CSV "\nyJ\n"; PrintHash (3,%yJ);
154	print CSV "\nxP\n"; PrintHash (89,%xP);
155	print CSV "\nyP\n"; PrintHash (89,%yP);
156	close (CSV);
157	print 'Wrote values of arrays to file "Macros4/Hashes4.csv".', "\n";
158}	# End of Skip printing hashes in Spreadsheet readable csv format
159
160\$Skip = "Yes";
161unless (\$Skip) {	# Option to skip printing, hashes by their keys
162	print "\n* * * Keys of hashes\n\n";
163	print "* * * \%xJ\n";
164	foreach \$key (keys (%xJ)) {
165		printf (\$HashKeyValue, \$key, \$xJ{\$key});
166	}
167	print "* * * \%yJ\n";
168	foreach \$key (keys (%yJ)) {
169		printf (\$HashKeyValue, \$key, \$yJ{\$key});
170	}
171	print "* * * \%L\n";
172	foreach \$key (keys (%L)) {
173		printf (\$HashKeyValue, \$key, \$L{\$key});
174	}
175	print "\n* * * COORDINATES\n";
176	print "* * * \%xP\n";
177	foreach \$key (keys (%xP)) {
178		printf (\$HashKeyValue, \$key, \$xP{\$key});
179	}
180	print "* * * \%yP\n";
181	foreach \$key (keys (%yP)) {
182		printf (\$HashKeyValue, \$key, \$yP{\$key});
183	}
184}	# End of Skip of printing hashes by their keys
185
186\$Skip = "Yes";
187unless (\$Skip) {	# Option to skip printing points and lengths
188print "* * * Major points (point, x, y);\n";
189printf (\$PointFormat, "M", \$xM, \$yM);
190printf (\$PointFormat, "N", \$xN, \$yN);
191printf (\$PointFormat, "A", \$xA, \$yA);
192printf (\$PointFormat, "B", \$xB, \$yB);
193printf (\$PointFormat, "D", \$xD, \$yD);
194printf (\$PointFormat, "E", \$xE, \$yE);
195printf (\$PointFormat, "F", \$xF, \$yF);
196printf (\$PointFormat, "G", \$xG, \$yG);
197# Main lengths
198print "\n* * * Main lengths\n";
199printf (\$LengthFormat,'MG',\$xG);
200printf (\$LengthFormat,'MN',\$MN);
201printf (\$LengthFormat,'MA',\$MA);
202printf (\$LengthFormat,'MB',\$MB);
203printf (\$LengthFormat,'AG',\$AG);
204printf (\$LengthFormat,'AB',\$AB);
205printf (\$LengthFormat,'BD',\$BD);
206printf (\$LengthFormat,'EF',\$EF);
207printf (\$LengthFormat,'FG',\$FG);
208printf (\$LengthFormat,'EFG',\$EFG);
209printf (\$LengthFormat,'ABDE',2 * \$EFG);
210print "\n";
211printf (\$LengthFormat,'Frigid parallel spacing',\$DeltaFrigid);
212printf (\$LengthFormat,'Equatorial meridian spacing',\$DeltaMEq);
213}	# End of Skip of printing major points and lengths
214
216# Scaffold triangle macro
217
218\$Skip = "Yes";
219unless (\$Skip) {	# Option to skip making macro files of half octant or one octant
220StartMacroFile (MACRO, "Macros4\/Triangle.txt");
221print MACRO ("Pen 9\r");
222print MACRO ("Line ", Coordinate(\$xM, \$yM), Coordinate (\$xN, \$yN), "\\\r");
223print MACRO (Coordinate(\$xG, \$yG), Coordinate (\$xM, \$yM), '[`;`]', "\r");
224EndMacroFile (MACRO);
225print 'Wrote file "Macros4/Triangle.txt".', "\n";
226
227# And the same thing in Gene's coordinate system
228StartMacroFile (MACRO, "Macros4\/TriangleG.txt");
229print MACRO ("Pen 9\r");
230print MACRO ("Line ", Coordinate(MJtoG(\$xM,\$yM)), Coordinate(MJtoG(\$xN,\$yN)), "\\\r");
231print MACRO (Coordinate(MJtoG(\$xG,\$yG)), Coordinate(MJtoG(\$xM,\$yM)), '[`;`]', "\r");
232EndMacroFile (MACRO);
233print 'Wrote file "Macros4/TriangleG.txt."', "\n";
234
235# Half-octant outline macro
236StartMacroFile (MACRO, "Macros4\/Half-Octant.txt");
237print MACRO ("Pen 12\r");
238print MACRO ("Line ", Coordinate(\$xA, \$yA), Coordinate (\$xB, \$yB), "\\\r");
239print MACRO (Coordinate(\$xD, \$yD), Coordinate (\$xE, \$yE), "\\\r");
240print MACRO (Coordinate(\$xF, \$yF), Coordinate (\$xG, \$yG), '[`;`]', "\r");
241EndMacroFile (MACRO);
242print 'Wrote file "Macros4/Half-Octant.txt".', "\n";
243
244# And the same thing in Gene's coordinate system
245StartMacroFile (MACRO, "Macros4\/Half-OctantG.txt");
246print MACRO ("Pen 12\r");
247print MACRO ("Line ", Coordinate(MJtoG(\$xA,\$yA)), Coordinate(MJtoG(\$xB, \$yB)), "\\\r");
248print MACRO (Coordinate(MJtoG(\$xD,\$yD)), Coordinate(MJtoG(\$xE, \$yE)), "\\\r");
249print MACRO (Coordinate(MJtoG(\$xF,\$yF)), Coordinate(MJtoG(\$xG, \$yG)), '[`;`]', "\r");
250EndMacroFile (MACRO);
251print 'Wrote file "Macros4/Half-OctantG.txt".', "\n";
252
253# Macro of 8 equilateral triangles
254StartMacroFile (MACRO, "Macros4\/MScaffold.txt");
255print MACRO ("Pen 9\r");
256MLinesMacro (MACRO,0,\$xN,-\$yN,1,\$xM,\$yM,1,\$xN,\$yN,1,\$xN,-\$yN,1,\$xM,\$yM,6,
257	\$xN,\$yN,1,\$xN,\$yN,2,\$xM,\$yM,1);
258MLinesMacro (MACRO,0,\$xN,\$yN,1,\$xM,\$yM,7,\$xN,\$yN,2,\$xN,\$yN,3,\$xM,\$yM,7);
259MLinesMacro (MACRO,0,\$xN,\$yN,2,\$xM,\$yM,3,\$xN,\$yN,3,\$xN,\$yN,4,\$xM,\$yM,3);
260MLinesMacro (MACRO,0,\$xN,\$yN,3,\$xM,\$yM,5,\$xN,\$yN,4);
261EndMacroFile (MACRO);
262print 'Wrote file "Macros4/MScaffold.txt".', "\n";
263
264# Macro of 8 Octants
265StartMacroFile (MACRO, "Macros4\/MOctants.txt");
266print MACRO ("Pen 12\r");
267MLinesMacro (MACRO,6,	\$xF,\$yF,\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD,\$xB,-\$yB,
268				\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,\$xF,\$yF);
269MLinesMacro (MACRO,1,	\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD,\$xB,-\$yB,
270				\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,\$xF,\$yF);
271MLinesMacro (MACRO,2,	\$xB,-\$yB,\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,
272				\$xF,\$yF,\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD);
273MLinesMacro (MACRO,7,	\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD,\$xB,-\$yB,
274				\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,\$xF,\$yF);
275MLinesMacro (MACRO,8,	\$xB,-\$yB,\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,
276				\$xF,\$yF,\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD);
277MLinesMacro (MACRO,3,	\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD,\$xB,-\$yB,
278				\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,\$xF,\$yF);
279MLinesMacro (MACRO,5,	\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD,\$xB,-\$yB,
280				\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,\$xF,\$yF);
281MLinesMacro (MACRO,4,	\$xB,-\$yB,\$xA,\$yA,\$xB,\$yB,\$xD,\$yD,\$xE,\$yE,
282				\$xF,\$yF,\$xF,-\$yF,\$xE,-\$yE,\$xD,-\$yD);
283EndMacroFile (MACRO);
284print 'Wrote file "Macros4/MOctants.txt".', "\n";
285
286# Meridians macros
287
288StartMacroFile (MAJOR, "Macros4\/MajorMeridians.txt");
289print MAJOR ("Pen 12\r");
290StartMacroFile (MINOR, "Macros4\/MinorMeridians.txt");
291print MINOR ("Pen 13\r");
292MeridianMacro ();
293EndMacroFile (MAJOR);
294EndMacroFile (MINOR);
295print 'Wrote files "Macros4/MajorMeridians.txt" and "Macros4/MinorMeridians.txt".', "\n";
296
297# Parallels macros
298StartMacroFile (MAJOR, "Macros4\/ParallelsMajor.txt");
299print MAJOR ("Pen 12\r");
300StartMacroFile (MINOR, "Macros4\/ParallelsMinor.txt");
301print MINOR ("Pen 13\r");
302ParallelMacro (1,89,0,45);
303EndMacroFile (MAJOR);
304EndMacroFile (MINOR);
305print 'Wrote files "Macros4/ParallelsMajor.txt" and "Macros4/ParallelsMinor.txt".', "\n";
306}	# End of Skip of making macro files of half octant or one octant
307
308\$Skip = "Yes";
309unless (\$Skip) {	# Option to skip making macro files of parallel points and of joints points
310# Macro of points for the intersections of the parallels with the meridians; this will be
311# a large file which will take quite awhile to plot; pens to be used are specified in the
312# subroutine. Meant only to check if the data seems right.
313StartMacroFile (MACRO, "Macros4\/Points.txt");
314PointsMacro ();
315EndMacroFile (MACRO);
316print 'Wrote file "Macros4/Points.txt".', "\n";
317
318# Macro of meridian joints and end points (4 points per meridian)
319StartMacroFile (MACRO, "Macros4\/Joints.txt");
320JointsMacro ();
321EndMacroFile (MACRO);
322print 'Wrote file "Macros4/Joints.txt".', "\n";
323}	# End of Skip of making macro files
324
325# Macros for all eight whole octants
326
327\$Skip = "Yes";
328unless (\$Skip) {	# Option to skip making macro files of eight octants
329# Make arrays for Octant; I should have done this much earlier in the program
330@xOct = (\$xA,\$xB,\$xD,\$xE,\$xF,\$xG);
331@yOct = (\$yA,\$yB,\$yD,\$yE,\$yF,\$yG);
332
333# Make a copy of arrays and hashes that I will want to convert; the other ones
334# will be changed to coordinates for different octant halves.
335@xOctC = @xOct;	@yOctC = @yOct;
336%xJC = %xJ;		%yJC = %yJ;
337%xPC = %xP;	%yPC = %yP;
338
339# Make grid files
340StartMacroFile (MINOR, "Macros4\/M-Grid-Minor.txt");
341print MINOR ("Pen 15\r");
342GridMacro (MINOR, -20000,20000,-10000,10000,100,1000);
343EndMacroFile (MINOR);
344
345StartMacroFile (MAJOR, "Macros4\/M-Grid-Major.txt");
346print MAJOR ("Pen 14\r");
347GridMacro (MAJOR, -20000,20000,-10000,10000,1000,0);
348EndMacroFile (MAJOR);
349
350StartMacroFile (MAJOR, "Macros4\/M-Grid-5Kred.txt");
351print MAJOR ("Pen 10\r");
352GridMacro (MAJOR, -20000,20000,-10000,10000,5000,0);
353EndMacroFile (MAJOR);
354
355print 'Wrote files "M-Grid-Minor.txt", "M-Grid-Major.txt" and "M-Grid-5Kred.txt" in "Macros4/".',
356	"\n";
357# Draw octants, half at a time
358StartMacroFile (MACRO, "Macros4\/M-Octants.txt");
359print MACRO ("Pen 11\r");
360foreach \$Oct (1 .. 8) {		# For octant 1 and octant 2
361	foreach \$i (0..5) {		# For each point (A,B,D,E,F,G)
362		(\$xOct[\$i],\$yOct[\$i]) = MJtoG (\$xOctC[\$i],-\$yOctC[\$i],\$Oct);	# West half
363	}
364	OctantMacro ();
365	foreach \$i (0..5) {
366		(\$xOct[\$i],\$yOct[\$i]) = MJtoG (\$xOctC[\$i],\$yOctC[\$i],\$Oct);	# East half
367	}
368	OctantMacro ();
369}
370EndMacroFile (MACRO);
371print 'Wrote file "Macros4/M-Octants.txt".', "\n";
372
373# Draw Meridians, half an octant at a time
374StartMacroFile (MINOR, "Macros4\/M-Meridians-Minor.txt");
375StartMacroFile (MAJOR, "Macros4\/M-Meridians-Major.txt");
376print MINOR ("Pen 13\r");
377print MAJOR ("Pen 12\r");
378
379foreach \$Oct (1 .. 8) {		# For each octant
380	# West half of octant
381	foreach \$m (0..45) {
382		\$M = \$m . ',' ;
383		foreach \$i (0..3) {
384			(\$xJ{\$M.\$i},\$yJ{\$M.\$i}) = MJtoG (\$xJC{\$M.\$i},-\$yJC{\$M.\$i},\$Oct);
385		}
386	}
387	MeridianMacro ();
388
389	# East half of octant
390	foreach \$m (0..45) {
391		\$M = \$m . ',' ;
392		foreach \$i (0..3) {
393			(\$xJ{\$M.\$i},\$yJ{\$M.\$i}) = MJtoG (\$xJC{\$M.\$i},\$yJC{\$M.\$i},\$Oct);
394		}
395	}
396	MeridianMacro ();
397}
398EndMacroFile (MINOR);
399EndMacroFile (MAJOR);
400print 'Wrote files "Macros4/M-Meridians-Minor.txt" and "Macros4/M-Meridians-Major.txt".', "\n";
401
402# Draw Parallels, half an octant at a time
403StartMacroFile (MINOR, "Macros4\/M-Parallels-Minor.txt");
404StartMacroFile (MAJOR, "Macros4\/M-Parallels-Major.txt");
405print MINOR ("Pen 13\r");
406print MAJOR ("Pen 12\r");
407
408foreach \$Oct (1 .. 8) {		# For each octant
409	# West half of octant
410	foreach \$m (0..45) {	# Convert coordinates
411		\$M = \$m . ',' ;
412		foreach \$p (0..89) {
413			(\$xP{\$M.\$p},\$yP{\$M.\$p}) = MJtoG (\$xPC{\$M.\$p},-\$yPC{\$M.\$p},\$Oct);
414		}
415	}
416	# All parallels from 1� to 89� of latitude; skip parallel 0� &#x2013; octant boundary
417	ParallelMacro (1, 89, 0, 45, ", Octant \$Oct West");
418
419	# East half of octant
420	foreach \$m (0..45) {	# Convert coordinates
421		\$M = \$m . ',' ;
422		foreach \$p (0..89) {
423			(\$xP{\$M.\$p},\$yP{\$M.\$p}) = MJtoG (\$xPC{\$M.\$p},\$yPC{\$M.\$p},\$Oct);
424		}
425	}
426	# All parallels from 1� to 89� of latitude; skip parallel 0� &#x2013; octant boundary
427	ParallelMacro (1, 89, 0, 45, ", Octant \$Oct East");
428}
429EndMacroFile (MINOR);
430EndMacroFile (MAJOR);
431print 'Wrote file "Macros4/M-Parallels-Minor.txt" and "Macros4/M-Parallels-Major.txt".', "\n";
432
433# Draw Parallels, half an octant at a time, but connecting only nodes of major meridians;
434# that is, not connecting nodes of minor meridians
435StartMacroFile (MINOR, "Macros4\/M-Parallels-Minor5.txt");
436StartMacroFile (MAJOR, "Macros4\/M-Parallels-Major5.txt");
437print MINOR ("Pen 13\r");
438print MAJOR ("Pen 12\r");
439
440foreach \$Oct (1 .. 8) {		# For each octant
441	# West half of octant
442	foreach \$m (0..45) {	# Convert coordinates
443		\$M = \$m . ',' ;
444		foreach \$p (0..89) {
445			(\$xP{\$M.\$p},\$yP{\$M.\$p}) = MJtoG (\$xPC{\$M.\$p},-\$yPC{\$M.\$p},\$Oct);
446		}
447	}
448	# All parallels from 1� to 89� of latitude; skip parallel 0� &#x2013; octant boundary
449	ParallelMacro5 (1, 89, 0, 45, ", Octant \$Oct West");
450
451	# East half of octant
452	foreach \$m (0..45) {	# Convert coordinates
453		\$M = \$m . ',' ;
454		foreach \$p (0..89) {
455			(\$xP{\$M.\$p},\$yP{\$M.\$p}) = MJtoG (\$xPC{\$M.\$p},\$yPC{\$M.\$p},\$Oct);
456		}
457	}
458	# All parallels from 1� to 89� of latitude; skip parallel 0� &#x2013; octant boundary
459	ParallelMacro5 (1, 89, 0, 45, ", Octant \$Oct East");
460}
461EndMacroFile (MINOR);
462EndMacroFile (MAJOR);
463print 'Wrote file "Macros4/M-Parallels-Major5.txt" and "Macros4/M-Parallels-Minor5.txt".', "\n";
464
465}	# End of skip of making macro files of eight octants
466
467# Read some Coastal Data, convert to M-map coordinates, and output a TurboCAD macro.
468# File read is of MAPGEN data format: two columns ASCII flat file with: longitude tab latitude;
469# at the start of each segment there is a line containing only "# -b".
470# This particular file covers longitude from -67� to -59� and latitude from 43� to 48�, that is,
471# Nova Scotia's area. File is 5781.dat.gz, is zipped, and has 92,227 lines.
472
473\$Skip = "Yes";
474unless (\$Skip) {	# Option to skip making macro files of Coastal Data
475print "Going to read land data.\n";
476\$maxX=-99999;	\$maxY=-99999;	\$minX=999999;	\$minY=999999;
477StartMacroFile (LAND, "Macros4\/NSmacro.txt");
478print LAND ("Pen 0\r");
479open (DATA, "zcat 5781.dat.gz | ");
480\$nData = 0;
481while (<DATA>) {
482	\$Line = \$_ ;
483	chomp(\$Line);
484	\$nData ++;
485	if (\$nData % 1000 == 0) { print "Read \$nData lines of land data so far.\n"; }
486	if (\$Line ne "\# -b") {
487		(\$Long,\$Lat) = split(/\t/,\$Line);
488		# Convert to M-map x and y coordinates
489		(\$x,\$y) = Real2Map (\$Long,\$Lat);
490		# Keep track of maximum and minimum values
491		if (\$x < \$minX) {\$minX = \$x;}
492		if (\$x > \$maxX) {\$maxX = \$x;}
493		if (\$y < \$minY) {\$minY = \$y;}
494		if (\$y > \$maxY) {\$maxY = \$y}
495		# Make arrays of points for this segment to be used by LandMacro
496		push (@Xs,\$x);
497		push (@Ys,\$y);
498	} elsif (defined(@Xs) ) {	# Do only if data points have already been read
499		LandMacro ();	# Draw the segment of boundary with arrays @Xs and @Ys
500		@Xs = ();	@Ys = ();	# Empty arrays, ready for next segment
501	}
502}
503close (DATA);
504if (\$Line ne "\# -b") { # Draw last segment, if it hasn't been drawn
505	LandMacro ();
506}
507EndMacroFile (LAND);
508print "Read a total of \$nData lines, and wrote file",' "Macros4/NSmacro.txt".', "\n";
509print \$minX, ' <= x <= ',\$maxX," and ",\$minY, ' <= y <= ', \$maxY, "\n";
510
511}	# End of skipping making macro for Coastal Data
512
513\$Skip = "";
514unless (\$Skip) {	# Option to skip making macro file of NS 1:2M of Coastal Data
515print "Going to read land data.\n";
516\$maxX=-99999;	\$maxY=-99999;	\$minX=999999;	\$minY=999999;
517StartMacroFile (LAND, "Macros4\/NS-2M-macro4.txt");
518print LAND ("Pen 0\r");
519open (DATA, "zcat NS-2M.dat.gz | ");
520\$nData = 0;
521while (<DATA>) {
522	\$Line = \$_ ;
523	chomp(\$Line);
524	\$nData ++;
525	if (\$nData % 1000 == 0) { print "Read \$nData lines of land data so far.\n"; }
526	if (\$Line ne "\# -b") {
527		(\$Long,\$Lat) = split(/\t/,\$Line);
528		# Convert to M-map x and y coordinates
529		(\$x,\$y) = Real2Map (\$Long,\$Lat);
530		# Keep track of maximum and minimum values
531		if (\$x < \$minX) {\$minX = \$x;}
532		if (\$x > \$maxX) {\$maxX = \$x;}
533		if (\$y < \$minY) {\$minY = \$y;}
534		if (\$y > \$maxY) {\$maxY = \$y}
535		# Make arrays of points for this segment to be used by LandMacro
536		push (@Xs,\$x);
537		push (@Ys,\$y);
538	} elsif (defined(@Xs) ) {	# Do only if data points have already been read
539		# Shift coordinates to see if TurboCAD plots them better
540		\$dx = -7600;
541		\$dy = 3700;
542		\$n = @Xs -1;
543		foreach \$i (0..\$n) {
544			\$Xs[\$i] = \$Xs[\$i] - \$dx;
545			\$Ys[\$i] = \$Ys[\$i] - \$dy;
546		}
547		LandMacro ();	# Draw the segment of boundary with arrays @Xs and @Ys
548		@Xs = ();	@Ys = ();	# Empty arrays, ready for next segment
549	}
550}
551close (DATA);
552if (\$Line ne "\# -b") { # Draw last segment, if it hasn't been drawn
553	LandMacro ();
554}
555EndMacroFile (LAND);
556print "Read a total of \$nData lines, and wrote file",' "Macros4/NS-2M4-macro .txt".', "\n";
557print \$minX, ' <= x <= ',\$maxX," and ",\$minY, ' <= y <= ', \$maxY, "\n";
558print "Need plot area ", \$maxX-\$dx, " by ", \$maxY-\$dy,".\n";
559
560}	# End of skipping making macro for NS, 2:M Coastal Data
561
562
563# -  -  -  -  -  -  -  -  -    SUBROUTINES         -  -  -  -  -  -  -  -  -  -  -  -  -
564sub PrintHash {
565	# Prints one of my two-dimension arrays in the form of a hash; it is assumed that
566	# the indices are separated by a comma, the first indece varies from 0 to 45 and the
567	# second indice varies from 0 to the first input parameter.
568	# Values are printed to filehandle CSV.
569	my (\$nI, %Hash) = @_;
570	my (\$i,\$m);
571	foreach \$i (0 .. 45) { printf CSV "\t%10d",\$i; } print CSV "\n";
572	foreach \$i (0 .. \$nI) {
573		print CSV \$i;
574		foreach \$m (0 .. 45) {
575			if (defined (\$Hash{\$m.",".\$i})) {printf CSV "\t%11.4f",\$Hash{\$m.",".\$i};
576			} else { printf CSV "\t undef"; }
577		}
578		print CSV "\n";
579	}
580}	# End of subroutine Hash
581
582sub LineIntersection {	# Written 2010-02-28
583# Subroutine/function to calculate coordinates of point of intersection of two lines which
584# are given by a point on the line and the line's slope angle in degrees.
585#
586# Return is an array of either one or three values:
587#	- if the lines don't intersect or are the same, the return is "Parallel";
588#	- if the lines intersect, the returns are: "Intersect" to mean that it worked, and
589#	the x and y coordinates of the point of intersection.
590#	- if not enough arguments were given or the angles were too large, then return is
591#	 "Error".
592#
593# Arguments should be 6, in this order:
594#	x and y coordinates of point of first line; slope of first line in degrees;
595#	x and y coordinates of point of second line; slope of second line in degrees;
596#
597# Equations used are from: slope of line = tangent angle = delta-y / delta-x, and the fact
598#	 that intersection point x,y is on both lines.
599
600use Math::Trig;
601my (\$nArguments,\$xp,\$yp,\$m1,\$m2);
602\$nArguments=@_ ;
603my (\$x1,\$y1,\$angle1,\$x2,\$y2,\$angle2) = @_ ;
604
605if (\$nArguments != 6) {	# Not the exact number of arguments
606	print "Sub LineIntersection requires 6 arguments but received \$nArguments \n";
607	return "Error";
608} elsif (\$angle1<-180 || \$angle2<-180 || \$angle1>180 || \$angle2>180) { # Angles too large
609	print "Sub LineIntersection: Angles should be between -180 and 180. \n";
610	return "Error";
611} elsif (\$angle1==\$angle2 || \$angle1 == (\$angle2+180) || \$angle1 == (\$angle2-180) ) {
612	return "Parallel";		# Lines are either parallel or the same line
613} elsif (	not defined(\$x1) || not defined(\$y1) || not defined(\$angle1) ||
614		not defined(\$x2) || not defined(\$y2) || not defined(\$angle2) ) {
615	print ("Undefined in LineInterpolation: \$x1,\$y1,\$angle1,\$x2,\$y2,\$angle2\n");
616	return "Error";
617} else {	# Lines intersect. Calculate point of intersection.
618	if (\$angle1 == 0 || \$angle1 == 180 || \$angle1 == -180) {		# Line 1 horizontal
619		if (\$angle2 == 90 || \$angle2 == -90) {	# Line 1 horizontal, line 2 vertical
620			return ("Intersect",\$y1,\$x2);
621		} else {					# Line 1 horizontal, line 2 slanted
622			\$xp = \$x2 - (\$y2 - \$y1) / tan(deg2rad(\$angle2));
623			return ("Intersect",\$xp,\$y1);
624		}
625	} elsif (\$angle1 == 90 || \$angle1 == -90) {	# Calculate when line 1 is vertical
626		if (\$angle2 == 0 || \$angle2 == 180 || \$angle2 == -180) {
627								# Line 1 vertical, line 2 horizontal
628			return ("Intersect",\$x1,\$y2);
629		} else {					# Line 1 vertical, line 2 slanted
630			\$yp = \$y2 - (\$x2 - \$x1) * tan(deg2rad(\$angle2));
631			return ("Intersect",\$x1,\$yp);
632		}
633	} elsif (\$angle2 == 0 || \$angle2 == 180 || \$angle2 == -180) {
634								# Line 1 slanted, line 2 horizontal
635			\$xp = \$x1 - (\$y1 - \$y2) / tan(deg2rad(\$angle1));
636			return ("Intersect",\$xp,\$y2);
637	} elsif (\$angle2 == 90 or \$angle2 == -90) {	# Line 1 slanted, line 2 vertical
638			\$yp = \$y1 - (\$x1-\$x2) * tan(deg2rad(\$angle1));
639			return ("Intersect",\$x2,\$yp);
640	} else {						# Line 1 slanted, line 2 slanted
643			\$xp = (\$m1 * \$x1 - \$m2 * \$x2 - \$y1 + \$y2) / (\$m1 - \$m2);
644			\$yp = \$m1 * \$xp - \$m1 * \$x1 + \$y1;
645			return ("Intersect",\$xp,\$yp);
646	}	#End of all the vertical, horizontal, slanted combinations
647}	# End of if-else &#x2013; lines intersect
648} # End of Subroutine LineIntersection
649
650sub Length {
651	# Input are x1,y1,x2,y2
652	my (\$x1,\$y1,\$x2,\$y2) = @_;
653	return sqrt( (\$x1-\$x2)**2 + (\$y1-\$y2)**2);
654}
655
656sub Interpolate {
657# Inputs are 4: length wanted, of total-segment-length, start, end;
658# Total-segment-length is different from (end - start); end and start may be x-coordinates,
659# or y-coordinates, while length takes into account the other coordinates.
660# (Start - End ) / Length = (Wanted - Start) / NewLength
661# Returns single value: Wanted.
662	my (\$NewLength, \$Length, \$Start, \$End) = @_;
663	my (\$Wanted);
664	\$Wanted = \$Start + (\$End - \$Start) * \$NewLength / \$Length;
665\$Skip = "Yes";
666unless (\$Skip) {	# Skip of not the printing of arguments and value being returned
667	my (\$nIn, \$i);
668	\$nIn = @_ - 1;
669#	print "Interpolate inputs are ", \$nIn+1, ": ";
670	foreach \$i (0 .. \$nIn) {
671		if (defined (\$_[\$i]) ) {
672			printf ("%11.4f, ", \$_[\$i]);
673		} else {
674			printf ("%11s, ", 'UNDEFINED!');
675		}
676	}
677	if (defined (\$Wanted) ) {
678#		print "Returning: ";
679		printf ("\>%11.4f",\$Wanted);
680	} else {
681		printf ("\>%11s", 'UNDEFINED!');
682	}
683	print "\n";
684}	# End of Skip printing arguments and value being returned
685	return \$Wanted;
686}	# End of subroutine Interpolate
687
688# -  -  -  -  -  -  -  -  -    SUBROUTINES relating to TurboCAD    -  -  -  -  -  -  -  -  -  -  -  -  -
689
690sub StartMacroFile {
691	my (\$FileHandle, \$ToOpen);
692	\$FileHandle = \$_ ;
693	\$ToOpen = '>' . \$_;
694	open (\$FileHandle, \$ToOpen);
695	print \$FileHandle ('#PenWas = %WT.P', "\r");
696}
697
698sub EndMacroFile {
699	my (\$FileHandle);
700	\$FileHandle = \$_ ;
701	print \$FileHandle ('Pen #PenWas', "\r");
702	print \$FileHandle ("Bell\rBell\r");
703	close (\$FileHandle);
704}
705
706sub Coordinate {
707	# returns string with coordinate format for TurboCAD;
708	# expects two arguments: x and y.
709	return '[a,' . sprintf ("%11.4f,%11.4f",\$_,\$_) . ',0,` `] '
710}
711
712sub PointsMacro {
713	# I will be lazy and use global variables instead of trying to pass hashes as
714	# parameters.
715	# Plotting all parallel intersections, that have already been calculated: arrays
716	# %xP, %yP.
717	# Subroutine will prepare TurboCad commands to plot each point. This will be written
718	# to filehandle MACRO
719	my (\$m, \$p, \$i, \$s, \$c);
720	\$s = 'PT Period Pen = #';
721	print MACRO ("\#MiMi = 13\r\#MjMi = 12\r\#MjMj = 9\r");	# Pen variables
722	print MACRO (\$s, 'MjMj ', Coordinate(\$xA, \$yA), '[`;`]', "\r" );	# Point at the pole
723	foreach \$m (0 .. 45) {
724		foreach \$p (0 .. 89) {
725			\$i = \$m . "," . \$p ;
726			if (defined (\$xP {\$i}) && defined (\$yP {\$i})) {
727				if (\$m % 5 == 0 && \$p % 5 == 0) {
728					\$c = 'MjMj';		# major meridian and major parallel
729				} elsif (\$m % 5 == 0 || \$p % 5 == 0) {
730					\$c = 'MjMi';	# major parallel or major meridian
731				} else {
732					\$c = 'MiMi';	# minor parallel with minor meridian
733				}
734				print MACRO (\$s,\$c,' ', Coordinate(\$xP{\$i},\$yP{\$i}), '[`;`]', "\r" );
735			}
736		}
737	}
738}		# End of subroutine PointsMacro.
739
740sub JointsMacro {
741	# I will be lazy and use global variables instead of trying to pass hashes as
742	# parameters.
743	# Plotting all joint intersections, that have already been calculated: arrays
744	# %xJ, %yJ.
745	# Subroutine will prepare TurboCad commands to plot each point. This will be written
746	# to filehandle MACRO.
747	# This subroutine is almost the same as PointsMacro; differences are:  different
748	# arrays and array indeces; not plotting "Point at the pole"; only two pens used;
749	# and there is no check for whether the variable is defined.
750
751	my (\$m, \$p, \$i, \$s, \$c);
752	\$s = 'PT Period Pen = #';
753	print MACRO ("\#Major = 9\r\#Minor = 12\r");	# Pen variables
754	foreach \$m (0 .. 45) {
755		if (\$m % 5 == 0) {
756			\$c = 'Major';		# Major meridian
757		} else {
758			\$c = 'Minor';		# Minor meridian
759		}
760		foreach \$p (0 .. 3) {
761			\$i = \$m . "," . \$p ;
762			print MACRO (\$s,\$c,' ', Coordinate(\$xJ{\$i},\$yJ{\$i}), '[`;`]', "\r" );
763		}
764	}
765}		# End of subroutine JointsMacro.
766
767
768sub MeridianMacro {
769	# I will be lazy and use global variables instead of trying to pass hashes as
770	# parameters.
771	# It is assumed that every meridian multiple of 5 will be written to file MAJOR;
772	# other meridians are written to file MINOR; Meridian numbers are from 0 to 45;
773	# will use hashes %xJ and %yJ, each of which is indexed as "i,j", where i= 0 .. 45,
774	# and j = 0 .. 4, for start polar point, two joints and equator point.
775	my (\$m, \$M, \$s, \$FileHandle);
776	foreach \$m (0 .. 45) {
777		\$M = \$m . ',' ;
778		if ( (\$m % 5) == 0 ) {		# Major meridian
779			\$FileHandle = MAJOR;
780		} else {				# Minor meridian
781			\$FileHandle = MINOR;
782		}
783		print \$FileHandle ('LI ', Coordinate(\$xJ {\$M . 0}, \$yJ {\$M . 0})); # First point
784		for (\$s = 1; \$s <= 3; \$s++ ) {	# for each line segment
785			print \$FileHandle ("\\\r",Coordinate(\$xJ {\$M . \$s},\$yJ {\$M . \$s}) );
786		}	# End of doing one meridian with 3 segments
787		print \$FileHandle ('[`;`]', "\r");		# End drawing of this multi-line
788	}	# End of for(m)-loop
789}	# End of Subroutine MeridianMacro
790
791sub ParallelArcs {
792	# This subroutine expects first and last parallels to do, and start and end meridians
793	# for the circular arc. If there is a 5th argument, it is a string to append to
794	# the comment line at the start of each parallel (e.g. ", Octant 1 West").
795	# It is assumed that centre of circle is point A, with global coordinates \$xA, \$yA.
796	# I will be lazy and use global variables instead of trying to pass hashes as
797	# parameters. Will use hashes %xP and %yP, each of which is indexed as "m,p".
798	# It is assumed that every parallel multiple of 5 will be written to file MAJOR;
799	# other meridians are written to file MINOR.
800	my (\$m, \$p, \$P, \$FileHandle, \$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo);
801	(\$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo) = @_;
802	if (not defined(\$Memo) ) {\$Memo = "";}
803
804	#	Make sure start is less than end, to be able to use foreach
805	if (\$pStart > \$pEnd) {
806		\$pStart = \$pEnd;	\$pEnd = \$_;
807	}
808	foreach \$p (\$pStart .. \$pEnd) {
809		\$P = ',' . \$p ;
810		if ( (\$p % 5) == 0 ) {		# Major parallel
811			\$FileHandle = MAJOR;
812		} else {				# Minor parallel
813			\$FileHandle = MINOR;
814		}
815		print \$FileHandle ('{ Parallel ', \$p, \$Memo,' }', "\r");
816		print \$FileHandle ('AR ', Coordinate (\$xA, \$yA), "\\\r");
817		print \$FileHandle (Coordinate (\$xP {\$mStart . \$P}, \$yP {\$mStart . \$P}) );
818		print \$FileHandle (Coordinate (\$xP {\$mEnd . \$P}, \$yP {\$mEnd . \$P} ) );
819		print \$FileHandle (' [`;`]', "\r");
820	}
821}	# End of subroutine ParallelArcs
822
823
824
825sub ParallelMacro {
826	# This subroutine expects first and last parallels to do, and start and end meridians
827	# for the segmented parallel. If there is a 5th argument, it is a string to append to
828	# the comment line at the start of each parallel (e.g. ", Octant 1 West").
829	# I will be lazy and use global variables instead of trying to pass hashes as
830	# parameters. Will use hashes %xP and %yP, each of which is indexed as "m,p".
831	# It is assumed that every parallel multiple of 5 will be written to file MAJOR;
832	# other parallels are written to file MINOR.
833	my (\$m, \$p, \$P, \$FileHandle, \$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo, \$mS, \$mE);
834	(\$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo) = @_;
835	if (not defined(\$Memo) ) {\$Memo = "";}
836
837	#	Make sure starts are less than ends, to be able to use foreach
838	if (\$pStart > \$pEnd) {
839		\$pStart = \$pEnd;	\$pEnd = \$_;
840	}
841	if (\$mStart > \$mEnd) {
842		\$mStart = \$mEnd;	\$mEnd = \$_;
843	}
844	\$mS = \$mStart + 1;	\$mE = \$mEnd - 1;
845
846	# For each parallel to be done
847	foreach \$p (\$pStart .. \$pEnd) {
848		\$P = ',' . \$p ;
849		if ( (\$p % 5) == 0 ) {		# Major parallel
850			\$FileHandle = MAJOR;
851		} else {				# Minor parallel
852			\$FileHandle = MINOR;
853		}
854		print \$FileHandle ('{ Parallel ', \$p, \$Memo, ' }', "\r");
855		# Start point of line
856		print \$FileHandle ('LI ', Coordinate (\$xP {\$mStart . \$P}, \$yP {\$mStart . \$P} ) );
857		foreach \$m (\$mS .. \$mE) {	# For each meridian
858			print \$FileHandle ("\\\r", Coordinate (\$xP {\$m . \$P}, \$yP {\$m . \$P} ) );
859			if ((\$m % 15) == 0) {
860				# Every 15 segments or so, terminate and restart line command
861				print \$FileHandle ('[`;`]', "\rLI ",
862					Coordinate (\$xP {\$m . \$P}, \$yP {\$m . \$P} ) );
863			}
864		} 	# End of foreach meridian
865		# Last point of parallel
866		print \$FileHandle (Coordinate (\$xP{\$mEnd.\$P}, \$yP{\$mEnd.\$P}), '[`;`]', "\r");
867	}	#End of foreach parallel
868}	# End of subroutine ParallelMacro
869
870sub ParallelMacro5 {
871	# This subroutine expects first and last parallels to do, and start and end meridians
872	# for the segmented parallel. It draws line segments for parallels by connecting nodes
873	# of major meridians, that is, every 5th meridian, skipping in between nodes.
874	# If there is a 5th argument, it is a string to append to the comment line at the start
875	# of each parallel (e.g. ", Octant 1 West").
876	# I will be lazy and use global variables instead of trying to pass hashes as
877	# parameters. Will use hashes %xP and %yP, each of which is indexed as "m,p".
878	# It is assumed that every parallel multiple of 5 will be written to file MAJOR;
879	# other meridians are written to file MINOR.
880	# This subroutine writes one parallel per LINE command (5� times 9 segments=45�).
881
882	my (\$m, \$p, \$P, \$FileHandle, \$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo, \$mS, \$mE);
883	(\$pStart, \$pEnd, \$mStart, \$mEnd, \$Memo) = @_;
884	if (not defined(\$Memo) ) {\$Memo = "";}
885
886	#	Make sure starts are less than ends, to be able to use foreach
887	if (\$pStart > \$pEnd) {
888		\$pStart = \$pEnd;	\$pEnd = \$_;
889	}
890	if (\$mStart > \$mEnd) {
891		\$mStart = \$mEnd;	\$mEnd = \$_;
892	}
893	\$mS = \$mStart + 1;	\$mE = \$mEnd - 1;
894
895	# For each parallel to be done
896	foreach \$p (\$pStart .. \$pEnd) {
897		\$P = ',' . \$p ;
898		if ( (\$p % 5) == 0 ) {		# Major parallel
899			\$FileHandle = MAJOR;
900		} else {				# Minor parallel
901			\$FileHandle = MINOR;
902		}
903		print \$FileHandle ('{ Parallel ', \$p, \$Memo, ' }', "\r");
904		# Start point of line
905		print \$FileHandle ('LI ', Coordinate (\$xP {\$mStart . \$P}, \$yP {\$mStart . \$P} ) );
906		foreach \$m (\$mS .. \$mE) {	# For each meridian
907			if ((\$m % 5) == 0) {	# only do for major meridians
908				print \$FileHandle ("\\\r", Coordinate(\$xP{\$m.\$P},\$yP{\$m.\$P} ) );
909			}
910		} 	# End of foreach meridian
911		# Last point of parallel
912		print \$FileHandle (Coordinate (\$xP{\$mEnd.\$P}, \$yP{\$mEnd.\$P}), '[`;`]', "\r");
913	}	#End of foreach parallel
914}	# End of subroutine ParallelMacro5
915
916sub MLinesMacro {
917	# Macro to draw concatenated lines, with coordinates from different octants.
918	# First argument is file handle. The second argument must be a number:
919	# If it is 0, the following arguments are sets of x,y, and octant number.
920	# If it is a number between 1 and 8, that is the number of the octant for all
921	# subsequent pairs of x,y coordinates.
922	# All coordinates are converted to M-map system.
923	my (\$FileHandle,\$Octant, @a) = @_ ;
924	my (\$nArg, \$i,\$x,\$y,\$x2,\$y2,\$n);
925	\$nArg = @a - 1;
926	\$n = 0;	# Flag for starting point
927	if (\$Octant == 0) {
928		for (\$i = 0; \$i < \$nArg; \$i=\$i+3) {
929			(\$x2,\$y2) = MJtoG (\$a[\$i], \$a[\$i+1], \$a[\$i+2]);
930			if (\$n == 0) {	# Starting point; update flag
931				print \$FileHandle ('LI ', Coordinate (\$x2,\$y2) );
932				\$n ++ ;
933			} else {
934				# Do line continuation, and this point
935				print \$FileHandle ("\\\r",Coordinate (\$x2,\$y2) );
936			}
937		}
938	} else {
939		for (\$i = 0; \$i < \$nArg; \$i=\$i+2) {
940			(\$x2,\$y2) = MJtoG (\$a[\$i], \$a[\$i+1], \$Octant);
941			if (\$n == 0) {	# Starting point; and update flag
942				print \$FileHandle ('LI ', Coordinate (\$x2,\$y2) );
943				\$n ++ ;
944			} else {
945				# Do line continuation, and this point
946				print \$FileHandle ("\\\r",Coordinate (\$x2,\$y2) );
947			}
948		}
949	}
950	print \$FileHandle (' [`;`]', "\r");
951}	# End of subroutine MLinesMacro
952
953sub OctantMacro {
954	# Macro to draw concatenated lines. Being lazy and using global arrays
955	# @xOct and @yOct, and writing to file handle MACRO.
956	my (\$FileHandle,@x, @y, \$n, \$i);
957	\$FileHandle = MACRO;
958	@x = @xOct;	@y = @yOct;
959	\$n = @xOct - 1;
960	print \$FileHandle ('LI ', Coordinate (\$x,\$y) );
961	foreach \$i (1 .. \$n) {
962		# Do line continuation, and this point
963		print \$FileHandle ("\\\r",Coordinate (\$x[\$i],\$y[\$i]) );
964	}
965	print \$FileHandle (' [`;`]', "\r");
966}	# End of subroutine LinesMacro
967
968sub GridMacro {
969	# Arguments are: filehandle, left, right, bottom, top, spacing, spacing of next
970	# larger grid or 0 to not skip major lines.
971	# Notice: every \$skip line is not plotted; it will be part of the grid of next larger size.
972	my (\$filehandle, \$left, \$right, \$bottom, \$top, \$delta, \$skip) = @_ ;
973	my (\$i, \$a, \$b, \$c);
974	\$a = 'LI [a,';
975	\$b = ',0,` `] [a,';
976	\$c = ',0,` `] [`;`]' . "\r";
977	\$i = \$bottom;
978	# If bottom is not a multiple of delta, start plotting at the first multiple of delta
979	if (\$i % \$delta !=0) { \$i = \$i + \$delta - (\$i % \$delta); }
980	print \$filehandle '{ Horizontal lines }', "\r" ;
981	while (\$i <= \$top) {	# horizontal lines
982		if (\$skip != 0 && \$i % \$skip == 0) {
983			\$i = \$i + \$delta;
984			 next;
985		}
986		print \$filehandle \$a,\$left,',',\$i,\$b,\$right,',',\$i,\$c;
987		\$i = \$i + \$delta;
988	}
989	\$i = \$left;
990	# If left is not a multiple of delta, start plotting at the first multiple of delta
991	if (\$i % \$delta !=0) { \$i = \$i + \$delta - (\$i % \$delta); }
992	print \$filehandle '{ Vertical lines }', "\r" ;
993	while (\$i <= \$right) {	# vertical lines
994		if (\$skip != 0 && \$i % \$skip == 0) {
995			\$i = \$i + \$delta;
996			 next;
997		}
998		print \$filehandle \$a,\$i,',',\$top,\$b,\$i,',',\$bottom,\$c;
999		\$i = \$i + \$delta;
1000	}
1001}	# End of subroutine GridMacro
1002
1003sub LandMacro {
1004	# Will be lazy and read global arrays @Xs and @Ys, and write to filehandle LAND.
1005	my (\$nPoints, \$FileHandle, \$n);
1006	\$nPoints = @Xs - 2;
1007	\$FileHandle = LAND;
1008	print \$FileHandle ('LI ', Coordinate (\$Xs , \$Ys ) );
1009	foreach \$n (1 .. \$nPoints) {
1010		print \$FileHandle ("\\\r", Coordinate (\$Xs[\$n], \$Ys [\$n] ) );
1011		if ((\$n % 15) == 0) {
1012			# Every 15 segments or so, terminate and restart line command
1013			print \$FileHandle ('[`;`]', "\rLI ", Coordinate (\$Xs[\$n], \$Ys [\$n] ) );
1014		}
1015	}
1016	print \$FileHandle (Coordinate (\$Xs[\$nPoints+1], \$Ys[\$nPoints+1]), '[`;`]', "\r");
1017}		# End of subroutine LandMacro
1018
1019# -  -  -  -  -  -  -  -  -    SUBROUTINES For Coordinate conversion   -  -  -  -  -  -  -  -  -  -  -  -
1020
1021sub MJtoG {
1022# Subroutine to convert (that is, do coordinate transformation of) x and y coordinates
1023# from Mary Jo's half-octant on its side to Gene's leaning, single octant coordinates, or
1024# to Gene's M-map (eight-octants) coordinates.
1025#
1026# Subroutine returns converted x and y coordinates.
1027#
1028# Arguments are:
1029# - x and y coordinates of point to convert.
1030# - Third argument, if given, is the Octant to convert to:
1031#	- no argument or 0 &#x2013; Gene's single-octant system, with y-axis on its left;
1032#	- 1, 2, 3 or 4 &#x2013; convert to M-map coordinates, respectively to first, second, third or
1033#		 fourth northern octant, from the left;
1034#	- 5, 6, 7 or 8 &#x2013; convert to M-map coordinates, respectively to fourth, first, second or
1035#		 third southern octant from the left.
1036#
1037# - In MJ's coordinates, point M is the origin, at (0,0), points M, A and G are on the positive
1038# x-axis, and point G is at (10000, 0).
1039# - In G's system, point M, L, J and P are on the positive y-axis, and point G is on the positive
1040# x-axis; in this system, point G is at coordinates (5000, 0).
1041# - From MJ's system to G's, there is a +60� rotation, and also a translation.
1042# - The M-map coordinate system is like G's system, except that the y-axis is 10000mm
1043# to the right, that is, the x-coordinates for the start octant are 10000mm smaller.
1044#
1045# I got the equations for rotation and translation from my pocketbook "The Universal
1046# Encyclopedia of Mathematics, with a Foreword by James R. Newman", �1964 by
1047# George Allen and Unwin, Ltd.; translated from original German language edition,
1048# pages 152, 153.
1049#
1050# This subroutine, in order to not have to recompute them each time the routine is called,
1051# uses a few global variables from the main program:
1054# \$yTranslate = 10000 * \$sin60;
1055
1056	my (\$nArgs, \$xnew, \$ynew);
1057	my (\$x, \$y, \$Octant) = @_ ;
1058	if (not defined (\$Octant) ) { \$Octant = 0; }
1059	if (\$Octant == 0) {
1060		(\$xnew, \$ynew) = Rotate (\$x, \$y, 60);
1061	} elsif (\$Octant == 1) {
1062		(\$xnew, \$ynew) = Rotate (\$x, \$y, 120);
1063		\$xnew = \$xnew - 10000;
1064	} elsif (\$Octant == 2) {
1065		(\$xnew, \$ynew) = Rotate (\$x, \$y, 60);
1066		\$xnew = \$xnew - 10000;
1067	} elsif (\$Octant == 3) {
1068		(\$xnew, \$ynew) = Rotate (\$x, \$y, 120);
1069		\$xnew = \$xnew + 10000;
1070	} elsif (\$Octant == 4) {
1071		(\$xnew, \$ynew) = Rotate (\$x, \$y, 60);
1072		\$xnew = \$xnew + 10000;
1073	} elsif (\$Octant == 5) {
1074		(\$xnew, \$ynew) = Rotate ((20000-\$x), \$y, 60);
1075		\$xnew = \$xnew + 10000;
1076	} elsif (\$Octant == 6) {
1077		(\$xnew, \$ynew) = Rotate ((20000-\$x), \$y, 120);
1078		\$xnew = \$xnew - 10000;
1079	} elsif (\$Octant == 7) {
1080		(\$xnew, \$ynew) = Rotate ((20000-\$x), \$y, 60);
1081		\$xnew = \$xnew - 10000;
1082	} elsif (\$Octant == 8) {
1083		(\$xnew, \$ynew) = Rotate ((20000-\$x), \$y, 120);
1084		\$xnew = \$xnew + 10000;
1085	} else {
1086		print "Error converting to M-map coordinates; there is no \$Octant octant!\n";
1087		return (\$x,\$y);
1088	}
1089	\$ynew = \$ynew + \$yTranslate;
1090	return (\$xnew, \$ynew);
1091
1092}	# End of subroutine MJtoM, which converts coordinates to octants on M-map
1093
1094sub Rotate {
1095	# Receives 3 arguments: x, y, and angle by which to rotate the coordinate system;
1096	# Expects that the axes will be rotated either 60� or 120�.
1097	# Returns new x and y values.
1098	# Uses \$cos60, \$sin60 global variables from main program, to minimize computations.
1099	my (\$x, \$y, \$angle) = @_ ;
1100	my (\$xnew, \$ynew);
1101	if ( \$angle == 60 ) {
1102		\$xnew = \$x * \$cos60 + \$y * \$sin60;
1103		\$ynew = -\$x * \$sin60 + \$y * \$cos60;
1104	} elsif ( \$angle == 120 ) {
1105		\$xnew = -\$x * \$cos60 + \$y * \$sin60;
1106		\$ynew = -\$x * \$sin60 - \$y * \$cos60;
1107	}
1108	return \$xnew, \$ynew;
1109}	# End of subroutine Rotate
1110
1111sub Real2Map {
1112	# Arguments are real world longitude and latitude for one point, in decimal degrees.
1113	# West longitudes and south latitudes have negative values.
1114	# Returns x and y values in M-Map coordinates.
1115	#
1116	# Subroutine uses the following global variables, without modifying them:
1117	# \$DeltaMEq, \$EF, \$AG, \$xA, \$yA, \$xE, \$yE, \$xF, \$yF, \$xG, \$xM, \$yM.
1118	# Also, it calls subroutines LineIntersection, Length, Interpolate and MJtoG.
1119
1120	# @South are southern octants corresponding to northern octants 1, 2, 3 and 4; the 0
1121	# is just a place holder to facilitate correspondence.
1122	my (@South) = (0,6,7,8,5);
1123	# \$m and \$p are the meridian and parallel numbers in scaffold half-octant setting;
1124	# \$Octant is the M-map octant of the real point; \$Sign is for east or west side of
1125	# scaffold octant; \$x and \$y are coordinates in scaffold half-octant setting; \$xNew
1126	# and \$yNew are the coordinates converted to M-Map system, which are returned.
1127	my (\$Octant, \$m, \$p, \$Sign, \$x, \$y, \$xNew, \$yNew);
1128	# Local arrays @xJ, @yJ and @L are equivalent to the global hashes %xJ, %yJ, %L,
1129	# but just for this particular meridian. The other variables are temporary.
1130	my (@xJ, @yJ, @L, \$distance, \$extra, \$flag);
1131	# Input values
1132	my (\$Long, \$Lat) = @_ ;
1133	# Determine the correct octant; Octant 1 is +160� to -110�; octant 4 is 70� to 160�
1134	\$Octant = int ( (\$Long + 200) / 90 ) + 1;
1135	# Make longitude fit within scaffold/template triangle, and determine if y value should
1136	# be positive or negative.
1137	\$m = ( (\$Long +200) - (90*(\$Octant - 1))) - 45;
1138	if (\$m < 0) {
1139		\$Sign = -1;	\$m = -\$m;
1140	} else {
1141		\$Sign = 1;
1142	}
1143	# Fix the octant number, if necessary
1144	if (\$Octant == 5) { \$Octant = 1; }
1145	if (\$Lat < 0) {
1146		\$Octant = \$South [\$Octant];
1147		\$p = -\$Lat;
1148	} else {
1149		\$p = \$Lat;
1150	}
1151#	print "\$Long, \$Lat convert to \$m, \$p, \$Sign, \$Octant.\n";
1152	# Now need to calculate x and y coordinates for \$m,\$p in scaffold triangle.
1153	if (\$m == 0) {	# Centre line meridian is straight, and easy to calculate
1154		\$y = 0;
1155		\$x = \$xG - \$p * \$AG / 90;
1156	} else {	# Complicated jointed meridian
1157		# Calculate coordinates of joints
1158		\$xJ  = \$xA;	\$yJ  = \$yA;	# Polar end of meridian
1159		\$distance = \$m * \$DeltaMEq;
1160		# Equatorial end of meridian
1161		if (\$distance <= \$yF) {
1162			\$xJ  = \$xG;
1163			\$yJ  = \$distance;
1164		} else {	# Past joint at point F
1165			\$extra = \$distance - \$yF;
1166			\$xJ  = Interpolate (\$extra, \$EF, \$xF, \$xE);
1167			\$yJ  = Interpolate (\$extra, \$EF, \$yF, \$yE);
1168		}
1169		# Calculate joints by line intersection from A and M
1170		(\$flag, \$xJ , \$yJ ) = LineIntersection (\$xM, \$yM, 2*\$m/3, \$xA, \$yA, \$m);
1171		if (\$flag ne "Intersect") { die "\$flag while calculating 1st joint at \$m.\n"; }
1172		(\$flag, \$xJ , \$yJ ) =
1173			LineIntersection (\$xM, \$yM, 2*\$m/3, \$xJ , \$yJ , \$m/3);
1174		if (\$flag ne "Intersect") { die "\$flag while calculating 2nd joint at \$m.\n"; }
1175		# Calculate length of meridian segments.
1176		\$L  = Length (\$xA, \$yA, \$xJ , \$yJ );
1177		\$L  = Length (\$xJ , \$yJ , \$xJ , \$yJ );
1178		\$L  = Length (\$xJ , \$yJ , \$xJ , \$yJ );
1179		\$L  = \$L  + \$L ;		# Torrid + Temperate
1180		\$L  = \$L  + \$L ;		# Total meridian length
1181		# Calculate parallel intersection &#x2013; location of this latitude on this meridian
1182		\$distance = \$p *  \$L  / 90;	# total length corresponds to 90�
1183		if (\$distance <= \$L  ) {		# Torrid segment
1184			\$x = Interpolate (\$distance, \$L , \$xJ , \$xJ );
1185			\$y = Interpolate (\$distance, \$L , \$yJ , \$yJ );
1186		} elsif (\$distance <= \$L  ) {	# Temperate seg., past one joint
1187			\$extra = \$distance - \$L ;
1188			\$x = Interpolate (\$extra, \$L , \$xJ , \$xJ );
1189			\$y = Interpolate (\$extra, \$L , \$yJ , \$yJ );
1190		} else {				# Frigid segment, past two joints
1191			\$extra = \$distance - \$L ;
1192			\$x = Interpolate (\$extra, \$L , \$xJ , \$xJ );
1193			\$y = Interpolate (\$extra, \$L , \$yJ , \$yJ );
1194		}
1195	}	# Now have \$x and \$y coordinates in scaffold half-octant coordinates
1196	# Convert to M-map coordinates
1197	(\$xNew, \$yNew) = MJtoG (\$x,\$Sign*\$y, \$Octant);
1198#	print "and x,y: \$x, \$y, converted to \$xNew, \$yNew.\n";
1199	return \$xNew, \$yNew;
1200}	# End of subroutine Real2Map
1201
1202# -  -  -  -  -  -  -  -  -    SUBROUTINES         -  -  -  -  -  -  -  -  -  -  -  -  -
```