从经纬度转换到笛卡尔坐标。

时间:2021-09-16 16:03:12

I have some earth-centered coordinate points given as latitude and longitude (WGS-84).

我有一些以地球为中心的坐标点作为纬度和经度(WGS-84)。

How can i convert them to Cartesian coordinates (x,y,z) with the origin at the center of the earth?

如何将它们转换为笛卡尔坐标(x,y,z)与地球中心的原点?

7 个解决方案

#1


34  

I have recently done something similar to this using the "Haversine Formula" on WGS-84 data, which is a derivative of the "Law of Haversines" with very satisfying results.

我最近在WGS-84数据上使用“Haversine公式”做了类似的事情,这是“Haversines定律”的导数,结果非常令人满意。

Yes, WGS-84 assumes the Earth is an ellipsoid, but I believe you only get about a 0.5% average error using an approach like the "Haversine Formula", which may be an acceptable amount of error in your case. You will always have some amount of error unless you're talking about a distance of a few feet and even then there is theoretically curvature of the Earth... If you require more rigidly WGS-84 compatible approach checkout the "Vincenty Formula."

是的,WGS-84假定地球是一个椭球体,但我相信,你只会得到0.5%的平均误差,使用的方法类似于“Haversine公式”,这可能是你的情况中一个可接受的错误量。你总是会有一些误差,除非你说的是几英尺的距离,即使是在理论上地球的曲率。如果您需要更严格的WGS-84兼容方法,请签出“Vincenty公式”。

I understand where starblue is coming from, but good software engineering is often about trade offs, so it all depends on the accuracy you require for what you are doing. For example, the result calculated from "Manhattan Distance Formula" versus the result from the "Distance Formula" can be better for certain situations as it is computationally less expensive. Think "which point is closest?" scenarios where you don't need a precise distance measurement.

我理解starblue的来源,但是好的软件工程经常是关于权衡的,所以这完全取决于你所做的事情的准确性。例如,从“曼哈顿距离公式”计算得出的结果与“距离公式”的结果相比,在某些情况下会更好,因为它在计算上更便宜。想想“哪个点最接近”,你不需要精确的距离测量。

Regarding, the "Haversine Formula" it is easy to implement and is nice because it is uses "Spherical Trigonometry" instead of a "Law of Cosines" based approach which is based on two-dimensional trigonometry, therefore you get a nice balance of accuracy over complexity.

关于“Haversine公式”,它很容易实现,而且很好,因为它使用的是“球面三角学”,而不是基于二维三角学的“余弦定理”,因此你可以在复杂度上得到精确的平衡。

A gentlemen by the name of Chris Veness has a great website at http://www.movable-type.co.uk/scripts/latlong.html that explains some the concepts you are interested in and demonstrates various programmatic implementations; this should answer your x/y conversion question as well.

一个以Chris的名字命名的先生有一个很好的网站:http://www.movabletype.co.uk/scripts/latlong.html,它解释了您感兴趣的一些概念,并展示了各种编程实现;这也可以回答你的x/y转换问题。

#2


93  

Here's the answer I found:

下面是我找到的答案:

Just to make the definition complete, in the Cartesian coordinate system:

为了使定义完整,在笛卡尔坐标系中

  • the x-axis goes through long,lat (0,0), so longitude 0 meets the equator;
  • x轴经过长,lat(0,0),所以经度0与赤道相交;
  • the y-axis goes through (0,90);
  • y轴穿过(0,90);
  • and the z-axis goes through the poles.
  • z轴穿过极点。

The conversion is:

转换:

x = R * cos(lat) * cos(lon)

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

z = R *sin(lat)

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

其中R为地球的近似半径(如6371KM)。

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

如果你的三角函数期望的是弧度(他们可能会这样做),你需要将经度和纬度转换成弧度。显然,您需要一个十进制的表示,而不是度数\分钟\秒(参见这里的转换)。

The formula for back conversion:

反向转换公式:

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

asin is of course arc sine. read about atan2 in wikipedia. Don’t forget to convert back from radians to degrees.

asin当然是反正弦。阅读*上的atan2。别忘了把弧度转换成角度。

This page gives c# code for this (note that it is very different from the formulas), and also some explanation and nice diagram of why this is correct,

这个页面给出了c#代码(注意它与公式非常不同),还有一些解释和很好的图表说明了为什么这是正确的,

#3


6  

Theory for convert GPS(WGS84) to Cartesian coordinates https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_geodetic_to_ECEF_coordinates

理论转换GPS(WGS84)笛卡尔坐标https://en.wikipedia.org/wiki/Geographic_coordinate_conversion From_geodetic_to_ECEF_coordinates

The following is what I am using:

以下是我所使用的:

  • Longitude in GPS(WGS84) and Cartesian coordinates are the same.
  • 经度GPS(WGS84)和笛卡尔坐标是相同的。
  • Latitude need be converted by WGS 84 ellipsoid parameters semi-major axis is 6378137 m, and
  • 纬度需要转换为WGS 84椭球参数半主轴为6378137 m,且。
  • Reciprocal of flattening is 298.257223563.
  • 扁平化的倒数是298.257223563。

