Geo Stars Library  0.9.3
Geodetic and Astrometry library
geoPoint.c File Reference

This file contains the geo library routines for: More...

#include <stdio.h>
#include <math.h>
#include <string.h>
#include <stdlib.h>
#include <ctype.h>
#include "geoStars.h"
Include dependency graph for geoPoint.c:

Go to the source code of this file.

Functions

double geoVersion (void)
 This function returns the version of the GeoStarsLib. More...
 
int geoEfg2XyzDiff (GEO_LOCATION *src_desc, GEO_LOCATION *tgt_desc, double xyz_disp[])
 This function returns the XYZ offset of the target point with respect to the source point, given the Earth Fixed Geodetic coordinates of the points. The EFG coordinates for the source must appear in the GEO_LOCATION record, which must be built previous to the call to this procedure. More...
 
int geoEfg2XyzDiff_packed (GEO_LOCATION *src_desc, int count, double efg_xyz[])
 This function returns the XYZ offset of the target point with respect to the source point, given the Earth Fixed Geodetic coordinates of the points. The EFG coordinates for the source must appear in the GEO_LOCATION record, which must be built previous to the call to this procedure. More...
 
double geoLlh2DiffX (double lat1, double lon1, double hgt1, int datum1, double lat2, double lon2, double hgt2, int datum2)
 
double geoLlh2DiffY (double lat1, double lon1, double hgt1, int datum1, double lat2, double lon2, double hgt2, int datum2)
 
double geoLlh2DiffZ (double lat1, double lon1, double hgt1, int datum1, double lat2, double lon2, double hgt2, int datum2)
 
void geoLlh2Efg (double lat, double lon, double height, int datum, double *e, double *f, double *g)
 This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). More...
 
void geoLlh2Efg_packed (GEO_DATUM *datum, int count, double llh_efg[])
 This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). More...
 
double geoLlh2E (double lat, double lon, double hgt, int datum)
 This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the E component. More...
 
double geoLlh2F (double lat, double lon, double hgt, int datum)
 This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the F component. More...
 
double geoLlh2G (double lat, double lon, double hgt, int datum)
 This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the G component. More...
 
void geoXyz2Rae (double xyz_in[], double rae_out[])
 Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2Rae calculates the range, azimuth, and elevation to that point. More...
 
void geoXyz2Rae_packed (int count, double xyz_rae[])
 Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2Rae calculates the range, azimuth, and elevation to that point. More...
 
double geoXyz2R (double x, double y, double z)
 Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the range to that point. More...
 
double geoXyz2A (double x, double y, double z)
 Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the azimuth to that point. More...
 
double geoXyz2E (double x, double y, double z)
 Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the elevation angle to that point. More...
 
void geoRae2Xyz (double rae_in[], double xyz_out[])
 This routine converts from Range, Azimuth, and Elevation into Cartesian coordinates X,Y,Z. More...
 
void geoRae2Xyz_packed (int count, double rae_xyz[])
 This routine iterates over count vertices in the array rae_xyz converting it from Range, Azimuth, and Elevation into Cartesian coordinates X,Y,Z. More...
 
void geoXyz2Efg (GEO_LOCATION *loc, double xyz_in[], double efg_out[])
 Ingests X, Y, Z, and site info and returns the EFG coordinates. More...
 
void geoXyz2Efg_packed (GEO_LOCATION *loc, int count, double xyz_efg[])
 This routine takes site info and iterates over count vertices in the tightly packed array xyz_efg converting it from cartesian XYZ into geocentric EFG coordinates. More...
 
void geoRae2Efg (GEO_LOCATION *loc, double aer_in[], double efg_out[])
 Ingests Range, Azimuth, Elevation and site info and returns the EFG coordinates that the RAE points to. More...
 
double geoDms2Rads (double deg, double min, double sec, char *sign)
 Converts degrees minutes seconds to radians. More...
 
double geoDms2DD (double deg, double min, double sec, char *sign)
 Converts degrees minutes seconds to Decimal Degrees. More...
 
double geoDecdms2Rads (double in)
 Convert decimal degrees, minutes, and seconds ("dmmss.s") to radians. More...
 
void geoRads2Dms (double rads, double *deg, double *min, double *sec, double *dir)
 Converts radians to degrees minutes seconds. More...
 
double geoRads2Decdms (double rads)
 Convert radians to decimal degrees, minutes, and seconds ("dddmmss.s"). More...
 
double geoRads2DD (double rads)
 Converts radians to Decimal Degrees. More...
 
double geoDD2Rads (double dd)
 Converts Decimal Degrees to radians. More...
 
void geoDD2Dms (double dd, double *deg, double *min, double *sec, double *dir)
 Converts Decimal Degrees to degrees minutes seconds. More...
 
double geoDD2Deg (double dd)
 
double geoDD2Min (double dd)
 
double geoDD2Sec (double dd)
 

