001/*
002 * To change this template, choose Tools | Templates
003 * and open the template in the editor.
004 */
005
006package armyc2.c5isr.JavaTacticalRenderer;
007import armyc2.c5isr.JavaLineArray.POINT2;
008import armyc2.c5isr.JavaLineArray.ref;
009import java.util.ArrayList;
010import armyc2.c5isr.renderer.utilities.ErrorLogger;
011import armyc2.c5isr.renderer.utilities.RendererException;
012import armyc2.c5isr.graphics2d.Rectangle2D;
013
014/**
015 * Class to calculate the geodesic based shapes for the Fire Support Areas
016*
017 */
018public final class mdlGeodesic {
019    private static final String _className = "mdlGeodesic";
020    private static final double sm_a    = 6378137;
021
022    private static double DegToRad(double deg) {
023        return deg / 180.0 * Math.PI;
024    }
025
026    private static double RadToDeg(double rad) {
027        return rad / Math.PI * 180.0;
028    }
029/**
030 * Returns the azimuth from true north between two points
031 * @param c1
032 * @param c2
033 * @return the azimuth from c1 to c2
034 */
035    public static double GetAzimuth(POINT2 c1, 
036            POINT2 c2) {//was private
037        double theta = 0;
038        try {
039            double lat1 = DegToRad(c1.y);
040            double lon1 = DegToRad(c1.x);
041            double lat2 = DegToRad(c2.y);
042            double lon2 = DegToRad(c2.x);
043            //formula
044            //θ = atan2( sin(Δlong).cos(lat2),
045            //cos(lat1).sin(lat2) − sin(lat1).cos(lat2).cos(Δlong) )
046            //var theta:Number = Math.atan2( Math.sin(lon2-lon1)*Math.cos(lat2),
047            //Math.cos(lat1)*Math.sin(lat2) − Math.sin(lat1)*Math.cos(lat2)*Math.cos(lon2-lon1) );
048            double y = Math.sin(lon2 - lon1);
049            y *= Math.cos(lat2);
050            double x = Math.cos(lat1);
051            x *= Math.sin(lat2);
052            double z = Math.sin(lat1);
053            z *= Math.cos(lat2);
054            z *= Math.cos(lon2 - lon1);
055            x = x - z;
056            theta = Math.atan2(y, x);
057            theta = RadToDeg(theta);
058        }
059        catch (Exception exc) {
060            //System.out.println(e.getMessage());
061            //clsUtility.WriteFile("Error in mdlGeodesic.GetAzimuth");
062               ErrorLogger.LogException(_className ,"GetAzimuth",
063                    new RendererException("Failed inside GetAzimuth", exc));
064        }
065        return theta;//RadToDeg(k);
066    }
067    /**
068     * Calculates the distance in meters between two geodesic points.
069     * Also calculates the azimuth from c1 to c2 and from c2 to c1.
070     *
071     * @param c1 the first point
072     * @param c2 the last point
073     * @param a12 OUT - an object with a member to hold the calculated azimuth in degrees from c1 to c2
074     * @param a21 OUT - an object with a member to hold the calculated azimuth in degrees from c2 to c1
075     * @return the distance in meters between c1 and c2
076     */
077    public static double geodesic_distance(POINT2 c1,
078            POINT2 c2,
079            ref<double[]> a12,
080            ref<double[]> a21) {
081        double h = 0;
082        try {
083            //formula
084            //R = earth’s radius (mean radius = 6,371km)
085            //Δlat = lat2− lat1
086            //Δlong = long2− long1
087            //a = sin²(Δlat/2) + cos(lat1).cos(lat2).sin²(Δlong/2)
088            //c = 2.atan2(√a, √(1−a))
089            //d = R.c
090            if(a12 != null && a21 !=null)
091            {
092                a12.value = new double[1];
093                a21.value = new double[1];
094                //set the azimuth
095                a12.value[0] = GetAzimuth(c1, c2);
096                a21.value[0] = GetAzimuth(c2, c1);
097            }
098            //c1.x+=360;
099            double dLat = DegToRad(c2.y - c1.y);
100            double dLon = DegToRad(c2.x - c1.x);
101
102            double b = 0, lat1 = 0, lat2 = 0, e = 0, f = 0, g = 0, k = 0;
103            b = Math.sin(dLat / 2);
104            lat1 = DegToRad(c1.y);
105            lat2 = DegToRad(c2.y);
106            e = Math.sin(dLon / 2);
107            f = Math.cos(lat1);
108            g = Math.cos(lat2);
109            //uncomment this to test calculation
110            //var a:Number = Math.sin(dLat / 2) * Math.sin(dLat / 2) + Math.cos(DegToRad(c1.y)) * Math.cos(DegToRad(c2.y)) * Math.sin(dLon / 2) * Math.sin(dLon / 2);
111            double a = b * b + f * g * e * e;
112            h = Math.sqrt(a);
113            k = Math.sqrt(1 - a);
114            h = 2 * Math.atan2(h, k);
115        }
116        catch (Exception exc) {
117            //System.out.println(e.getMessage());
118            //clsUtility.WriteFile("Error in mdlGeodesic.geodesic_distance");
119               ErrorLogger.LogException(_className ,"geodesic_distance",
120                    new RendererException("Failed inside geodesic_distance", exc));
121        }
122        return sm_a * h;
123    }
124    /**
125     * Calculates a geodesic point and given distance and azimuth from the srating geodesic point
126     *
127     * @param start the starting point
128     * @param distance the distance in meters
129     * @param azimuth the azimuth or bearing in degrees
130     *
131     * @return the calculated point
132     */
133    public static POINT2 geodesic_coordinate(POINT2 start,
134            double distance,
135            double azimuth) {
136        POINT2 pt = null;
137        try
138        {
139        //formula
140        //lat2 = asin(sin(lat1)*cos(d/R) + cos(lat1)*sin(d/R)*cos(θ))
141        //lon2 = lon1 + atan2(sin(θ)*sin(d/R)*cos(lat1), cos(d/R)−sin(lat1)*sin(lat2))
142
143        double a = 0, b = 0, c = 0, d = 0, e = 0, f = 0, g = 0, h = 0,
144                j = 0, k = 0, l = 0, m = 0, n = 0, p = 0, q = 0;
145
146        a = DegToRad(start.y);
147        b = Math.cos(a);
148        c = DegToRad(azimuth);
149        d = Math.sin(a);
150        e = Math.cos(distance / sm_a);
151        f = Math.sin(distance / sm_a);
152        g = Math.cos(c);
153        //uncomment to test calculation
154        //var lat2:Number = RadToDeg(Math.asin(Math.sin(DegToRad(start.y)) * Math.cos(DegToRad(distance / sm_a)) + Math.cos(DegToRad(start.y)) * Math.sin(DegToRad(distance / sm_a)) * Math.cos(DegToRad(azimuth))));
155        //lat2 = asin(sin(lat1)*cos(d/R) + cos(lat1)*sin(d/R)*cos(θ))
156        //var lat2:Number = RadToDeg(Math.asin(Math.sin(DegToRad(start.y)) * Math.cos(distance / sm_a) + Math.cos(DegToRad(start.y)) * Math.sin(distance / sm_a) * Math.cos(DegToRad(azimuth))));
157        //double lat2 = RadToDeg(Math.asin(Math.sin(DegToRad(start.y)) * Math.cos(distance / sm_a) + Math.cos(DegToRad(start.y)) * Math.sin(distance / sm_a) * Math.cos(DegToRad(azimuth))));
158        double lat = RadToDeg(Math.asin(d * e + b * f * g));
159        h = Math.sin(c);
160        k = Math.sin(h);
161        l = Math.cos(a);
162        m = DegToRad(lat);
163        n = Math.sin(m);
164        p = Math.atan2(h * f * b, e - d * n);
165        //uncomment to test calculation
166        //var lon2:Number = start.x + DegToRad(Math.atan2(Math.sin(DegToRad(azimuth)) * Math.sin(DegToRad(distance / sm_a)) * Math.cos(DegToRad(start.y)), Math.cos(DegToRad(distance / sm_a)) - Math.sin(DegToRad(start.y)) * Math.sin(DegToRad(lat))));
167        //lon2 = lon1 + atan2(sin(θ)*sin(d/R)*cos(lat1), cos(d/R)−sin(lat1)*sin(lat2))
168        //var lon2:Number = start.x + RadToDeg(Math.atan2(Math.sin(DegToRad(azimuth)) * Math.sin(distance / sm_a) * Math.cos(DegToRad(start.y)), Math.cos(distance / sm_a) - Math.sin(DegToRad(start.y)) * Math.sin(DegToRad(lat2))));
169        double lon = start.x + RadToDeg(p);
170        pt = new POINT2(lon, lat);
171        }
172        catch (Exception exc) {
173            //clsUtility.WriteFile("Error in mdlGeodesic.geodesic_distance");
174               ErrorLogger.LogException(_className ,"geodesic_coordinate",
175                    new RendererException("Failed inside geodesic_coordinate", exc));
176        }
177        return pt;
178    }
179    /**
180     * Calculates an arc from geodesic point and uses them for the change 1 circular symbols
181     *
182     * @param pPoints array of 3 points, currently the last 2 points are the same. The first point
183     * is the center and the next point defines the radius.
184     *
185     * @return points for the geodesic circle
186     */
187    public static ArrayList<POINT2> GetGeodesicArc(POINT2[] pPoints) {
188        ArrayList<POINT2> pPoints2 = new ArrayList();
189        try {
190            if (pPoints == null) {
191                return null;
192            }
193            if (pPoints.length < 3) {
194                return null;
195            }
196
197            POINT2 ptCenter = new POINT2(pPoints[0]);
198            POINT2 pt1 = new POINT2(pPoints[1]);
199            POINT2 pt2 = new POINT2(pPoints[2]);
200            POINT2 ptTemp = null;
201            ref<double[]> a12b = new ref();
202            double dist2 = 0.0;
203            double dist1 = 0.0;
204            ref<double[]> a12 = new ref();
205            ref<double[]> a21 = new ref();
206            //distance and azimuth from the center to the 1st point
207            dist1 = geodesic_distance(ptCenter, pt1, a12, a21);
208            double saveAzimuth = a21.value[0];
209            //distance and azimuth from the center to the 2nd point
210            dist2 = geodesic_distance(ptCenter, pt2, a12b, a21);
211            //if the points are nearly the same we want 360 degree range fan
212            if (Math.abs(a21.value[0] - saveAzimuth) <= 1) {
213                if (a12.value[0] < 360) {
214                    a12.value[0] += 360;
215                }
216
217                a12b.value[0] = a12.value[0] + 360;
218            }
219
220            ref<double[]> a12c = new ref();
221            int j = 0;
222            if (a12b.value[0] < 0) {
223                a12b.value[0] = 360 + a12b.value[0];
224            }
225            if (a12.value[0] < 0) {
226                a12.value[0] = 360 + a12.value[0];
227            }
228            if (a12b.value[0] < a12.value[0]) {
229                a12b.value[0] = a12b.value[0] + 360;
230            }
231            a12c.value=new double[1];
232            for (j = 0; j <= 100; j++) {
233
234                a12c.value[0] = a12.value[0] + ((double) j / 100.0) * (a12b.value[0] - a12.value[0]);
235                ptTemp = geodesic_coordinate(ptCenter, dist1, a12c.value[0]);
236                pPoints2.add(ptTemp);
237            }
238
239            //if the points are nearly the same we want 360 degree range fan
240            //with no line from the center
241            if (Math.abs(a21.value[0] - saveAzimuth) > 1) {
242                pPoints2.add(ptCenter);
243            }
244
245            if (a12.value[0] < a12b.value[0]) {
246                pPoints2.add(pt1);
247            } else {
248                pPoints2.add(pt2);
249            }
250        } catch (Exception exc) {
251            //clsUtility.WriteFile("Error in mdlGeodesic.GetGeodesicArc");
252               ErrorLogger.LogException(_className ,"GetGeodesicArc",
253                    new RendererException("Failed inside GetGeodesicArc", exc));
254        }
255        return pPoints2;
256    }
257    /**
258     * Calculates the sector points for a sector range fan.
259     *
260     * @param pPoints array of 3 points. The first point
261     * is the center and the next two points define either side of the sector
262     * @param pPoints2 OUT - the calculated geodesic sector points
263     *
264     * @return true if the sector is a circle
265     */
266    public static boolean GetGeodesicArc2(ArrayList<POINT2> pPoints,
267            ArrayList<POINT2> pPoints2) {
268        boolean circle = false;
269        try {
270            POINT2 ptCenter = new POINT2(pPoints.get(0)), pt1 = new POINT2(pPoints.get(1)), pt2 = new POINT2(pPoints.get(2));
271
272            ref<double[]> a12b = new ref();
273            //double dist2 = 0d;
274            double dist1 = 0d;
275            ref<double[]> a12 = new ref();
276            ref<double[]> a21 = new ref();
277            //double lat2c = 0.0;
278            //distance and azimuth from the center to the 1st point
279            //geodesic_distance(lonCenter, latCenter, lon1, lat1, ref dist1, ref a12, ref a21);
280            dist1 = geodesic_distance(ptCenter, pt1, a12, a21);
281            double saveAzimuth = a21.value[0];
282            //distance and azimuth from the center to the 2nd point
283            //geodesic_distance(lonCenter, latCenter, lon2, lat2, ref dist2, ref a12b, ref a21);
284            double dist2 = geodesic_distance(ptCenter, pt2, a12b, a21);
285            //if the points are nearly the same we want 360 degree range fan
286            if (Math.abs(a21.value[0] - saveAzimuth) <= 1) {
287                if (a12.value[0] < 360) {
288                    a12.value[0] += 360;
289                }
290                a12b.value[0] = a12.value[0] + 360;
291                circle = true;
292            }
293
294            //assume caller has set pPoints2 as new Array
295
296            ref<double[]> a12c = new ref();
297            a12c.value = new double[1];
298            int j = 0;
299            POINT2 pPoint = new POINT2();
300            if (a12b.value[0] < 0) {
301                a12b.value[0] = 360 + a12b.value[0];
302            }
303            if (a12.value[0] < 0) {
304                a12.value[0] = 360 + a12.value[0];
305            }
306            if (a12b.value[0] < a12.value[0]) {
307                a12b.value[0] = a12b.value[0] + 360;
308            }
309            for (j = 0; j <= 100; j++) {
310
311                a12c.value[0] = a12.value[0] + ((double) j / 100) * (a12b.value[0] - a12.value[0]);
312                pPoint = geodesic_coordinate(ptCenter, dist1, a12c.value[0]);
313                pPoints2.add(pPoint);
314            }
315        }
316        catch (Exception exc) {
317            //System.out.println(e.getMessage());
318            //clsUtility.WriteFile("Error in mdlGeodesic.GetGeodesicArc2");
319               ErrorLogger.LogException(_className ,"GetGeodesicArc2",
320                    new RendererException("Failed inside GetGeodesicArc2", exc));
321        }
322        return circle;
323    }
324    /**
325     * returns intersection of two lines, each defined by a point and a bearing
326     * <a href="http://creativecommons.org/licenses/by/3.0/"><img alt="Creative Commons License" style="border-width:0" src="http://i.creativecommons.org/l/by/3.0/88x31.png"></a><br>This work is licensed under a <a href="http://creativecommons.org/licenses/by/3.0/">Creative Commons Attribution 3.0 Unported License</a>.
327     * @param p1 1st point
328     * @param brng1 first line bearing in degrees from true north
329     * @param p2 2nd point
330     * @param brng2 2nd point bearing in degrees from true north
331     * @return
332     * @deprecated
333     */
334    public static POINT2 IntersectLines(POINT2 p1,
335            double brng1, 
336            POINT2 p2,
337            double brng2) {
338        POINT2 ptResult = null;
339        try {
340            double lat1 = DegToRad(p1.y);//p1._lat.toRad();
341            double lon1 = DegToRad(p1.x);//p1._lon.toRad();
342            double lat2 = DegToRad(p2.y);//p2._lat.toRad();
343            double lon2 = DegToRad(p2.x);//p2._lon.toRad();
344            double brng13 = DegToRad(brng1);//brng1.toRad();
345            double brng23 = DegToRad(brng2);//brng2.toRad();
346            double dLat = lat2 - lat1;
347            double dLon = lon2 - lon1;
348
349
350            double dist12 = 2 * Math.asin(Math.sqrt(Math.sin(dLat / 2) * Math.sin(dLat / 2) +
351                    Math.cos(lat1) * Math.cos(lat2) * Math.sin(dLon / 2) * Math.sin(dLon / 2)));
352
353            if (dist12 == 0) {
354                return null;
355            }
356
357            double brngA = Math.acos((Math.sin(lat2) - Math.sin(lat1) * Math.cos(dist12)) /
358                    (Math.sin(dist12) * Math.cos(lat1)));
359
360            if (Double.isNaN(brngA)) {
361                brngA = 0;  // protect against rounding
362            }
363            double brngB = Math.acos((Math.sin(lat1) - Math.sin(lat2) * Math.cos(dist12)) /
364                    (Math.sin(dist12) * Math.cos(lat2)));
365
366            double brng12 = 0, brng21 = 0;
367            if (Math.sin(lon2 - lon1) > 0) {
368                brng12 = brngA;
369                brng21 = 2 * Math.PI - brngB;
370            } else {
371                brng12 = 2 * Math.PI - brngA;
372                brng21 = brngB;
373            }
374
375            double alpha1 = (brng13 - brng12 + Math.PI) % (2 * Math.PI) - Math.PI;  // angle 2-1-3
376            double alpha2 = (brng21 - brng23 + Math.PI) % (2 * Math.PI) - Math.PI;  // angle 1-2-3
377
378            if (Math.sin(alpha1) == 0 && Math.sin(alpha2) == 0) {
379                return null;  // infinite intersections
380            }
381            if (Math.sin(alpha1) * Math.sin(alpha2) < 0) {
382                return null;       // ambiguous intersection
383            }
384            //alpha1 = Math.abs(alpha1);
385            //alpha2 = Math.abs(alpha2);  // ... Ed Williams takes abs of alpha1/alpha2, but seems to break calculation?
386            double alpha3 = Math.acos(-Math.cos(alpha1) * Math.cos(alpha2) +
387                    Math.sin(alpha1) * Math.sin(alpha2) * Math.cos(dist12));
388
389            double dist13 = Math.atan2(Math.sin(dist12) * Math.sin(alpha1) * Math.sin(alpha2),
390                    Math.cos(alpha2) + Math.cos(alpha1) * Math.cos(alpha3));
391
392            double lat3 = Math.asin(Math.sin(lat1) * Math.cos(dist13) +
393                    Math.cos(lat1) * Math.sin(dist13) * Math.cos(brng13));
394            double dLon13 = Math.atan2(Math.sin(brng13) * Math.sin(dist13) * Math.cos(lat1),
395                    Math.cos(dist13) - Math.sin(lat1) * Math.sin(lat3));
396            double lon3 = lon1 + dLon13;
397            lon3 = (lon3 + Math.PI) % (2 * Math.PI) - Math.PI;  // normalise to -180..180º
398
399            //return new POINT2(lat3.toDeg(), lon3.toDeg());
400            ptResult = new POINT2(RadToDeg(lon3), RadToDeg(lat3));
401
402        } catch (Exception exc) {
403            ErrorLogger.LogException(_className, "IntersectLines",
404                    new RendererException("Failed inside IntersectLines", exc));
405        }
406        return ptResult;
407    }
408    /**
409     * Normalizes geo points for arrays which span the IDL
410     *
411     * @param geoPoints
412     * @return
413     */
414    public static ArrayList<POINT2> normalize_points(ArrayList<POINT2> geoPoints) {
415        ArrayList<POINT2> normalizedPts = null;
416        try {
417            if (geoPoints == null || geoPoints.isEmpty()) {
418                return normalizedPts;
419            }
420
421            int j = 0;
422            double minx = geoPoints.get(0).x;
423            double maxx = minx;
424            boolean spansIDL = false;
425            POINT2 pt = null;
426            int n=geoPoints.size();
427            //for (j = 1; j < geoPoints.size(); j++) 
428            for (j = 1; j < n; j++) 
429            {
430                pt = geoPoints.get(j);
431                if (pt.x < minx) {
432                    minx = pt.x;
433                }
434                if (pt.x > maxx) {
435                    maxx = pt.x;
436                }
437            }
438            if (maxx - minx > 180) {
439                spansIDL = true;
440            }
441
442            if (!spansIDL) {
443                return geoPoints;
444            }
445
446            normalizedPts = new ArrayList();
447            n=geoPoints.size();
448            //for (j = 0; j < geoPoints.size(); j++) 
449            for (j = 0; j < n; j++) 
450            {
451                pt = geoPoints.get(j);
452                if (pt.x < 0) {
453                    pt.x += 360;
454                }
455                normalizedPts.add(pt);
456            }
457        } catch (Exception exc) {
458            ErrorLogger.LogException(_className, "normalize_pts",
459                    new RendererException("Failed inside normalize_pts", exc));
460        }
461        return normalizedPts;
462    }
463
464    /**
465     * calculates the geodesic MBR, intended for regular shaped areas
466     *
467     * @param geoPoints
468     * @return
469     */
470    public static Rectangle2D.Double geodesic_mbr(ArrayList<POINT2> geoPoints) {
471        Rectangle2D.Double rect2d = null;
472        try {
473            if (geoPoints == null || geoPoints.isEmpty()) {
474                return rect2d;
475            }
476            
477            ArrayList<POINT2>normalizedPts=normalize_points(geoPoints);
478            double ulx=normalizedPts.get(0).x;
479            double lrx=ulx;
480            double uly=normalizedPts.get(0).y;
481            double lry=uly;
482            int j=0;
483            POINT2 pt=null;
484            int n=normalizedPts.size();
485            //for(j=1;j<normalizedPts.size();j++)
486            for(j=1;j<n;j++)
487            {
488                pt=normalizedPts.get(j);
489                if(pt.x<ulx)
490                    ulx=pt.x;
491                if(pt.x>lrx)
492                    lrx=pt.x;
493            
494                if(pt.y>uly)
495                    uly=pt.y;
496                if(pt.y<lry)
497                    lry=pt.y;
498            }
499            POINT2 ul=new POINT2(ulx,uly);
500            POINT2 ur=new POINT2(lrx,uly);
501            POINT2 lr=new POINT2(lrx,lry);
502            double width=geodesic_distance(ul,ur,null,null);
503            double height=geodesic_distance(ur,lr,null,null);
504            rect2d=new Rectangle2D.Double(ulx,uly,width,height);
505        } catch (Exception exc) {
506            ErrorLogger.LogException(_className, "geodesic_mbr",
507                    new RendererException("Failed inside geodesic_mbr", exc));
508        }
509        return rect2d;
510    }
511
512    /**
513     * Currently used by AddModifiers for greater accuracy on center labels
514     *
515     * @param geoPoints
516     * @return
517     */
518    public static POINT2 geodesic_center(ArrayList<POINT2> geoPoints) {
519        POINT2 pt = null;
520        try {
521            if(geoPoints==null || geoPoints.isEmpty())
522                return pt;
523            
524            Rectangle2D.Double rect2d=geodesic_mbr(geoPoints);
525            double deltax=rect2d.getWidth()/2;
526            double deltay=rect2d.getHeight()/2;
527            POINT2 ul=new POINT2(rect2d.x,rect2d.y);
528            //first walk east by deltax
529            POINT2 ptEast=geodesic_coordinate(ul,deltax,90);
530            //next walk south by deltay;
531            pt=geodesic_coordinate(ptEast,deltay,180);
532            
533        } catch (Exception exc) {
534            ErrorLogger.LogException(_className, "geodesic_center",
535                    new RendererException("Failed inside geodesic_center", exc));
536        }
537        return pt;
538    }
539    /**
540     * rotates a point from a center point in degrees
541     * @param ptCenter center point to rotate about
542     * @param ptRotate point to rotate
543     * @param rotation rotation angle in degrees
544     * @return 
545     */
546    private static POINT2 geoRotatePoint(POINT2 ptCenter, POINT2 ptRotate, double rotation)
547    {
548        try
549        {
550            double bearing=GetAzimuth(ptCenter, ptRotate);
551            double dist=geodesic_distance(ptCenter,ptRotate,null,null);
552            return geodesic_coordinate(ptCenter,dist,bearing+rotation);
553        }
554        catch (Exception exc) {
555            ErrorLogger.LogException(_className, "geoRotatePoint",
556                    new RendererException("Failed inside geoRotatePoint", exc));
557        }
558        return null;
559    }
560    /**
561     * Calculates points for a geodesic ellipse and rotates the points by rotation
562     * @param ptCenter
563     * @param majorRadius
564     * @param minorRadius
565     * @param rotation  rotation angle in degrees
566     * @return 
567     */
568    public static POINT2[] getGeoEllipse(POINT2 ptCenter, double majorRadius, double minorRadius, double rotation)
569    {        
570        POINT2[]pEllipsePoints=null;
571        try
572        {
573            pEllipsePoints=new POINT2[37];
574            //int l=0;
575            POINT2 pt=null;            
576            double dFactor, azimuth=0,a=0,b=0,dist=0,bearing=0;
577            POINT2 ptLongitude=null,ptLatitude=null;
578            for (int l = 1; l < 37; l++)
579            {
580                dFactor = (10.0 * l) * Math.PI / 180.0;                
581                a=majorRadius * Math.cos(dFactor);
582                b=minorRadius * Math.sin(dFactor);
583                //dist=Math.sqrt(a*a+b*b);
584                //azimuth = (10.0 * l);// * Math.PI / 180.0;  
585                //azimuth=90-azimuth;
586                //pt = geodesic_coordinate(ptCenter,dist,azimuth);                
587                //pt = geodesic_coordinate(ptCenter,dist,azimuth);                
588                ptLongitude=geodesic_coordinate(ptCenter,a,90);
589                ptLatitude=geodesic_coordinate(ptCenter,b,0);
590                //pt=new POINT2(ptLatitude.x,ptLongitude.y);
591                pt=new POINT2(ptLongitude.x,ptLatitude.y);
592                //pEllipsePoints[l-1]=pt;
593                pEllipsePoints[l-1]=geoRotatePoint(ptCenter,pt,-rotation);
594            }            
595            pEllipsePoints[36]=new POINT2(pEllipsePoints[0]);
596        }
597        catch(Exception exc)
598        {
599            ErrorLogger.LogException(_className, "GetGeoEllipse",
600                    new RendererException("GetGeoEllipse", exc));
601        }
602        return pEllipsePoints;
603    }
604}