I attached a VB code I wrote:

我附上了我写的VB代码:

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

Please notice that the h is altitude above the WGS 84 ellipsoid.

请注意,h是WGS 84椭球面的高度。

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

通常GPS会给我们高于MSL高度的H。利用重力势模型EGM96 (Lemoine et al, 1998),将MSL高度转换为WGS 84椭球面的高度h。这是通过插入geoid高度文件网格来完成的,它的空间分辨率为15弧分钟。

Or if you have some level professional GPS has Altitude H (msl,heigh above mean sea level) and UNDULATION,the relationship between the geoid and the ellipsoid (m) of the chosen datum output from internal table. you can get h = H(msl) + undulation

或者,如果你有某种级别的专业GPS,它有高度H (msl,高于平均海平面)和波动,geoid和椭球体(m)之间的关系从内部表中选择的数据输出。你可以得到h = h (msl) +波动。

To XYZ by Cartesian coordinates:

对于XYZ,笛卡尔坐标:

x = R * cos(lat) * cos(lon)

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

z = R *sin(lat)

#4


6  

Why implement something which has already been implemented and test-proven?

为什么要实现已经实现和测试的东西?

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

其中一个是NetTopologySuite,它是JTS拓扑套件的。net端口。

Specifically, you have a severe flaw in your calculation. The earth is not a perfect sphere, and the approximation of the earth's radius might not cut it for precise measurements.

具体来说,你的计算有一个严重的缺陷。地球不是一个完美的球体,地球半径的近似值可能无法精确测量。

If in some cases it's acceptable to use homebrew functions, GIS is a good example of a field in which it is much preferred to use a reliable, test-proven library.

如果在某些情况下使用homebrew函数是可以接受的,GIS就是一个很好的例子,在这个领域中,它更倾向于使用可靠的、经过测试的库。

#5


3  

If you care about getting coordinates based on an ellipsoid rather than a sphere, take a look at http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF - it gives the formulae as well as the WGS84 constants you need for the conversion.

如果您关注的是基于椭球体而不是球体的坐标,请查看http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF—它给出了转换所需的公式和WGS84常量。

The formulae there also take into account the altitude relative to the reference ellipsoid surface (useful if you are getting altitude data from a GPS device).

公式中还考虑到相对于参考椭球面的高度(如果你从GPS设备获得高度数据,这是有用的)。

#6


2  

The proj.4 software provides a command line program that can do the conversion, e.g.

项目。软件提供了一个可以进行转换的命令行程序。

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 function pj_geodetic_to_geocentric will do the conversion without having to set up a projection object first.

它还提供了一个C API。特别是,函数pj_geodetic_to_geocentric将完成转换,而不必先设置投影对象。

#7


1  

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);

#1


34  

I have recently done something similar to this using the "Haversine Formula" on WGS-84 data, which is a derivative of the "Law of Haversines" with very satisfying results.

我最近在WGS-84数据上使用“Haversine公式”做了类似的事情,这是“Haversines定律”的导数,结果非常令人满意。

Yes, WGS-84 assumes the Earth is an ellipsoid, but I believe you only get about a 0.5% average error using an approach like the "Haversine Formula", which may be an acceptable amount of error in your case. You will always have some amount of error unless you're talking about a distance of a few feet and even then there is theoretically curvature of the Earth... If you require more rigidly WGS-84 compatible approach checkout the "Vincenty Formula."

是的,WGS-84假定地球是一个椭球体,但我相信,你只会得到0.5%的平均误差,使用的方法类似于“Haversine公式”,这可能是你的情况中一个可接受的错误量。你总是会有一些误差,除非你说的是几英尺的距离,即使是在理论上地球的曲率。如果您需要更严格的WGS-84兼容方法,请签出“Vincenty公式”。

I understand where starblue is coming from, but good software engineering is often about trade offs, so it all depends on the accuracy you require for what you are doing. For example, the result calculated from "Manhattan Distance Formula" versus the result from the "Distance Formula" can be better for certain situations as it is computationally less expensive. Think "which point is closest?" scenarios where you don't need a precise distance measurement.

我理解starblue的来源,但是好的软件工程经常是关于权衡的,所以这完全取决于你所做的事情的准确性。例如,从“曼哈顿距离公式”计算得出的结果与“距离公式”的结果相比,在某些情况下会更好,因为它在计算上更便宜。想想“哪个点最接近”,你不需要精确的距离测量。

Regarding, the "Haversine Formula" it is easy to implement and is nice because it is uses "Spherical Trigonometry" instead of a "Law of Cosines" based approach which is based on two-dimensional trigonometry, therefore you get a nice balance of accuracy over complexity.

关于“Haversine公式”,它很容易实现,而且很好,因为它使用的是“球面三角学”,而不是基于二维三角学的“余弦定理”,因此你可以在复杂度上得到精确的平衡。

A gentlemen by the name of Chris Veness has a great website at http://www.movable-type.co.uk/scripts/latlong.html that explains some the concepts you are interested in and demonstrates various programmatic implementations; this should answer your x/y conversion question as well.

