Convert from Longitude / Latitude to Cartesian Coordinates

Question

I have coordinate coordinates with the earth set as latitude and longitude ( WGS-84 ).

How can I convert them to Cartesian coordinates (x, y, z) with the origin in the center of the Earth?

Answer 1

I recently did something similar to this using the Haversin Formula according to WGS-84, which is derived from the Haversines Law with very satisfactory results.

Yes, WGS-84 suggests that the Earth is an ellipsoid, but I believe that you get an approximately 0.5% average error using an approach such as the Haversin Formula, which may be an acceptable error in your case. You will always have some mistake if you are not talking about a distance of several feet, and even then theoretically the curvature of the Earth ... If you need a more stringent verification of compliance with WGS-84 on the "Vincenty Formula".

I understand where starblue comes from , but good software development is often a compromise, so it all depends on the accuracy you require from what you do. For example, a result calculated using the Manhattan Distance Formula formula compared to a Distance Formula result may be better for certain situations, as it is estimated to be less expensive. Think "which moment is closer?" Scenarios in which you do not need accurate distance measurement.

As for the Haversin Formula, it is easy to implement and good, because it uses Spherical Trigonometry instead of the Cosine Law approach based on two-dimensional trigonometry, so you get a nice balance of accuracy in complexity.

Gentlemen named Chris Veness have a great website at http://www.movable-type.co.uk/scripts/latlong.html that explains some of the concepts you are interested in and demonstrate various software implementations; this should also answer your question about changing x / y.

Answer 2

Here is the answer I found:

To make the definition complete, in a Cartesian coordinate system:

• the x axis goes through long, lat (0,0), so longitude 0 meets the equator;
• the y axis passes (0.90);
• and the z axis goes through the poles.

Conversion:

x = R * cos(lat) * cos(lon) y = R * cos(lat) * sin(lon) z = R *sin(lat) 

Where R is the approximate radius of the earth (e.g. 6371KM).

If your trigonometric functions expect radians (which they probably do), you will need to first convert your longitude and latitude to radians. You obviously need a decimal representation, not degrees \ minutes \ seconds (see conversion here ).

The formula for the inverse transformation is:

  lat = asin(z / R) lon = atan2(y, x) 

asin is, of course, the sine. reads about atan2 on wikipedia . Remember to translate back from radians to degrees.

This page gives C # code for this (note that it is very different from the formulas), as well as some explanation and a good diagram of the reasons for this being correct.

Answer 3

Theory of Converting GPS（WGS84) to Cartesian Coordinates https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_geodetic_to_ECEF_coordinates

Below I use the following:

• GPS longitude (WGS84) and Cartesian coordinates are the same.
• The latitude should be converted according to the parameters of the ellipsoid WGS 84, the half-axis - 6378137 m, and
• Mutual alignment is 298.257223563.

I attached the VB code that I wrote:

 Imports System.Math 'Input GPSLatitude is WGS84 Latitude,h is altitude above the WGS 84 ellipsoid Public Function GetSphericalLatitude(ByVal GPSLatitude As Double, ByVal h As Double) As Double Dim A As Double = 6378137 'semi-major axis Dim f As Double = 1 / 298.257223563 '1/f Reciprocal of flattening Dim e2 As Double = f * (2 - f) Dim Rc As Double = A / (Sqrt(1 - e2 * (Sin(GPSLatitude * PI / 180) ^ 2))) Dim p As Double = (Rc + h) * Cos(GPSLatitude * PI / 180) Dim z As Double = (Rc * (1 - e2) + h) * Sin(GPSLatitude * PI / 180) Dim r As Double = Sqrt(p ^ 2 + z ^ 2) Dim SphericalLatitude As Double = Asin(z / r) * 180 / PI Return SphericalLatitude End Function 

Note that h is the height above WGS 84 ellipsoid .

Usually GPS will give us h above MSL height. The MSL height must be converted to a height h above the WGS 84 ellipsoid using the EGM96 geopotential model (Lemoine et al, 1998).
This is done by interpolating the grid of the geoid height file with a spatial resolution of 15 minutes of arc.

Or, if you have a GPS professional level, there is a height h (msl, altitude above average sea level) and UNDULATION , the relationship between the geoid and ellipsoid (m) selected source point, derived from the internal table. you can get h = H(msl) + undulation

To XYZ Cartesian coordinates:

 x = R * cos(lat) * cos(lon) y = R * cos(lat) * sin(lon) z = R *sin(lat) 

Answer 4

Why introduce something that has already been implemented and tested for testing?

C #, for one, has NetTopologySuite , which is the .NET port of the JTS topology topology.

In particular, you have a serious flaw in your calculations. Earth is not an ideal sphere, and approximation of the radius of the earth may not cut it out for accurate measurements.

If in some cases it is acceptable to use the homebrew functions, GIS is a good example of a field in which it is much preferable to use a reliable, tested testing library.

Answer 5

If you want to get coordinates based on an ellipsoid, not a sphere, see http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF - it gives the formulas as well as the WGS84 constants needed for the conversion.

Formulas that also take into account the height relative to the surface of the reference ellipsoid (useful if you get height data from a GPS device).

Answer 6

proj.4 software provides a command line program that can perform conversion, for example

 LAT=40 LON=-110 echo $LON $LAT | cs2cs +proj=latlong +datum=WGS84 +to +proj=geocent +datum=WGS84 

It also provides a C API . In particular, the pj_geodetic_to_geocentric function pj_geodetic_to_geocentric perform the conversion without having to first create a projection object.

Answer 7

 Coordinate[] coordinates = new Coordinate[3]; coordinates[0] = new Coordinate(102, 26); coordinates[1] = new Coordinate(103, 25.12); coordinates[2] = new Coordinate(104, 16.11); CoordinateSequence coordinateSequence = new CoordinateArraySequence(coordinates); Geometry geo = new LineString(coordinateSequence, geometryFactory); CoordinateReferenceSystem wgs84 = DefaultGeographicCRS.WGS84; CoordinateReferenceSystem cartesinaCrs = DefaultGeocentricCRS.CARTESIAN; MathTransform mathTransform = CRS.findMathTransform(wgs84, cartesinaCrs, true); Geometry geo1 = JTS.transform(geo, mathTransform); 

Answer 8

You can do it in Java.

 public List<Double> convertGpsToECEF(double lat, double longi, float alt) { double a=6378.1; double b=6356.8; double N; double e= 1-(Math.pow(b, 2)/Math.pow(a, 2)); N= a/(Math.sqrt(1.0-(e*Math.pow(Math.sin(Math.toRadians(lat)), 2)))); double cosLatRad=Math.cos(Math.toRadians(lat)); double cosLongiRad=Math.cos(Math.toRadians(longi)); double sinLatRad=Math.sin(Math.toRadians(lat)); double sinLongiRad=Math.sin(Math.toRadians(longi)); double x =(N+0.001*alt)*cosLatRad*cosLongiRad; double y =(N+0.001*alt)*cosLatRad*sinLongiRad; double z =((Math.pow(b, 2)/Math.pow(a, 2))*N+0.001*alt)*sinLatRad; List<Double> ecef= new ArrayList<>(); ecef.add(x); ecef.add(y); ecef.add(z); return ecef; }