-
-
Save klucar/1536194 to your computer and use it in GitHub Desktop.
/* | |
* | |
* ECEF - Earth Centered Earth Fixed | |
* | |
* LLA - Lat Lon Alt | |
* | |
* ported from matlab code at | |
* https://gist.github.com/1536054 | |
* and | |
* https://gist.github.com/1536056 | |
*/ | |
// WGS84 ellipsoid constants | |
private final double a = 6378137; // radius | |
private final double e = 8.1819190842622e-2; // eccentricity | |
private final double asq = Math.pow(a,2); | |
private final double esq = Math.pow(e,2); | |
private double[] ecef2lla(double[] ecef){ | |
double x = ecef[0]; | |
double y = ecef[1]; | |
double z = ecef[2]; | |
double b = Math.sqrt( asq * (1-esq) ); | |
double bsq = Math.pow(b,2); | |
double ep = Math.sqrt( (asq - bsq)/bsq); | |
double p = Math.sqrt( Math.pow(x,2) + Math.pow(y,2) ); | |
double th = Math.atan2(a*z, b*p); | |
double lon = Math.atan2(y,x); | |
double lat = Math.atan2( (z + Math.pow(ep,2)*b*Math.pow(Math.sin(th),3) ), (p - esq*a*Math.pow(Math.cos(th),3)) ); | |
double N = a/( Math.sqrt(1-esq*Math.pow(Math.sin(lat),2)) ); | |
double alt = p / Math.cos(lat) - N; | |
// mod lat to 0-2pi | |
lon = lon % (2*Math.PI); | |
// correction for altitude near poles left out. | |
double[] ret = {lat, lon, alt}; | |
return ret; | |
} | |
private double[] lla2ecef(double[] lla){ | |
double lat = lla[0]; | |
double lon = lla[1]; | |
double alt = lla[2]; | |
double N = a / Math.sqrt(1 - esq * Math.pow(Math.sin(lat),2) ); | |
double x = (N+alt) * Math.cos(lat) * Math.cos(lon); | |
double y = (N+alt) * Math.cos(lat) * Math.sin(lon); | |
double z = ((1-esq) * N + alt) * Math.sin(lat); | |
double[] ret = {x, y, z}; | |
return ret; | |
} | |
Thanks for the catch. I updated line 28.
As for the radius, 6378137 meters is the equatorial radius of the earth ellipsoid as defined by the WGS-84 model.
http://en.wikipedia.org/wiki/WGS-84#A_new_World_Geodetic_System:_WGS_84
Hello,
Do you know and can you point some articles or literature where there is described this conversion?
I need it to my thesis and I have problem with finding something usefull :(
Note that lat and lon are in radians, not degrees! Use Math.toDegress() and Math.toRadians().
pthon code
def lla2ecef(lon,lat,alt):
a = 6378137
e = 8.1819190842622e-2
esq=pow(e,2)
N = a / sqrt(1 - esq * pow(sin(radians(lat)), 2))
x = (N + alt) * cos(radians(lat)) * cos(radians(lon))
y = (N + alt) * cos(radians(lat) )* sin(radians(lon))
z = ((1 - esq) * N + alt) * sin(radians(lat))
return x,y,z
The ecef2lla
method cannot possibly be correct, because there is no closed form solution for this, it has to be solved iteratively. See:
https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_ECEF_to_geodetic_coordinates
Hey! You're missing the Math class name for the sqrt operation at line 28.
Thanks for your code!
Edit, added question: what units are you considering? earth radius is 6371 in km. :/