1 /*
2  * Copyright (c) 2013-2021, The OSKAR Developers.
3  * See the LICENSE file at the top-level directory of this distribution.
4  */
5 
6 #include "convert/oskar_convert_geodetic_spherical_to_ecef.h"
7 #include <math.h>
8 
9 #ifdef __cplusplus
10 extern "C" {
11 #endif
12 
oskar_convert_geodetic_spherical_to_ecef(int num_points,const double * lon_rad,const double * lat_rad,const double * alt_metres,double * x,double * y,double * z)13 void oskar_convert_geodetic_spherical_to_ecef(int num_points,
14         const double* lon_rad, const double* lat_rad, const double* alt_metres,
15         double* x, double* y, double* z)
16 {
17     const double a = 6378137.000; /* Equatorial radius (semi-major axis). */
18     const double b = 6356752.314; /* Polar radius (semi-minor axis). */
19     const double e2 = 1.0 - (b*b) / (a*a);
20     int i = 0;
21     for (i = 0; i < num_points; ++i)
22     {
23         const double lat_ = lat_rad[i];
24         const double lon_ = lon_rad[i];
25         const double alt_ = alt_metres[i];
26         const double sin_lat = sin(lat_);
27         const double cos_lat = cos(lat_);
28         const double sin_lon = sin(lon_);
29         const double cos_lon = cos(lon_);
30         const double n_phi = a / sqrt(1.0 - e2 * sin_lat * sin_lat);
31         x[i] = (n_phi + alt_) * cos_lat * cos_lon;
32         y[i] = (n_phi + alt_) * cos_lat * sin_lon;
33         z[i] = ((1.0 - e2) * n_phi + alt_) * sin_lat;
34     }
35 }
36 
37 #ifdef __cplusplus
38 }
39 #endif
40