Detailed Description

This file contains the geo library routines for:

  • Conversion between coordinate systems
  • Angle conversions
  • Tangential plane calculations

Definition in file geoPoint.c.

Function Documentation

double geoDD2Deg ( double  dd)

Definition at line 888 of file geoPoint.c.

References geoDD2Dms().

889 {
890  double deg, min, sec,dir;
891  geoDD2Dms(dd,&deg, &min, &sec,&dir);
892  return(deg*dir);
893 }

Here is the call graph for this function:

void geoDD2Dms ( double  dd,
double *  deg,
double *  min,
double *  sec,
double *  dir 
)

Converts Decimal Degrees to degrees minutes seconds.

Parameters
doubledecimal degrees
doubledeg, min, sec
chardir : -1.0 or 1.0
Returns
nothing

Definition at line 865 of file geoPoint.c.

Referenced by geoDD2Deg(), geoDD2Min(), and geoDD2Sec().

867 {
868  double temp;
869  double fraction;
870 
871  if (dd < 0.0)
872  *dir = -1.0;
873  else
874  *dir = 1.0;
875 
876  dd = fabs (dd);
877 
878  //temp = RAD_TO_DEG * rads;
879  fraction = modf(dd,deg);
880 
881  temp = fraction * 60.0;
882  fraction = modf(temp,min);
883 
884  *sec = fraction * 60.0;
885 }
double geoDD2Min ( double  dd)

Definition at line 894 of file geoPoint.c.

References geoDD2Dms().

895 {
896  double deg, min, sec,dir;
897  geoDD2Dms(dd,&deg, &min, &sec,&dir);
898  return(min);
899 }

Here is the call graph for this function:

double geoDD2Rads ( double  dd)

Converts Decimal Degrees to radians.

Parameters
doubledecimal degrees
Returns
radians

Definition at line 849 of file geoPoint.c.

References DEG_TO_RAD.

850 {
851  return(dd * DEG_TO_RAD);
852 }
double geoDD2Sec ( double  dd)

Definition at line 900 of file geoPoint.c.

References geoDD2Dms().

901 {
902  double deg, min, sec,dir;
903  geoDD2Dms(dd,&deg, &min, &sec,&dir);
904  return(sec);
905 }

Here is the call graph for this function:

double geoDecdms2Rads ( double  in)

Convert decimal degrees, minutes, and seconds ("dmmss.s") to radians.

Parameters
doublein; // decimal deg,min,sec
Returns
radians

Definition at line 756 of file geoPoint.c.

References DEG_TO_RAD.

757 {
758  double t1,m,d,s,sign;
759 
760  if (in < 0.0)
761  sign = -1.0;
762  else
763  sign = 1.0;
764 
765  s = modf(fabs(in)/100.0, &t1) * 100.0;
766  m = modf(t1/100.0, &d) * 100.0;
767 
768  /* Return radians to calling function. */
769  return( (d + (m/60.0) + (s/3600.0)) * sign * DEG_TO_RAD);
770 }
double geoDms2DD ( double  deg,
double  min,
double  sec,
char *  sign 
)

Converts degrees minutes seconds to Decimal Degrees.

Parameters
doubledeg, min, sec
charsign : N,E,S,W,-,+
Returns
decimal degrees

Definition at line 723 of file geoPoint.c.

References MIN_TO_DEG, and SEC_TO_DEG.

724 {
725  double direction;
726 
727  switch (toupper(*sign))
728  {
729  case '-':
730  case 'W': /* If coordinate is West or South, returned */
731  case 'S':
732  direction = -1.0; /* radian will have negative value. */
733  break;
734  case '+':
735  case 'E': /* If coordinate is East or North, returned */
736  case 'N':
737  direction = 1.0; /* radian will have positive value. */
738  break;
739  default:
740  direction = 1.0; /* If no compass direction entered, returned*/
741  } /* radian will be assumed positive. */
742 
743  /* Return radians to calling function. */
744  return( direction *
745  ( fabs(deg) + (min * MIN_TO_DEG) + (sec * SEC_TO_DEG)));
746 }
double geoDms2Rads ( double  deg,
double  min,
double  sec,
char *  sign 
)

Converts degrees minutes seconds to radians.

Parameters
doubledeg, min, sec
charsign : N,E,S,W,-,+
Returns
radians

Definition at line 692 of file geoPoint.c.

References DEG_TO_RAD, MIN_TO_DEG, and SEC_TO_DEG.

