如何计算给定lat / lng位置的边界框？

15 回答

• 10

  我建议将地球表面局部近似为一个球体，其半径由给定纬度的WGS84椭球给出 . 我怀疑latMin和latMax的精确计算需要椭圆函数，并且不会产生明显的精度增加（WGS84本身就是近似值） .

  我的实现如下（它是用Python编写的;我还没有测试过）：

  # degrees to radians
  def deg2rad(degrees):
      return math.pi*degrees/180.0
  # radians to degrees
  def rad2deg(radians):
      return 180.0*radians/math.pi
  
  # Semi-axes of WGS-84 geoidal reference
  WGS84_a = 6378137.0  # Major semiaxis [m]
  WGS84_b = 6356752.3  # Minor semiaxis [m]
  
  # Earth radius at a given latitude, according to the WGS-84 ellipsoid [m]
  def WGS84EarthRadius(lat):
      # http://en.wikipedia.org/wiki/Earth_radius
      An = WGS84_a*WGS84_a * math.cos(lat)
      Bn = WGS84_b*WGS84_b * math.sin(lat)
      Ad = WGS84_a * math.cos(lat)
      Bd = WGS84_b * math.sin(lat)
      return math.sqrt( (An*An + Bn*Bn)/(Ad*Ad + Bd*Bd) )
  
  # Bounding box surrounding the point at given coordinates,
  # assuming local approximation of Earth surface as a sphere
  # of radius given by WGS84
  def boundingBox(latitudeInDegrees, longitudeInDegrees, halfSideInKm):
      lat = deg2rad(latitudeInDegrees)
      lon = deg2rad(longitudeInDegrees)
      halfSide = 1000*halfSideInKm
  
      # Radius of Earth at given latitude
      radius = WGS84EarthRadius(lat)
      # Radius of the parallel at given latitude
      pradius = radius*math.cos(lat)
  
      latMin = lat - halfSide/radius
      latMax = lat + halfSide/radius
      lonMin = lon - halfSide/pradius
      lonMax = lon + halfSide/pradius
  
      return (rad2deg(latMin), rad2deg(lonMin), rad2deg(latMax), rad2deg(lonMax))
  

  编辑：以下代码将（度，素数，秒）转换为度数的度数分数，反之亦然（未经测试）：

  def dps2deg(degrees, primes, seconds):
      return degrees + primes/60.0 + seconds/3600.0
  
  def deg2dps(degrees):
      intdeg = math.floor(degrees)
      primes = (degrees - intdeg)*60.0
      intpri = math.floor(primes)
      seconds = (primes - intpri)*60.0
      intsec = round(seconds)
      return (int(intdeg), int(intpri), int(intsec))
  

  回复于 2024-05-04T03:02:28+08:00
• 6

  我写了一篇关于找到边界坐标的文章：

  http://JanMatuschek.de/LatitudeLongitudeBoundingCoordinates

  本文解释了这些公式，并提供了Java实现 . （这也说明了为什么费德里科的最小/最大经度公式不准确 . ）

  回复于 2024-05-04T03:02:28+08:00
