Re: Earth map transformation in Ada ?

[Gary Scott <scottg@flash.net>]

If not too much trouble (or rude of me to ask), I would be overjoyed to
receive such a package.  I've been working on the same thing and have
searched high and low for one so that I didn't have to reinvent the
wheel.


Geoff Bull wrote:
> 
> Reinert Korsnes wrote:
> >
> > Does someone know about available (Earth) map stransformations in Ada ?
> >
> > like:
> >
> > lat,long  <-> x,y (Polar Stereographic)
> >
> > reinert
> 
> I have a package for doing wgs84 calculations, which you may find useful
> as a starting point. I can't remember wehther it produces
> correct answers. Here is an excerpt, I can mail you the whole
> package if you like:
> 
> generic
>    type Float_Type is digits <>;
> package WGS84_Package is
> 
>      --  first  numerical eccentricity
>       e1sqr : constant := (a ** 2 - b ** 2) / (a ** 2) ;
> 
>       --  second numerical eccentricity
>       e2sqr : constant := (a ** 2 - b ** 2) / (b ** 2);
> 
> ...
>    subtype Metres is Float_Type;
>    subtype Radians is Float_Type;
>    subtype Degrees is Float_Type;
> 
>    type ECEF is -- earth centreed earth fixed rectangular coordinate
> system
>       record
>          X : Metres;  --  though equator at prime meridian
>          Y : Metres;  --  though equator 90 degrees east of prime
> meridian
>          Z : Metres;  --  through north pole
>       end record;
> 
> 
>    type Polar is
>       record
>          Latitude : Radians;
>          Longitude : Radians;
>          Altitude : Metres;
>       end record;
> 
>    type Azimuth_Elevation is
>       record
>          Azimuth : Radians;
>          Elevation : Radians;
>          Distance : Metres;
>       end record;
> 
> ...
> 
> package body WGS84_Package is
> 
>     --  converts Latitude, Longitude and Altitude (above ellipsoid) to
> XYZ
>    function To_ECEF (From : in Polar) return ECEF is
>       N : Float_Type := a / sqrt (1.0 - e1sqr * sin (From.Latitude) **
> 2);
>    begin -- To_ECEF
>       return
>         (X => (N + From.Altitude) * cos (From.Latitude) * cos
> (From.Longitude),
> 
>          Y => (N + From.Altitude) * cos (From.Latitude) * sin
> (From.Longitude),
> 
>          Z => (N * (1.0 - e1sqr) + From.Altitude) * sin (From.Latitude)
>          );
>    end To_ECEF;
> 
>    --  converts XYZ to Latitude, Longitude, Altitude (above ellipsoid)
>    function To_Polar (From : in ECEF) return Polar is
>       -- length projected onto equatorial plane;
>       P : Metres := sqrt (From.X ** 2 + From.Y ** 2);
> 
>       T : Float_Type := arctan ((From.Z * a) / (p * b));
> 
>       N : Float_Type;
> 
>       Latitude, Longitude : Radians;
>       Altitude : Metres;
>    begin -- To_Polar
>       Latitude := arctan ((From.Z + e2sqr * b * Sin (T) ** 3)
>                       / (P - e1sqr * a * Cos (T) ** 3));
> 
>       if From.X = 0.0 then
>          Longitude := Pi / 2.0;
>       else
>          Longitude := arctan (From.Y / From.X);
>       end if;
> 
>       N :=  A / Sqrt (1.0 - e1sqr * sin (Latitude) ** 2);
>       Altitude := P / cos (Latitude) - N;
> 
>       return (Latitude, Longitude, Altitude);
> 
>    end To_Polar;
> 
> ...


-- 

Gary Scott
mailto:scottg@flash.net

mailto:webmaster@fortranlib.com
http://www.fortranlib.com

Support the GNU Fortran G95 Project:  http://xena.eas.asu.edu/~andy/