dRonin  adbada4
dRonin firmware
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
tranmerc.c
Go to the documentation of this file.
1 /***************************************************************************/
2 /* RSC IDENTIFIER: TRANSVERSE MERCATOR
3  *
4  * ABSTRACT
5  *
6  * This component provides conversions between Geodetic coordinates
7  * (latitude and longitude) and Transverse Mercator projection coordinates
8  * (easting and northing).
9  *
10  * ERROR HANDLING
11  *
12  * This component checks parameters for valid values. If an invalid value
13  * is found the error code is combined with the current error code using
14  * the bitwise or. This combining allows multiple error codes to be
15  * returned. The possible error codes are:
16  *
17  * TRANMERC_NO_ERROR : No errors occurred in function
18  * TRANMERC_LAT_ERROR : Latitude outside of valid range
19  * (-90 to 90 degrees)
20  * TRANMERC_LON_ERROR : longitude outside of valid range
21  * (-180 to 360 degrees, and within
22  * +/-90 of Central Meridian)
23  * TRANMERC_EASTING_ERROR : Easting outside of valid range
24  * (depending on ellipsoid and
25  * projection parameters)
26  * TRANMERC_NORTHING_ERROR : Northing outside of valid range
27  * (depending on ellipsoid and
28  * projection parameters)
29  * TRANMERC_ORIGIN_LAT_ERROR : Origin latitude outside of valid range
30  * (-90 to 90 degrees)
31  * TRANMERC_CENT_MER_ERROR : Central meridian outside of valid range
32  * (-180 to 360 degrees)
33  * TRANMERC_A_ERROR : Semi-major axis less than or equal to zero
34  * TRANMERC_INV_F_ERROR : Inverse flattening outside of valid range
35  * (250 to 350)
36  * TRANMERC_SCALE_FACTOR_ERROR : Scale factor outside of valid
37  * range (0.3 to 3.0)
38  * TM_LON_WARNING : Distortion will result if longitude is more
39  * than 9 degrees from the Central Meridian
40  *
41  * REUSE NOTES
42  *
43  * TRANSVERSE MERCATOR is intended for reuse by any application that
44  * performs a Transverse Mercator projection or its inverse.
45  *
46  * REFERENCES
47  *
48  * Further information on TRANSVERSE MERCATOR can be found in the
49  * Reuse Manual.
50  *
51  * TRANSVERSE MERCATOR originated from :
52  * U.S. Army Topographic Engineering Center
53  * Geospatial Information Division
54  * 7701 Telegraph Road
55  * Alexandria, VA 22310-3864
56  *
57  * LICENSES
58  *
59  * None apply to this component.
60  *
61  * RESTRICTIONS
62  *
63  * TRANSVERSE MERCATOR has no restrictions.
64  *
65  * ENVIRONMENT
66  *
67  * TRANSVERSE MERCATOR was tested and certified in the following
68  * environments:
69  *
70  * 1. Solaris 2.5 with GCC, version 2.8.1
71  * 2. Windows 95 with MS Visual C++, version 6
72  *
73  * MODIFICATIONS
74  *
75  * Date Description
76  * ---- -----------
77  * 10-02-97 Original Code
78  * 03-02-97 Re-engineered Code
79  *
80  */
81 
82 
83 /***************************************************************************/
84 /*
85  * INCLUDES
86  */
87 
88 #include <math.h>
89 #include "tranmerc.h"
90 
91 /*
92  * math.h - Standard C math library
93  * tranmerc.h - Is for prototype error checking
94  */
95 
96 
97 /***************************************************************************/
98 /* DEFINES
99  *
100  */
101 
102 #define PI 3.14159265358979323e0 /* PI */
103 #define PI_OVER_2 (PI/2.0e0) /* PI over 2 */
104 #define MAX_LAT ((PI * 89.99)/180.0) /* 89.99 degrees in radians */
105 #define MAX_DELTA_int ((PI * 90)/180.0) /* 90 degrees in radians */
106 #define MIN_SCALE_FACTOR 0.3
107 #define MAX_SCALE_FACTOR 3.0
108 
109 #define SPHTMD(Latitude) ((double) (TranMerc_ap * Latitude \
110  - TranMerc_bp * sin(2.e0 * Latitude) + TranMerc_cp * sin(4.e0 * Latitude) \
111  - TranMerc_dp * sin(6.e0 * Latitude) + TranMerc_ep * sin(8.e0 * Latitude) ) )
112 
113 #define SPHSN(Latitude) ((double) (TranMerc_a / sqrt( 1.e0 - TranMerc_es * \
114  pow(sin(Latitude), 2))))
115 
116 #define SPHSR(Latitude) ((double) (TranMerc_a * (1.e0 - TranMerc_es) / \
117  pow(DENOM(Latitude), 3)))
118 
119 #define DENOM(Latitude) ((double) (sqrt(1.e0 - TranMerc_es * pow(sin(Latitude),2))))
120 
121 
122 /**************************************************************************/
123 /* GLOBAL DECLARATIONS
124  *
125  */
126 
127 /* Ellipsoid Parameters, default to WGS 84 */
128 static double TranMerc_a = 6378137.0; /* Semi-major axis of ellipsoid in meters */
129 static double TranMerc_f = 1 / 298.257223563; /* Flattening of ellipsoid */
130 static double TranMerc_es = 0.0066943799901413800; /* Eccentricity (0.08181919084262188000) squared */
131 static double TranMerc_ebs = 0.0067394967565869; /* Second Eccentricity squared */
132 
133 /* Transverse_Mercator projection Parameters */
134 static double TranMerc_Origin_Lat = 0.0; /* Latitude of origin in radians */
135 static double TranMerc_Origin_int = 0.0; /* longitude of origin in radians */
136 static double TranMerc_False_Northing = 0.0; /* False northing in meters */
137 static double TranMerc_False_Easting = 0.0; /* False easting in meters */
138 static double TranMerc_Scale_Factor = 1.0; /* Scale factor */
139 
140 /* Isometeric to geodetic latitude parameters, default to WGS 84 */
141 static double TranMerc_ap = 6367449.1458008;
142 static double TranMerc_bp = 16038.508696861;
143 static double TranMerc_cp = 16.832613334334;
144 static double TranMerc_dp = 0.021984404273757;
145 static double TranMerc_ep = 3.1148371319283e-005;
146 
147 /* Maximum variance for easting and northing values for WGS 84. */
148 static double TranMerc_Delta_Easting = 40000000.0;
149 static double TranMerc_Delta_Northing = 40000000.0;
150 
151 /* These state variables are for optimization purposes. The only function
152  * that should modify them is Set_Tranverse_Mercator_Parameters. */
153 
154 
155 /************************************************************************/
156 /* FUNCTIONS
157  *
158  */
159 
160 
161 int Set_Transverse_Mercator_Parameters(double a,
162  double f,
163  double Origin_Latitude,
164  double Central_Meridian,
165  double False_Easting,
166  double False_Northing,
167  double Scale_Factor)
168 
169 { /* BEGIN Set_Tranverse_Mercator_Parameters */
170  /*
171  * The function Set_Tranverse_Mercator_Parameters receives the ellipsoid
172  * parameters and Tranverse Mercator projection parameters as inputs, and
173  * sets the corresponding state variables. If any errors occur, the error
174  * code(s) are returned by the function, otherwise TRANMERC_NO_ERROR is
175  * returned.
176  *
177  * a : Semi-major axis of ellipsoid, in meters (input)
178  * f : Flattening of ellipsoid (input)
179  * Origin_Latitude : Latitude in radians at the origin of the (input)
180  * projection
181  * Central_Meridian : longitude in radians at the center of the (input)
182  * projection
183  * False_Easting : Easting/X at the center of the projection (input)
184  * False_Northing : Northing/Y at the center of the projection (input)
185  * Scale_Factor : Projection scale factor (input)
186  */
187 
188  double tn; /* True Meridianal distance constant */
189  double tn2;
190  double tn3;
191  double tn4;
192  double tn5;
193  double dummy_northing;
194  double TranMerc_b; /* Semi-minor axis of ellipsoid, in meters */
195  double inv_f = 1 / f;
196  int Error_Code = TRANMERC_NO_ERROR;
197 
198  if (a <= 0.0)
199  { /* Semi-major axis must be greater than zero */
200  Error_Code |= TRANMERC_A_ERROR;
201  }
202  if ((inv_f < 250) || (inv_f > 350))
203  { /* Inverse flattening must be between 250 and 350 */
204  Error_Code |= TRANMERC_INV_F_ERROR;
205  }
206  if ((Origin_Latitude < -PI_OVER_2) || (Origin_Latitude > PI_OVER_2))
207  { /* origin latitude out of range */
208  Error_Code |= TRANMERC_ORIGIN_LAT_ERROR;
209  }
210  if ((Central_Meridian < -PI) || (Central_Meridian > (2*PI)))
211  { /* origin longitude out of range */
212  Error_Code |= TRANMERC_CENT_MER_ERROR;
213  }
214  if ((Scale_Factor < MIN_SCALE_FACTOR) || (Scale_Factor > MAX_SCALE_FACTOR))
215  {
216  Error_Code |= TRANMERC_SCALE_FACTOR_ERROR;
217  }
218  if (!Error_Code)
219  { /* no errors */
220  TranMerc_a = a;
221  TranMerc_f = f;
222  TranMerc_Origin_Lat = Origin_Latitude;
223  if (Central_Meridian > PI)
224  Central_Meridian -= (2*PI);
225  TranMerc_Origin_int = Central_Meridian;
226  TranMerc_False_Northing = False_Northing;
227  TranMerc_False_Easting = False_Easting;
228  TranMerc_Scale_Factor = Scale_Factor;
229 
230  /* Eccentricity Squared */
231  TranMerc_es = 2 * TranMerc_f - TranMerc_f * TranMerc_f;
232  /* Second Eccentricity Squared */
233  TranMerc_ebs = (1 / (1 - TranMerc_es)) - 1;
234 
235  TranMerc_b = TranMerc_a * (1 - TranMerc_f);
236  /*True meridianal constants */
237  tn = (TranMerc_a - TranMerc_b) / (TranMerc_a + TranMerc_b);
238  tn2 = tn * tn;
239  tn3 = tn2 * tn;
240  tn4 = tn3 * tn;
241  tn5 = tn4 * tn;
242 
243  TranMerc_ap = TranMerc_a * (1.e0 - tn + 5.e0 * (tn2 - tn3)/4.e0
244  + 81.e0 * (tn4 - tn5)/64.e0 );
245  TranMerc_bp = 3.e0 * TranMerc_a * (tn - tn2 + 7.e0 * (tn3 - tn4)
246  /8.e0 + 55.e0 * tn5/64.e0 )/2.e0;
247  TranMerc_cp = 15.e0 * TranMerc_a * (tn2 - tn3 + 3.e0 * (tn4 - tn5 )/4.e0) /16.0;
248  TranMerc_dp = 35.e0 * TranMerc_a * (tn3 - tn4 + 11.e0 * tn5 / 16.e0) / 48.e0;
249  TranMerc_ep = 315.e0 * TranMerc_a * (tn4 - tn5) / 512.e0;
250  Convert_Geodetic_To_Transverse_Mercator(MAX_LAT,
251  MAX_DELTA_int + Central_Meridian,
252  &TranMerc_Delta_Easting,
253  &TranMerc_Delta_Northing);
254  Convert_Geodetic_To_Transverse_Mercator(0,
255  MAX_DELTA_int + Central_Meridian,
256  &TranMerc_Delta_Easting,
257  &dummy_northing);
258  TranMerc_Delta_Northing++;
259  TranMerc_Delta_Easting++;
260 
261  } /* END OF if(!Error_Code) */
262  return (Error_Code);
263 } /* END of Set_Transverse_Mercator_Parameters */
264 
265 
266 void Get_Transverse_Mercator_Parameters(double *a,
267  double *f,
268  double *Origin_Latitude,
269  double *Central_Meridian,
270  double *False_Easting,
271  double *False_Northing,
272  double *Scale_Factor)
273 
274 { /* BEGIN Get_Tranverse_Mercator_Parameters */
275  /*
276  * The function Get_Transverse_Mercator_Parameters returns the current
277  * ellipsoid and Transverse Mercator projection parameters.
278  *
279  * a : Semi-major axis of ellipsoid, in meters (output)
280  * f : Flattening of ellipsoid (output)
281  * Origin_Latitude : Latitude in radians at the origin of the (output)
282  * projection
283  * Central_Meridian : longitude in radians at the center of the (output)
284  * projection
285  * False_Easting : Easting/X at the center of the projection (output)
286  * False_Northing : Northing/Y at the center of the projection (output)
287  * Scale_Factor : Projection scale factor (output)
288  */
289 
290  *a = TranMerc_a;
291  *f = TranMerc_f;
292  *Origin_Latitude = TranMerc_Origin_Lat;
293  *Central_Meridian = TranMerc_Origin_int;
294  *False_Easting = TranMerc_False_Easting;
295  *False_Northing = TranMerc_False_Northing;
296  *Scale_Factor = TranMerc_Scale_Factor;
297  return;
298 } /* END OF Get_Tranverse_Mercator_Parameters */
299 
300 
301 
302 int Convert_Geodetic_To_Transverse_Mercator (double Latitude,
303  double longitude,
304  double *Easting,
305  double *Northing)
306 
307 { /* BEGIN Convert_Geodetic_To_Transverse_Mercator */
308 
309  /*
310  * The function Convert_Geodetic_To_Transverse_Mercator converts geodetic
311  * (latitude and longitude) coordinates to Transverse Mercator projection
312  * (easting and northing) coordinates, according to the current ellipsoid
313  * and Transverse Mercator projection coordinates. If any errors occur, the
314  * error code(s) are returned by the function, otherwise TRANMERC_NO_ERROR is
315  * returned.
316  *
317  * Latitude : Latitude in radians (input)
318  * longitude : longitude in radians (input)
319  * Easting : Easting/X in meters (output)
320  * Northing : Northing/Y in meters (output)
321  */
322 
323  double c; /* Cosine of latitude */
324  double c2;
325  double c3;
326  double c5;
327  double c7;
328  double dlam; /* Delta longitude - Difference in longitude */
329  double eta; /* constant - TranMerc_ebs *c *c */
330  double eta2;
331  double eta3;
332  double eta4;
333  double s; /* Sine of latitude */
334  double sn; /* Radius of curvature in the prime vertical */
335  double t; /* Tangent of latitude */
336  double tan2;
337  double tan3;
338  double tan4;
339  double tan5;
340  double tan6;
341  double t1; /* Term in coordinate conversion formula - GP to Y */
342  double t2; /* Term in coordinate conversion formula - GP to Y */
343  double t3; /* Term in coordinate conversion formula - GP to Y */
344  double t4; /* Term in coordinate conversion formula - GP to Y */
345  double t5; /* Term in coordinate conversion formula - GP to Y */
346  double t6; /* Term in coordinate conversion formula - GP to Y */
347  double t7; /* Term in coordinate conversion formula - GP to Y */
348  double t8; /* Term in coordinate conversion formula - GP to Y */
349  double t9; /* Term in coordinate conversion formula - GP to Y */
350  double tmd; /* True Meridional distance */
351  double tmdo; /* True Meridional distance for latitude of origin */
352  int Error_Code = TRANMERC_NO_ERROR;
353  double temp_Origin;
354  double temp_int;
355 
356  if ((Latitude < -MAX_LAT) || (Latitude > MAX_LAT))
357  { /* Latitude out of range */
358  Error_Code|= TRANMERC_LAT_ERROR;
359  }
360  if (longitude > PI)
361  longitude -= (2 * PI);
362  if ((longitude < (TranMerc_Origin_int - MAX_DELTA_int))
363  || (longitude > (TranMerc_Origin_int + MAX_DELTA_int)))
364  {
365  if (longitude < 0)
366  temp_int = longitude + 2 * PI;
367  else
368  temp_int = longitude;
369  if (TranMerc_Origin_int < 0)
370  temp_Origin = TranMerc_Origin_int + 2 * PI;
371  else
372  temp_Origin = TranMerc_Origin_int;
373  if ((temp_int < (temp_Origin - MAX_DELTA_int))
374  || (temp_int > (temp_Origin + MAX_DELTA_int)))
375  Error_Code|= TRANMERC_LON_ERROR;
376  }
377  if (!Error_Code)
378  { /* no errors */
379 
380  /*
381  * Delta longitude
382  */
383  dlam = longitude - TranMerc_Origin_int;
384 
385  if (fabs(dlam) > (9.0 * PI / 180))
386  { /* Distortion will result if longitude is more than 9 degrees from the Central Meridian */
387  Error_Code |= TRANMERC_LON_WARNING;
388  }
389 
390  if (dlam > PI)
391  dlam -= (2 * PI);
392  if (dlam < -PI)
393  dlam += (2 * PI);
394  if (fabs(dlam) < 2.e-10)
395  dlam = 0.0;
396 
397  s = sin(Latitude);
398  c = cos(Latitude);
399  c2 = c * c;
400  c3 = c2 * c;
401  c5 = c3 * c2;
402  c7 = c5 * c2;
403  t = tan (Latitude);
404  tan2 = t * t;
405  tan3 = tan2 * t;
406  tan4 = tan3 * t;
407  tan5 = tan4 * t;
408  tan6 = tan5 * t;
409  eta = TranMerc_ebs * c2;
410  eta2 = eta * eta;
411  eta3 = eta2 * eta;
412  eta4 = eta3 * eta;
413 
414  /* radius of curvature in prime vertical */
415  sn = SPHSN(Latitude);
416 
417  /* True Meridianal Distances */
418  tmd = SPHTMD(Latitude);
419 
420  /* Origin */
421  tmdo = SPHTMD (TranMerc_Origin_Lat);
422 
423  /* northing */
424  t1 = (tmd - tmdo) * TranMerc_Scale_Factor;
425  t2 = sn * s * c * TranMerc_Scale_Factor/ 2.e0;
426  t3 = sn * s * c3 * TranMerc_Scale_Factor * (5.e0 - tan2 + 9.e0 * eta
427  + 4.e0 * eta2) /24.e0;
428 
429  t4 = sn * s * c5 * TranMerc_Scale_Factor * (61.e0 - 58.e0 * tan2
430  + tan4 + 270.e0 * eta - 330.e0 * tan2 * eta + 445.e0 * eta2
431  + 324.e0 * eta3 -680.e0 * tan2 * eta2 + 88.e0 * eta4
432  -600.e0 * tan2 * eta3 - 192.e0 * tan2 * eta4) / 720.e0;
433 
434  t5 = sn * s * c7 * TranMerc_Scale_Factor * (1385.e0 - 3111.e0 *
435  tan2 + 543.e0 * tan4 - tan6) / 40320.e0;
436 
437  *Northing = TranMerc_False_Northing + t1 + pow(dlam,2.e0) * t2
438  + pow(dlam,4.e0) * t3 + pow(dlam,6.e0) * t4
439  + pow(dlam,8.e0) * t5;
440 
441  /* Easting */
442  t6 = sn * c * TranMerc_Scale_Factor;
443  t7 = sn * c3 * TranMerc_Scale_Factor * (1.e0 - tan2 + eta ) /6.e0;
444  t8 = sn * c5 * TranMerc_Scale_Factor * (5.e0 - 18.e0 * tan2 + tan4
445  + 14.e0 * eta - 58.e0 * tan2 * eta + 13.e0 * eta2 + 4.e0 * eta3
446  - 64.e0 * tan2 * eta2 - 24.e0 * tan2 * eta3 )/ 120.e0;
447  t9 = sn * c7 * TranMerc_Scale_Factor * ( 61.e0 - 479.e0 * tan2
448  + 179.e0 * tan4 - tan6 ) /5040.e0;
449 
450  *Easting = TranMerc_False_Easting + dlam * t6 + pow(dlam,3.e0) * t7
451  + pow(dlam,5.e0) * t8 + pow(dlam,7.e0) * t9;
452  }
453  return (Error_Code);
454 } /* END OF Convert_Geodetic_To_Transverse_Mercator */
455 
456 
457 int Convert_Transverse_Mercator_To_Geodetic (
458  double Easting,
459  double Northing,
460  double *Latitude,
461  double *longitude)
462 { /* BEGIN Convert_Transverse_Mercator_To_Geodetic */
463 
464  /*
465  * The function Convert_Transverse_Mercator_To_Geodetic converts Transverse
466  * Mercator projection (easting and northing) coordinates to geodetic
467  * (latitude and longitude) coordinates, according to the current ellipsoid
468  * and Transverse Mercator projection parameters. If any errors occur, the
469  * error code(s) are returned by the function, otherwise TRANMERC_NO_ERROR is
470  * returned.
471  *
472  * Easting : Easting/X in meters (input)
473  * Northing : Northing/Y in meters (input)
474  * Latitude : Latitude in radians (output)
475  * longitude : longitude in radians (output)
476  */
477 
478  double c; /* Cosine of latitude */
479  double de; /* Delta easting - Difference in Easting (Easting-Fe) */
480  double dlam; /* Delta longitude - Difference in longitude */
481  double eta; /* constant - TranMerc_ebs *c *c */
482  double eta2;
483  double eta3;
484  double eta4;
485  double ftphi; /* Footpoint latitude */
486  int i; /* Loop iterator */
487  double sn; /* Radius of curvature in the prime vertical */
488  double sr; /* Radius of curvature in the meridian */
489  double t; /* Tangent of latitude */
490  double tan2;
491  double tan4;
492  double t10; /* Term in coordinate conversion formula - GP to Y */
493  double t11; /* Term in coordinate conversion formula - GP to Y */
494  double t12; /* Term in coordinate conversion formula - GP to Y */
495  double t13; /* Term in coordinate conversion formula - GP to Y */
496  double t14; /* Term in coordinate conversion formula - GP to Y */
497  double t15; /* Term in coordinate conversion formula - GP to Y */
498  double t16; /* Term in coordinate conversion formula - GP to Y */
499  double t17; /* Term in coordinate conversion formula - GP to Y */
500  double tmd; /* True Meridional distance */
501  double tmdo; /* True Meridional distance for latitude of origin */
502  int Error_Code = TRANMERC_NO_ERROR;
503 
504  if ((Easting < (TranMerc_False_Easting - TranMerc_Delta_Easting))
505  ||(Easting > (TranMerc_False_Easting + TranMerc_Delta_Easting)))
506  { /* Easting out of range */
507  Error_Code |= TRANMERC_EASTING_ERROR;
508  }
509  if ((Northing < (TranMerc_False_Northing - TranMerc_Delta_Northing))
510  || (Northing > (TranMerc_False_Northing + TranMerc_Delta_Northing)))
511  { /* Northing out of range */
512  Error_Code |= TRANMERC_NORTHING_ERROR;
513  }
514 
515  if (!Error_Code)
516  {
517  /* True Meridional Distances for latitude of origin */
518  tmdo = SPHTMD(TranMerc_Origin_Lat);
519 
520  /* Origin */
521  tmd = tmdo + (Northing - TranMerc_False_Northing) / TranMerc_Scale_Factor;
522 
523  /* First Estimate */
524  sr = SPHSR(0.e0);
525  ftphi = tmd/sr;
526 
527  for (i = 0; i < 5 ; i++)
528  {
529  t10 = SPHTMD (ftphi);
530  sr = SPHSR(ftphi);
531  ftphi = ftphi + (tmd - t10) / sr;
532  }
533 
534  /* Radius of Curvature in the meridian */
535  sr = SPHSR(ftphi);
536 
537  /* Radius of Curvature in the meridian */
538  sn = SPHSN(ftphi);
539 
540  /* Sine Cosine terms */
541  c = cos(ftphi);
542 
543  /* Tangent Value */
544  t = tan(ftphi);
545  tan2 = t * t;
546  tan4 = tan2 * tan2;
547  eta = TranMerc_ebs * pow(c,2);
548  eta2 = eta * eta;
549  eta3 = eta2 * eta;
550  eta4 = eta3 * eta;
551  de = Easting - TranMerc_False_Easting;
552  if (fabs(de) < 0.0001)
553  de = 0.0;
554 
555  /* Latitude */
556  t10 = t / (2.e0 * sr * sn * pow(TranMerc_Scale_Factor, 2));
557  t11 = t * (5.e0 + 3.e0 * tan2 + eta - 4.e0 * pow(eta,2)
558  - 9.e0 * tan2 * eta) / (24.e0 * sr * pow(sn,3)
559  * pow(TranMerc_Scale_Factor,4));
560  t12 = t * (61.e0 + 90.e0 * tan2 + 46.e0 * eta + 45.E0 * tan4
561  - 252.e0 * tan2 * eta - 3.e0 * eta2 + 100.e0
562  * eta3 - 66.e0 * tan2 * eta2 - 90.e0 * tan4
563  * eta + 88.e0 * eta4 + 225.e0 * tan4 * eta2
564  + 84.e0 * tan2* eta3 - 192.e0 * tan2 * eta4)
565  / ( 720.e0 * sr * pow(sn,5) * pow(TranMerc_Scale_Factor, 6) );
566  t13 = t * ( 1385.e0 + 3633.e0 * tan2 + 4095.e0 * tan4 + 1575.e0
567  * pow(t,6))/ (40320.e0 * sr * pow(sn,7) * pow(TranMerc_Scale_Factor,8));
568  *Latitude = ftphi - pow(de,2) * t10 + pow(de,4) * t11 - pow(de,6) * t12
569  + pow(de,8) * t13;
570 
571  t14 = 1.e0 / (sn * c * TranMerc_Scale_Factor);
572 
573  t15 = (1.e0 + 2.e0 * tan2 + eta) / (6.e0 * pow(sn,3) * c *
574  pow(TranMerc_Scale_Factor,3));
575 
576  t16 = (5.e0 + 6.e0 * eta + 28.e0 * tan2 - 3.e0 * eta2
577  + 8.e0 * tan2 * eta + 24.e0 * tan4 - 4.e0
578  * eta3 + 4.e0 * tan2 * eta2 + 24.e0
579  * tan2 * eta3) / (120.e0 * pow(sn,5) * c
580  * pow(TranMerc_Scale_Factor,5));
581 
582  t17 = (61.e0 + 662.e0 * tan2 + 1320.e0 * tan4 + 720.e0
583  * pow(t,6)) / (5040.e0 * pow(sn,7) * c
584  * pow(TranMerc_Scale_Factor,7));
585 
586  /* Difference in longitude */
587  dlam = de * t14 - pow(de,3) * t15 + pow(de,5) * t16 - pow(de,7) * t17;
588 
589  /* longitude */
590  (*longitude) = TranMerc_Origin_int + dlam;
591 
592  if((fabs)(*Latitude) > (90.0 * PI / 180.0))
593  Error_Code |= TRANMERC_NORTHING_ERROR;
594 
595  if((*longitude) > (PI))
596  {
597  *longitude -= (2 * PI);
598  if((fabs)(*longitude) > PI)
599  Error_Code |= TRANMERC_EASTING_ERROR;
600  }
601  else if((*longitude) < (-PI))
602  {
603  *longitude += (2 * PI);
604  if((fabs)(*longitude) > PI)
605  Error_Code |= TRANMERC_EASTING_ERROR;
606  }
607 
608  if (fabs(dlam) > (9.0 * PI / 180) * cos(*Latitude))
609  { /* Distortion will result if longitude is more than 9 degrees from the Central Meridian at the equator */
610  /* and decreases to 0 degrees at the poles */
611  /* As you move towards the poles, distortion will become more significant */
612  Error_Code |= TRANMERC_LON_WARNING;
613  }
614  }
615  return (Error_Code);
616 } /* END OF Convert_Transverse_Mercator_To_Geodetic */
TRANMERC_ORIGIN_LAT_ERROR
#define TRANMERC_ORIGIN_LAT_ERROR
Definition: tranmerc.h:96
Convert_Transverse_Mercator_To_Geodetic
int Convert_Transverse_Mercator_To_Geodetic(double Easting, double Northing, double *Latitude, double *longitude)
Definition: tranmerc.c:457
TranMerc_False_Easting
static double TranMerc_False_Easting
Definition: tranmerc.c:137
TRANMERC_LON_WARNING
#define TRANMERC_LON_WARNING
Definition: tranmerc.h:101
TRANMERC_EASTING_ERROR
#define TRANMERC_EASTING_ERROR
Definition: tranmerc.h:94
TranMerc_False_Northing
static double TranMerc_False_Northing
Definition: tranmerc.c:136
PI_OVER_2
#define PI_OVER_2
Definition: tranmerc.c:103
MAX_LAT
#define MAX_LAT
Definition: tranmerc.c:104
TranMerc_Delta_Easting
static double TranMerc_Delta_Easting
Definition: tranmerc.c:148
TRANMERC_NORTHING_ERROR
#define TRANMERC_NORTHING_ERROR
Definition: tranmerc.h:95
SPHSN
#define SPHSN(Latitude)
Definition: tranmerc.c:113
Set_Transverse_Mercator_Parameters
int Set_Transverse_Mercator_Parameters(double a, double f, double Origin_Latitude, double Central_Meridian, double False_Easting, double False_Northing, double Scale_Factor)
Definition: tranmerc.c:161
PI
#define PI
Definition: tranmerc.c:102
TRANMERC_SCALE_FACTOR_ERROR
#define TRANMERC_SCALE_FACTOR_ERROR
Definition: tranmerc.h:100
MIN_SCALE_FACTOR
#define MIN_SCALE_FACTOR
Definition: tranmerc.c:106
TRANMERC_A_ERROR
#define TRANMERC_A_ERROR
Definition: tranmerc.h:98
TRANMERC_LON_ERROR
#define TRANMERC_LON_ERROR
Definition: tranmerc.h:93
TranMerc_es
static double TranMerc_es
Definition: tranmerc.c:130
SPHSR
#define SPHSR(Latitude)
Definition: tranmerc.c:116
TranMerc_cp
static double TranMerc_cp
Definition: tranmerc.c:143
longitude
int32_t longitude
Definition: pios_crossfire.h:42
TRANMERC_NO_ERROR
#define TRANMERC_NO_ERROR
Definition: tranmerc.h:91
TranMerc_Origin_int
static double TranMerc_Origin_int
Definition: tranmerc.c:135
TranMerc_ebs
static double TranMerc_ebs
Definition: tranmerc.c:131
TranMerc_Delta_Northing
static double TranMerc_Delta_Northing
Definition: tranmerc.c:149
TranMerc_f
static double TranMerc_f
Definition: tranmerc.c:129
TRANMERC_INV_F_ERROR
#define TRANMERC_INV_F_ERROR
Definition: tranmerc.h:99
TranMerc_a
static double TranMerc_a
Definition: tranmerc.c:128
Convert_Geodetic_To_Transverse_Mercator
int Convert_Geodetic_To_Transverse_Mercator(double Latitude, double longitude, double *Easting, double *Northing)
Definition: tranmerc.c:302
TRANMERC_LAT_ERROR
#define TRANMERC_LAT_ERROR
Definition: tranmerc.h:92
TranMerc_Origin_Lat
static double TranMerc_Origin_Lat
Definition: tranmerc.c:134
MAX_DELTA_int
#define MAX_DELTA_int
Definition: tranmerc.c:105
TranMerc_ep
static double TranMerc_ep
Definition: tranmerc.c:145
i
uint8_t i
Definition: msp_messages.h:97
px_mkfw.f
tuple f
Definition: px_mkfw.py:81
SPHTMD
#define SPHTMD(Latitude)
Definition: tranmerc.c:109
TranMerc_ap
static double TranMerc_ap
Definition: tranmerc.c:141
TranMerc_bp
static double TranMerc_bp
Definition: tranmerc.c:142
TRANMERC_CENT_MER_ERROR
#define TRANMERC_CENT_MER_ERROR
Definition: tranmerc.h:97
Get_Transverse_Mercator_Parameters
void Get_Transverse_Mercator_Parameters(double *a, double *f, double *Origin_Latitude, double *Central_Meridian, double *False_Easting, double *False_Northing, double *Scale_Factor)
Definition: tranmerc.c:266
MAX_SCALE_FACTOR
#define MAX_SCALE_FACTOR
Definition: tranmerc.c:107
TranMerc_Scale_Factor
static double TranMerc_Scale_Factor
Definition: tranmerc.c:138
TranMerc_dp
static double TranMerc_dp
Definition: tranmerc.c:144