• 0

  在这里，我已经将Federico A. Ramponi对C＃的回答转换为任何感兴趣的人：

  public class MapPoint
  {
      public double Longitude { get; set; } // In Degrees
      public double Latitude { get; set; } // In Degrees
  }
  
  public class BoundingBox
  {
      public MapPoint MinPoint { get; set; }
      public MapPoint MaxPoint { get; set; }
  }        
  
  // Semi-axes of WGS-84 geoidal reference
  private const double WGS84_a = 6378137.0; // Major semiaxis [m]
  private const double WGS84_b = 6356752.3; // Minor semiaxis [m]
  
  // 'halfSideInKm' is the half length of the bounding box you want in kilometers.
  public static BoundingBox GetBoundingBox(MapPoint point, double halfSideInKm)
  {            
      // Bounding box surrounding the point at given coordinates,
      // assuming local approximation of Earth surface as a sphere
      // of radius given by WGS84
      var lat = Deg2rad(point.Latitude);
      var lon = Deg2rad(point.Longitude);
      var halfSide = 1000 * halfSideInKm;
  
      // Radius of Earth at given latitude
      var radius = WGS84EarthRadius(lat);
      // Radius of the parallel at given latitude
      var pradius = radius * Math.Cos(lat);
  
      var latMin = lat - halfSide / radius;
      var latMax = lat + halfSide / radius;
      var lonMin = lon - halfSide / pradius;
      var lonMax = lon + halfSide / pradius;
  
      return new BoundingBox { 
          MinPoint = new MapPoint { Latitude = Rad2deg(latMin), Longitude = Rad2deg(lonMin) },
          MaxPoint = new MapPoint { Latitude = Rad2deg(latMax), Longitude = Rad2deg(lonMax) }
      };            
  }
  
  // degrees to radians
  private static double Deg2rad(double degrees)
  {
      return Math.PI * degrees / 180.0;
  }
  
  // radians to degrees
  private static double Rad2deg(double radians)
  {
      return 180.0 * radians / Math.PI;
  }
  
  // Earth radius at a given latitude, according to the WGS-84 ellipsoid [m]
  private static double WGS84EarthRadius(double lat)
  {
      // http://en.wikipedia.org/wiki/Earth_radius
      var An = WGS84_a * WGS84_a * Math.Cos(lat);
      var Bn = WGS84_b * WGS84_b * Math.Sin(lat);
      var Ad = WGS84_a * Math.Cos(lat);
      var Bd = WGS84_b * Math.Sin(lat);
      return Math.Sqrt((An*An + Bn*Bn) / (Ad*Ad + Bd*Bd));
  }
  

  回复于 2024-05-04T03:02:28+08:00
• 0

  我编写了一个JavaScript函数，它返回一个方形边界框的四个坐标，给定一个距离和一对坐标：

  'use strict';
  
  /**
   * @param {number} distance - distance (km) from the point represented by centerPoint
   * @param {array} centerPoint - two-dimensional array containing center coords [latitude, longitude]
   * @description
   *   Computes the bounding coordinates of all points on the surface of a sphere
   *   that has a great circle distance to the point represented by the centerPoint
   *   argument that is less or equal to the distance argument.
   *   Technique from: Jan Matuschek <http://JanMatuschek.de/LatitudeLongitudeBoundingCoordinates>
   * @author Alex Salisbury
  */
  
  getBoundingBox = function (centerPoint, distance) {
    var MIN_LAT, MAX_LAT, MIN_LON, MAX_LON, R, radDist, degLat, degLon, radLat, radLon, minLat, maxLat, minLon, maxLon, deltaLon;
    if (distance < 0) {
      return 'Illegal arguments';
    }
    // helper functions (degrees<–>radians)
    Number.prototype.degToRad = function () {
      return this * (Math.PI / 180);
    };
    Number.prototype.radToDeg = function () {
      return (180 * this) / Math.PI;
    };
    // coordinate limits
    MIN_LAT = (-90).degToRad();
    MAX_LAT = (90).degToRad();
    MIN_LON = (-180).degToRad();
    MAX_LON = (180).degToRad();
    // Earth's radius (km)
    R = 6378.1;
    // angular distance in radians on a great circle
    radDist = distance / R;
    // center point coordinates (deg)
    degLat = centerPoint[0];
    degLon = centerPoint[1];
    // center point coordinates (rad)
    radLat = degLat.degToRad();
    radLon = degLon.degToRad();
    // minimum and maximum latitudes for given distance
    minLat = radLat - radDist;
    maxLat = radLat + radDist;
    // minimum and maximum longitudes for given distance
    minLon = void 0;
    maxLon = void 0;
    // define deltaLon to help determine min and max longitudes
    deltaLon = Math.asin(Math.sin(radDist) / Math.cos(radLat));
    if (minLat > MIN_LAT && maxLat < MAX_LAT) {
      minLon = radLon - deltaLon;
      maxLon = radLon + deltaLon;
      if (minLon < MIN_LON) {
        minLon = minLon + 2 * Math.PI;
      }
      if (maxLon > MAX_LON) {
        maxLon = maxLon - 2 * Math.PI;
      }
    }
    // a pole is within the given distance
    else {
      minLat = Math.max(minLat, MIN_LAT);
      maxLat = Math.min(maxLat, MAX_LAT);
      minLon = MIN_LON;
      maxLon = MAX_LON;
    }
    return [
      minLon.radToDeg(),
      minLat.radToDeg(),
      maxLon.radToDeg(),
      maxLat.radToDeg()
    ];
  };
  

  回复于 2024-05-04T03:02:28+08:00
