1 //*******************************************************************
2 //
3 // License:  See top LICENSE.txt file.
4 //
5 // Author:  Garrett Potts
6 //
7 // Description:
8 //
9 // Calls Lamberts Conformal Conic projection code.
10 //*******************************************************************
11 //  $Id: ossimLambertConformalConicProjection.cpp 19640 2011-05-25 15:58:00Z oscarkramer $
12 
13 #include <iostream>
14 #include <iomanip>
15 #include <ossim/projection/ossimLambertConformalConicProjection.h>
16 #include <ossim/base/ossimKeywordNames.h>
17 
18 using namespace std;
19 
20 RTTI_DEF1(ossimLambertConformalConicProjection, "ossimLambertConformalConicProjection", ossimMapProjection)
21 
22 /***************************************************************************/
23 /*                               DEFINES
24  *
25  */
26 #ifndef PI_OVER_2
27 #  define PI_OVER_2  ( M_PI / 2.0)
28 #endif
29 #ifndef TWO_PI
30 #  define TWO_PI     (2.0 * M_PI)
31 #endif
32 #define MAX_LAT    (( M_PI *  89.99972222222222) / 180.0)  /* 89 59 59.0 degrees in radians */
33 #define LAMBERT_m(clat,essin)                  (clat / sqrt(1.0 - essin * essin))
34 #define LAMBERT_t(lat,essin)                   tan(PI_OVER_4 - lat / 2) /				\
35 										            pow((1.0 - essin) / (1.0 + essin), es_OVER_2)
36 #define ES_SIN(sinlat)                         (es * sinlat)
37 /**************************************************************************/
38 /*                               GLOBAL DECLARATIONS
39  *
40  */
41 
42 const double PI_OVER_4 = (M_PI / 4.0);
43 
44 #define LAMBERT_NO_ERROR           0x0000
45 #define LAMBERT_LAT_ERROR          0x0001
46 #define LAMBERT_LON_ERROR          0x0002
47 #define LAMBERT_EASTING_ERROR      0x0004
48 #define LAMBERT_NORTHING_ERROR     0x0008
49 #define LAMBERT_FIRST_STDP_ERROR   0x0010
50 #define LAMBERT_SECOND_STDP_ERROR  0x0020
51 #define LAMBERT_ORIGIN_LAT_ERROR   0x0040
52 #define LAMBERT_CENT_MER_ERROR     0x0080
53 #define LAMBERT_A_ERROR            0x0100
54 #define LAMBERT_B_ERROR            0x0200
55 #define LAMBERT_A_LESS_B_ERROR     0x0400
56 #define LAMBERT_HEMISPHERE_ERROR   0x0800
57 #define LAMBERT_FIRST_SECOND_ERROR 0x1000
58 
ossimLambertConformalConicProjection(const ossimEllipsoid & ellipsoid,const ossimGpt & origin)59 ossimLambertConformalConicProjection::ossimLambertConformalConicProjection(const ossimEllipsoid& ellipsoid,
60                                                                            const ossimGpt& origin)
61    :ossimMapProjection(ellipsoid, origin)
62 {
63    setDefaults();
64    update();
65 }
66 
ossimLambertConformalConicProjection(const ossimEllipsoid & ellipsoid,const ossimGpt & origin,double stdParallel1,double stdParallel2,double falseEasting,double falseNorthing)67 ossimLambertConformalConicProjection::ossimLambertConformalConicProjection(const ossimEllipsoid& ellipsoid,
68                                                                            const ossimGpt& origin,
69                                                                            double stdParallel1,
70                                                                            double stdParallel2,
71                                                                            double falseEasting,
72                                                                            double falseNorthing)
73    :ossimMapProjection(ellipsoid, origin)
74 {
75    Lambert_Std_Parallel_1 = stdParallel1*RAD_PER_DEG;
76    Lambert_Std_Parallel_2 = stdParallel2*RAD_PER_DEG;
77    Lambert_False_Easting  = falseEasting;
78    Lambert_False_Northing = falseNorthing;
79 
80    update();
81 }
82 
~ossimLambertConformalConicProjection()83 ossimLambertConformalConicProjection::~ossimLambertConformalConicProjection()
84 {
85 }
86 
dup() const87 ossimObject* ossimLambertConformalConicProjection::dup()const
88 {
89    return new ossimLambertConformalConicProjection(*this);
90 }
91 
update()92 void ossimLambertConformalConicProjection::update()
93 {
94    Set_Lambert_Parameters(theEllipsoid.getA(),
95                           theEllipsoid.getFlattening(),
96                           theOrigin.latr(),
97                           theOrigin.lonr(),
98                           Lambert_Std_Parallel_1,
99                           Lambert_Std_Parallel_2,
100                           Lambert_False_Easting,
101                           Lambert_False_Northing);
102 
103    theFalseEastingNorthing.x = Lambert_False_Easting;
104    theFalseEastingNorthing.y = Lambert_False_Northing;
105 
106    ossimMapProjection::update();
107 }
108 
setStandardParallel1(double degree)109 void ossimLambertConformalConicProjection::setStandardParallel1(double degree)
110 {
111    Lambert_Std_Parallel_1 = degree*RAD_PER_DEG;
112    update();
113 }
114 
setStandardParallel2(double degree)115 void ossimLambertConformalConicProjection::setStandardParallel2(double degree)
116 {
117    Lambert_Std_Parallel_2 = degree*RAD_PER_DEG;
118    update();
119 }
120 
setStandardParallels(double parallel1Degree,double parallel2Degree)121 void ossimLambertConformalConicProjection::setStandardParallels(double parallel1Degree,
122                                                                 double parallel2Degree)
123 {
124    Lambert_Std_Parallel_1 = parallel1Degree*RAD_PER_DEG;
125    Lambert_Std_Parallel_2 = parallel2Degree*RAD_PER_DEG;
126    update();
127 
128 }
129 
setFalseEasting(double falseEasting)130 void ossimLambertConformalConicProjection::setFalseEasting(double falseEasting)
131 {
132    Lambert_False_Easting = falseEasting;
133    update();
134 }
135 
setFalseNorthing(double falseNorthing)136 void ossimLambertConformalConicProjection::setFalseNorthing(double falseNorthing)
137 {
138    Lambert_False_Northing = falseNorthing;
139    update();
140 }
141 
setFalseEastingNorthing(double falseEasting,double falseNorthing)142 void ossimLambertConformalConicProjection::setFalseEastingNorthing(double falseEasting,
143                                                                    double falseNorthing)
144 {
145    Lambert_False_Easting = falseEasting;
146    Lambert_False_Northing = falseNorthing;
147    update();
148 }
149 
setParameters(double parallel1,double parallel2,double falseEasting,double falseNorthing)150 void ossimLambertConformalConicProjection::setParameters(double parallel1,
151                                                          double parallel2,
152                                                          double falseEasting,
153                                                          double falseNorthing)
154 {
155    Lambert_False_Easting = falseEasting;
156    Lambert_False_Northing = falseNorthing;
157    Lambert_Std_Parallel_1 = parallel1*RAD_PER_DEG;
158    Lambert_Std_Parallel_2 = parallel2*RAD_PER_DEG;
159    update();
160 }
161 
162 
setDefaults()163 void ossimLambertConformalConicProjection::setDefaults()
164 {
165    Lambert_Std_Parallel_1 = 40*RAD_PER_DEG;
166    Lambert_Std_Parallel_2 = 50*RAD_PER_DEG;
167    Lambert_False_Northing = 0.0;
168    Lambert_False_Easting  = 0.0;
169 }
170 
inverse(const ossimDpt & eastingNorthing) const171 ossimGpt ossimLambertConformalConicProjection::inverse(const ossimDpt &eastingNorthing)const
172 {
173    double lat = 0.0;
174    double lon = 0.0;
175 
176    Convert_Lambert_To_Geodetic(eastingNorthing.x,
177                                eastingNorthing.y,
178                                &lat,
179                                &lon);
180 
181    return ossimGpt(lat*DEG_PER_RAD, lon*DEG_PER_RAD, 0.0, theDatum);
182 }
183 
forward(const ossimGpt & latLon) const184 ossimDpt ossimLambertConformalConicProjection::forward(const ossimGpt &latLon)const
185 {
186    double easting  = 0.0;
187    double northing = 0.0;
188    ossimGpt gpt = latLon;
189 
190    if (theDatum)
191    {
192       if (theDatum->code() != latLon.datum()->code())
193       {
194          gpt.changeDatum(theDatum); // Shift to our datum.
195       }
196    }
197 
198    Convert_Geodetic_To_Lambert(gpt.latr(),
199                                gpt.lonr(),
200                                &easting,
201                                &northing);
202 
203    return ossimDpt(easting, northing);
204 }
205 
saveState(ossimKeywordlist & kwl,const char * prefix) const206 bool ossimLambertConformalConicProjection::saveState(ossimKeywordlist& kwl, const char* prefix) const
207 {
208    kwl.add(prefix,
209            ossimKeywordNames::STD_PARALLEL_1_KW,
210            Lambert_Std_Parallel_1*DEG_PER_RAD,
211            true);
212 
213    kwl.add(prefix,
214            ossimKeywordNames::STD_PARALLEL_2_KW,
215            Lambert_Std_Parallel_2*DEG_PER_RAD,
216            true);
217 
218    return ossimMapProjection::saveState(kwl, prefix);
219 }
220 
loadState(const ossimKeywordlist & kwl,const char * prefix)221 bool ossimLambertConformalConicProjection::loadState(const ossimKeywordlist& kwl, const char* prefix)
222 {
223    bool flag = ossimMapProjection::loadState(kwl, prefix);
224 
225    const char* type          = kwl.find(prefix, ossimKeywordNames::TYPE_KW);
226    const char* stdParallel1  = kwl.find(prefix, ossimKeywordNames::STD_PARALLEL_1_KW);
227    const char* stdParallel2  = kwl.find(prefix, ossimKeywordNames::STD_PARALLEL_2_KW);
228 
229    setDefaults();
230 
231    if(ossimString(type) == STATIC_TYPE_NAME(ossimLambertConformalConicProjection))
232    {
233       Lambert_False_Easting  = theFalseEastingNorthing.x;
234       Lambert_False_Northing = theFalseEastingNorthing.y;
235 
236       if(stdParallel1)
237       {
238          Lambert_Std_Parallel_1 = ossimString(stdParallel1).toDouble()*RAD_PER_DEG;
239       }
240       if(stdParallel2)
241       {
242          Lambert_Std_Parallel_2 = ossimString(stdParallel2).toDouble()*RAD_PER_DEG;
243       }
244    }
245    update();
246    return flag;
247 }
248 
print(std::ostream & out) const249 std::ostream& ossimLambertConformalConicProjection::print(
250    std::ostream& out) const
251 {
252    // Capture the original flags.
253    std::ios_base::fmtflags f = out.flags();
254 
255    out << setiosflags(ios::fixed) << setprecision(15);
256 
257    out << "// ossimLambertConformalConicProjection::print\n"
258        << ossimKeywordNames::STD_PARALLEL_1_KW << ": "
259        << Lambert_Std_Parallel_1*DEG_PER_RAD << "\n"
260        << ossimKeywordNames::STD_PARALLEL_2_KW << ": "
261       << Lambert_Std_Parallel_2*DEG_PER_RAD << std::endl;
262 
263    // Reset flags.
264    out.setf(f);
265 
266    return ossimMapProjection::print(out);
267 }
268 
269 /************************************************************************/
270 /*                              FUNCTIONS
271  *
272  */
273 
Set_Lambert_Parameters(double a,double f,double Origin_Latitude,double Central_Meridian,double Std_Parallel_1,double Std_Parallel_2,double False_Easting,double False_Northing)274 long ossimLambertConformalConicProjection::Set_Lambert_Parameters(double a,
275                                                                   double f,
276                                                                   double Origin_Latitude,
277                                                                   double Central_Meridian,
278                                                                   double Std_Parallel_1,
279                                                                   double Std_Parallel_2,
280                                                                   double False_Easting,
281                                                                   double False_Northing)
282 
283 { /* BEGIN Set_Lambert_Parameters */
284 /*
285  * The function Set_Lambert_Parameters receives the ellipsoid parameters and
286  * Lambert Conformal Conic projection parameters as inputs, and sets the
287  * corresponding state variables.  If any errors occur, the error code(s)
288  * are returned by the function, otherwise LAMBERT_NO_ERROR is returned.
289  *
290  *   a                   : Semi-major axis of ellipsoid, in meters   (input)
291  *   f                   : Flattening of ellipsoid						       (input)
292  *   Origin_Latitude     : Latitude of origin, in radians            (input)
293  *   Central_Meridian    : Longitude of origin, in radians           (input)
294  *   Std_Parallel_1      : First standard parallel, in radians       (input)
295  *   Std_Parallel_2      : Second standard parallel, in radians      (input)
296  *   False_Easting       : False easting, in meters                  (input)
297  *   False_Northing      : False northing, in meters                 (input)
298  */
299 
300   double slat, slat1, clat;
301   double es_sin;
302   double t0, t1, t2;
303   double m1, m2;
304 //  double inv_f = 1 / f;
305   long Error_Code = LAMBERT_NO_ERROR;
306 
307 //   if (a <= 0.0)
308 //   { /* Semi-major axis must be greater than zero */
309 //     Error_Code |= LAMBERT_A_ERROR;
310 //   }
311 //   if ((inv_f < 250) || (inv_f > 350))
312 //   { /* Inverse flattening must be between 250 and 350 */
313 //     Error_Code |= LAMBERT_INV_F_ERROR;
314 //   }
315 //   if ((Origin_Latitude < -MAX_LAT) || (Origin_Latitude > MAX_LAT))
316 //   { /* Origin Latitude out of range */
317 //     Error_Code |= LAMBERT_ORIGIN_LAT_ERROR;
318 //   }
319 //   if ((Std_Parallel_1 < -MAX_LAT) || (Std_Parallel_1 > MAX_LAT))
320 //   { /* First Standard Parallel out of range */
321 //     Error_Code |= LAMBERT_FIRST_STDP_ERROR;
322 //   }
323 //   if ((Std_Parallel_2 < -MAX_LAT) || (Std_Parallel_2 > MAX_LAT))
324 //   { /* Second Standard Parallel out of range */
325 //     Error_Code |= LAMBERT_SECOND_STDP_ERROR;
326 //   }
327 //   if ((Std_Parallel_1 == 0) && (Std_Parallel_2 == 0))
328 //   { /* First & Second Standard Parallels are both 0 */
329 //     Error_Code |= LAMBERT_FIRST_SECOND_ERROR;
330 //   }
331 //   if (Std_Parallel_1 == -Std_Parallel_2)
332 //   { /* Parallels are the negation of each other */
333 //     Error_Code |= LAMBERT_HEMISPHERE_ERROR;
334 //   }
335 //   if ((Central_Meridian < -PI) || (Central_Meridian > TWO_PI))
336 //   { /* Origin Longitude out of range */
337 //     Error_Code |= LAMBERT_CENT_MER_ERROR;
338 //   }
339 
340   if (!Error_Code)
341   { /* no errors */
342 
343     Lambert_a = a;
344     Lambert_f = f;
345     Lambert_Origin_Lat = Origin_Latitude;
346     Lambert_Std_Parallel_1 = Std_Parallel_1;
347     Lambert_Std_Parallel_2 = Std_Parallel_2;
348 //     if (Central_Meridian > PI)
349 //       Central_Meridian -= TWO_PI;
350     Lambert_Origin_Long = Central_Meridian;
351     Lambert_False_Easting = False_Easting;
352     Lambert_False_Northing = False_Northing;
353 
354     es2 = 2 * Lambert_f - Lambert_f * Lambert_f;
355     es = sqrt(es2);
356     es_OVER_2 = es / 2.0;
357 
358     slat = sin(Lambert_Origin_Lat);
359     es_sin = ES_SIN(slat);
360     t0 = LAMBERT_t(Lambert_Origin_Lat, es_sin);
361 
362     slat1 = sin(Lambert_Std_Parallel_1);
363     clat = cos(Lambert_Std_Parallel_1);
364     es_sin = ES_SIN(slat1);
365     m1 = LAMBERT_m(clat, es_sin);
366     t1 = LAMBERT_t(Lambert_Std_Parallel_1, es_sin);
367 
368 
369     if (fabs(Lambert_Std_Parallel_1 - Lambert_Std_Parallel_2) > 1.0e-10)
370     {
371       slat = sin(Lambert_Std_Parallel_2);
372       clat = cos(Lambert_Std_Parallel_2);
373       es_sin = ES_SIN(slat);
374       m2 = LAMBERT_m(clat, es_sin);
375       t2 = LAMBERT_t(Lambert_Std_Parallel_2, es_sin);
376       n = log(m1 / m2) / log(t1 / t2);
377     }
378     else
379       n = slat1;
380     F = m1 / (n * pow(t1, n));
381     Lambert_aF = Lambert_a * F;
382     if ((t0 == 0) && (n < 0))
383       rho0 = 0.0;
384     else
385       rho0 = Lambert_aF * pow(t0, n);
386 
387   }
388   return (Error_Code);
389 } /* END OF Set_Lambert_Parameters */
390 
391 
Get_Lambert_Parameters(double * a,double * f,double * Origin_Latitude,double * Central_Meridian,double * Std_Parallel_1,double * Std_Parallel_2,double * False_Easting,double * False_Northing) const392 void ossimLambertConformalConicProjection::Get_Lambert_Parameters(double *a,
393                                                                   double *f,
394                                                                   double *Origin_Latitude,
395                                                                   double *Central_Meridian,
396                                                                   double *Std_Parallel_1,
397                                                                   double *Std_Parallel_2,
398                                                                   double *False_Easting,
399                                                                   double *False_Northing)const
400 
401 { /* BEGIN Get_Lambert_Parameters */
402 /*
403  * The function Get_Lambert_Parameters returns the current ellipsoid
404  * parameters and Lambert Conformal Conic projection parameters.
405  *
406  *   a                   : Semi-major axis of ellipsoid, in meters   (output)
407  *   f                   : Flattening of ellipsoid					         (output)
408  *   Origin_Latitude     : Latitude of origin, in radians            (output)
409  *   Central_Meridian    : Longitude of origin, in radians           (output)
410  *   Std_Parallel_1      : First standard parallel, in radians       (output)
411  *   Std_Parallel_2      : Second standard parallel, in radians      (output)
412  *   False_Easting       : False easting, in meters                  (output)
413  *   False_Northing      : False northing, in meters                 (output)
414  */
415 
416 
417   *a = Lambert_a;
418   *f = Lambert_f;
419   *Std_Parallel_1 = Lambert_Std_Parallel_1;
420   *Std_Parallel_2 = Lambert_Std_Parallel_2;
421   *Origin_Latitude = Lambert_Origin_Lat;
422   *Central_Meridian = Lambert_Origin_Long;
423   *False_Easting = Lambert_False_Easting;
424   *False_Northing = Lambert_False_Northing;
425 
426   return;
427 } /* END OF Get_Lambert_Parameters */
428 
429 
Convert_Geodetic_To_Lambert(double Latitude,double Longitude,double * Easting,double * Northing) const430 long ossimLambertConformalConicProjection::Convert_Geodetic_To_Lambert (double Latitude,
431                                                                         double Longitude,
432                                                                         double *Easting,
433                                                                         double *Northing)const
434 
435 { /* BEGIN Convert_Geodetic_To_Lambert */
436 /*
437  * The function Convert_Geodetic_To_Lambert converts Geodetic (latitude and
438  * longitude) coordinates to Lambert Conformal Conic projection (easting
439  * and northing) coordinates, according to the current ellipsoid and
440  * Lambert Conformal Conic projection parameters.  If any errors occur, the
441  * error code(s) are returned by the function, otherwise LAMBERT_NO_ERROR is
442  * returned.
443  *
444  *    Latitude         : Latitude, in radians                         (input)
445  *    Longitude        : Longitude, in radians                        (input)
446  *    Easting          : Easting (X), in meters                       (output)
447  *    Northing         : Northing (Y), in meters                      (output)
448  */
449 
450   double slat;
451   double es_sin;
452   double t;
453   double rho;
454   double dlam;
455   double theta;
456   long  Error_Code = LAMBERT_NO_ERROR;
457 
458 //   if ((Latitude < -PI_OVER_2) || (Latitude > PI_OVER_2))
459 //   {  /* Latitude out of range */
460 //     Error_Code|= LAMBERT_LAT_ERROR;
461 //   }
462 //   if ((Longitude < -PI) || (Longitude > TWO_PI))
463 //   {  /* Longitude out of range */
464 //     Error_Code|= LAMBERT_LON_ERROR;
465 //   }
466 
467   if (!Error_Code)
468   { /* no errors */
469 
470     if (fabs(fabs(Latitude) - PI_OVER_2) > 1.0e-10)
471     {
472       slat = sin(Latitude);
473       es_sin = ES_SIN(slat);
474       t = LAMBERT_t(Latitude, es_sin);
475       rho = Lambert_aF * pow(t, n);
476     }
477     else
478     {
479       if ((Latitude * n) <= 0)
480       { /* Point can not be projected */
481         Error_Code |= LAMBERT_LAT_ERROR;
482         return (Error_Code);
483       }
484       rho = 0.0;
485     }
486 
487     dlam = Longitude - Lambert_Origin_Long;
488 
489 //     if (dlam > PI)
490 //     {
491 //       dlam -= TWO_PI;
492 //     }
493 //     if (dlam < -PI)
494 //     {
495 //       dlam += TWO_PI;
496 //     }
497 
498     theta = n * dlam;
499 
500     *Easting = rho * sin(theta) + Lambert_False_Easting;
501     *Northing = rho0 - rho * cos(theta) + Lambert_False_Northing;
502 
503   }
504   return (Error_Code);
505 } /* END OF Convert_Geodetic_To_Lambert */
506 
507 
508 
Convert_Lambert_To_Geodetic(double Easting,double Northing,double * Latitude,double * Longitude) const509 long ossimLambertConformalConicProjection::Convert_Lambert_To_Geodetic (double Easting,
510                                                                         double Northing,
511                                                                         double *Latitude,
512                                                                         double *Longitude)const
513 
514 { /* BEGIN Convert_Lambert_To_Geodetic */
515 /*
516  * The function Convert_Lambert_To_Geodetic converts Lambert Conformal
517  * Conic projection (easting and northing) coordinates to Geodetic
518  * (latitude and longitude) coordinates, according to the current ellipsoid
519  * and Lambert Conformal Conic projection parameters.  If any errors occur,
520  * the error code(s) are returned by the function, otherwise LAMBERT_NO_ERROR
521  * is returned.
522  *
523  *    Easting          : Easting (X), in meters                       (input)
524  *    Northing         : Northing (Y), in meters                      (input)
525  *    Latitude         : Latitude, in radians                         (output)
526  *    Longitude        : Longitude, in radians                        (output)
527  */
528 
529 
530   double dy, dx;
531   double rho, rho0_MINUS_dy;
532   double t;
533   double PHI;
534   double tempPHI = 0.0;
535   double sin_PHI;
536   double es_sin;
537   double theta = 0.0;
538   double tolerance = 4.85e-10;
539   long Error_Code = LAMBERT_NO_ERROR;
540 
541 //   if ((Easting < (Lambert_False_Easting - Lambert_Delta_Easting))
542 //       ||(Easting > (Lambert_False_Easting + Lambert_Delta_Easting)))
543 //   { /* Easting out of range  */
544 //     Error_Code |= LAMBERT_EASTING_ERROR;
545 //   }
546 //   if ((Northing < (Lambert_False_Northing - Lambert_Delta_Northing))
547 //       || (Northing > (Lambert_False_Northing + Lambert_Delta_Northing)))
548 //   { /* Northing out of range */
549 //     Error_Code |= LAMBERT_NORTHING_ERROR;
550 //   }
551 
552   if (!Error_Code)
553   { /* no errors */
554 
555     dy = Northing - Lambert_False_Northing;
556     dx = Easting - Lambert_False_Easting;
557     rho0_MINUS_dy = rho0 - dy;
558     rho = sqrt(dx * dx + (rho0_MINUS_dy) * (rho0_MINUS_dy));
559 
560     if (n < 0.0)
561     {
562       rho *= -1.0;
563       dy *= -1.0;
564       dx *= -1.0;
565       rho0_MINUS_dy *= -1.0;
566     }
567 
568     if (rho != 0.0)
569     {
570       theta = atan2(dx, rho0_MINUS_dy);
571       t = pow(rho / (Lambert_aF) , 1.0 / n);
572       PHI = PI_OVER_2 - 2.0 * atan(t);
573       while (fabs(PHI - tempPHI) > tolerance)
574       {
575         tempPHI = PHI;
576         sin_PHI = sin(PHI);
577         es_sin = ES_SIN(sin_PHI);
578         PHI = PI_OVER_2 - 2.0 * atan(t * pow((1.0 - es_sin) / (1.0 + es_sin), es_OVER_2));
579       }
580       *Latitude = PHI;
581       *Longitude = theta / n + Lambert_Origin_Long;
582 
583       if (fabs(*Latitude) < 2.0e-7)  /* force lat to 0 to avoid -0 degrees */
584         *Latitude = 0.0;
585       if (*Latitude > PI_OVER_2)  /* force distorted values to 90, -90 degrees */
586         *Latitude = PI_OVER_2;
587       else if (*Latitude < -PI_OVER_2)
588         *Latitude = -PI_OVER_2;
589 
590       if (*Longitude > M_PI)
591       {
592         if (*Longitude - M_PI < 3.5e-6)
593           *Longitude = M_PI;
594 //         else
595 //           *Longitude -= TWO_PI;
596       }
597       if (*Longitude < -M_PI)
598       {
599         if (fabs(*Longitude + M_PI) < 3.5e-6)
600           *Longitude = -M_PI;
601 //         else
602 //           *Longitude += TWO_PI;
603       }
604 
605       if (fabs(*Longitude) < 2.0e-7)  /* force lon to 0 to avoid -0 degrees */
606         *Longitude = 0.0;
607       if (*Longitude > M_PI)  /* force distorted values to 180, -180 degrees */
608         *Longitude = M_PI;
609       else if (*Longitude < -M_PI)
610         *Longitude = -M_PI;
611     }
612     else
613     {
614       if (n > 0.0)
615         *Latitude = PI_OVER_2;
616       else
617         *Latitude = -PI_OVER_2;
618       *Longitude = Lambert_Origin_Long;
619     }
620   }
621   return (Error_Code);
622 } /* END OF Convert_Lambert_To_Geodetic */
623 
getStandardParallel1() const624 double ossimLambertConformalConicProjection::getStandardParallel1()const
625 {
626    return  Lambert_Std_Parallel_1/RAD_PER_DEG;
627 }
628 
getStandardParallel2() const629 double ossimLambertConformalConicProjection::getStandardParallel2()const
630 {
631    return  Lambert_Std_Parallel_2/RAD_PER_DEG;
632 }
633 
634 //*************************************************************************************************
635 //! Returns TRUE if principal parameters are within epsilon tolerance.
636 //*************************************************************************************************
operator ==(const ossimProjection & proj) const637 bool ossimLambertConformalConicProjection::operator==(const ossimProjection& proj) const
638 {
639    if (!ossimMapProjection::operator==(proj))
640       return false;
641 
642    const ossimLambertConformalConicProjection* p = dynamic_cast<const ossimLambertConformalConicProjection*>(&proj);
643    if (!p) return false;
644 
645    if (!ossim::almostEqual(Lambert_Std_Parallel_1,p->Lambert_Std_Parallel_1)) return false;
646    if (!ossim::almostEqual(Lambert_Std_Parallel_2,p->Lambert_Std_Parallel_2)) return false;
647 
648    return true;
649 }
650