693 {
694  double direction;
695 
696  switch (toupper(*sign))
697  {
698  case '-':
699  case 'W': /* If coordinate is West or South, returned */
700  case 'S':
701  direction = -1.0; /* radian will have negative value. */
702  break;
703  case '+':
704  case 'E': /* If coordinate is East or North, returned */
705  case 'N':
706  direction = 1.0; /* radian will have positive value. */
707  break;
708  default:
709  direction = 1.0; /* If no compass direction entered, returned*/
710  } /* radian will be assumed positive. */
711 
712  /* Return radians to calling function. */
713  return( direction * DEG_TO_RAD *
714  ( fabs(deg) + (min * MIN_TO_DEG) + (sec * SEC_TO_DEG)));
715 }
int geoEfg2XyzDiff ( GEO_LOCATION src_desc,
GEO_LOCATION tgt_desc,
double  xyz_disp[] 
)

This function returns the XYZ offset of the target point with respect to the source point, given the Earth Fixed Geodetic coordinates of the points. The EFG coordinates for the source must appear in the GEO_LOCATION record, which must be built previous to the call to this procedure.

Parameters
GEO_LOCATION*src_desc
GEO_LOCATION*tgt_desc
doublexyz_disp[]
Return values
GEO_OKon success
GEO_ERRORon error
Note
This routine will allow you input site coordinates from two different datums. It is up to the caller to make sure the datums are the same (if it matters).

Definition at line 50 of file geoPoint.c.

References GEO_LOCATION::clat, GEO_LOCATION::clon, GEO_LOCATION::e, GEO_LOCATION::f, GEO_LOCATION::g, GEO_OK, GEO_X, GEO_Y, GEO_Z, GEO_LOCATION::slat, and GEO_LOCATION::slon.

Referenced by geoLlh2DiffX(), geoLlh2DiffY(), and geoLlh2DiffZ().

53 {
54  double delta_e, delta_f, delta_g, intermed_val;
55 
56  // if(src_desc->datum != tgt_desc->datum) return(GEO_ERROR);
57 
58  delta_e = tgt_desc->e - src_desc->e;
59  delta_f = tgt_desc->f - src_desc->f;
60  delta_g = tgt_desc->g - src_desc->g;
61 
62  intermed_val = (-src_desc->clon * delta_e) -
63  ( src_desc->slon * delta_f);
64 
65  xyz_disp[GEO_X] =
66  (-src_desc->slon * delta_e) +
67  ( src_desc->clon * delta_f);
68  xyz_disp[GEO_Y] =
69  (src_desc->slat * intermed_val) +
70  (src_desc->clat * delta_g);
71  xyz_disp[GEO_Z] =
72  (-src_desc->clat * intermed_val) +
73  ( src_desc->slat * delta_g);
74 
75  return(GEO_OK);
76 
77 } /* End procedure site_xyz_diff */
int geoEfg2XyzDiff_packed ( GEO_LOCATION src_desc,
int  count,
double  efg_xyz[] 
)

This function returns the XYZ offset of the target point with respect to the source point, given the Earth Fixed Geodetic coordinates of the points. The EFG coordinates for the source must appear in the GEO_LOCATION record, which must be built previous to the call to this procedure.

Parameters
GEO_LOCATION*src_desc
GEO_LOCATION*tgt_desc
doublexyz_disp[]
Return values
GEO_OKon success
GEO_ERRORon error
Note
This routine will allow you input site coordinates from two different datums. It is up to the caller to make sure the datums are the same (if it matters).

Definition at line 96 of file geoPoint.c.

References GEO_LOCATION::clat, GEO_LOCATION::clon, GEO_LOCATION::e, GEO_LOCATION::f, GEO_LOCATION::g, GEO_E, GEO_F, GEO_G, GEO_OK, GEO_X, GEO_Y, GEO_Z, GEO_LOCATION::slat, and GEO_LOCATION::slon.

97 {
98  double delta_e, delta_f, delta_g, intermed_val;
99 
100  double* vector;
101 
102  // calculate end
103  double* vend = efg_xyz + (3 * count);
104 
105  /* Convert all the sets. */
106  for (vector = efg_xyz; vector < vend; vector += 3)
107  {
108  delta_e = vector[GEO_E] - src_desc->e;
109  delta_f = vector[GEO_F] - src_desc->f;
110  delta_g = vector[GEO_G] - src_desc->g;
111 
112  intermed_val = (-src_desc->clon * delta_e) -
113  ( src_desc->slon * delta_f);
114 
115  vector[GEO_X] =
116  (-src_desc->slon * delta_e) +
117  ( src_desc->clon * delta_f);
118  vector[GEO_Y] =
119  (src_desc->slat * intermed_val) +
120  (src_desc->clat * delta_g);
121  vector[GEO_Z] =
122  (-src_desc->clat * intermed_val) +
123  ( src_desc->slat * delta_g);
124  }
125 
126  return(GEO_OK);
127 
128 } /* End procedure efg_xyz_diff */
double geoLlh2DiffX ( double  lat1,
double  lon1,
double  hgt1,
int  datum1,
double  lat2,
double  lon2,
double  hgt2,
int  datum2 
)

Definition at line 130 of file geoPoint.c.