• 4

  你正在寻找一个椭球公式 .

  我发现开始编码的最佳位置是基于CPAN的Geo :: Ellipsoid库 . 它为您提供了创建测试的基线，并将结果与结果进行比较 . 我在之前的雇主那里用它作为PHP类似库的基础 .

  Geo::Ellipsoid

  看一下 location 方法 . 叫它两次，你有你的bbox .

  你没有发布你正在使用的语言 . 可能已经有一个地理编码库可供您使用 .

  哦，如果你现在还没弄明白，谷歌 Map 使用的是WGS84椭圆体 .

  回复于 2024-05-04T03:02:28+08:00
• 56

  我改编了一个PHP脚本，我发现这样做 . 你可以用它来找到一个点周围的一个角落（例如，20公里外） . 我的具体示例适用于Google Maps API：

  http://www.richardpeacock.com/blog/2011/11/draw-box-around-coordinate-google-maps-based-miles-or-kilometers

  回复于 2024-05-04T03:02:28+08:00
• 0

  由于我需要非常粗略的估计，所以为了在弹性搜索查询中过滤掉一些不必要的文档，我采用了以下公式：

  Min.lat = Given.Lat - (0.009 x N)
  Max.lat = Given.Lat + (0.009 x N)
  Min.lon = Given.lon - (0.009 x N)
  Max.lon = Given.lon + (0.009 x N)
  

  从给定位置需要N = kms . 对于你的情况N = 10

  不准确但方便 .

  回复于 2024-05-04T03:02:28+08:00
• 1

  这是一个使用javascript的简单实现，它基于纬度到kms的转换 1 degree latitude ~ 111.2 km .

  我正在计算从给定纬度和经度宽度为10km的 Map 边界 .

  function getBoundsFromLatLng(lat, lng){
       var lat_change = 10/111.2;
       var lon_change = Math.abs(Math.cos(lat*(Math.PI/180)));
       var bounds = { 
           lat_min : lat - lat_change,
           lon_min : lng - lon_change,
           lat_max : lat + lat_change,
           lon_max : lng + lon_change
       };
       return bounds;
  }
  

  回复于 2024-05-04T03:02:28+08:00
• 0

  @Jan Philip Matuschek的插图非常好的解释 . （请把他的答案投票，而不是这个;我加上这个，因为我花了一点时间来理解原来的答案）

  寻找最近邻居的优化的边界框技术需要导出距离d处的点P的最小和最大纬度，经度对 . 所有落在这些点之外的点肯定距离点大于d . 这里需要注意的一点是Jan Philip Matuschek解释中突出显示的交叉纬度的计算 . 交叉纬度不在点P的纬度处，而是略微偏离它 . 这是确定距离d的点P的正确最小和最大边界经度时经常遗漏但重要的部分 . 这在验证中也是有用的 .

  P的（纬度，经度）与（纬度，经度）之间的半径距离等于距离d .

  Python gist here https://gist.github.com/alexcpn/f95ae83a7ee0293a5225

  

  回复于 2024-05-04T03:02:28+08:00
