1 /*
2  * tclgeomap.c --
3  *
4  *	This file defines the structures and functions that implement the
5  *	tclgeomap extension to Tcl, which adds the ability to store and
6  *	manipulate geographic data in Tcl.
7  *
8  * Copyright (c) 2004 Gordon D. Carrie.  All rights reserved.
9  *
10  * Licensed under the Open Software License version 2.1
11  *
12  * Please address questions and feedback to user0@tkgeomap.org
13  *
14  * @(#) $Id: tclgeomap.c,v 1.14 2007/04/20 15:08:58 tkgeomap Exp $
15  *
16  ********************************************
17  *
18  */
19 
20 #include "tclgeomap.h"
21 #include "tclgeomapInt.h"
22 
23 /*
24  * Forward declarations.
25  */
26 
27 static int		version _ANSI_ARGS_((ClientData clientData,
28                                 Tcl_Interp *interp, int objc,
29 				Tcl_Obj *CONST objv[]));
30 static int		setEarthRadius _ANSI_ARGS_((ClientData clientData,
31                                 Tcl_Interp *interp, int objc,
32 				Tcl_Obj *CONST objv[]));
33 static int		geoPtIsSomewhere _ANSI_ARGS_((ClientData clientData,
34                                 Tcl_Interp *interp, int objc,
35 				Tcl_Obj *CONST objv[]));
36 static int		mapPtIsSomewhere _ANSI_ARGS_((ClientData clientData,
37                                 Tcl_Interp *interp, int objc,
38 				Tcl_Obj *CONST objv[]));
39 static int		lonbtwn _ANSI_ARGS_((ClientData clientData,
40                                 Tcl_Interp *interp, int objc,
41 				Tcl_Obj *CONST objv[]));
42 static int		cross _ANSI_ARGS_((ClientData clientData,
43                                 Tcl_Interp *interp, int objc,
44 				Tcl_Obj *CONST objv[]));
45 static int		rotate _ANSI_ARGS_((ClientData clientData,
46                                 Tcl_Interp *interp, int objc,
47 				Tcl_Obj *CONST objv[]));
48 static int		scale _ANSI_ARGS_((ClientData clientData,
49                                 Tcl_Interp *interp, int objc,
50 				Tcl_Obj *CONST objv[]));
51 static int		domnlonpt _ANSI_ARGS_((ClientData clientData,
52                                 Tcl_Interp *interp, int objc,
53 				Tcl_Obj *CONST objv[]));
54 static int		domnlat _ANSI_ARGS_((ClientData clientData,
55                                 Tcl_Interp *interp, int objc,
56 				Tcl_Obj *CONST objv[]));
57 static int		domnlon _ANSI_ARGS_((ClientData clientData,
58                                 Tcl_Interp *interp, int objc,
59 				Tcl_Obj *CONST objv[]));
60 static int		gwchpt _ANSI_ARGS_((ClientData clientData,
61                                 Tcl_Interp *interp, int objc,
62 				Tcl_Obj *CONST objv[]));
63 static int		gwchlon _ANSI_ARGS_((ClientData clientData,
64                                 Tcl_Interp *interp, int objc,
65 				Tcl_Obj *CONST objv[]));
66 static int		dmsToDec _ANSI_ARGS_((ClientData clientData,
67                                 Tcl_Interp *interp, int objc,
68 				Tcl_Obj *CONST objv[]));
69 static int		decToDM _ANSI_ARGS_((ClientData clientData,
70                                 Tcl_Interp *interp, int objc,
71 				Tcl_Obj *CONST objv[]));
72 static int		decToDMS _ANSI_ARGS_((ClientData clientData,
73                                 Tcl_Interp *interp, int objc,
74 				Tcl_Obj *CONST objv[]));
75 static int		cartg _ANSI_ARGS_((ClientData clientData,
76                                 Tcl_Interp *interp, int objc,
77 				Tcl_Obj *CONST objv[]));
78 static int		centroid _ANSI_ARGS_((ClientData clientData,
79                                 Tcl_Interp *interp, int objc,
80 				Tcl_Obj *CONST objv[]));
81 static int		jul_tm _ANSI_ARGS_((ClientData clientData,
82                                 Tcl_Interp *interp, int objc,
83 				Tcl_Obj *CONST objv[]));
84 static int		cal_tm _ANSI_ARGS_((ClientData clientData,
85                                 Tcl_Interp *interp, int objc,
86 				Tcl_Obj *CONST objv[]));
87 
88 /*
89  * The following functions and data structure define a Tcl object
90  * that stores a {latitude longitude} pair.
91  */
92 
93 static Tcl_FreeInternalRepProc freeGeoPt;
94 static Tcl_DupInternalRepProc dupGeoPt;
95 static Tcl_UpdateStringProc updateGeoPtString;
96 static Tcl_SetFromAnyProc setGeoPtFromAny;
97 
98 static Tcl_ObjType GeoPtType = {
99     "GeoPt",
100     freeGeoPt,
101     dupGeoPt,
102     updateGeoPtString,
103     setGeoPtFromAny
104 };
105 
106 /*
107  * The following functions and data structure define a Tcl object
108  * that stores a point on a two dimensional map.
109  */
110 
111 static Tcl_FreeInternalRepProc freeMapPt;
112 static Tcl_DupInternalRepProc dupMapPt;
113 static Tcl_UpdateStringProc updateMapPtString;
114 static Tcl_SetFromAnyProc setMapPtFromAny;
115 
116 static Tcl_ObjType MapPtType = {
117     "MapPt",
118     freeMapPt,
119     dupMapPt,
120     updateMapPtString,
121     setMapPtFromAny
122 };
123 
124 /*
125  * The following functions and data structure define a Tcl object
126  * that stores a point on a two dimensional map.
127  */
128 
129 static Tcl_FreeInternalRepProc freeGeoTime;
130 static Tcl_DupInternalRepProc dupGeoTime;
131 static Tcl_UpdateStringProc updateGeoTimeString;
132 static Tcl_SetFromAnyProc setGeoTimeFromAny;
133 
134 static Tcl_ObjType GeoTimeType = {
135     "GeoTime",
136     freeGeoTime,
137     dupGeoTime,
138     updateGeoTimeString,
139     setGeoTimeFromAny
140 };
141 
142 /*
143  *------------------------------------------------------------------------
144  *
145  * Tclgeomap_Init --
146  *
147  *	This procedure initializes the Tclgeomap interface and provides
148  *	the tclgeomap package.
149  *
150  * Results:
151  *	A standard Tcl result.
152  *
153  * Side effects:
154  *	Initializes other interfaces.
155  *
156  *
157  *------------------------------------------------------------------------
158  */
159 
160 int
Tclgeomap_Init(interp)161 Tclgeomap_Init(interp)
162     Tcl_Interp *interp;		/* Current Tcl interpreter */
163 {
164     static int loaded;		/* Tell if package already loaded */
165 
166     if (loaded) {
167 	return TCL_OK;
168     }
169 
170     /*
171      * Initialize other interfaces.
172      */
173 
174     if (Tcl_InitStubs(interp, TCL_VERSION, 0) == NULL) {
175 	return TCL_ERROR;
176     }
177     if (Tcl_PkgRequire(interp, "Tcl", TCL_VERSION, 0) == NULL) {
178 	if (TCL_VERSION[0] == '7') {
179 	    if (Tcl_PkgRequire(interp, "Tcl", "8.0", 0) == NULL) {
180 		return TCL_ERROR;
181 	    }
182 	}
183     }
184     if (TclgeomapProjInit(interp) != TCL_OK) {
185 	return TCL_ERROR;
186     }
187     if (TclgeomapTimeInit(interp) != TCL_OK) {
188 	return TCL_ERROR;
189     }
190     if (TclgeomapPlaceInit(interp) != TCL_OK) {
191 	return TCL_ERROR;
192     }
193     if (TclgeomapLnArrInit(interp) != TCL_OK) {
194 	return TCL_ERROR;
195     }
196     if (Tcl_PkgProvide(interp, "tclgeomap", TCLGEOMAP_VERSION) != TCL_OK) {
197 	return TCL_ERROR;
198     }
199 
200     /*
201      * Initialize general geography objects and functions.
202      */
203 
204     Tcl_RegisterObjType(&GeoPtType);
205     Tcl_RegisterObjType(&MapPtType);
206     Tcl_CreateObjCommand(interp, "::geomap::version", version, NULL, NULL);
207     Tcl_CreateObjCommand(interp, "::geomap::georadius", setEarthRadius, NULL,
208 	    NULL);
209     Tcl_CreateObjCommand(interp, "::geomap::latlonok", geoPtIsSomewhere, NULL,
210 	    NULL);
211     Tcl_CreateObjCommand(interp, "::geomap::mapptok", mapPtIsSomewhere, NULL,
212 	    NULL);
213     Tcl_CreateObjCommand(interp, "::geomap::lonbtwn", lonbtwn, NULL, NULL);
214     Tcl_CreateObjCommand(interp, "::geomap::gclcross", cross, NULL, NULL);
215     Tcl_CreateObjCommand(interp, "::geomap::rotatpt", rotate, NULL, NULL);
216     Tcl_CreateObjCommand(interp, "::geomap::scalept", scale, NULL, NULL);
217     Tcl_CreateObjCommand(interp, "::geomap::domnlonpt", domnlonpt, NULL, NULL);
218     Tcl_CreateObjCommand(interp, "::geomap::domnlat", domnlat, NULL, NULL);
219     Tcl_CreateObjCommand(interp, "::geomap::domnlon", domnlon, NULL, NULL);
220     Tcl_CreateObjCommand(interp, "::geomap::gwchpt", gwchpt, NULL, NULL);
221     Tcl_CreateObjCommand(interp, "::geomap::gwchlon", gwchlon, NULL, NULL);
222     Tcl_CreateObjCommand(interp, "::geomap::dmstodec", dmsToDec, NULL, NULL);
223     Tcl_CreateObjCommand(interp, "::geomap::dectodm", decToDM, NULL, NULL);
224     Tcl_CreateObjCommand(interp, "::geomap::dectodms", decToDMS, NULL, NULL);
225     Tcl_CreateObjCommand(interp, "::geomap::cartg", cartg, NULL, NULL);
226     Tcl_CreateObjCommand(interp, "::geomap::centroid", centroid, NULL, NULL);
227     Tcl_CreateObjCommand(interp, "::geomap::jul_tm", jul_tm, NULL, NULL);
228     Tcl_CreateObjCommand(interp, "::geomap::cal_tm", cal_tm, NULL, NULL);
229 
230     loaded = 1;
231     return TCL_OK;
232 }
233 
234 /*
235  *------------------------------------------------------------------------
236  *
237  * Tclgeomap_NewGeoPtObj --
238  *
239  *	This procedure creates a new GeoPt object and initializes it from
240  *	the argument GeoPt value.
241  *
242  * Results:
243  *	The newly created object is returned. This object will have an
244  *	invalid string representation and a ref count 0.
245  *
246  * Side effects:
247  *	None.
248  *
249  *------------------------------------------------------------------------
250  */
251 
252 Tcl_Obj *
Tclgeomap_NewGeoPtObj(geoPt)253 Tclgeomap_NewGeoPtObj(geoPt)
254     GeoPt geoPt;	/* GeoPt value for the new object. */
255 {
256     GeoPt *geoPtPtr;	/* GeoPt to store in the new object */
257     Tcl_Obj *objPtr;	/* The new object */
258 
259     objPtr = Tcl_NewObj();
260     geoPtPtr = (GeoPt *)CKALLOC(sizeof(GeoPt));
261     objPtr->internalRep.otherValuePtr = geoPtPtr;
262     *(GeoPt *)objPtr->internalRep.otherValuePtr = geoPt;
263     objPtr->typePtr = &GeoPtType;
264     objPtr->bytes = NULL;
265     return objPtr;
266 }
267 
268 /*
269  *------------------------------------------------------------------------
270  *
271  * freeGeoPt --
272  *
273  *	This procedure frees the memory associated with the argument
274  *	GeoPt object.
275  *
276  * Results:
277  *	None.
278  *
279  * Side effects:
280  *	Deallocates the internalRep of a geoPt object.
281  *
282  *------------------------------------------------------------------------
283  */
284 
285 void
freeGeoPt(objPtr)286 freeGeoPt(objPtr)
287     Tcl_Obj *objPtr;
288 {
289     CKFREE(objPtr->internalRep.otherValuePtr);
290     objPtr->typePtr = NULL;
291 }
292 
293 /*
294  *------------------------------------------------------------------------
295  *
296  * dupGeoPt --
297  *
298  *	This procedure copies a GeoPt object into another.
299  *
300  * Results:
301  *	The internal representation of the destination object is freed and
302  *	its string representation is invalidated.  A new GeoPt structure
303  *	is allocated and stored in the destination object.
304  *
305  * Side effects:
306  *	None.
307  *
308  *------------------------------------------------------------------------
309  */
310 
311 void
dupGeoPt(srcObjPtr,dupObjPtr)312 dupGeoPt(srcObjPtr, dupObjPtr)
313     Tcl_Obj *srcObjPtr;		/* Source object */
314     Tcl_Obj *dupObjPtr;		/* Destination object */
315 {
316     GeoPt *srcGeoPtPtr;	/* Value stored in source object */
317     GeoPt *dupGeoPtPtr;	/* Value to store in destination object */
318 
319     if (dupObjPtr
320 	    && dupObjPtr->typePtr
321 	    && dupObjPtr->typePtr->freeIntRepProc) {
322 	dupObjPtr->typePtr->freeIntRepProc(dupObjPtr);
323     }
324     Tcl_InvalidateStringRep(dupObjPtr);
325     srcGeoPtPtr = (GeoPt *)srcObjPtr->internalRep.otherValuePtr;
326     dupGeoPtPtr = (GeoPt *)CKALLOC(sizeof(GeoPt));
327     memcpy(dupGeoPtPtr, srcGeoPtPtr, sizeof(GeoPt));
328     dupObjPtr->internalRep.otherValuePtr = dupGeoPtPtr;
329     dupObjPtr->typePtr = &GeoPtType;
330 }
331 
332 /*
333  *------------------------------------------------------------------------
334  *
335  * updateGeoPtString --
336  *
337  *	This procedure updates the string representation of a GeoPt object.
338  *
339  * Results:
340  *	The string member of the argument object is modified.
341  *
342  * Side effects:
343  *	None.
344  *
345  *------------------------------------------------------------------------
346  */
347 
348 void
updateGeoPtString(objPtr)349 updateGeoPtString(objPtr)
350     Tcl_Obj *objPtr;		/* Input GeoPt object */
351 {
352     double dLat, dLon;		/* Latitude and longitude of the geopoint */
353     char *latLon[2];		/* Strings where "lat lon" will be printed */
354     char lat[TCL_DOUBLE_SPACE];	/* latLon[0] */
355     char lon[TCL_DOUBLE_SPACE];	/* latLon[1] */
356     GeoPt *geoPtPtr;		/* GeoPt value from objPtr */
357 
358     latLon[0] = lat;
359     latLon[1] = lon;
360     geoPtPtr = (GeoPt *)objPtr->internalRep.otherValuePtr;
361     GeoPtGetDeg(*geoPtPtr, &dLat, &dLon);
362     Tcl_PrintDouble(NULL, dLat, lat);
363     Tcl_PrintDouble(NULL, dLon, lon);
364     objPtr->bytes = Tcl_Merge(2, latLon);
365     objPtr->length = strlen(objPtr->bytes);
366 }
367 
368 /*
369  *------------------------------------------------------------------------
370  *
371  * setGeoPtFromAny --
372  *
373  *	This procedure attempts to  convert an object into a GeoPt object.
374  *
375  * Results:
376  *	A standard Tcl result.
377  *
378  * Side effects:
379  *	A new GeoPt value is stored in objPtr's internal representation.
380  *
381  *------------------------------------------------------------------------
382  */
383 
384 int
setGeoPtFromAny(interp,objPtr)385 setGeoPtFromAny(interp, objPtr)
386     Tcl_Interp *interp;		/* Current interpreter */
387     Tcl_Obj *objPtr;		/* Object to convert to GeoPt object */
388 {
389     Tcl_Obj **listPtr;		/* Two element list obtained from objPtr's
390 				 * original value */
391     int n;			/* Number of elements in list at listPtr.
392 				 * Should be 2. */
393     double lat, lon;		/* Coordinates obtained from objPtr's original
394 				 * value */
395     GeoPt *geoPtPtr;		/* New value to assign to objPtr's internal
396 				 * representation */
397 
398     if (Tcl_ListObjGetElements(interp, objPtr, &n, &listPtr) != TCL_OK
399 	    || n != 2
400 	    || Tcl_GetDoubleFromObj(interp, listPtr[0], &lat) != TCL_OK
401 	    || Tcl_GetDoubleFromObj(interp, listPtr[1], &lon) != TCL_OK) {
402 	if (interp) {
403 	    Tcl_AppendResult(interp, "Could not read ", objPtr->bytes,
404 		    " as a lat-lon pair.  ", NULL);
405 	}
406 	return TCL_ERROR;
407     }
408     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
409 	objPtr->typePtr->freeIntRepProc(objPtr);
410     }
411     geoPtPtr = (GeoPt *)CKALLOC(sizeof(GeoPt));
412     *geoPtPtr = GeoPtFmDeg(lat, lon);
413     objPtr->internalRep.otherValuePtr = geoPtPtr;
414     objPtr->typePtr = &GeoPtType;
415 
416     return TCL_OK;
417 }
418 
419 /*
420  *------------------------------------------------------------------------
421  *
422  * Tclgeomap_SetGeoPtObj --
423  *
424  *	This procedure converts an object to a GeoPt object with a given
425  *	value.
426  *
427  * Results:
428  *	The previous contents of objPtr are invalidated and freed.  A new
429  *	GeoPt is allocated and stored in objPtr.
430  *
431  * Side effects:
432  *	None.
433  *------------------------------------------------------------------------
434  */
435 
436 void
Tclgeomap_SetGeoPtObj(objPtr,geoPt)437 Tclgeomap_SetGeoPtObj(objPtr, geoPt)
438     Tcl_Obj *objPtr;	/* Object to become a GeoPt object */
439     GeoPt geoPt;	/* Value to store in objPtr */
440 {
441     GeoPt *geoPtPtr;	/* New value for objPtr's internalRep */
442 
443     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
444 	objPtr->typePtr->freeIntRepProc(objPtr);
445     }
446     Tcl_InvalidateStringRep(objPtr);
447     geoPtPtr = (GeoPt *)CKALLOC(sizeof(GeoPt));
448     objPtr->internalRep.otherValuePtr = geoPtPtr;
449     *(GeoPt *)objPtr->internalRep.otherValuePtr = geoPt;
450     objPtr->typePtr = &GeoPtType;
451 }
452 
453 /*
454  *------------------------------------------------------------------------
455  *
456  * Tclgeomap_GetGeoPtFromObj --
457  *
458  *	This procedure retrieves the GeoPt stored in a GeoPt object.
459  *
460  * Results:
461  *	If successful, the GeoPt value stored in objPtr is copied into
462  *	geoPtPtr and TCL_OK is returned.  Otherwise, TCL_ERROR is returned.
463  *
464  * Side effects:
465  *	None.
466  *
467  *------------------------------------------------------------------------
468  */
469 
470 int
Tclgeomap_GetGeoPtFromObj(interp,objPtr,geoPtPtr)471 Tclgeomap_GetGeoPtFromObj(interp, objPtr, geoPtPtr)
472     Tcl_Interp *interp;		/* Current interpreter */
473     Tcl_Obj *objPtr;		/* GeoPt object with GeoPt value */
474     GeoPt *geoPtPtr;		/* Pointer to received value from objPtr */
475 {
476     if (   objPtr->typePtr != &GeoPtType
477 	    && setGeoPtFromAny(interp, objPtr) != TCL_OK) {
478 	/*
479 	 * Conversion error.  setGeoPtFromAny provides error message
480 	 */
481 
482 	return TCL_ERROR;
483     }
484     *geoPtPtr = *(GeoPt *)objPtr->internalRep.otherValuePtr;
485     return TCL_OK;
486 }
487 
488 /*
489  *------------------------------------------------------------------------
490  *
491  * Tclgeomap_NewMapPtObj --
492  *
493  *	This procedure creates a new MapPt object and initializes it from
494  *	the argument MapPt value.
495  *
496  * Results:
497  *	The newly created object is returned. This object will have an
498  *	invalid string representation and a ref count 0.
499  *
500  * Side effects:
501  *	None.
502  *
503  *------------------------------------------------------------------------
504  */
505 Tcl_Obj *
Tclgeomap_NewMapPtObj(mapPt)506 Tclgeomap_NewMapPtObj(mapPt)
507     MapPt mapPt;	/* MapPt value for the new object */
508 {
509     MapPt *mapPtPtr;	/* MapPt to store in the new object */
510     Tcl_Obj *objPtr;	/* Return value */
511 
512     objPtr = Tcl_NewObj();
513     mapPtPtr = (MapPt *)CKALLOC(sizeof(MapPt));
514     objPtr->internalRep.otherValuePtr = mapPtPtr;
515     *(MapPt *)objPtr->internalRep.otherValuePtr = mapPt;
516     objPtr->typePtr = &MapPtType;
517     objPtr->bytes = NULL;
518     return objPtr;
519 }
520 
521 /*
522  *------------------------------------------------------------------------
523  *
524  * freeMapPt --
525  *
526  *	This procedure frees the memory associated with the argument
527  *	MapPt object.
528  *
529  * Results:
530  *	None.
531  *
532  * Side effects:
533  *	Deallocates the internalRep of a mapPt object.
534  *
535  *------------------------------------------------------------------------
536  */
537 
538 void
freeMapPt(objPtr)539 freeMapPt(objPtr)
540     Tcl_Obj *objPtr;
541 {
542     CKFREE(objPtr->internalRep.otherValuePtr);
543     objPtr->typePtr = NULL;
544 }
545 
546 /*
547  *------------------------------------------------------------------------
548  *
549  * dupMapPt --
550  *
551  *	This procedure copies a MapPt object into another.
552  *
553  * Results:
554  *	The internal representation of the destination object is freed and
555  *	its string representation is invalidated.  A new MapPt structure
556  *	is allocated and stored in the destination object.
557  *
558  * Side effects:
559  *	None.
560  *
561  *------------------------------------------------------------------------
562  */
563 
564 void
dupMapPt(srcObjPtr,dupObjPtr)565 dupMapPt(srcObjPtr, dupObjPtr)
566     Tcl_Obj *srcObjPtr;		/* Source object */
567     Tcl_Obj *dupObjPtr;		/* Destination object */
568 {
569     MapPt *srcMapPtPtr;		/* Value stored in source object */
570     MapPt *dupMapPtPtr;		/* Value to store in destination object */
571 
572     if (dupObjPtr
573 	    && dupObjPtr->typePtr
574 	    && dupObjPtr->typePtr->freeIntRepProc) {
575 	dupObjPtr->typePtr->freeIntRepProc(dupObjPtr);
576     }
577     Tcl_InvalidateStringRep(dupObjPtr);
578     srcMapPtPtr = (MapPt *)srcObjPtr->internalRep.otherValuePtr,
579     dupMapPtPtr = (MapPt *)CKALLOC(sizeof(MapPt));
580     dupMapPtPtr->ord = srcMapPtPtr->ord;
581     dupMapPtPtr->abs = srcMapPtPtr->abs;
582     dupObjPtr->internalRep.otherValuePtr = dupMapPtPtr;
583     dupObjPtr->typePtr = &MapPtType;
584 }
585 
586 /*
587  *------------------------------------------------------------------------
588  *
589  * updateMapPtString --
590  *
591  *	This procedure updates the string representation of a MapPt object.
592  *
593  * Results:
594  *	The string member of the argument object is modified.
595  *
596  * Side effects:
597  *	None.
598  *
599  *------------------------------------------------------------------------
600  */
601 
602 void
updateMapPtString(objPtr)603 updateMapPtString(objPtr)
604     Tcl_Obj *objPtr;		/* Input MapPt object */
605 {
606     char *coords[2];		/* Strings where "abs ord" will be printed */
607     char abs[TCL_DOUBLE_SPACE];	/* coords[0] */
608     char ord[TCL_DOUBLE_SPACE];	/* coords[1] */
609     MapPt *mapPtPtr;		/* MapPt value from objPtr */
610 
611     coords[0] = abs;
612     coords[1] = ord;
613     mapPtPtr = (MapPt *)objPtr->internalRep.otherValuePtr;
614     Tcl_PrintDouble(NULL, mapPtPtr->abs, abs);
615     Tcl_PrintDouble(NULL, mapPtPtr->ord, ord);
616     objPtr->bytes = Tcl_Merge(2, coords);
617     objPtr->length = strlen(objPtr->bytes);
618 }
619 
620 /*
621  *------------------------------------------------------------------------
622  *
623  * setMapPtFromAny --
624  *
625  *	This procedure attempts to  convert an object into a MapPt object.
626  *
627  * Results:
628  *	A standard Tcl result.
629  *
630  * Side effects:
631  *	A new MapPt value is stored in objPtr's internal representation.
632  *
633  *------------------------------------------------------------------------
634  */
635 
636 int
setMapPtFromAny(interp,objPtr)637 setMapPtFromAny(interp, objPtr)
638     Tcl_Interp *interp;		/* Current interpreter */
639     Tcl_Obj *objPtr;		/* Object to convert to MapPt object */
640 {
641     Tcl_Obj **listPtr;		/* Two element list obtained from objPtr's
642 				 * original value */
643     int n;			/* Number of elements in list at listPtr.
644 				 * Should be 2. */
645     double abs, ord;		/* Coordinates obtained from objPtr's original
646 				 * value */
647     MapPt *mapPtPtr;		/* New value to assign to objPtr's internal
648 				 * representation */
649 
650     if (Tcl_ListObjGetElements(interp, objPtr, &n, &listPtr) != TCL_OK
651 	    || n != 2
652 	    || Tcl_GetDoubleFromObj(interp, listPtr[0], &abs) != TCL_OK
653 	    || Tcl_GetDoubleFromObj(interp, listPtr[1], &ord) != TCL_OK) {
654 	if (interp) {
655 	    Tcl_AppendResult(interp, "Could not read ", objPtr->bytes,
656 		    " as a map point.  ", NULL);
657 	}
658 	return TCL_ERROR;
659     }
660     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
661 	objPtr->typePtr->freeIntRepProc(objPtr);
662     }
663     mapPtPtr = (MapPt *)CKALLOC(sizeof(MapPt));
664     mapPtPtr->abs = abs;
665     mapPtPtr->ord = ord;
666     objPtr->internalRep.otherValuePtr = mapPtPtr;
667     objPtr->typePtr = &MapPtType;
668     return TCL_OK;
669 }
670 
671 /*
672  *------------------------------------------------------------------------
673  *
674  * Tclgeomap_SetMapPtObj --
675  *
676  *	This procedure converts an object to a MapPt object with a given
677  *	value.
678  *
679  * Results:
680  *	The previous contents of objPtr are invalidated and freed.  A new
681  *	MapPt is allocated and stored in objPtr.
682  *
683  * Side effects:
684  *	None.
685  *------------------------------------------------------------------------
686  */
687 
688 void
Tclgeomap_SetMapPtObj(objPtr,mapPt)689 Tclgeomap_SetMapPtObj(objPtr, mapPt)
690     Tcl_Obj *objPtr;	/* Object to become a MapPt object */
691     MapPt mapPt;	/* Value to store in objPtr */
692 {
693     MapPt *mapPtPtr;	/* New value for objPtr's internalRep */
694 
695     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
696 	objPtr->typePtr->freeIntRepProc(objPtr);
697     }
698     Tcl_InvalidateStringRep(objPtr);
699     mapPtPtr = (MapPt *)CKALLOC(sizeof(MapPt));
700     objPtr->internalRep.otherValuePtr = mapPtPtr;
701     *(MapPt *)objPtr->internalRep.otherValuePtr = mapPt;
702     objPtr->typePtr = &MapPtType;
703 }
704 
705 /*
706  *------------------------------------------------------------------------
707  *
708  * Tclgeomap_GetMapPtFromObj --
709  *
710  *	This procedure retrieves the MapPt stored in a MapPt object.
711  *
712  * Results:
713  *	If successful, the MapPt value stored in objPtr is copied into
714  *	geoPtPtr and TCL_OK is returned.  Otherwise, TCL_ERROR is returned.
715  *
716  * Side effects:
717  *	None.
718  *
719  *------------------------------------------------------------------------
720  */
721 
722 int
Tclgeomap_GetMapPtFromObj(interp,objPtr,mapPtPtr)723 Tclgeomap_GetMapPtFromObj(interp, objPtr, mapPtPtr)
724     Tcl_Interp *interp;
725     Tcl_Obj *objPtr;
726     MapPt *mapPtPtr;
727 {
728     /* Put MapPt value from objPtr into mapPtPtr */
729     if (   objPtr->typePtr != &MapPtType
730 	    && setMapPtFromAny(interp, objPtr) != TCL_OK) {
731 	/* Conversion error.  setMapPtFromAny provides error message */
732 	return TCL_ERROR;
733     }
734     *mapPtPtr = *(MapPt *)objPtr->internalRep.otherValuePtr;
735     return TCL_OK;
736 }
737 
738 /*
739  *------------------------------------------------------------------------
740  *
741  * Tclgeomap_NewGeoTimeObj --
742  *
743  *	This procedure creates a new GeoTime object and initializes it from
744  *	the argument GeoTime value.
745  *
746  * Results:
747  *	The newly created object is returned. This object will have an
748  *	invalid string representation and a ref count 0.
749  *
750  * Side effects:
751  *	None.
752  *
753  *------------------------------------------------------------------------
754  */
755 
756 Tcl_Obj *
Tclgeomap_NewGeoTimeObj(jul)757 Tclgeomap_NewGeoTimeObj(jul)
758     struct GeoTime_Jul jul;	/* GeoTime value for the new object */
759 {
760     struct GeoTime_Jul *julPtr;	/* GeoTime to store in the new object */
761     Tcl_Obj *objPtr;		/* Return value */
762 
763     objPtr = Tcl_NewObj();
764     julPtr = (struct GeoTime_Jul *)CKALLOC(sizeof(struct GeoTime_Jul));
765     objPtr->internalRep.otherValuePtr = julPtr;
766     *(struct GeoTime_Jul *)objPtr->internalRep.otherValuePtr = jul;
767     objPtr->typePtr = &GeoTimeType;
768     objPtr->bytes = NULL;
769     return objPtr;
770 }
771 
772 /*
773  *------------------------------------------------------------------------
774  *
775  * freeGeoTime --
776  *
777  *	This procedure frees the memory associated with the argument
778  *	GeoTime object.
779  *
780  * Results:
781  *	None.
782  *
783  * Side effects:
784  *	Deallocates the internalRep of a geoTime object.
785  *
786  *------------------------------------------------------------------------
787  */
788 
789 void
freeGeoTime(objPtr)790 freeGeoTime(objPtr)
791     Tcl_Obj *objPtr;
792 {
793     CKFREE(objPtr->internalRep.otherValuePtr);
794     objPtr->typePtr = NULL;
795 }
796 
797 /*
798  *------------------------------------------------------------------------
799  *
800  * dupGeoTime --
801  *
802  *	This procedure copies a GeoTime object into another.
803  *
804  * Results:
805  *	The internal representation of the destination object is freed and
806  *	its string representation is invalidated.  A new GeoTime structure
807  *	is allocated and stored in the destination object.
808  *
809  * Side effects:
810  *	None.
811  *
812  *------------------------------------------------------------------------
813  */
814 
815 void
dupGeoTime(srcObjPtr,dupObjPtr)816 dupGeoTime(srcObjPtr, dupObjPtr)
817     Tcl_Obj *srcObjPtr;			/* Source object */
818     Tcl_Obj *dupObjPtr;			/* Destination object */
819 {
820     struct GeoTime_Jul *srcGeoTimePtr;	/* Value stored in source object */
821     struct GeoTime_Jul *dupGeoTimePtr;	/* Value for destination object */
822 
823     if (dupObjPtr && dupObjPtr->typePtr && dupObjPtr->typePtr->freeIntRepProc) {
824 	dupObjPtr->typePtr->freeIntRepProc(dupObjPtr);
825     }
826     Tcl_InvalidateStringRep(dupObjPtr);
827     srcGeoTimePtr = (struct GeoTime_Jul *)srcObjPtr->internalRep.otherValuePtr,
828     dupGeoTimePtr = (struct GeoTime_Jul *)CKALLOC(sizeof(struct GeoTime_Jul));
829     dupGeoTimePtr->day = srcGeoTimePtr->day;
830     dupGeoTimePtr->second = srcGeoTimePtr->second;
831     dupObjPtr->internalRep.otherValuePtr = dupGeoTimePtr;
832     dupObjPtr->typePtr = &GeoTimeType;
833 }
834 
835 /*
836  *------------------------------------------------------------------------
837  *
838  * updateGeoTimeString --
839  *
840  *	This procedure updates the string representation of a GeoTime object.
841  *
842  * Results:
843  *	The string member of the argument object is modified.
844  *
845  * Side effects:
846  *	None.
847  *
848  *------------------------------------------------------------------------
849  */
850 
851 void
updateGeoTimeString(objPtr)852 updateGeoTimeString(objPtr)
853     Tcl_Obj *objPtr;			/* Input GeoTime object */
854 {
855     struct GeoTime_Jul *julPtr;		/* GeoTime value from objPtr */
856     char *elems[6];			/* Strings to print to */
857     char year[TCL_INTEGER_SPACE];	/* coords[0] */
858     char month[TCL_INTEGER_SPACE];	/* coords[1] */
859     char day[TCL_INTEGER_SPACE];	/* coords[2] */
860     char hour[TCL_INTEGER_SPACE];	/* coords[3] */
861     char minute[TCL_INTEGER_SPACE];	/* coords[4] */
862     char second[TCL_DOUBLE_SPACE];	/* coords[5] */
863     struct GeoTime_Cal cal;		/* Calendar version of julPtr */
864 
865     julPtr = (struct GeoTime_Jul *)objPtr->internalRep.otherValuePtr;
866     cal = GeoTime_JulToCal(*julPtr);
867     sprintf(year, "%d", cal.year);
868     sprintf(month, "%d", cal.month);
869     sprintf(day, "%d", cal.day);
870     sprintf(hour, "%d", cal.hour);
871     sprintf(minute, "%d", cal.minute);
872     Tcl_PrintDouble(NULL, cal.second, second);
873     elems[0] = year;
874     elems[1] = month;
875     elems[2] = day;
876     elems[3] = hour;
877     elems[4] = minute;
878     elems[5] = second;
879     objPtr->bytes = Tcl_Merge(6, elems);
880     objPtr->length = strlen(objPtr->bytes);
881 }
882 
883 /*
884  *------------------------------------------------------------------------
885  *
886  * setGeoTimeFromAny --
887  *
888  *	This procedure attempts to  convert an object into a GeoTime object.
889  *
890  * Results:
891  *	A standard Tcl result.
892  *
893  * Side effects:
894  *	A new GeoTime value is stored in objPtr's internal representation.
895  *
896  *------------------------------------------------------------------------
897  */
898 
899 int
setGeoTimeFromAny(interp,objPtr)900 setGeoTimeFromAny(interp, objPtr)
901     Tcl_Interp *interp;		/* Current interpreter */
902     Tcl_Obj *objPtr;		/* Object to convert to GeoTime object */
903 {
904     Tcl_Obj **listPtr;		/* Six element list with
905 				 * {year month day hour minute second} obtained
906 				 * from objPtr's original value */
907     int n;			/* Number of elements in list at listPtr.
908 				 * Should be 6. */
909     struct GeoTime_Cal cal;	/* Holder for time */
910     struct GeoTime_Jul *julPtr; /* New value to assign to objPtr's internal
911 				 * representation */
912     int year, month, day, hour, minute;
913     double second;
914 
915     if (Tcl_ListObjGetElements(interp, objPtr, &n, &listPtr) != TCL_OK
916 	    || n != 6
917 	    || Tcl_GetIntFromObj(interp, listPtr[0], &year) != TCL_OK
918 	    || Tcl_GetIntFromObj(interp, listPtr[1], &month) != TCL_OK
919 	    || Tcl_GetIntFromObj(interp, listPtr[2], &day) != TCL_OK
920 	    || Tcl_GetIntFromObj(interp, listPtr[3], &hour) != TCL_OK
921 	    || Tcl_GetIntFromObj(interp, listPtr[4], &minute) != TCL_OK
922 	    || Tcl_GetDoubleFromObj(interp, listPtr[5], &second) != TCL_OK) {
923 	if (interp) {
924 	    Tcl_AppendResult(interp, "Expected {year month day hour minute "
925 		    "second}.  Got:", objPtr->bytes, NULL);
926 	}
927 	return TCL_ERROR;
928     }
929     cal = GeoTime_CalSet(year, month, day, hour, minute, second);
930     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
931 	objPtr->typePtr->freeIntRepProc(objPtr);
932     }
933     julPtr = (struct GeoTime_Jul *)CKALLOC(sizeof(struct GeoTime_Jul));
934     *julPtr = GeoTime_CalToJul(cal);
935     objPtr->internalRep.otherValuePtr = julPtr;
936     objPtr->typePtr = &GeoTimeType;
937     return TCL_OK;
938 }
939 
940 /*
941  *------------------------------------------------------------------------
942  *
943  * Tclgeomap_SetGeoTimeObj --
944  *
945  *	This procedure converts an object to a GeoTime object with a given
946  *	value.
947  *
948  * Results:
949  *	The previous contents of objPtr are invalidated and freed.  A new
950  *	GeoTime is allocated and stored in objPtr.
951  *
952  * Side effects:
953  *	None.
954  *------------------------------------------------------------------------
955  */
956 
957 void
Tclgeomap_SetGeoTimeObj(objPtr,jul)958 Tclgeomap_SetGeoTimeObj(objPtr, jul)
959     Tcl_Obj *objPtr;			/* Object to become a GeoTime object */
960     struct GeoTime_Jul jul;		/* Value to store in objPtr */
961 {
962     struct GeoTime_Jul *julPtr;		/* New value for objPtr's internalRep */
963 
964     if (objPtr && objPtr->typePtr && objPtr->typePtr->freeIntRepProc) {
965 	objPtr->typePtr->freeIntRepProc(objPtr);
966     }
967     Tcl_InvalidateStringRep(objPtr);
968     julPtr = (struct GeoTime_Jul *)CKALLOC(sizeof(struct GeoTime_Jul));
969     objPtr->internalRep.otherValuePtr = julPtr;
970     *(struct GeoTime_Jul *)objPtr->internalRep.otherValuePtr = jul;
971     objPtr->typePtr = &GeoTimeType;
972 }
973 
974 /*
975  *------------------------------------------------------------------------
976  *
977  * Tclgeomap_GetGeoTimeFromObj --
978  *
979  *	This procedure retrieves the GeoTime stored in a GeoTime object.
980  *
981  * Results:
982  *	If successful, the GeoTime value stored in objPtr is copied into
983  *	geoPtPtr and TCL_OK is returned.  Otherwise, TCL_ERROR is returned.
984  *
985  * Side effects:
986  *	None.
987  *
988  *------------------------------------------------------------------------
989  */
990 
991 int
Tclgeomap_GetGeoTimeFromObj(interp,objPtr,julPtr)992 Tclgeomap_GetGeoTimeFromObj(interp, objPtr, julPtr)
993     Tcl_Interp *interp;
994     Tcl_Obj *objPtr;
995     struct GeoTime_Jul *julPtr;
996 {
997     if (objPtr->typePtr != &GeoTimeType
998 	    && setGeoTimeFromAny(interp, objPtr) != TCL_OK) {
999 	/*
1000 	 * Conversion error.  setGeoTimeFromAny provides error message
1001 	 */
1002 
1003 	return TCL_ERROR;
1004     }
1005     *julPtr = *(struct GeoTime_Jul *)objPtr->internalRep.otherValuePtr;
1006     return TCL_OK;
1007 }
1008 
1009 /*
1010  *------------------------------------------------------------------------
1011  *
1012  * version --
1013  *
1014  *	This procedure prints version and copyright information.
1015  *
1016  * Results:
1017  *	A standard Tcl result.
1018  *
1019  * Side effects:
1020  * 	Set result to a copyright notice.
1021  *
1022  *------------------------------------------------------------------------
1023  */
1024 
1025 int
version(clientData,interp,objc,objv)1026 version(clientData, interp, objc, objv)
1027     ClientData clientData;	/* Not used */
1028     Tcl_Interp *interp;		/* Current interpreter */
1029     int objc;			/* Number of arguments */
1030     Tcl_Obj *CONST objv[];	/* Argument objects */
1031 {
1032     static char vsn[] =
1033 	"tclgeomap/tkgeomap, version " TCLGEOMAP_VERSION "\n"
1034 	"Copyright (C) 2003 Gordon D. Carrie\n"
1035 	"This software comes with ABSOLUTELY NO WARRANTY.\n"
1036 	"Some or all of this program is free software.  You should have\n"
1037 	"received or been offered source code with this program.\n"
1038 	"See the source code and user documentation for warranty and \n"
1039 	"distribution terms.";
1040     Tcl_SetResult(interp, vsn, TCL_STATIC);
1041     return TCL_OK;
1042 }
1043 
1044 /*
1045  *------------------------------------------------------------------------
1046  *
1047  * setEarthRadius --
1048  *
1049  *	This procedure gets or sets the Earth's radius.
1050  *
1051  * Results:
1052  *	A standard Tcl result.
1053  *
1054  * Side effects:
1055  *	If a new value is given, this procedure modifies the value of Earth's
1056  *	radius used in all subsequent geography calculations.
1057  *
1058  *------------------------------------------------------------------------
1059  */
1060 
1061 int
setEarthRadius(clientData,interp,objc,objv)1062 setEarthRadius(clientData, interp, objc, objv)
1063     ClientData clientData;	/* Not used */
1064     Tcl_Interp *interp;		/* Current interpreter */
1065     int objc;			/* Number of arguments */
1066     Tcl_Obj *CONST objv[];	/* Argument objects */
1067 {
1068     char *units[] =  {
1069 	"nmiles", "smiles", "km", "meters", "cm", NULL
1070     };				/* Unit specifiers */
1071     enum unitIdx {
1072 	NMILES,   SMILES,   KM,   METERS,   CM
1073     };				/* Corresponding indices */
1074     int idx;			/* Index for unit specifier on command line */
1075 
1076     if (objc == 2) {
1077 	/*
1078 	 * Command is of form "georadius unit".  Set result to radius of Earth
1079 	 * in that unit.
1080 	 */
1081 
1082 	if (Tcl_GetIndexFromObj(interp, objv[1], units, "unit", 0, &idx)
1083 		!= TCL_OK) {
1084 	    return TCL_ERROR;
1085 	}
1086 	switch ((enum unitIdx)idx) {
1087 	    case NMILES:
1088 		Tcl_SetObjResult(interp, Tcl_NewDoubleObj(REarth() / NMI));
1089 		break;
1090 	    case SMILES:
1091 		Tcl_SetObjResult(interp, Tcl_NewDoubleObj(REarth() / SMI));
1092 		break;
1093 	    case KM:
1094 		Tcl_SetObjResult(interp, Tcl_NewDoubleObj(REarth() * 0.001));
1095 		break;
1096 	    case METERS:
1097 		Tcl_SetObjResult(interp, Tcl_NewDoubleObj(REarth()));
1098 		break;
1099 	    case CM:
1100 		Tcl_SetObjResult(interp, Tcl_NewDoubleObj(REarth() * 100.0));
1101 		break;
1102 	}
1103     } else if (objc == 3) {
1104 	/*
1105 	 * Command is of form "georadius value unit".  Set Earth radius to
1106 	 * given value, adjusting for unit specification.
1107 	 */
1108 
1109 	double r;		/* New value for Earth radius */
1110 	if (Tcl_GetDoubleFromObj(interp, objv[1], &r) != TCL_OK) {
1111 	    return TCL_ERROR;
1112 	}
1113 	if (Tcl_GetIndexFromObj(interp, objv[2], units, "unit", 0, &idx)
1114 		!= TCL_OK) {
1115 	    return TCL_ERROR;
1116 	}
1117 	switch ((enum unitIdx)idx) {
1118 	    case NMILES:
1119 		SetREarth(r * NMI);
1120 		break;
1121 	    case SMILES:
1122 		SetREarth(r * SMI);
1123 		break;
1124 	    case KM:
1125 		SetREarth(r * 1000.0);
1126 		break;
1127 	    case METERS:
1128 		SetREarth(r);
1129 		break;
1130 	    case CM:
1131 		SetREarth(r * 0.01);
1132 		break;
1133 	}
1134     } else {
1135 	Tcl_WrongNumArgs(interp, objc, objv, "?radius? unit");
1136 	return TCL_ERROR;
1137     }
1138     return TCL_OK;
1139 }
1140 
1141 /*
1142  *------------------------------------------------------------------------
1143  *
1144  * geoPtIsSomewhere --
1145  *
1146  *	This is the callback for the latlonok command.
1147  *
1148  * Results:
1149  *	A standard Tcl result.
1150  *
1151  * Side effects:
1152  * 	See the user documentation.
1153  *
1154  *------------------------------------------------------------------------
1155  */
1156 
1157 int
geoPtIsSomewhere(clientData,interp,objc,objv)1158 geoPtIsSomewhere(clientData, interp, objc, objv)
1159     ClientData clientData;	/* Not used */
1160     Tcl_Interp *interp;		/* Current interpreter */
1161     int objc;			/* Number of arguments */
1162     Tcl_Obj *CONST objv[];	/* Argument objects */
1163 {
1164     GeoPt geoPt;		/* GeoPt to evaluate */
1165 
1166     if (objc != 2) {
1167 	Tcl_WrongNumArgs(interp, 1, objv, " point");
1168 	return TCL_ERROR;
1169     }
1170     if (Tclgeomap_GetGeoPtFromObj(NULL, objv[1], &geoPt) != TCL_OK) {
1171 	Tcl_SetObjResult(interp, Tcl_NewBooleanObj(0));
1172 	return TCL_OK;
1173     }
1174     Tcl_SetObjResult(interp, Tcl_NewBooleanObj(GeoPtIsSomewhere(geoPt)));
1175     return TCL_OK;
1176 }
1177 
1178 /*
1179  *------------------------------------------------------------------------
1180  *
1181  * mapPtIsSomewhere --
1182  *
1183  *	This is the callback for the mapptok command.
1184  *
1185  * Results:
1186  *	A standard Tcl result.
1187  *
1188  * Side effects:
1189  * 	See the user documentation.
1190  *
1191  *------------------------------------------------------------------------
1192  */
1193 
1194 int
mapPtIsSomewhere(clientData,interp,objc,objv)1195 mapPtIsSomewhere(clientData, interp, objc, objv)
1196     ClientData clientData;	/* Not used */
1197     Tcl_Interp *interp;		/* Current interpreter */
1198     int objc;			/* Number of arguments */
1199     Tcl_Obj *CONST objv[];	/* Argument objects */
1200 {
1201     MapPt mapPt;		/* MapPt to evaluate */
1202 
1203     if (objc != 2) {
1204 	Tcl_WrongNumArgs(interp, 1, objv, " point");
1205 	return TCL_ERROR;
1206     }
1207     if (Tclgeomap_GetMapPtFromObj(NULL, objv[1], &mapPt) != TCL_OK) {
1208 	Tcl_SetObjResult(interp, Tcl_NewBooleanObj(0));
1209 	return TCL_OK;
1210     }
1211     Tcl_SetObjResult(interp, Tcl_NewBooleanObj(MapPtIsSomewhere(mapPt)));
1212     return TCL_OK;
1213 }
1214 
1215 /*
1216  *------------------------------------------------------------------------
1217  *
1218  * lonbtwn --
1219  *
1220  *	This is the callback for the lonbtwn command.
1221  *
1222  * Results:
1223  *	A standard Tcl result.
1224  *
1225  * Side effects:
1226  * 	See the user documentation.
1227  *
1228  *------------------------------------------------------------------------
1229  */
1230 
1231 int
lonbtwn(clientData,interp,objc,objv)1232 lonbtwn(clientData, interp, objc, objv)
1233     ClientData clientData;	/* Not used */
1234     Tcl_Interp *interp;		/* Current interpreter */
1235     int objc;			/* Number of arguments */
1236     Tcl_Obj *CONST objv[];	/* Argument objects */
1237 {
1238     double lon;			/* Longitude to consider */
1239     double lon0, lon1;		/* Longitudes that lon might be between */
1240     Tcl_Obj **pairPtr;		/* List from command line containing
1241 				 * "lon0 lon1" */
1242     int pairCnt;		/* Number of elements in pairPtr */
1243 
1244     if (objc != 3
1245 	    || Tcl_GetDoubleFromObj(interp, objv[1], &lon) != TCL_OK
1246 	    || Tcl_ListObjGetElements(interp, objv[2], &pairCnt, &pairPtr)
1247 	    != TCL_OK
1248 	    || pairCnt != 2
1249 	    || Tcl_GetDoubleFromObj(interp, pairPtr[0], &lon0) != TCL_OK
1250 	    || Tcl_GetDoubleFromObj(interp, pairPtr[1], &lon1) != TCL_OK) {
1251 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1252 		" lon {lon0 lon1}", NULL);
1253 	return TCL_ERROR;
1254     }
1255     Tcl_SetObjResult(interp, Tcl_NewBooleanObj(LonBtwn(AngleFmDeg(lon),
1256 		    AngleFmDeg(lon0), AngleFmDeg(lon1))));
1257     return TCL_OK;
1258 }
1259 
1260 /*
1261  *------------------------------------------------------------------------
1262  *
1263  * cross --
1264  *
1265  *	This is the callback for the gclcross command.
1266  *
1267  * Results:
1268  *	A standard Tcl result.
1269  *
1270  * Side effects:
1271  * 	See the user documentation.
1272  *
1273  *------------------------------------------------------------------------
1274  */
1275 
1276 int
cross(clientData,interp,objc,objv)1277 cross(clientData, interp, objc, objv)
1278     ClientData clientData;	/* Not used */
1279     Tcl_Interp *interp;		/* Current interpreter */
1280     int objc;			/* Number of arguments */
1281     Tcl_Obj *CONST objv[];	/* Argument objects */
1282 {
1283     GeoPt ln1Pt1, ln1Pt2;		/* Points on 1st line */
1284     GeoPt ln2Pt1, ln2Pt2;		/* Points on 2nd line */
1285     GeoPt geoPt;			/* Return value */
1286 
1287     if (objc != 5) {
1288 	Tcl_WrongNumArgs(interp, 1, objv, "{line1_lat1 line1_lon1} "
1289 		"{line1_lat2 line1_lon2} {line2_lat1 line2_lon1} "
1290 		"{line2_lat2 line2_lon2}");
1291 	return TCL_ERROR;
1292     }
1293     if (Tclgeomap_GetGeoPtFromObj(interp, objv[1], &ln1Pt1) != TCL_OK
1294 	    || Tclgeomap_GetGeoPtFromObj(interp, objv[2], &ln1Pt2) != TCL_OK
1295 	    || Tclgeomap_GetGeoPtFromObj(interp, objv[3], &ln2Pt1) != TCL_OK
1296 	    || Tclgeomap_GetGeoPtFromObj(interp, objv[4], &ln2Pt2) != TCL_OK
1297 	    ) {
1298 	return TCL_ERROR;
1299     }
1300     geoPt = GCircleX(ln1Pt1, ln1Pt2, ln2Pt1, ln2Pt2);
1301     if (GeoPtIsNowhere(geoPt)) {
1302 	Tcl_AppendResult(interp, "Undefined intersection", NULL);
1303 	return TCL_ERROR;
1304     } else {
1305 	Tcl_SetObjResult(interp, Tclgeomap_NewGeoPtObj(geoPt));
1306 	return TCL_OK;
1307     }
1308 }
1309 
1310 /*
1311  *------------------------------------------------------------------------
1312  *
1313  * rotate --
1314  *
1315  *	This is the callback for the rotatpt command.
1316  *
1317  * Results:
1318  *	A standard Tcl result.
1319  *
1320  * Side effects:
1321  * 	See the user documentation.
1322  *
1323  *------------------------------------------------------------------------
1324  */
1325 
1326 int
rotate(clientData,interp,objc,objv)1327 rotate(clientData, interp, objc, objv)
1328     ClientData clientData;	/* Not used */
1329     Tcl_Interp *interp;		/* Current interpreter */
1330     int objc;			/* Number of arguments */
1331     Tcl_Obj *CONST objv[];	/* Argument objects */
1332 {
1333     MapPt
1334 	mapPt,			/* Input value */
1335 	rPt;			/* Return value */
1336     double
1337 	theta,			/* Angle to rotate point by */
1338 	costh, sinth;		/* Computational constants */
1339 
1340     if (objc != 3
1341 	    || Tclgeomap_GetMapPtFromObj(interp, objv[1], &mapPt) != TCL_OK
1342 	    || Tcl_GetDoubleFromObj(interp, objv[2], &theta) != TCL_OK) {
1343 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1344 		"mapPt theta", NULL);
1345 	return TCL_ERROR;
1346     }
1347     if (MapPtIsNowhere(mapPt)) {
1348 	Tcl_SetObjResult(interp, Tclgeomap_NewMapPtObj(mapPt));
1349     }
1350     costh = cos(theta * RADPERDEG);
1351     sinth = sin(theta * RADPERDEG);
1352     rPt.abs = mapPt.abs * costh + mapPt.ord * sinth;
1353     rPt.ord = mapPt.ord * costh - mapPt.abs * sinth;
1354     Tcl_SetObjResult(interp, Tclgeomap_NewMapPtObj(rPt));
1355     return TCL_OK;
1356 }
1357 
1358 /*
1359  *------------------------------------------------------------------------
1360  *
1361  * scale --
1362  *
1363  *	This is the callback for the scalept command.
1364  *
1365  * Results:
1366  *	A standard Tcl result.
1367  *
1368  * Side effects:
1369  * 	See the user documentation.
1370  *
1371  *------------------------------------------------------------------------
1372  */
1373 
1374 int
scale(clientData,interp,objc,objv)1375 scale(clientData, interp, objc, objv)
1376     ClientData clientData;	/* Not used */
1377     Tcl_Interp *interp;		/* Current interpreter */
1378     int objc;			/* Number of arguments */
1379     Tcl_Obj *CONST objv[];	/* Argument objects */
1380 {
1381     MapPt mapPt;		/* Input value */
1382     double fac;			/* Scale factor */
1383 
1384     if (   objc != 3
1385 	    || Tclgeomap_GetMapPtFromObj(interp, objv[1], &mapPt) != TCL_OK
1386 	    || Tcl_GetDoubleFromObj(interp, objv[2], &fac)) {
1387 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1388 		"mapPt fac", NULL);
1389 	return TCL_ERROR;
1390     }
1391     if (MapPtIsNowhere(mapPt)) {
1392 	Tcl_SetObjResult(interp, Tclgeomap_NewMapPtObj(mapPt));
1393     }
1394     Tcl_SetObjResult(interp, Tclgeomap_NewMapPtObj(ScaleMapPt(mapPt, fac)));
1395     return TCL_OK;
1396 }
1397 
1398 /*
1399  *------------------------------------------------------------------------
1400  *
1401  * domnlat --
1402  *
1403  *	This is the callback for the domnlat command.
1404  *
1405  * Results:
1406  *	A standard Tcl result.
1407  *
1408  * Side effects:
1409  * 	See the user documentation.
1410  *
1411  *------------------------------------------------------------------------
1412  */
1413 
1414 int
domnlat(clientData,interp,objc,objv)1415 domnlat(clientData, interp, objc, objv)
1416     ClientData clientData;	/* Not used */
1417     Tcl_Interp *interp;		/* Current interpreter */
1418     int objc;			/* Number of arguments */
1419     Tcl_Obj *CONST objv[];	/* Argument objects */
1420 {
1421     double latDeg;		/* Input values */
1422     Angle lat;
1423 
1424     if (objc != 2
1425 	    || Tcl_GetDoubleFromObj(interp, objv[1], &latDeg) != TCL_OK) {
1426 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1427 		" latitude", NULL);
1428 	return TCL_ERROR;
1429     }
1430     lat = AngleFmDeg(latDeg);
1431     Tcl_SetObjResult(interp, Tcl_NewDoubleObj(AngleToDeg(DomainLat(lat))));
1432     return TCL_OK;
1433 }
1434 
1435 /*
1436  *------------------------------------------------------------------------
1437  *
1438  * domnlon --
1439  *
1440  *	This is the callback for the domnlon command.
1441  *
1442  * Results:
1443  *	A standard Tcl result.
1444  *
1445  * Side effects:
1446  * 	See the user documentation.
1447  *
1448  *------------------------------------------------------------------------
1449  */
1450 
1451 int
domnlon(clientData,interp,objc,objv)1452 domnlon(clientData, interp, objc, objv)
1453     ClientData clientData;	/* Not used */
1454     Tcl_Interp *interp;		/* Current interpreter */
1455     int objc;			/* Number of arguments */
1456     Tcl_Obj *CONST objv[];	/* Argument objects */
1457 {
1458     double lonDeg, rLonDeg;	/* Input values */
1459     Angle lon, rLon;
1460 
1461     if (objc != 3
1462 	    || Tcl_GetDoubleFromObj(interp, objv[1], &lonDeg) != TCL_OK
1463 	    || Tcl_GetDoubleFromObj(interp, objv[2], &rLonDeg) != TCL_OK) {
1464 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1465 		" lonDeg reflon", NULL);
1466 	return TCL_ERROR;
1467     }
1468     lon = AngleFmDeg(lonDeg);
1469     rLon = AngleFmDeg(rLonDeg);
1470     Tcl_SetObjResult(interp,
1471 	    Tcl_NewDoubleObj(AngleToDeg(DomainLon(lon, rLon))));
1472     return TCL_OK;
1473 }
1474 
1475 /*
1476  *------------------------------------------------------------------------
1477  *
1478  * domnlonpt --
1479  *
1480  *	This is the callback for the domnlonpt command.
1481  *
1482  * Results:
1483  *	A standard Tcl result.
1484  *
1485  * Side effects:
1486  * 	See the user documentation.
1487  *
1488  *------------------------------------------------------------------------
1489  */
1490 
1491 int
domnlonpt(clientData,interp,objc,objv)1492 domnlonpt(clientData, interp, objc, objv)
1493     ClientData clientData;	/* Not used */
1494     Tcl_Interp *interp;		/* Current interpreter */
1495     int objc;			/* Number of arguments */
1496     Tcl_Obj *CONST objv[];	/* Argument objects */
1497 {
1498     GeoPt geoPt;		/* Input value */
1499     double rlon;		/* Reference longitude */
1500 
1501     if (objc != 3
1502 	    || Tclgeomap_GetGeoPtFromObj(interp, objv[1], &geoPt) != TCL_OK
1503 	    || Tcl_GetDoubleFromObj(interp, objv[2], &rlon) != TCL_OK) {
1504 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1505 		" geoPt rlon", NULL);
1506 	return TCL_ERROR;
1507     }
1508     if (GeoPtIsNowhere(geoPt)) {
1509 	Tcl_SetObjResult(interp, Tclgeomap_NewGeoPtObj(geoPt));
1510 	return TCL_OK;
1511     }
1512     Tcl_SetObjResult(interp,
1513 	    Tclgeomap_NewGeoPtObj(DomainLonPt(geoPt, AngleFmDeg(rlon))));
1514     return TCL_OK;
1515 }
1516 
1517 /*
1518  *------------------------------------------------------------------------
1519  *
1520  * gwchlon --
1521  *
1522  *	This is the callback for the gwchlon command.
1523  *
1524  * Results:
1525  *	A standard Tcl result.
1526  *
1527  * Side effects:
1528  * 	See the user documentation.
1529  *
1530  *------------------------------------------------------------------------
1531  */
1532 
1533 int
gwchlon(clientData,interp,objc,objv)1534 gwchlon(clientData, interp, objc, objv)
1535     ClientData clientData;	/* Not used */
1536     Tcl_Interp *interp;		/* Current interpreter */
1537     int objc;			/* Number of arguments */
1538     Tcl_Obj *CONST objv[];	/* Argument objects */
1539 {
1540     double lon;			/* Input value */
1541 
1542     if (objc != 2
1543 	    || Tcl_GetDoubleFromObj(interp, objv[1], &lon) != TCL_OK) {
1544 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1545 		" lon", NULL);
1546 	return TCL_ERROR;
1547     }
1548     Tcl_SetObjResult(interp,
1549 	    Tcl_NewDoubleObj(AngleToDeg(GwchLon(AngleFmDeg(lon)))));
1550     return TCL_OK;
1551 }
1552 
1553 /*
1554  *------------------------------------------------------------------------
1555  *
1556  * gwchpt --
1557  *
1558  *	This is the callback for the gwchpt command.
1559  *
1560  * Results:
1561  *	A standard Tcl result.
1562  *
1563  * Side effects:
1564  * 	See the user documentation.
1565  *
1566  *------------------------------------------------------------------------
1567  */
1568 
1569 int
gwchpt(clientData,interp,objc,objv)1570 gwchpt(clientData, interp, objc, objv)
1571     ClientData clientData;	/* Not used */
1572     Tcl_Interp *interp;		/* Current interpreter */
1573     int objc;			/* Number of arguments */
1574     Tcl_Obj *CONST objv[];	/* Argument objects */
1575 {
1576     GeoPt geoPt;		/* Input value */
1577 
1578     if (objc != 2 || Tclgeomap_GetGeoPtFromObj(interp, objv[1], &geoPt)
1579 	    != TCL_OK) {
1580 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
1581 		" geoPt", NULL);
1582 	return TCL_ERROR;
1583     }
1584     if (GeoPtIsNowhere(geoPt)) {
1585 	Tcl_SetObjResult(interp, Tclgeomap_NewGeoPtObj(geoPt));
1586 	return TCL_OK;
1587     }
1588     Tcl_SetObjResult(interp, Tclgeomap_NewGeoPtObj(GwchLonPt(geoPt)));
1589     return TCL_OK;
1590 }
1591 
1592 /*
1593  *------------------------------------------------------------------------
1594  *
1595  * dmsToDec --
1596  *
1597  *	This is the callback for the dmstodec command.
1598  *
1599  * Results:
1600  *	A standard Tcl result.
1601  *
1602  * Side effects:
1603  * 	See the user documentation.
1604  *
1605  *------------------------------------------------------------------------
1606  */
1607 
1608 int
dmsToDec(clientData,interp,objc,objv)1609 dmsToDec(clientData, interp, objc, objv)
1610     ClientData clientData;	/* Not used */
1611     Tcl_Interp *interp;		/* Current interpreter */
1612     int objc;			/* Number of arguments */
1613     Tcl_Obj *CONST objv[];	/* Argument objects */
1614 {
1615     char
1616 	*cmdNm,			/* Command name */
1617 	*usmg, *w1, *w2;	/* Usage message */
1618     Tcl_Obj **dmsPtr;		/* Elements from command line */
1619     int n;			/* Element count from command line */
1620     double latd, latm, lats;	/* Latitude degrees, minutes, seconds */
1621     double lond, lonm, lons;	/* Longitude degrees, minutes, seconds */
1622     double d, m, s;		/* Degrees, minutes, seconds */
1623     char *ns, *we;		/* Hemisphere: "N" or "S" or "W" or "E"
1624 				 * or lower case */
1625     double lat, lon;		/* Latitude and longitude for point */
1626     double deg;			/* Return value if defining an angle */
1627 
1628     /*
1629      * Initialize variables
1630      */
1631 
1632     cmdNm = Tcl_GetStringFromObj(objv[0], NULL);
1633     w1 = "{deg ?min ?sec?? NorS deg ?min ?sec?? WorE}";
1634     w2 = "{deg ?min ?sec??}";
1635     usmg = CKALLOC(2 * strlen(cmdNm) + strlen(w1) + strlen(w2) + 18);
1636     sprintf(usmg, "Usage: \"%s %s\" or \"%s %s\"", cmdNm, w1, cmdNm, w2);
1637     latd = latm = lats = lond = lonm = lons = d = m = s = 0.0;
1638 
1639     if (objc != 2) {
1640 	Tcl_WrongNumArgs(interp, 1, objv, usmg);
1641 	return TCL_ERROR;
1642     }
1643     if (Tcl_ListObjGetElements(interp, objv[1], &n, &dmsPtr) != TCL_OK) {
1644 	Tcl_AppendResult(interp, usmg, NULL);
1645 	return TCL_ERROR;
1646     }
1647     switch (n) {
1648 	case 8:
1649 	    /*
1650 	     * Form is {deg min sec NorS deg min sec WorE}
1651 	     */
1652 
1653 	    if (Tcl_GetDoubleFromObj(interp, dmsPtr[0], &latd) != TCL_OK
1654 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[1], &latm) != TCL_OK
1655 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[2], &lats) != TCL_OK
1656 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[4], &lond) != TCL_OK
1657 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[5], &lonm) != TCL_OK
1658 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[6], &lons) != TCL_OK)
1659 		return TCL_ERROR;
1660 	    ns = Tcl_GetStringFromObj(dmsPtr[3], NULL);
1661 	    if ( !strchr("NnSs", *ns) ) {
1662 		Tcl_AppendResult(interp, usmg, NULL);
1663 		return TCL_ERROR;
1664 	    }
1665 	    we = Tcl_GetStringFromObj(dmsPtr[7], NULL);
1666 	    if ( !strchr("WwEe", *we) ) {
1667 		Tcl_AppendResult(interp, usmg, NULL);
1668 		return TCL_ERROR;
1669 	    }
1670 	    lat = (latd + latm / 60.0 + lats / 3600.0)
1671 		* (*ns == 'N' || *ns == 'n' ? 1.0 : -1.0);
1672 	    lon = (lond + lonm / 60.0 + lons / 3600.0)
1673 		* (*we == 'E' || *we == 'e' ? 1.0 : -1.0);
1674 	    Tcl_SetObjResult(interp,
1675 		    Tclgeomap_NewGeoPtObj(GeoPtFmDeg(lat, lon)));
1676 	    break;
1677 	case 6:
1678 	    /*
1679 	     * Form is {deg min NorS deg min WorE}
1680 	     */
1681 
1682 	    if (Tcl_GetDoubleFromObj(interp, dmsPtr[0], &latd) != TCL_OK
1683 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[1], &latm) != TCL_OK
1684 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[3], &lond) != TCL_OK
1685 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[4], &lonm)
1686 		    != TCL_OK) {
1687 		return TCL_ERROR;
1688 	    }
1689 	    ns = Tcl_GetStringFromObj(dmsPtr[2], NULL);
1690 	    if ( !strchr("NnSs", *ns) ) {
1691 		Tcl_AppendResult(interp, usmg, NULL);
1692 		return TCL_ERROR;
1693 	    }
1694 	    we = Tcl_GetStringFromObj(dmsPtr[5], NULL);
1695 	    if ( !strchr("WeEe", *we) ) {
1696 		Tcl_AppendResult(interp, usmg, NULL);
1697 		return TCL_ERROR;
1698 	    }
1699 	    lat = (latd + latm / 60.0)
1700 		* (*ns == 'N' || *ns == 'n' ? 1.0 : -1.0);
1701 	    lon = (lond + lonm / 60.0)
1702 		* (*we == 'E' || *we == 'e' ? 1.0 : -1.0);
1703 	    Tcl_SetObjResult(interp,
1704 		    Tclgeomap_NewGeoPtObj(GeoPtFmDeg(lat, lon)));
1705 	    break;
1706 	case 4:
1707 	    /*
1708 	     * Form is {deg NorS deg WorE}
1709 	     */
1710 
1711 	    if (Tcl_GetDoubleFromObj(interp, dmsPtr[0], &latd) != TCL_OK
1712 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[2], &lond)
1713 		    != TCL_OK) {
1714 		return TCL_ERROR;
1715 	    }
1716 	    ns = Tcl_GetStringFromObj(dmsPtr[1], NULL);
1717 	    if ( !strchr("NnSs", *ns) ) {
1718 		Tcl_AppendResult(interp, usmg, NULL);
1719 		return TCL_ERROR;
1720 	    }
1721 	    we = Tcl_GetStringFromObj(dmsPtr[3], NULL);
1722 	    if ( !strchr("WwEe", *ns) ) {
1723 		Tcl_AppendResult(interp, usmg, NULL);
1724 		return TCL_ERROR;
1725 	    }
1726 	    lat = latd * (*ns == 'N' || *ns == 'n' ? 1.0 : -1.0);
1727 	    lon = lond * (*we == 'E' || *we == 'e' ? 1.0 : -1.0);
1728 	    Tcl_SetObjResult(interp,
1729 		    Tclgeomap_NewGeoPtObj(GeoPtFmDeg(lat, lon)));
1730 	    break;
1731 	case 3:
1732 	    /*
1733 	     * Form is {deg min sec}
1734 	     */
1735 
1736 	    if (Tcl_GetDoubleFromObj(interp, dmsPtr[0], &d) != TCL_OK
1737 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[1], &m) != TCL_OK
1738 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[2], &s) != TCL_OK) {
1739 		return TCL_ERROR;
1740 	    }
1741 	    deg = d + m / 60.0 + s / 3600.0;
1742 	    Tcl_SetObjResult(interp, Tcl_NewDoubleObj(deg));
1743 	    break;
1744 	case 2:
1745 	    /*
1746 	     * Form is {deg min}
1747 	     */
1748 
1749 	    if (   Tcl_GetDoubleFromObj(interp, dmsPtr[0], &d) != TCL_OK
1750 		    || Tcl_GetDoubleFromObj(interp, dmsPtr[1], &m) != TCL_OK) {
1751 		return TCL_ERROR;
1752 	    }
1753 	    deg = d + m / 60.0;
1754 	    Tcl_SetObjResult(interp, Tcl_NewDoubleObj(deg));
1755 	    break;
1756 	default:
1757 	    Tcl_AppendResult(interp, usmg, NULL);
1758 	    return TCL_ERROR;
1759     }
1760     CKFREE(usmg);
1761     return TCL_OK;
1762 }
1763 
1764 /*
1765  *------------------------------------------------------------------------
1766  *
1767  * decToDM --
1768  *
1769  *	This is the callback for the dectodm command.
1770  *
1771  * Results:
1772  *	A standard Tcl result.
1773  *
1774  * Side effects:
1775  * 	See the user documentation.
1776  *
1777  *------------------------------------------------------------------------
1778  */
1779 
1780 int
decToDM(clientData,interp,objc,objv)1781 decToDM(clientData, interp, objc, objv)
1782     ClientData clientData;	/* Not used */
1783     Tcl_Interp *interp;		/* Current interpreter */
1784     int objc;			/* Number of arguments */
1785     Tcl_Obj *CONST objv[];	/* Argument objects */
1786 {
1787     GeoPt geoPt;		/* Input point */
1788     double lat, lon;		/* Latitude and longitude from geoPt */
1789     double deg;			/* Input angle */
1790     char *ns, *we;		/* Hemisphere specifiers */
1791     int latDeg, latMin;		/* Latitude of input point */
1792     int lonDeg, lonMin;		/* Longitude of input point */
1793     int iDeg, iMin;		/* Input angle */
1794     Tcl_Obj *result;		/* Result */
1795 
1796     result = Tcl_NewObj();
1797     if (objc != 2) {
1798 	Tcl_WrongNumArgs(interp, 1, objv, "degOR{lat lon}");
1799 	return TCL_ERROR;
1800     }
1801     if (Tclgeomap_GetGeoPtFromObj(NULL, objv[1], &geoPt) == TCL_OK) {
1802 	/*
1803 	 * Input is a geographic point.
1804 	 */
1805 
1806 	GeoPtGetDeg(geoPt, &lat, &lon);
1807 	ns = lat < 0.0 ? "S" : "N";
1808 	lat = fabs(lat);
1809 	latDeg = (int)lat;
1810 	latMin = (int)(fmod(lat * 60.0, 60.0) + 0.5);
1811 	if (latMin == 60) {
1812 	    latDeg++;
1813 	    latMin = 0;
1814 	}
1815 	we = lon < 0.0 ? "W" : "E";
1816 	lon = fabs(lon);
1817 	lonDeg = (int)lon;
1818 	lonMin = (int)(fmod(lon * 60.0, 60.0) + 0.5);
1819 	if (lonMin == 60) {
1820 	    lonDeg++;
1821 	    lonMin = 0;
1822 	}
1823 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(latDeg));
1824 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(latMin));
1825 	Tcl_ListObjAppendElement(interp, result, Tcl_NewStringObj(ns, -1));
1826 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(lonDeg));
1827 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(lonMin));
1828 	Tcl_ListObjAppendElement(interp, result, Tcl_NewStringObj(we, -1));
1829 
1830     } else if (Tcl_GetDoubleFromObj(NULL, objv[1], &deg) == TCL_OK) {
1831 	/*
1832 	 * Input is a number of degrees.
1833 	 */
1834 
1835 	deg = fabs(deg);
1836 	iDeg = (int)deg;
1837 	iMin = (int)(fmod(deg * 60.0, 60.0) + 0.5);
1838 	if (iMin == 60) {
1839 	    iDeg++;
1840 	    iMin = 0;
1841 	}
1842 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(iDeg));
1843 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(iMin));
1844 
1845     } else {
1846 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetString(objv[0]),
1847 		" degOR{lat lon}", NULL);
1848 	return TCL_ERROR;
1849     }
1850     Tcl_SetObjResult(interp, result);
1851     return TCL_OK;
1852 }
1853 
1854 /*
1855  *------------------------------------------------------------------------
1856  *
1857  * decToDMS --
1858  *
1859  *	This is the callback for the dectodms command.
1860  *
1861  * Results:
1862  *	A standard Tcl result.
1863  *
1864  * Side effects:
1865  * 	See the user documentation.
1866  *
1867  *------------------------------------------------------------------------
1868  */
1869 
1870 int
decToDMS(clientData,interp,objc,objv)1871 decToDMS(clientData, interp, objc, objv)
1872     ClientData clientData;	/* Not used */
1873     Tcl_Interp *interp;		/* Current interpreter */
1874     int objc;			/* Number of arguments */
1875     Tcl_Obj *CONST objv[];	/* Argument objects */
1876 {
1877     GeoPt geoPt;		/* Input point */
1878     double lat, lon;		/* Latitude and longitude from geoPt */
1879     double deg, min;		/* Input angle */
1880     char *ns, *we;		/* Hemisphere specifiers */
1881     int latDeg, latMin, latSec;	/* Latitude of input point */
1882     int lonDeg, lonMin, lonSec;	/* Longitude of input point */
1883     int iDeg, iMin, iSec;	/* Input angle */
1884     Tcl_Obj *result;		/* Result */
1885 
1886     result = Tcl_NewObj();
1887     if (objc != 2) {
1888 	Tcl_WrongNumArgs(interp, 1, objv, " degOR{lat lon}");
1889 	return TCL_ERROR;
1890     }
1891     if (Tclgeomap_GetGeoPtFromObj(NULL, objv[1], &geoPt) == TCL_OK) {
1892 	/*
1893 	 * Input is a geographic point.
1894 	 */
1895 
1896 	GeoPtGetDeg(geoPt, &lat, &lon);
1897 	ns = lat < 0.0 ? "S" : "N";
1898 	lat = fabs(lat);
1899 	latDeg = (int)lat;
1900 	min = fmod(lat * 60.0, 60.0);
1901 	latSec = (int)(fmod(min * 60.0, 60.0) + 0.5);
1902 	latMin = (int)min;
1903 	if (latSec == 60) {
1904 	    latMin++;
1905 	    latSec = 0;
1906 	}
1907 	if (latMin == 60) {
1908 	    latDeg++;
1909 	    latMin = 0;
1910 	}
1911 	we = lon < 0.0 ? "W" : "E";
1912 	lon = fabs(lon);
1913 	lonDeg = (int)lon;
1914 	min = fmod(lon * 60.0, 60.0);
1915 	lonSec = (int)(fmod(min * 60.0, 60.0) + 0.5);
1916 	lonMin = (int)min;
1917 	if (lonSec == 60) {
1918 	    lonMin++;
1919 	    lonSec = 0;
1920 	}
1921 	if (lonMin == 60) {
1922 	    lonDeg++;
1923 	    lonMin = 0;
1924 	}
1925 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(latDeg));
1926 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(latMin));
1927 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(latSec));
1928 	Tcl_ListObjAppendElement(interp, result, Tcl_NewStringObj(ns, -1));
1929 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(lonDeg));
1930 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(lonMin));
1931 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(lonSec));
1932 	Tcl_ListObjAppendElement(interp, result, Tcl_NewStringObj(we, -1));
1933 
1934     } else if (Tcl_GetDoubleFromObj(NULL, objv[1], &deg) == TCL_OK) {
1935 	/*
1936 	 * Input is a number of degrees.
1937 	 */
1938 
1939 	deg = fabs(deg);
1940 	iDeg = (int)deg;
1941 	min = fmod(deg * 60.0, 60.0);
1942 	iSec = (int)(fmod(min * 60.0, 60.0) + 0.5);
1943 	iMin = (int)min;
1944 	if (iSec == 60) {
1945 	    iMin++;
1946 	    iSec = 0;
1947 	}
1948 	if (iMin == 60) {
1949 	    iDeg++;
1950 	    iMin = 0;
1951 	}
1952 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(iDeg));
1953 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(iMin));
1954 	Tcl_ListObjAppendElement(interp, result, Tcl_NewIntObj(iSec));
1955 
1956     } else {
1957 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetString(objv[0]),
1958 		" degOR{lat lon}", NULL);
1959 	return TCL_ERROR;
1960     }
1961     Tcl_SetObjResult(interp, result);
1962     return TCL_OK;
1963 }
1964 
1965 /*
1966  *------------------------------------------------------------------------
1967  *
1968  * cartg --
1969  *
1970  *	This is the callback for the cartg command.
1971  *
1972  * Results:
1973  *	A standard Tcl result.
1974  *
1975  * Side effects:
1976  * 	See the user documentation.
1977  *
1978  *------------------------------------------------------------------------
1979  */
1980 
1981 int
cartg(clientData,interp,objc,objv)1982 cartg(clientData, interp, objc, objv)
1983     ClientData clientData;	/* Not used */
1984     Tcl_Interp *interp;		/* Current interpreter */
1985     int objc;			/* Number of arguments */
1986     Tcl_Obj *CONST objv[];	/* Argument objects */
1987 {
1988     double scale;		/* Scale as a double */
1989     char cscale[13];		/* Scale as a string in cartographic form */
1990     int numer, denom;		/* Integers from cscale */
1991 
1992     if (objc != 2) {
1993 	Tcl_WrongNumArgs(interp, 1, objv, "{1:XXXXXXX}ORdoubleValue");
1994 	return TCL_ERROR;
1995     }
1996     if (Tcl_GetDoubleFromObj(interp, objv[1], &scale) == TCL_OK) {
1997 	/*
1998 	 * Input is a double value.
1999 	 */
2000 
2001 	sprintf(cscale, "1:%-9d", (int)(1.0 / scale + 0.5));
2002 	Tcl_SetObjResult(interp, Tcl_NewStringObj(cscale, -1));
2003 	return TCL_OK;
2004     } else if (sscanf(Tcl_GetStringFromObj(objv[1], NULL), "%d:%d",
2005 		&numer, &denom)) {
2006 	/*
2007 	 * Input is a string of form "int:int"
2008 	 */
2009 
2010 	scale = (double)numer / denom;
2011 	Tcl_SetObjResult(interp, Tcl_NewDoubleObj(scale));
2012 	return TCL_OK;
2013     } else {
2014 	Tcl_AppendResult(interp, "Usage: ", Tcl_GetStringFromObj(objv[0], NULL),
2015 		" {1:XXXXXXX}ORdoubleValue", NULL);
2016 	return TCL_ERROR;
2017     }
2018 }
2019 
2020 /*
2021  *------------------------------------------------------------------------
2022  *
2023  * centroid --
2024  *
2025  *	This is the callback for the centroid command.
2026  *
2027  * Results:
2028  *	A standard Tcl result.
2029  *
2030  * Side effects:
2031  * 	See the user documentation.
2032  *
2033  *------------------------------------------------------------------------
2034  */
2035 
2036 int
centroid(clientData,interp,objc,objv)2037 centroid(clientData, interp, objc, objv)
2038     ClientData clientData;	/* Not used */
2039     Tcl_Interp *interp;		/* Current interpreter */
2040     int objc;			/* Number of arguments */
2041     Tcl_Obj *CONST objv[];	/* Argument objects */
2042 {
2043     Tcl_Obj **ptsPtr;		/* List of lat-lon values from command line */
2044     int n;			/* Number of lat-lon values */
2045     Tcl_Obj **p, **pe;
2046     GeoPt geoPt;		/* lat-lon from command line */
2047     double lat, lon;		/* Latitude and longitude of geoPt, radians */
2048     double cos_lat;		/* Computational constant */
2049     double x, y, z;		/* Geocentric cartesian coordinates of geoPt */
2050     GeoPt centroid;		/* Result */
2051 
2052     if (objc != 2) {
2053 	Tcl_WrongNumArgs(interp, 1, objv, "point_list");
2054 	return TCL_ERROR;
2055     }
2056     if (Tcl_ListObjGetElements(interp, objv[1], &n, &ptsPtr) != TCL_OK) {
2057 	Tcl_AppendResult(interp, "Could not get points list", NULL);
2058 	return TCL_ERROR;
2059     }
2060 
2061     for (p = ptsPtr, pe = p + n, x = y = z = 0.0; p < pe; p++) {
2062 	if (Tclgeomap_GetGeoPtFromObj(interp, *p, &geoPt) != TCL_OK) {
2063 	    return TCL_ERROR;
2064 	}
2065 	GeoPtGetRad(geoPt, &lat, &lon);
2066 	cos_lat = cos(lat);
2067 	x += cos_lat * sin(lon);
2068 	y += cos_lat * cos(lon);
2069 	z += sin(lat);
2070     }
2071     x /= n;
2072     y /= n;
2073     z /= n;
2074     lon = atan2(x, y);
2075     lat = asin(z);
2076     centroid = GeoPtFmRad(lat, lon);
2077     Tcl_SetObjResult(interp, Tclgeomap_NewGeoPtObj(centroid));
2078     return TCL_OK;
2079 }
2080 
2081 /*
2082  *-------------------------------------------------------------------------
2083  *
2084  * jul_tm --
2085  *
2086  *	This is the callback for the jul_tm command.
2087  *
2088  * Results:
2089  * 	A standard Tcl result.
2090  *
2091  * Side effects:
2092  * 	The interpreter result is set to a list of form {day seconds}
2093  * 	corresponding to year, month, day, hour, minute, second on
2094  * 	the command line.
2095  *
2096  *-------------------------------------------------------------------------
2097  */
2098 
2099 static int
jul_tm(clientData,interp,objc,objv)2100 jul_tm(clientData, interp, objc, objv)
2101     ClientData clientData;		/* Not used */
2102     Tcl_Interp *interp;			/* Current interpreter */
2103     int objc;				/* Number of arguments */
2104     Tcl_Obj *CONST objv[];		/* Argument objects */
2105 {
2106     struct GeoTime_Jul jul;		/* Days and seconds in cal */
2107     Tcl_Obj *result;			/* Return value {day sec} */
2108     int year, month, day, hour, minute;
2109     double second;
2110 
2111     if (objc != 7) {
2112 	Tcl_WrongNumArgs(interp, 1, objv, "year month day hour minute second");
2113 	return TCL_ERROR;
2114     }
2115     if (Tcl_GetIntFromObj(interp, objv[1], &year)) {
2116 	return TCL_ERROR;
2117     }
2118     if (Tcl_GetIntFromObj(interp, objv[2], &month)) {
2119 	return TCL_ERROR;
2120     }
2121     if (Tcl_GetIntFromObj(interp, objv[3], &day)) {
2122 	return TCL_ERROR;
2123     }
2124     if (Tcl_GetIntFromObj(interp, objv[4], &hour)) {
2125 	return TCL_ERROR;
2126     }
2127     if (Tcl_GetIntFromObj(interp, objv[5], &minute)) {
2128 	return TCL_ERROR;
2129     }
2130     if (Tcl_GetDoubleFromObj(interp, objv[6], &second)) {
2131 	return TCL_ERROR;
2132     }
2133     jul = GeoTime_CalToJul(GeoTime_CalSet(year, month, day,
2134 		hour, minute, second));
2135     result = Tcl_NewObj();
2136     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(jul.day));
2137     Tcl_ListObjAppendElement(NULL, result, Tcl_NewDoubleObj(jul.second));
2138     Tcl_SetObjResult(interp, result);
2139     return TCL_OK;
2140 }
2141 
2142 /*
2143  *-------------------------------------------------------------------------
2144  *
2145  * cal_tm --
2146  *
2147  *	This is the callback for the cal_tm command.
2148  *
2149  * Results:
2150  * 	A standard Tcl result.
2151  *
2152  * Side effects:
2153  * 	The interpreter result is set to a list of form
2154  * 	{year month day hour minute second}
2155  * 	corresponding to day and seconds on the command line.
2156  *
2157  *-------------------------------------------------------------------------
2158  */
2159 
2160 static int
cal_tm(clientData,interp,objc,objv)2161 cal_tm(clientData, interp, objc, objv)
2162     ClientData clientData;		/* Not used */
2163     Tcl_Interp *interp;			/* Current interpreter */
2164     int objc;				/* Number of arguments */
2165     Tcl_Obj *CONST objv[];		/* Argument objects */
2166 {
2167     int day;
2168     double sec;
2169 
2170     if (objc != 3) {
2171 	Tcl_WrongNumArgs(interp, 1, objv, "days seconds");
2172 	return TCL_ERROR;
2173     }
2174     if (Tcl_GetIntFromObj(interp, objv[1], &day)) {
2175 	return TCL_ERROR;
2176     }
2177     if (Tcl_GetDoubleFromObj(interp, objv[2], &sec)) {
2178 	return TCL_ERROR;
2179     }
2180     Tcl_SetObjResult(interp,
2181 	    Tclgeomap_JulDayToList(GeoTime_JulSet(day, sec)));
2182     return TCL_OK;
2183 }
2184 
2185 /*
2186  *-------------------------------------------------------------------------
2187  *
2188  * Tclgeomap_JulDayToList --
2189  *
2190  * 	The function returns the calendar representation of a Julian day
2191  * 	in list format.
2192  *
2193  * Results:
2194  * 	See the user documentation.
2195  *
2196  * Side effects:
2197  * 	See the user documentation.
2198  *
2199  *-------------------------------------------------------------------------
2200  */
2201 
2202 Tcl_Obj *
Tclgeomap_JulDayToList(jul)2203 Tclgeomap_JulDayToList(jul)
2204     struct GeoTime_Jul jul;
2205 {
2206     struct GeoTime_Cal cal;	/* Cal representation of jul */
2207     Tcl_Obj *result;		/* Return value {year mon day hr min sec} */
2208 
2209     cal = GeoTime_JulToCal(jul);
2210     result = Tcl_NewObj();
2211     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(cal.year));
2212     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(cal.month));
2213     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(cal.day));
2214     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(cal.hour));
2215     Tcl_ListObjAppendElement(NULL, result, Tcl_NewIntObj(cal.minute));
2216     Tcl_ListObjAppendElement(NULL, result, Tcl_NewDoubleObj(cal.second));
2217     return result;
2218 }
2219 
2220 /*
2221  *-------------------------------------------------------------------------
2222  *
2223  * Tclgeomap_AppendTime --
2224  *
2225  *	This convenience function converts a Julian date representation
2226  *	to {year month day hour minute second}.
2227  *
2228  * Results:
2229  * 	None.
2230  *
2231  * Side effects:
2232  * 	See the user documentation.
2233  *
2234  *-------------------------------------------------------------------------
2235  */
2236 
2237 void
Tclgeomap_AppendTime(list,jul)2238 Tclgeomap_AppendTime(list, jul)
2239     Tcl_Obj *list;
2240     struct GeoTime_Jul jul;
2241 {
2242     Tcl_ListObjAppendElement(NULL, list, Tclgeomap_JulDayToList(jul));
2243 }
2244