一个以Chris的名字命名的先生有一个很好的网站:http://www.movabletype.co.uk/scripts/latlong.html,它解释了您感兴趣的一些概念,并展示了各种编程实现;这也可以回答你的x/y转换问题。

#2


93  

Here's the answer I found:

下面是我找到的答案:

Just to make the definition complete, in the Cartesian coordinate system:

为了使定义完整,在笛卡尔坐标系中

  • the x-axis goes through long,lat (0,0), so longitude 0 meets the equator;
  • x轴经过长,lat(0,0),所以经度0与赤道相交;
  • the y-axis goes through (0,90);
  • y轴穿过(0,90);
  • and the z-axis goes through the poles.
  • z轴穿过极点。

The conversion is:

转换:

x = R * cos(lat) * cos(lon)

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

z = R *sin(lat)

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

其中R为地球的近似半径(如6371KM)。

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

如果你的三角函数期望的是弧度(他们可能会这样做),你需要将经度和纬度转换成弧度。显然,您需要一个十进制的表示,而不是度数\分钟\秒(参见这里的转换)。

The formula for back conversion:

反向转换公式:

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

asin is of course arc sine. read about atan2 in wikipedia. Don’t forget to convert back from radians to degrees.

asin当然是反正弦。阅读*上的atan2。别忘了把弧度转换成角度。

This page gives c# code for this (note that it is very different from the formulas), and also some explanation and nice diagram of why this is correct,

这个页面给出了c#代码(注意它与公式非常不同),还有一些解释和很好的图表说明了为什么这是正确的,

#3


6  

Theory for convert GPS(WGS84) to Cartesian coordinates https://en.wikipedia.org/wiki/Geographic_coordinate_conversion#From_geodetic_to_ECEF_coordinates

理论转换GPS(WGS84)笛卡尔坐标https://en.wikipedia.org/wiki/Geographic_coordinate_conversion From_geodetic_to_ECEF_coordinates

The following is what I am using:

以下是我所使用的:

  • Longitude in GPS(WGS84) and Cartesian coordinates are the same.
  • 经度GPS(WGS84)和笛卡尔坐标是相同的。
  • Latitude need be converted by WGS 84 ellipsoid parameters semi-major axis is 6378137 m, and
  • 纬度需要转换为WGS 84椭球参数半主轴为6378137 m,且。
  • Reciprocal of flattening is 298.257223563.
  • 扁平化的倒数是298.257223563。

I attached a VB code I wrote:

我附上了我写的VB代码:

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

Please notice that the h is altitude above the WGS 84 ellipsoid.

请注意,h是WGS 84椭球面的高度。

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

通常GPS会给我们高于MSL高度的H。利用重力势模型EGM96 (Lemoine et al, 1998),将MSL高度转换为WGS 84椭球面的高度h。这是通过插入geoid高度文件网格来完成的,它的空间分辨率为15弧分钟。

Or if you have some level professional GPS has Altitude H (msl,heigh above mean sea level) and UNDULATION,the relationship between the geoid and the ellipsoid (m) of the chosen datum output from internal table. you can get h = H(msl) + undulation

或者,如果你有某种级别的专业GPS,它有高度H (msl,高于平均海平面)和波动,geoid和椭球体(m)之间的关系从内部表中选择的数据输出。你可以得到h = h (msl) +波动。

To XYZ by Cartesian coordinates:

对于XYZ,笛卡尔坐标:

x = R * cos(lat) * cos(lon)

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

z = R *sin(lat)

#4


6  

Why implement something which has already been implemented and test-proven?

为什么要实现已经实现和测试的东西?

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

其中一个是NetTopologySuite,它是JTS拓扑套件的。net端口。

Specifically, you have a severe flaw in your calculation. The earth is not a perfect sphere, and the approximation of the earth's radius might not cut it for precise measurements.

具体来说,你的计算有一个严重的缺陷。地球不是一个完美的球体,地球半径的近似值可能无法精确测量。

If in some cases it's acceptable to use homebrew functions, GIS is a good example of a field in which it is much preferred to use a reliable, test-proven library.

如果在某些情况下使用homebrew函数是可以接受的,GIS就是一个很好的例子,在这个领域中,它更倾向于使用可靠的、经过测试的库。

#5


3  

If you care about getting coordinates based on an ellipsoid rather than a sphere, take a look at http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF - it gives the formulae as well as the WGS84 constants you need for the conversion.

如果您关注的是基于椭球体而不是球体的坐标,请查看http://en.wikipedia.org/wiki/Geodetic_system#From_geodetic_to_ECEF—它给出了转换所需的公式和WGS84常量。

The formulae there also take into account the altitude relative to the reference ellipsoid surface (useful if you are getting altitude data from a GPS device).

公式中还考虑到相对于参考椭球面的高度(如果你从GPS设备获得高度数据,这是有用的)。

#6


2  

The proj.4 software provides a command line program that can do the conversion, e.g.

项目。软件提供了一个可以进行转换的命令行程序。

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 function pj_geodetic_to_geocentric will do the conversion without having to set up a projection object first.

它还提供了一个C API。特别是,函数pj_geodetic_to_geocentric将完成转换,而不必先设置投影对象。

#7


1  

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);