给定当前点、距离和方位角，获取纬度/经度

Question

给定当前点、距离和方位角，获取纬度/经度

114

给定一个经纬度坐标点，以及距离（以公里为单位）和方位角（以弧度表示的角度），我想要计算新的经纬度坐标。这个网站一次又一次地出现，但是我就是无法使公式对我起作用。

上述链接中提供的公式如下：

lat2 = asin(sin(lat1)*cos(d/R) + cos(lat1)*sin(d/R)*cos(θ))

lon2 = lon1 + atan2(sin(θ)*sin(d/R)*cos(lat1), cos(d/R)−sin(lat1)*sin(lat2))

上述公式适用于MSExcel，其中-

asin          = arc sin()   
d             = distance (in any unit)   
R             = Radius of the earth (in the same unit as above)  
and hence d/r = is the angular distance (in radians)  
atan2(a,b)    = arc tan(b/a)  
θ is the bearing (in radians, clockwise from north);

这是我在Python中得到的代码。

import math

R = 6378.1 #Radius of the Earth
brng = 1.57 #Bearing is 90 degrees converted to radians.
d = 15 #Distance in km

#lat2  52.20444 - the lat result I'm hoping for
#lon2  0.36056 - the long result I'm hoping for.

lat1 = 52.20472 * (math.pi * 180) #Current lat point converted to radians
lon1 = 0.14056 * (math.pi * 180) #Current long point converted to radians

lat2 = math.asin( math.sin(lat1)*math.cos(d/R) +
             math.cos(lat1)*math.sin(d/R)*math.cos(brng))

lon2 = lon1 + math.atan2(math.sin(brng)*math.sin(d/R)*math.cos(lat1),
                     math.cos(d/R)-math.sin(lat1)*math.sin(lat2))

print(lat2)
print(lon2)

我得到

lat2 = 0.472492248844 
lon2 = 79.4821662373

- David M

1

@GWW，我得到了一个不合理的答案。这个不合理的原因是因为我没有将答案转换回度数。代码已更改并包含在原帖中作为编辑。 - David M

6

你应该将编辑后的内容直接提交为答案，并接受该答案，以便更清楚地表明您解决了自己的问题。否则，Stack Overflow会对您未解决的问题进行惩罚，这会使未来的用户更有可能不愿回答您的问题。 - Cerin

如果您使用numpy对象，您将获得更好的精度和结果。 - Mike Pennington

应该是“lat1 = 52.20472 * (math.pi */180)”吗？ - Tushar Seth

如果方位角为90度，为什么纬度会改变呢？这不只是沿着经线移动吗？ - danjampro

可能可以忽略不计，但如果您需要绝对精度，您可以将地球半径设置为3个小数位，R = 6378.137 - 来自Google Maps的几何库文档“默认半径是地球半径6378137米”- https://developers.google.com/maps/documentation/javascript/reference/geometry - Esteban Ortega

15个回答

41

这个geopy库支持此操作：

import geopy
from geopy.distance import VincentyDistance

# given: lat1, lon1, b = bearing in degrees, d = distance in kilometers

origin = geopy.Point(lat1, lon1)
destination = VincentyDistance(kilometers=d).destination(origin, b)

lat2, lon2 = destination.latitude, destination.longitude

通过https://dev59.com/b1LTa4cB1Zd3GeqPeuoI#4531227找到。

- elliot42

这个库存在一些距离问题需要解决：https://github.com/geopy/geopy/pull/144 - jpoehnelt

2

请注意，自v2.0.0以来API已更改。请改用geopy.distance.geodesic： https://dev59.com/hVIG5IYBdhLWcg3wzEzT#62866744 - danjampro

33

这个问题被称为大地测量学中的正问题。

这是一个非常普遍的问题，也是不断引起困惑的原因。原因是大多数人都在寻找简单明了的答案。但事实上并没有简单的答案，因为大多数提出这个问题的人都没有提供足够的信息，而这是因为他们并不知道：

地球不是一个完美的球体，因为它的极点会被压扁/压缩
由于 (1)，地球的半径R不是恒定的。参见此处。
地球不是完全光滑的（存在高度差异）等等。
由于板块运动，一个地理点的纬度/经度位置可能每年变化几毫米(至少)。

因此，在各种几何模型中使用许多不同的假设，这些假设根据您所需的精度以不同的方式应用。因此，要回答这个问题，您需要考虑您想要的精度是什么。

一些例子：

我只是想找到最近几公里的大约位置，针对在 0-70 deg N|S 的纬度间的小距离 (< 100 km)。（地球是平坦模型。）
我需要一个可行的答案，可以适用于全球的任何地方，但只精确到几米。
我需要一个超级精确的定位，这个定位可以有效准确到纳米[nm]的原子尺度。
我需要快速简便易行且计算量不大的答案。