• 50

  我正在研究边界框问题，作为查找静态LAT，LONG点的SrcRad半径内的所有点的副问题 . 已经有很多计算使用

  maxLon = $lon + rad2deg($rad/$R/cos(deg2rad($lat)));
  minLon = $lon - rad2deg($rad/$R/cos(deg2rad($lat)));
  

  计算经度界限，但我发现这并没有给出所需的所有答案 . 因为你真正想做的是

  (SrcRad/RadEarth)/cos(deg2rad(lat))
  

  我知道，我知道答案应该是一样的，但我发现事实并非如此 . 它似乎通过不确定我正在做（SRCrad / RadEarth）首先然后除去Cos部分我留下了一些位置点 .

  获得所有边界框点后，如果您有一个计算给定纬度的点到点距离的函数，那么很容易只获得距离固定点一定距离半径的那些点 . 这就是我做的 . 我知道它需要一些额外的步骤，但它帮助了我

  -- GLOBAL Constants
  gc_pi CONSTANT REAL := 3.14159265359;  -- Pi
  
  -- Conversion Factor Constants
  gc_rad_to_degs          CONSTANT NUMBER := 180/gc_pi; -- Conversion for Radians to Degrees 180/pi
  gc_deg_to_rads          CONSTANT NUMBER := gc_pi/180; --Conversion of Degrees to Radians
  
  lv_stat_lat    -- The static latitude point that I am searching from 
  lv_stat_long   -- The static longitude point that I am searching from 
  
  -- Angular radius ratio in radians
  lv_ang_radius := lv_search_radius / lv_earth_radius;
  lv_bb_maxlat := lv_stat_lat + (gc_rad_to_deg * lv_ang_radius);
  lv_bb_minlat := lv_stat_lat - (gc_rad_to_deg * lv_ang_radius);
  
  --Here's the tricky part, accounting for the Longitude getting smaller as we move up the latitiude scale
  -- I seperated the parts of the equation to make it easier to debug and understand
  -- I may not be a smart man but I know what the right answer is... :-)
  
  lv_int_calc := gc_deg_to_rads * lv_stat_lat;
  lv_int_calc := COS(lv_int_calc);
  lv_int_calc := lv_ang_radius/lv_int_calc;
  lv_int_calc := gc_rad_to_degs*lv_int_calc;
  
  lv_bb_maxlong := lv_stat_long + lv_int_calc;
  lv_bb_minlong := lv_stat_long - lv_int_calc;
  
  -- Now select the values from your location datatable 
  SELECT *  FROM (
  SELECT cityaliasname, city, state, zipcode, latitude, longitude, 
  -- The actual distance in miles
  spherecos_pnttopntdist(lv_stat_lat, lv_stat_long, latitude, longitude, 'M') as miles_dist    
  FROM Location_Table 
  WHERE latitude between lv_bb_minlat AND lv_bb_maxlat
  AND   longitude between lv_bb_minlong and lv_bb_maxlong)
  WHERE miles_dist <= lv_limit_distance_miles
  order by miles_dist
  ;
  

  回复于 2024-05-04T03:02:28+08:00
• 2

  这很简单，只需访问panoramio网站，然后从panoramio网站打开世界 Map . 然后前往所需的纬度和经度的指定位置 .

  然后，您在地址栏中找到了纬度和经度，例如在此地址中 .

  http://www.panoramio.com/map#lt=32.739485&ln=70.491211&z=9&k=1&a=1&tab=1&pl=all

  lt = 32.739485 =>纬度ln = 70.491211 =>经度

  这个Panoramio JavaScript API小部件在纬度/经度对周围创建一个边界框，然后在这些边界内返回所有照片 .

  另一种类型的Panoramio JavaScript API小部件，您还可以使用example and code is here更改背景颜色 .

  它没有在作曲情绪中表现出来 . 它在出版后表现出来 .

  <div dir="ltr" style="text-align: center;" trbidi="on">
  <script src="https://ssl.panoramio.com/wapi/wapi.js?v=1&amp;hl=en"></script>
  <div id="wapiblock" style="float: right; margin: 10px 15px"></div>
  <script type="text/javascript">
  var myRequest = {
    'tag': 'kahna',
    'rect': {'sw': {'lat': -30, 'lng': 10.5}, 'ne': {'lat': 50.5, 'lng': 30}}
  };
    var myOptions = {
    'width': 300,
    'height': 200
  };
  var wapiblock = document.getElementById('wapiblock');
  var photo_widget = new panoramio.PhotoWidget('wapiblock', myRequest, myOptions);
  photo_widget.setPosition(0);
  </script>
  </div>
  

  回复于 2024-05-04T03:02:28+08:00