References GEO_X, geoEfg2XyzDiff(), and geoInitLocation().

131 {
132  double xyz[3];
133  GEO_LOCATION site1, site2;
134 
135  geoInitLocation(&site1, lat1, lon1, hgt1, datum1, "site1");
136  geoInitLocation(&site2, lat2, lon2, hgt2, datum2, "site2");
137  geoEfg2XyzDiff(&site1, &site2,xyz);
138  return(xyz[GEO_X]);
139 
140 }

Here is the call graph for this function:

double geoLlh2DiffY ( double  lat1,
double  lon1,
double  hgt1,
int  datum1,
double  lat2,
double  lon2,
double  hgt2,
int  datum2 
)

Definition at line 142 of file geoPoint.c.

References GEO_Y, geoEfg2XyzDiff(), and geoInitLocation().

143 {
144  double xyz[3];
145  GEO_LOCATION site1, site2;
146 
147  geoInitLocation(&site1, lat1, lon1, hgt1, datum1, "site1");
148  geoInitLocation(&site2, lat2, lon2, hgt2, datum2, "site2");
149  geoEfg2XyzDiff(&site1, &site2,xyz);
150  return(xyz[GEO_Y]);
151 
152 }

Here is the call graph for this function:

double geoLlh2DiffZ ( double  lat1,
double  lon1,
double  hgt1,
int  datum1,
double  lat2,
double  lon2,
double  hgt2,
int  datum2 
)

Definition at line 153 of file geoPoint.c.

References GEO_Z, geoEfg2XyzDiff(), and geoInitLocation().

154 {
155  double xyz[3];
156  GEO_LOCATION site1, site2;
157 
158  geoInitLocation(&site1, lat1, lon1, hgt1, datum1, "site1");
159  geoInitLocation(&site2, lat2, lon2, hgt2, datum2, "site2");
160  geoEfg2XyzDiff(&site1, &site2,xyz);
161  return(xyz[GEO_Z]);
162 
163 }

Here is the call graph for this function:

double geoLlh2E ( double  lat,
double  lon,
double  hgt,
int  datum 
)

This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the E component.

Parameters
doublelat : Latitude in RADIANS
doublelon : Longitude in RADIANS
doubleheight : Height in METERS
intdatum1
Returns
double *e : E(x) in METERS

Definition at line 262 of file geoPoint.c.

References GEO_LOCATION::e, and geoInitLocation().

264 {
265  GEO_LOCATION site;
266 
267  geoInitLocation(&site, lat, lon, hgt, datum, "site");
268  return(site.e);
269 }

Here is the call graph for this function:

void geoLlh2Efg ( double  lat,
double  lon,
double  height,
int  datum,
double *  e,
double *  f,
double *  g 
)

This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G).

Parameters
doublelat : Latitude in RADIANS
doublelon : Longitude in RADIANS
doubleheight : Height in METERS
intdatum1
double*e : E(x) in METERS
double*f : F(y) in METERS
double*g : G(z) in METERS
Returns
nothing

Definition at line 181 of file geoPoint.c.

References geoGetEllipsoid().

184 {
185  double N,a,e2,ee2, b, flat;
186  double slat, clat;
187 
188  clat = cos(lat);
189  slat = sin(lat);
190 
191  /* Get the ellipsoid parameters */
192  geoGetEllipsoid(&a,&b,&e2,&ee2,&flat,datum);
193 
194  /* Compute the radius of curvature */
195  N = a / (sqrt(1.0 - e2 * pow(slat,2.0)));
196 
197  /* Compute the EFG(XYZ) coordinates (earth centered) */
198  *e = (N + height) * clat * cos(lon);
199  *f = (N + height) * clat * sin(lon);
200  *g = (N * (1.0 - e2) + height) * sin(lat);
201 }

Here is the call graph for this function:

void geoLlh2Efg_packed ( GEO_DATUM datum,
int  count,
double  llh_efg[] 
)

This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G).

Parameters
GEO_DATUM*datum Datum of the lat/lon/height
intcount Number of vertices
doublellh_efg[] Pointer to the first coordinate of the first vertex of the array.
Returns
nothing

Definition at line 214 of file geoPoint.c.

References GEO_DATUM::a, GEO_DATUM::e2, GEO_E, GEO_F, GEO_G, GEO_HGT, GEO_LAT, and GEO_LON.