所以你可以在算法选择上有很多选择。此外，每种编程语言都有自己的实现或“包”，加上模型数量和模型开发人员的特定需求不计其数。在这里，为了实际目的，最好忽略除javascript之外的任何其他语言，因为它非常接近伪代码的性质。因此，它可以轻松地转换为任何其他语言，只需进行最少的更改。

然后，主要的模型是：

欧几里得/平地模型：适用于距离非常短的情况下（约10公里以下）
球面模型：适用于经度距离很大，但纬度差异很小的情况。流行的模型是：
- Haversine：在[km]尺度上具有米级精度，代码非常简单。
椭球模型：在任何纬度/经度和距离上最精确，但仍然是数字逼近，取决于您需要什么精度。一些流行的模型包括：
- Lambert：在1000多km的范围内具有约10米的精度。
- Paul D.Thomas：Andoyer-Lambert逼近
- Vincenty：毫米级精度和计算效率
- Kerney：纳米级精度

参考资料：

- not2qubit

16

可能有点晚回答，但在测试了其他答案后，发现它们不正确。这是我们系统中使用的PHP代码。可以在所有方向上工作。

PHP代码：

lat1 = 起点纬度（以度为单位）

long1 = 起点经度（以度为单位）

d = 距离（以公里为单位）

angle = 方位角（以度为单位）

function get_gps_distance($lat1,$long1,$d,$angle)
{
    # Earth Radious in KM
    $R = 6378.14;

    # Degree to Radian
    $latitude1 = $lat1 * (M_PI/180);
    $longitude1 = $long1 * (M_PI/180);
    $brng = $angle * (M_PI/180);

    $latitude2 = asin(sin($latitude1)*cos($d/$R) + cos($latitude1)*sin($d/$R)*cos($brng));
    $longitude2 = $longitude1 + atan2(sin($brng)*sin($d/$R)*cos($latitude1),cos($d/$R)-sin($latitude1)*sin($latitude2));

    # back to degrees
    $latitude2 = $latitude2 * (180/M_PI);
    $longitude2 = $longitude2 * (180/M_PI);

    # 6 decimal for Leaflet and other system compatibility
   $lat2 = round ($latitude2,6);
   $long2 = round ($longitude2,6);

   // Push in array and get back
   $tab[0] = $lat2;
   $tab[1] = $long2;
   return $tab;
 }

- Peter

2

看起来不错，但我认为请求者希望得到Python的东西。错了吗？ - user2360915

可能更好地命名为get_gps_coord或类似的名称。你不是在获取距离，而是将其提供给函数。但感谢您提供这个，这正是我要找的。许多搜索结果返回计算坐标之间距离的内容（误报）。谢谢！ - Dev Null

太棒了！感谢您的贡献！ - jkamizato

3

“6,378.14公里”似乎是地球的最大半径。平均值约为“6,371.0公里”，这可能会使计算更加准确。 - JVE999

感谢你为我节省了一些时间。 - Robert Moskal

大多数答案无法正常工作是因为最后一项出现了小错误；它们没有使用旧纬度和新纬度，而是两次使用已经更新的纬度。 - afp

13

我将Brad的答案移植到纯JS答案，没有使用Bing maps。

https://jsfiddle.net/kodisha/8a3hcjtd/

    // ----------------------------------------
    // Calculate new Lat/Lng from original points
    // on a distance and bearing (angle)
    // ----------------------------------------
    let llFromDistance = function(latitude, longitude, distance, bearing) {
      // taken from: https://dev59.com/42w05IYBdhLWcg3weBru#46410871 
      // distance in KM, bearing in degrees
    
      const R = 6378.1; // Radius of the Earth
      const brng = bearing * Math.PI / 180; // Convert bearing to radian
      let lat = latitude * Math.PI / 180; // Current coords to radians
      let lon = longitude * Math.PI / 180;
    
      // Do the math magic
      lat = Math.asin(Math.sin(lat) * Math.cos(distance / R) + Math.cos(lat) * Math.sin(distance / R) * Math.cos(brng));
      lon += Math.atan2(Math.sin(brng) * Math.sin(distance / R) * Math.cos(lat), Math.cos(distance / R) - Math.sin(lat) * Math.sin(lat));
    
      // Coords back to degrees and return
      return [(lat * 180 / Math.PI), (lon * 180 / Math.PI)];
    
    }
    
    let pointsOnMapCircle = function(latitude, longitude, distance, numPoints) {
      const points = [];
      for (let i = 0; i <= numPoints - 1; i++) {
        const bearing = Math.round((360 / numPoints) * i);
        console.log(bearing, i);
        const newPoints = llFromDistance(latitude, longitude, distance, bearing);
        points.push(newPoints);
      }
      return points;
    }
    
    const points = pointsOnMapCircle(41.890242042122836, 12.492358982563019, 0.2, 8);
    let geoJSON = {
      "type": "FeatureCollection",
      "features": []
    };
    points.forEach((p) => {
      geoJSON.features.push({
        "type": "Feature",
        "properties": {},
        "geometry": {
          "type": "Point",
          "coordinates": [
            p[1],
            p[0]
          ]
        }
      });
    });
    
    document.getElementById('res').innerHTML = JSON.stringify(geoJSON, true, 2);