• 4

  如果有人感兴趣，我在这里转换了Federico A. Ramponi对PHP的回答：

  <?php
  # deg2rad and rad2deg are already within PHP
  
  # Semi-axes of WGS-84 geoidal reference
  $WGS84_a = 6378137.0;  # Major semiaxis [m]
  $WGS84_b = 6356752.3;  # Minor semiaxis [m]
  
  # Earth radius at a given latitude, according to the WGS-84 ellipsoid [m]
  function WGS84EarthRadius($lat)
  {
      global $WGS84_a, $WGS84_b;
  
      $an = $WGS84_a * $WGS84_a * cos($lat);
      $bn = $WGS84_b * $WGS84_b * sin($lat);
      $ad = $WGS84_a * cos($lat);
      $bd = $WGS84_b * sin($lat);
  
      return sqrt(($an*$an + $bn*$bn)/($ad*$ad + $bd*$bd));
  }
  
  # Bounding box surrounding the point at given coordinates,
  # assuming local approximation of Earth surface as a sphere
  # of radius given by WGS84
  function boundingBox($latitudeInDegrees, $longitudeInDegrees, $halfSideInKm)
  {
      $lat = deg2rad($latitudeInDegrees);
      $lon = deg2rad($longitudeInDegrees);
      $halfSide = 1000 * $halfSideInKm;
  
      # Radius of Earth at given latitude
      $radius = WGS84EarthRadius($lat);
      # Radius of the parallel at given latitude
      $pradius = $radius*cos($lat);
  
      $latMin = $lat - $halfSide / $radius;
      $latMax = $lat + $halfSide / $radius;
      $lonMin = $lon - $halfSide / $pradius;
      $lonMax = $lon + $halfSide / $pradius;
  
      return array(rad2deg($latMin), rad2deg($lonMin), rad2deg($latMax), rad2deg($lonMax));
  }
  ?>
  

  回复于 2024-05-04T03:02:28+08:00
• 0

  感谢@Fedrico A.为Phyton实现，我已将其移植到Objective C类别类中 . 这是：

  #import "LocationService+Bounds.h"
  
  //Semi-axes of WGS-84 geoidal reference
  const double WGS84_a = 6378137.0; //Major semiaxis [m]
  const double WGS84_b = 6356752.3; //Minor semiaxis [m]
  
  @implementation LocationService (Bounds)
  
  struct BoundsLocation {
      double maxLatitude;
      double minLatitude;
      double maxLongitude;
      double minLongitude;
  };
  
  + (struct BoundsLocation)locationBoundsWithLatitude:(double)aLatitude longitude:(double)aLongitude maxDistanceKm:(NSInteger)aMaxKmDistance {
      return [self boundingBoxWithLatitude:aLatitude longitude:aLongitude halfDistanceKm:aMaxKmDistance/2];
  }
  
  #pragma mark - Algorithm 
  
  + (struct BoundsLocation)boundingBoxWithLatitude:(double)aLatitude longitude:(double)aLongitude halfDistanceKm:(double)aDistanceKm {
      double radianLatitude = [self degreesToRadians:aLatitude];
      double radianLongitude = [self degreesToRadians:aLongitude];
      double halfDistanceMeters = aDistanceKm*1000;
  
  
      double earthRadius = [self earthRadiusAtLatitude:radianLatitude];
      double parallelRadius = earthRadius*cosl(radianLatitude);
  
      double radianMinLatitude = radianLatitude - halfDistanceMeters/earthRadius;
      double radianMaxLatitude = radianLatitude + halfDistanceMeters/earthRadius;
      double radianMinLongitude = radianLongitude - halfDistanceMeters/parallelRadius;
      double radianMaxLongitude = radianLongitude + halfDistanceMeters/parallelRadius;
  
      struct BoundsLocation bounds;
      bounds.minLatitude = [self radiansToDegrees:radianMinLatitude];
      bounds.maxLatitude = [self radiansToDegrees:radianMaxLatitude];
      bounds.minLongitude = [self radiansToDegrees:radianMinLongitude];
      bounds.maxLongitude = [self radiansToDegrees:radianMaxLongitude];
  
      return bounds;
  }
  
  + (double)earthRadiusAtLatitude:(double)aRadianLatitude {
      double An = WGS84_a * WGS84_a * cosl(aRadianLatitude);
      double Bn = WGS84_b * WGS84_b * sinl(aRadianLatitude);
      double Ad = WGS84_a * cosl(aRadianLatitude);
      double Bd = WGS84_b * sinl(aRadianLatitude);
      return sqrtl( ((An * An) + (Bn * Bn))/((Ad * Ad) + (Bd * Bd)) );
  }
  
  + (double)degreesToRadians:(double)aDegrees {
      return M_PI*aDegrees/180.0;
  }
  
  + (double)radiansToDegrees:(double)aRadians {
      return 180.0*aRadians/M_PI;
  }
  
  
  
  @end
  

  我测试了它，似乎工作得很好 . Struct BoundsLocation应该被一个类替换，我只是在这里分享它 .

  回复于 2024-05-04T03:02:28+08:00