215 {
216  double N, slat, clat, slon, clon;
217 
218  double* vector;
219 
220  // calculate end
221  double* vend = llh_efg + (3 * count);
222 
223  /* Convert all the sets. */
224  for (vector = llh_efg; vector < vend; vector += 3)
225  {
226 #ifdef HAVE_SINCOS
227  sincos(vector[GEO_LAT], &slat, &clat);
228  sincos(vector[GEO_LON], &slon, &clon);
229 
230 #else
231  slat = sin(vector[GEO_LAT]);
232  clat = cos(vector[GEO_LAT]);
233  slon = sin(vector[GEO_LON]);
234  clon = cos(vector[GEO_LON]);
235 #endif
236 
237  /* Compute the radius of curvature */
238  N = datum->a / (sqrt(1.0 - datum->e2 * slat * slat));
239 
240  /* Compute the EFG(XYZ) coordinates (earth centered) */
241  vector[GEO_E] = (N + vector[GEO_HGT]) * clat * clon;
242  vector[GEO_F] = (N + vector[GEO_HGT]) * clat * slon;
243  vector[GEO_G] = (N * (1.0 - datum->e2) + vector[GEO_HGT]) * slat;
244  }
245 }
double geoLlh2F ( double  lat,
double  lon,
double  hgt,
int  datum 
)

This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the F component.

Parameters
doublelat : Latitude in RADIANS
doublelon : Longitude in RADIANS
doubleheight : Height in METERS
intdatum1
Returns
double *f : F(x) in METERS

Definition at line 286 of file geoPoint.c.

References GEO_LOCATION::f, and geoInitLocation().

288 {
289  GEO_LOCATION site;
290 
291  geoInitLocation(&site, lat, lon, hgt, datum, "site");
292  return(site.f);
293 }

Here is the call graph for this function:

double geoLlh2G ( double  lat,
double  lon,
double  hgt,
int  datum 
)

This routine will convert geodetic coordinates (latitude $\phi$, longitude $\lambda$, and ellipsoid height $h$) into earth centered Cartesian coordinates (E,F,G). This returns the G component.

Parameters
doublelat : Latitude in RADIANS
doublelon : Longitude in RADIANS
doubleheight : Height in METERS
intdatum1
Returns
double *g : G(x) in METERS

Definition at line 310 of file geoPoint.c.

References GEO_LOCATION::g, and geoInitLocation().

312 {
313  GEO_LOCATION site;
314 
315  geoInitLocation(&site, lat, lon, hgt, datum, "site");
316  return(site.g);
317 }

Here is the call graph for this function:

double geoRads2DD ( double  rads)

Converts radians to Decimal Degrees.

Parameters
doubleradians
Returns
double decimal degrees

Definition at line 836 of file geoPoint.c.

References RAD_TO_DEG.

837 {
838  return(rads * RAD_TO_DEG);
839 }
double geoRads2Decdms ( double  rads)

Convert radians to decimal degrees, minutes, and seconds ("dddmmss.s").

Parameters
doublerads
Returns
double Decimal Deg/min/sec

Definition at line 810 of file geoPoint.c.

References RAD_TO_DEG.

811 {
812  double d,m,s,sign;
813  double frac;
814 
815  if (rads < 0.0)
816  sign = -1.0;
817  else
818  sign = 1.0;
819 
820  frac = modf(fabs(rads * RAD_TO_DEG),&d);
821  frac = modf(frac * 60.0,&m);
822  s = frac * 60.0;
823 
824  /* Return dddmmss.s to calling function. */
825  return(((d*10000.0)+(m*100)+s)*sign);
826 }
void geoRads2Dms ( double  rads,
double *  deg,
double *  min,
double *  sec,
double *  dir 
)

Converts radians to degrees minutes seconds.

Parameters
doublerads
doubledeg, min, sec
chardir : -1.0 or 1.0
Returns
nothing

Definition at line 782 of file geoPoint.c.

References RAD_TO_DEG.

784 {
785  double temp;
786  double fraction;
787 
788  if (rads < 0.0)
789  *dir = -1.0;
790  else
791  *dir = 1.0;
792 
793  rads = fabs (rads);
794 
795  temp = RAD_TO_DEG * rads;
796  fraction = modf(temp,deg);
797 
798  temp = fraction * 60.0;
799  fraction = modf(temp,min);
800 
801  *sec = fraction * 60.0;
802 }
void geoRae2Efg ( GEO_LOCATION loc,
double  aer_in[],
double  efg_out[] 
)

Ingests Range, Azimuth, Elevation and site info and returns the EFG coordinates that the RAE points to.

Parameters
GEO_LOCATION*loc
doubleaer_in[]
doubleefg_out[]
Returns
nothing

Definition at line 654 of file geoPoint.c.

References GEO_LOCATION::clat, GEO_LOCATION::clon, GEO_LOCATION::e, GEO_LOCATION::f, GEO_LOCATION::g, GEO_E, GEO_F, GEO_G, GEO_X, GEO_Y, GEO_Z, geoRae2Xyz(), GEO_LOCATION::slat, and GEO_LOCATION::slon.