此外，我添加了 geoJSON 导出功能，因此您可以将导出的 geoJSON 简单地粘贴到 http://geojson.io/#map=17/41.89017/12.49171 来立即查看结果。

结果：

- kodisha

GeoJSON 的地图对我在地图上定位位置非常有帮助。 - cloudscomputes

谢谢@kodisha，你的fiddle帮了我很多！ - Douglas Schmidt

和我在上一个答案中的评论一样，我认为经度计算的最后一部分可能是错误的，因为在计算 lon 之前已经更新了变量 lat，即术语 Math.sin(lat) * Math.sin(lat) 实际上并没有分别使用旧的和新的纬度。 - afp

可以用十进制度数表示距离吗？ - Cristián Vargas Acevedo

8

使用geopy的快速方法

from geopy import distance
#distance.distance(unit=15).destination((lat,lon),bering) 
#Exemples
distance.distance(nautical=15).destination((-24,-42),90) 
distance.distance(miles=15).destination((-24,-42),90)
distance.distance(kilometers=15).destination((-24,-42),90)

- Plinio Bueno Andrade Silva

3

不说明你所用的计算方法，答案基本上是毫无意义的。 - not2qubit

2

@not2qubit 不管 @plinio-bueno-andrade-silva 是否知道，geopy.distance.distance 目前使用的是测地线距离。geopy 更具体地说，默认使用的椭球模型是 WGS-84 椭球体，“这是全球最准确的”。 - hlongmore

3

我将Python移植到了Javascript。它返回一个Bing Maps Location对象，你可以自由更改它。

getLocationXDistanceFromLocation: function(latitude, longitude, distance, bearing) {
    // distance in KM, bearing in degrees

    var R = 6378.1,                         // Radius of the Earth
        brng = Math.radians(bearing)       // Convert bearing to radian
        lat = Math.radians(latitude),       // Current coords to radians
        lon = Math.radians(longitude);

    // Do the math magic
    lat = Math.asin(Math.sin(lat) * Math.cos(distance / R) + Math.cos(lat) * Math.sin(distance / R) * Math.cos(brng));
    lon += Math.atan2(Math.sin(brng) * Math.sin(distance / R) * Math.cos(lat), Math.cos(distance/R)-Math.sin(lat)*Math.sin(lat));

    // Coords back to degrees and return
    return new Microsoft.Maps.Location(Math.degrees(lat), Math.degrees(lon));

},

- Brad Mathews

请发布包含运行所需内容的可用代码，例如此代码似乎依赖于 Microsoft.Maps。在哪里可以找到或者如何安装它？ - not2qubit

如果你的程序使用Bing地图，那么你只能使用Bing（Microsoft）地图。只需获取Math.degrees(lat)和Math.degrees(lon)值，并根据你的应用程序需要进行处理即可。 - Brad Mathews

3

lon1和lat1是以度为单位的坐标

brng = 弧度制的方位角

d = 距离，单位为公里

R = 地球半径，单位为公里

lat2 = math.degrees((d/R) * math.cos(brng)) + lat1
long2 = math.degrees((d/(R*math.sin(math.radians(lat2)))) * math.sin(brng)) + long1

我用PHP实现了你的算法和我的算法，并进行了基准测试。这个版本运行时间大约只有原来的50%。生成的结果完全相同，因此它在数学上是等效的。

我没有测试上面的Python代码，所以可能会有语法错误。

- Chris

无法工作。从南到北，结果是正确的，但在“东西”方向上是错误的。 - Peter

2

感谢@kodisha，这里提供了Swift版本，但是对于地球半径的计算进行了改进和更精确的计算：

extension CLLocationCoordinate2D {
  
