Geo Stars Library  0.9.3
Geodetic and Astrometry library
geoEllips.c
Go to the documentation of this file.
1 /*!
2 \file geoEllips.c
3 \brief This file contains geo library location setup and initialization
4  of ellipsoid values. There are 23 ellipsoids to choose from, of which
5  the WGS 84 ellipsoid is the default ellipsoid.
6 */
7 
8 
9 #include <stdio.h>
10 #include <math.h>
11 #include <string.h>
12 #include <time.h>
13 #include "geoStars.h"
14 
15 #ifdef WIN32
16 #include <windows.h>
17 #endif
18 
19 #ifdef WIN32
20 #ifndef GEO_LIB
21 BOOL APIENTRY
22 DllMain(HANDLE hModule, DWORD ul_reason_for_call, LPVOID lpReserved) {
23  switch (ul_reason_for_call) {
24  case DLL_PROCESS_ATTACH :
25  // FreeImage_Initialise(FALSE);
26  break;
27 
28  case DLL_PROCESS_DETACH :
29  // FreeImage_DeInitialise();
30  break;
31 
32  case DLL_THREAD_ATTACH :
33  case DLL_THREAD_DETACH :
34  break;
35  }
36 
37  return TRUE;
38 }
39 #endif // GEO_LIB
40 #endif // WIN32
41 //-----------------------------------------------------------------
42 /*!
43 \brief This structure array contains the all of the ellipsoids used in this library
44 
45 */
47 {
48 /* Ellipsoid Name ID Major Axis Inverse Flattening */
49 /* ----------------------------------------------------------------------------------------- */
50 { "Default Ellipsoid (WGS 1984) ","00", 6378137.0 ,298.257223563 }, /*!< Default datum (WGS84) */
51 { "Airy 1830 ","AA", 6377563.396 ,299.3249646 },
52 { "Australian National ","AN", 6378160.0 ,298.25 },
53 { "Bessel 1841 ","BR", 6377397.155 ,299.1528128 },
54 { "Bessel 1841 (Namibia) ","BN", 6377483.865 ,299.1528128 },
55 { "Clarke 1866 ","CC", 6378206.4 ,294.9786982 },
56 { "Clarke 1880 ","CD", 6378249.145 ,293.465 },
57 { "Everest (Brunei, E. Malaysia) ","EB", 6377298.556 ,300.8017 },
58 { "Everest 1830 ","EA", 6377276.345 ,300.8017 },
59 { "Everest 1956 (India and Nepal) ","EC", 6377301.243 ,300.8017 },
60 { "Everest (Pakistan) ","EF", 6377309.613 ,300.8017 },
61 { "Everest 1948 (W. Malaysia and Singapore)","EE", 6377304.063 ,300.8017 },
62 { "Everest 1969 (W. Malaysia) ","ED", 6377295.664 ,300.8017 },
63 { "Geodetic Reference System 1980 ","RF", 6378137.0 ,298.257222101 },
64 { "Helmert 1906 ","HE", 6378200.0 ,298.3 },
65 { "Hough 1960 ","HO", 6378270.0 ,297.0 },
66 { "Indonesian 1974 ","ID", 6378160.0 ,298.247 },
67 { "International 1924 ","IN", 6378388.0 ,297.0 },
68 { "Krassovsky 1940 ","KA", 6378245.0 ,298.3 },
69 { "Modified Airy ","AM", 6377340.189 ,299.3249646 },
70 { "Modified Fischer 1960 ","FA", 6378155.0 ,298.3 },
71 { "South American 1969 ","SA", 6378160.0 ,298.25 },
72 { "WGS 1972 ","WD", 6378135.0 ,298.26 },
73 { "WGS 1984 ","WE", 6378137.0 ,298.257223563 }
74 
75 };
76 
77 
78 
79 
80 //-----------------------------------------------------------------
81 /*!
82 \brief This routine needs to be called when a site (or location) is initialized.
83 Several of the routines use the information in the structure that this routine fills
84 
85 \param GEO_LOCATION *l
86 \param double lat
87 \param double lon
88 \param double hgt
89 \param int datum
90 \param char *name //Name of the site location
91 \retval GEO_OK on success
92 \retval GEO_ERROR on error
93 */
94 
95 
96 int DLL_CALLCONV geoInitLocation(GEO_LOCATION *l, double lat, double lon, double hgt, int datum, char *name)
97 {
98  double N,a,e2,ee2, b, flat,Nh;
99  double dec;
100  // struct tm *newtime;
101  // time_t aclock;
102 
103 
104  /* Validate and get datum values */
105  if(datum > GEO_DATUM_MAX) return(GEO_ERROR);
106  geoGetEllipsoid(&a,&b,&e2,&ee2,&flat,datum); /* Calls to read in "to" ellipsoid data */
107 
108  /* Initialize ellipsoid values in the location descriptor */
109  l->datum.a = a;
110  l->datum.b = b;
111  l->datum.e2 = e2;
112  l->datum.flat = flat ;
113  l->datum.ee2 = ee2;
114  l->datum.m1e2 = 1.0 - e2;
115 
116 
117  /* Initialize lat/lon values */
118  if(lon>180.0)lon=lon-360.0; // normalize to -180 to +180
119  l->lat = lat;
120  l->lon = lon;
121  l->hgt = hgt;
122 
123  l->rlat = lat * DEG_TO_RAD; /* make radians */
124  l->rlon = lon * DEG_TO_RAD;
125 
126  /* Precompute sin/cos values */
127 #ifdef HAVE_SINCOS
128  sincos(l->rlat,&l->slat,&l->clat);
129  sincos(l->rlon,&l->slon,&l->clon);
130 #else
131  l->slat = sin(l->rlat);
132  l->slon = sin(l->rlon);
133  l->clat = cos(l->rlat);
134  l->clon = cos(l->rlon);
135 #endif
136  l->tlat = tan(l->rlat);
137  l->clonclat = l->clon * l->clat;
138  l->slonslat = l->slon * l->slat;
139  l->clonslat = l->clon * l->slat;
140  l->slonclat = l->slon * l->clat;
141 
142  /* Compute geocentric coordinates */
143 
144  /* Compute the radius of curvature */
145  N = a / (sqrt(1.0 - e2 * pow(l->slat,2.0)));
146 
147  /* Compute the EFG(XYZ) coordinates (earth centered) */
148  Nh = N + hgt;
149  l->e = Nh * l->clonclat; //l->clat * l->clon;
150  l->f = Nh * l->slonclat; //l->clat * l->slon;
151  l->g = (N * (l->datum.m1e2) + hgt) * l->slat;
152  l->efg[GEO_E] = l->e;
153  l->efg[GEO_F] = l->f;
154  l->efg[GEO_G] = l->g;
155 
156  /* Compute the Magnetic Declination */
157  geoMagFillDec(l,&dec);
158 
159  /* Time Zone */
160  l->timezone = 0.0; //use geoSetTimeZone to set this value
161 
162  /* Save datum and site info */
163  l->datum.datum_num = datum;
164  strcpy(l->name, name);
165 
166  return(GEO_OK);
167 }
168 
169 void DLL_CALLCONV geoSetTimeZone(GEO_LOCATION *l, double tz, int dst)
170 {
171  l->timezone = tz;
172  l->dst = dst;
173 }
174 
175 //-----------------------------------------------------------------
176 /*!
177 \brief This routine computes essential datum values from basic parameters
178 obtained from the \e ellips structure.
179 
180 \param double *a : Major axis ( Meters )
181 \param double *b : Minor axis ( Meters )
182 \param double *e2 : Eccentricity
183 \param double *ee2 : Eccentricity prime
184 \param double *f : Flattening
185 \param int datum : Datum to use
186 \return nothing
187 */
188 
189 void DLL_CALLCONV geoGetEllipsoid(double *a,double *b,double *e2,double *ee2,double *f,int datum)
190 {
191  /* - Major Axis. */
192  *a = ellips[datum].a;
193 
194  /* - Earth Flattening. */
195  *f = GEO_FL(ellips[datum].f1);
196 
197  /* - Minor axis value. */
198  *b = GEO_B(ellips[datum].a, (ellips[datum].f1));
199 
200  /* - Eccentricity squared. */
201  *e2 = GEO_E2(ellips[datum].a, (ellips[datum].f1));
202 
203  /* - Eccentricity squared prime. */
204  *ee2 = GEO_E2P(ellips[datum].a, (ellips[datum].f1));
205 
206 
207 #ifdef GEO_DEBUG
208  printf("DEBUG: a=%f, b=%f, f=%f, e2=%f, e2p=%f \n",*a, *b, *f, *e2, *ee2);
209 #endif
210 }
211 
212 
213 /*
214  Deflection of the vertical
215  latitude longitude Xi Eta Hor Lap
216  Station Name ddd mm ss.sssss ddd mm ss.sssss arc-sec arc-sec arc-sec
217  USER LOCATION 31 17 32.00000 85 51 3.00000 -4.22 0.30 -0.18
218 
219 
220 
221 */
222