657 {
658  double c1, c2, c3;
659  double xyz_val[3];
660 
661  /* Convert the RAE value to XYZ: */
662  geoRae2Xyz(aer_in, xyz_val);
663 
664  /* Do the matrix multiplication: */
665 
666  c1 = -loc->slon * xyz_val[GEO_X] +
667  -loc->clon * loc->slat * xyz_val[GEO_Y] +
668  loc->clon * loc->clat * xyz_val[GEO_Z];
669 
670  c2 = loc->clon * xyz_val[GEO_X] +
671  -loc->slon * loc->slat * xyz_val[GEO_Y] +
672  loc->slon * loc->clat * xyz_val[GEO_Z];
673 
674  c3 = loc->clat * xyz_val[GEO_Y] +
675  loc->slat * xyz_val[GEO_Z];
676 
677  /* Add resultant matrix to local EFG to get remote EFG: */
678  efg_out[GEO_E] = c1 + loc->e;
679  efg_out[GEO_F] = c2 + loc->f;
680  efg_out[GEO_G] = c3 + loc->g;
681 
682 } /* End procedure aer_to_efg */

Here is the call graph for this function:

void geoRae2Xyz ( double  rae_in[],
double  xyz_out[] 
)

This routine converts from Range, Azimuth, and Elevation into Cartesian coordinates X,Y,Z.

Parameters
doublerae_in[]
doublexyz_out[]
Returns
nothing

Definition at line 501 of file geoPoint.c.

References GEO_AZ, GEO_EL, GEO_RNG, GEO_X, GEO_Y, and GEO_Z.

Referenced by geoRae2Efg().

503 {
504  double r_cos_e;
505 
506  r_cos_e = rae_in[GEO_RNG] * cos(rae_in[GEO_EL]);
507 
508  xyz_out[GEO_X] = sin(rae_in[GEO_AZ]) * r_cos_e;
509  xyz_out[GEO_Y] = cos(rae_in[GEO_AZ]) * r_cos_e;
510  xyz_out[GEO_Z] = rae_in[GEO_RNG] * sin(rae_in[GEO_EL]);
511 }
void geoRae2Xyz_packed ( int  count,
double  rae_xyz[] 
)

This routine iterates over count vertices in the array rae_xyz converting it from Range, Azimuth, and Elevation into Cartesian coordinates X,Y,Z.

Parameters
countNumber of vertices
rae_xyz[]Pointer to the first coordinate of the first vertex of the array.
Returns
nothing

Definition at line 522 of file geoPoint.c.

References GEO_AZ, GEO_EL, GEO_RNG, GEO_X, GEO_Y, and GEO_Z.

524 {
525  double r_cos_e;
526  double r_sin_e;
527 #ifdef HAVE_SINCOS
528  double sin_a;
529 #endif
530  double cos_a;
531  double* vector;
532 
533  // calculate end
534  double* vend = rae_xyz + (3 * count);
535 
536  /* Convert all the sets. */
537  for (vector = rae_xyz; vector < vend; vector += 3)
538  {
539 #ifdef HAVE_SINCOS
540  sincos(vector[GEO_EL], &r_sin_e, &r_cos_e);
541  sincos(vector[GEO_AZ], &sin_a, &cos_a);
542  r_cos_e *= vector[GEO_RNG];
543  r_sin_e *= vector[GEO_RNG];
544 
545  vector[GEO_X] = sin_a * r_cos_e;
546  vector[GEO_Y] = cos_a * r_cos_e;
547  vector[GEO_Z] = r_sin_e;
548 #else
549  r_cos_e = vector[GEO_RNG] * cos(vector[GEO_EL]);
550  r_sin_e *= vector[GEO_RNG] * sin(vector[GEO_EL]);
551  cos_a = cos(vector[GEO_AZ]);
552 
553  vector[GEO_X] = sin(vector[GEO_AZ]) * r_cos_e;
554  vector[GEO_Y] = cos_a * r_cos_e;
555  vector[GEO_Z] = vector[GEO_RNG] ;
556 #endif
557  }
558 }
double geoVersion ( void  )

This function returns the version of the GeoStarsLib.

Return values
GEOSTARSLIB_VERSIONas a double

Definition at line 25 of file geoPoint.c.

References GEOSTARSLIB_VERSION.

26 {
27  return(GEOSTARSLIB_VERSION);
28 }
double geoXyz2A ( double  x,
double  y,
double  z 
)

Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the azimuth to that point.

Notes:

  • X is understood to be the east-west displacement of the point, with east being the positive direction.
  • Y is the north-south displacement, with north being positive.
  • Z is the vertical displacement of the point.
  • Azimuth is in decimal degrees
Parameters
doublexyz_in[]
Returns
double azimuth

Definition at line 455 of file geoPoint.c.

References GEO_AZ, geoXyz2Rae(), and RAD_TO_DEG.

456 {
457  double rae[3],xyz[3];
458  xyz[0]=x;
459  xyz[1]=y;
460  xyz[2]=z;
461  geoXyz2Rae(xyz,rae);
462  return(rae[GEO_AZ]*RAD_TO_DEG);
463 }