  func earthRadius() -> CLLocationDistance {
    let earthRadiusInMetersAtSeaLevel = 6378137.0
    let earthRadiusInMetersAtPole = 6356752.314
    
    let r1 = earthRadiusInMetersAtSeaLevel
    let r2 = earthRadiusInMetersAtPole
    let beta = latitude
    
    let earthRadiuseAtGivenLatitude = (
      ( pow(pow(r1, 2) * cos(beta), 2) + pow(pow(r2, 2) * sin(beta), 2) ) /
        ( pow(r1 * cos(beta), 2) + pow(r2 * sin(beta), 2) )
    )
    .squareRoot()
    
    return earthRadiuseAtGivenLatitude
  }
  
  func locationByAdding(
    distance: CLLocationDistance,
    bearing: CLLocationDegrees
  ) -> CLLocationCoordinate2D {
    let latitude = self.latitude
    let longitude = self.longitude
    
    let earthRadiusInMeters = self.earthRadius()
    let brng = bearing.degreesToRadians
    var lat = latitude.degreesToRadians
    var lon = longitude.degreesToRadians
    
    lat = asin(
      sin(lat) * cos(distance / earthRadiusInMeters) +
        cos(lat) * sin(distance / earthRadiusInMeters) * cos(brng)
    )
    lon += atan2(
      sin(brng) * sin(distance / earthRadiusInMeters) * cos(lat),
      cos(distance / earthRadiusInMeters) - sin(lat) * sin(lat)
    )
    
    let newCoordinate = CLLocationCoordinate2D(
      latitude: lat.radiansToDegrees,
      longitude: lon.radiansToDegrees
    )
    
    return newCoordinate
  }
}

extension FloatingPoint {
  var degreesToRadians: Self { self * .pi / 180 }
  var radiansToDegrees: Self { self * 180 / .pi }
}

- hbk

1

我认为经度计算的最后一部分可能是错误的，因为在计算lon之前变量lat已经被更新，即术语sin(lat) * sin(lat)实际上并没有分别使用旧和新的纬度。 - afp

@afp 我同意，那肯定是个错误，但很容易修复。 - Eric

1

如果有人需要，我将@David M的答案移植到了Java......我得到了一个稍微不同的结果为52.20462299620793，0.360433887489931。

    double R = 6378.1;  //Radius of the Earth
    double brng = 1.57;  //Bearing is 90 degrees converted to radians.
    double d = 15;  //Distance in km

    double lat2 = 52.20444; // - the lat result I'm hoping for
    double lon2 = 0.36056; // - the long result I'm hoping for.

    double lat1 = Math.toRadians(52.20472); //Current lat point converted to radians
    double lon1 = Math.toRadians(0.14056); //Current long point converted to radians

    lat2 = Math.asin( Math.sin(lat1)*Math.cos(d/R) +
            Math.cos(lat1)*Math.sin(d/R)*Math.cos(brng));

    lon2 = lon1 + Math.atan2(Math.sin(brng)*Math.sin(d/R)*Math.cos(lat1),
            Math.cos(d/R)-Math.sin(lat1)*Math.sin(lat2));

    lat2 = Math.toDegrees(lat2);
    lon2 = Math.toDegrees(lon2);

    System.out.println(lat2 + ", " + lon2);

- Alessandro

可能这是最正确的答案，因为在计算 lon2 表达式的最后一项 Math.sin(lat1)*Math.sin(lat2) 时，它分别正确使用了旧和新纬度。因此结果略有不同。 - afp

网页内容由stack overflow 提供, 点击上面的

可以查看英文原文，
原文链接

- David M · Accepted Answer

需要将弧度制的答案转换为角度制。以下是可行的代码：

from math import asin, atan2, cos, degrees, radians, sin

def get_point_at_distance(lat1, lon1, d, bearing, R=6371):
    """
    lat: initial latitude, in degrees
    lon: initial longitude, in degrees
    d: target distance from initial
    bearing: (true) heading in degrees
    R: optional radius of sphere, defaults to mean radius of earth

    Returns new lat/lon coordinate {d}km from initial, in degrees
    """
    lat1 = radians(lat1)
    lon1 = radians(lon1)
    a = radians(bearing)
    lat2 = asin(sin(lat1) * cos(d/R) + cos(lat1) * sin(d/R) * cos(a))
    lon2 = lon1 + atan2(
        sin(a) * sin(d/R) * cos(lat1),
        cos(d/R) - sin(lat1) * sin(lat2)
    )
    return (degrees(lat2), degrees(lon2),)

lat = 52.20472
lon = 0.14056
distance = 15
bearing = 90
lat2, lon2 = get_point_at_distance(lat, lon, distance, bearing)

# lat2  52.20444 - the lat result I'm hoping for
# lon2  0.36056 - the long result I'm hoping for.

print(lat2, lon2)
# prints "52.20451523755824 0.36067845713550956"