• 3

  All of the above answer are only partially correct . 特别是在像澳大利亚这样的地区，他们总是包括极点并计算一个非常大的矩形，甚至10kms .

  特别是Jan Philip Matuschek在http://janmatuschek.de/LatitudeLongitudeBoundingCoordinates#UsingIndex的算法包括一个非常大的矩形（-37，-90，-180,180），几乎在澳大利亚的每个点 . 这会影响数据库中的大量用户，并且必须为几乎一半国家的所有用户计算距离 .

  我发现 Drupal API Earth Algorithm by Rochester Institute of Technology 在极点和其他地方都能更好地工作，并且更容易实现 .

  https://www.rit.edu/drupal/api/drupal/sites%21all%21modules%21location%21earth.inc/7.54

  使用上述算法中的 earth_latitude_range 和 earth_longitude_range 来计算边界矩形

  并使用 distance calculation formula documented by google maps 来计算距离

  https://developers.google.com/maps/solutions/store-locator/clothing-store-locator#outputting-data-as-xml-using-php

  要按公里而不是里程搜索，请将3959替换为6371.对于（Lat，Lng）=（37，-122）和带有列lat和lng的Markers表，公式为：

  SELECT id, ( 3959 * acos( cos( radians(37) ) * cos( radians( lat ) ) * cos( radians( lng ) - radians(-122) ) + sin( radians(37) ) * sin( radians( lat ) ) ) ) AS distance FROM markers HAVING distance < 25 ORDER BY distance LIMIT 0 , 20;
  

  在https://stackoverflow.com/a/45950426/5076414阅读我的详细答案

  回复于 2024-05-04T03:02:28+08:00
• 28

  以下是Federico Ramponi在Go中的回答 . 注意：没有错误检查:(

  import (
      "math"
  )
  
  // Semi-axes of WGS-84 geoidal reference
  const (
      // Major semiaxis (meters)
      WGS84A = 6378137.0
      // Minor semiaxis (meters)
      WGS84B = 6356752.3
  )
  
  // BoundingBox represents the geo-polygon that encompasses the given point and radius
  type BoundingBox struct {
      LatMin float64
      LatMax float64
      LonMin float64
      LonMax float64
  }
  
  // Convert a degree value to radians
  func deg2Rad(deg float64) float64 {
      return math.Pi * deg / 180.0
  }
  
  // Convert a radian value to degrees
  func rad2Deg(rad float64) float64 {
      return 180.0 * rad / math.Pi
  }
  
  // Get the Earth's radius in meters at a given latitude based on the WGS84 ellipsoid
  func getWgs84EarthRadius(lat float64) float64 {
      an := WGS84A * WGS84A * math.Cos(lat)
      bn := WGS84B * WGS84B * math.Sin(lat)
  
      ad := WGS84A * math.Cos(lat)
      bd := WGS84B * math.Sin(lat)
  
      return math.Sqrt((an*an + bn*bn) / (ad*ad + bd*bd))
  }
  
  // GetBoundingBox returns a BoundingBox encompassing the given lat/long point and radius
  func GetBoundingBox(latDeg float64, longDeg float64, radiusKm float64) BoundingBox {
      lat := deg2Rad(latDeg)
      lon := deg2Rad(longDeg)
      halfSide := 1000 * radiusKm
  
      // Radius of Earth at given latitude
      radius := getWgs84EarthRadius(lat)
  
      pradius := radius * math.Cos(lat)
  
      latMin := lat - halfSide/radius
      latMax := lat + halfSide/radius
      lonMin := lon - halfSide/pradius
      lonMax := lon + halfSide/pradius
  
      return BoundingBox{
          LatMin: rad2Deg(latMin),
          LatMax: rad2Deg(latMax),
          LonMin: rad2Deg(lonMin),
          LonMax: rad2Deg(lonMax),
      }
  }
  

  回复于 2024-05-04T03:02:28+08:00