Here is the call graph for this function:

double geoXyz2E ( double  x,
double  y,
double  z 
)

Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the elevation angle to that point.

Notes:

  • X is understood to be the east-west displacement of the point, with east being the positive direction.
  • Y is the north-south displacement, with north being positive.
  • Z is the vertical displacement of the point.
  • Elevation is in decimal degrees
Parameters
doublexyz_in[]
Returns
double elevation

Definition at line 482 of file geoPoint.c.

References GEO_EL, geoXyz2Rae(), and RAD_TO_DEG.

483 {
484  double rae[3],xyz[3];
485  xyz[0]=x;
486  xyz[1]=y;
487  xyz[2]=z;
488  geoXyz2Rae(xyz,rae);
489  return(rae[GEO_EL]*RAD_TO_DEG);
490 }

Here is the call graph for this function:

void geoXyz2Efg ( GEO_LOCATION loc,
double  xyz_in[],
double  efg_out[] 
)

Ingests X, Y, Z, and site info and returns the EFG coordinates.

Parameters
GEO_LOCATION*loc
doublexyz_in[]
doubleefg_out[]
Returns
nothing

Definition at line 570 of file geoPoint.c.

References GEO_LOCATION::clat, GEO_LOCATION::clon, GEO_LOCATION::e, GEO_LOCATION::f, GEO_LOCATION::g, GEO_E, GEO_F, GEO_G, GEO_X, GEO_Y, GEO_Z, GEO_LOCATION::slat, and GEO_LOCATION::slon.

573 {
574  double c1, c2, c3;
575 
576  /* Do the matrix multiplication: */
577 
578  c1 = -loc->slon * xyz_in[GEO_X] +
579  -loc->clon * loc->slat * xyz_in[GEO_Y] +
580  loc->clon * loc->clat * xyz_in[GEO_Z];
581 
582  c2 = loc->clon * xyz_in[GEO_X] +
583  -loc->slon * loc->slat * xyz_in[GEO_Y] +
584  loc->slon * loc->clat * xyz_in[GEO_Z];
585 
586  c3 = loc->clat * xyz_in[GEO_Y] +
587  loc->slat * xyz_in[GEO_Z];
588 
589  /* Add resultant matrix to local EFG to get remote EFG: */
590  efg_out[GEO_E] = c1 + loc->e;
591  efg_out[GEO_F] = c2 + loc->f;
592  efg_out[GEO_G] = c3 + loc->g;
593 
594 } /* End procedure aer_to_efg */
void geoXyz2Efg_packed ( GEO_LOCATION loc,
int  count,
double  xyz_efg[] 
)

This routine takes site info and iterates over count vertices in the tightly packed array xyz_efg converting it from cartesian XYZ into geocentric EFG coordinates.

Parameters
locThe location that the XYZ grid is based around
countThe number of vertices in the array
xyz_efg[]Pointer to the first vertex of the array
Returns
nothing

Definition at line 606 of file geoPoint.c.

References GEO_LOCATION::clat, GEO_LOCATION::clon, GEO_LOCATION::e, GEO_LOCATION::f, GEO_LOCATION::g, GEO_E, GEO_F, GEO_G, GEO_X, GEO_Y, GEO_Z, GEO_LOCATION::slat, and GEO_LOCATION::slon.

609 {
610  double c1, c2, c3;
611  double nclonslat = -loc->clon * loc->slat;
612  double clonclat = loc->clon * loc->clat;
613  double nslonslat = -loc->slon * loc->slat;
614  double slonclat = loc->slon * loc->clat;
615  double* vector;
616 
617  // calculate end
618  double* vend = xyz_efg + (3 * count);
619 
620  /* Convert all the sets. */
621  for (vector = xyz_efg; vector < vend; vector += 3)
622  {
623  /* Do the matrix multiplication: */
624  c1 = -loc->slon * vector[GEO_X] +
625  nclonslat * vector[GEO_Y] +
626  clonclat * vector[GEO_Z];
627 
628  c2 = loc->clon * vector[GEO_X] +
629  nslonslat * vector[GEO_Y] +
630  slonclat * vector[GEO_Z];
631 
632  c3 = loc->clat * vector[GEO_Y] +
633  loc->slat * vector[GEO_Z];
634 
635  /* Add resultant matrix to local EFG to get remote EFG: */
636  vector[GEO_E] = c1 + loc->e;
637  vector[GEO_F] = c2 + loc->f;
638  vector[GEO_G] = c3 + loc->g;
639  }
640 
641 } /* End procedure xyz_to_efg */
double geoXyz2R ( double  x,
double  y,
double  z 
)

Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2R calculates the range to that point.

Notes:

  • X is understood to be the east-west displacement of the point, with east being the positive direction.
  • Y is the north-south displacement, with north being positive.
  • Z is the vertical displacement of the point.
  • Range is in meters.
Parameters
doublexyz_in[]
Returns
double range

Definition at line 428 of file geoPoint.c.

References GEO_RNG, and geoXyz2Rae().

429 {
430  double rae[3],xyz[3];
431  xyz[0]=x;
432  xyz[1]=y;
433  xyz[2]=z;
434  geoXyz2Rae(xyz,rae);
435  return(rae[GEO_RNG]);
436 }

Here is the call graph for this function:

void geoXyz2Rae ( double  xyz_in[],
double  rae_out[] 
)

Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2Rae calculates the range, azimuth, and elevation to that point.

Notes:

  • X is understood to be the east-west displacement of the point, with east being the positive direction.
  • Y is the north-south displacement, with north being positive.
  • Z is the vertical displacement of the point.
  • Range is in meters. Azimuth and elevation are in radians
Parameters
doublexyz_in[]
doublerae_out[]
Returns
nothing

Definition at line 339 of file geoPoint.c.

References GEO_AZ, GEO_EL, GEO_RNG, GEO_X, GEO_Y, GEO_Z, and M_PI.

Referenced by geoXyz2A(), geoXyz2E(), and geoXyz2R().

341 {
342  double horz_dist;
343 
344  /* Determine the range: */
345  rae_out[GEO_RNG] =
346  sqrt(pow(xyz_in[GEO_X],2.0) + pow(xyz_in[GEO_Y],2.0) + pow(xyz_in[GEO_Z],2.0));
347 
348  /* Determine the azimuth: */
349  rae_out[GEO_AZ] = atan2(xyz_in[GEO_X], xyz_in[GEO_Y]);
350  //if((xyz_in[GEO_X] >= 0.0) && (xyz_in[GEO_Y] < 0.0)) rae_out[GEO_AZ] += M_PI;
351  if((xyz_in[GEO_X] < 0.0) && (xyz_in[GEO_Y] < 0.0)) rae_out[GEO_AZ] += (2.0 * M_PI);
352  if((xyz_in[GEO_X] < 0.0) && (xyz_in[GEO_Y] > 0.0)) rae_out[GEO_AZ] += (2.0 * M_PI);
353 
354  /* Determine the elevation: */
355  horz_dist = sqrt(pow(xyz_in[GEO_X],2.0) + pow(xyz_in[GEO_Y],2.0));
356  rae_out[GEO_EL] = atan2(xyz_in[GEO_Z], horz_dist);
357  if(rae_out[GEO_EL] < 0.0) rae_out[GEO_EL] += (2.0 * M_PI);
358 
359 } /* End procedure xyz_to_rae */
void geoXyz2Rae_packed ( int  count,
double  xyz_rae[] 
)

Given the X, Y, and Z coordinates (in meters) of a point in space, the procedure geoXyz2Rae calculates the range, azimuth, and elevation to that point.

Notes:

  • X is understood to be the east-west displacement of the point, with east being the positive direction.
  • Y is the north-south displacement, with north being positive.
  • Z is the vertical displacement of the point.
  • Range is in meters. Azimuth and elevation are in radians
    Parameters
    countNumber of vertices
    xyz_rae[]Pointer to the first coordinate of the first vertex of the array.
    Returns
    nothing

Definition at line 378 of file geoPoint.c.

References GEO_AZ, GEO_EL, GEO_RNG, GEO_X, GEO_Y, GEO_Z, and M_PI.

380 {
381  double horz_dist;
382  double rng, az;
383  double* vector;
384 
385  // calculate end
386  double* vend = xyz_rae + (3 * count);
387 
388  /* Convert all the sets. */
389  for (vector = xyz_rae; vector < vend; vector += 3)
390  {
391  /* Determine the range: */
392  rng =
393  sqrt(vector[GEO_X]*vector[GEO_X] + vector[GEO_Y]*vector[GEO_Y] + vector[GEO_Z]*vector[GEO_Z]);
394 
395  /* Determine the azimuth: */
396  az = atan2(vector[GEO_X], vector[GEO_Y]);
397 
398  //if((vector[GEO_X] < 0.0) && (vector[GEO_Y] < 0.0)) az += (2.0 * M_PI);
399  //if((vector[GEO_X] < 0.0) && (vector[GEO_Y] > 0.0)) az += (2.0 * M_PI);
400  if (az < -0.0f) az += (2.0 * M_PI);
401 
402  /* Determine the elevation: */
403  horz_dist = sqrt(vector[GEO_X]*vector[GEO_X] + vector[GEO_Y]*vector[GEO_Y]);
404  vector[GEO_RNG] = rng;
405  vector[GEO_AZ] = az;
406  vector[GEO_EL] = atan2(vector[GEO_Z], horz_dist);
407  if(vector[GEO_EL] < -0.0f) vector[GEO_EL] += (2.0 * M_PI);
408  }
409 } /* End procedure xyz_to_rae */