【问题标题】:How to calculate the bounding box for a given lat/lng location?如何计算给定纬度/经度位置的边界框?
【发布时间】:2010-09-19 06:58:54
【问题描述】:

我已经给出了一个由纬度和经度定义的位置。 现在我想计算一个边界框,例如距离该点 10 公里。

边界框应定义为 latmin, lngmin 和 latmax, lngmax。

我需要这些东西才能使用panoramio API

有人知道如何获得积分的公式吗?

编辑: 伙计们,我正在寻找一个公式/函数,它以 lat & lng 作为输入并返回一个边界框作为 latmin & lngmin 和 latmax & latmin。 Mysql, php, c#, javascript都可以,伪代码也可以。

编辑:我不是在寻找显示 2 点距离的解决方案

【问题讨论】:

  • 如果您在某处使用地理数据库,它们肯定集成了边界框计算。例如,您甚至可以查看 PostGIS/GEOS 的来源。

标签: geocoding panoramio


【解决方案1】:

您正在寻找一个椭球公式。

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

Geo::Ellipsoid

看看location 方法。调用它两次,你就得到了你的 bbox。

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

哦,如果您现在还没有弄清楚,Google 地图使用 WGS84 椭球。

【讨论】:

    【解决方案2】:

    我建议将地球表面局部近似为一个球体,其半径由给定纬度的 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))
    

    【讨论】:

    • 正如建议的 CPAN 库的文档中所指出的,这仅对 halfSide
    • 这在两极附近有效吗?似乎没有,因为它看起来以 latMin pi (对于北极)结束?我认为,当您在杆的 halfSide 内时,您需要返回一个边界框,其中包括所有经度和通常为远离杆的一侧和靠近杆的一侧计算的纬度。
    • 这里是来自 JanMatuschek.de 规范的 PHP 实现:github.com/anthonymartin/GeoLocation.class.php
    • 我在下面添加了这个答案的 C# 实现。
    • @FedericoA.Ramponi 这里的 haldSideinKm 是什么?不明白...我必须在这个参数中传递什么,地图中两点之间的半径还是什么?
    【解决方案3】:

    我写了一篇关于寻找边界坐标的文章:

    http://JanMatuschek.de/LatitudeLongitudeBoundingCoordinates

    文章解释了这些公式,还提供了一个 Java 实现。 (这也说明了为什么 Federico 的最小/最大经度公式不准确。)

    【讨论】:

    • 我已经为您的 GeoLocation 类创建了一个 PHP 端口。可以在这里找到:pastie.org/5416584
    • 我已经上传到github了:github.com/anthonymartin/GeoLocation.class.php
    • 这还能回答问题吗?如果我们只有 1 个起点,我们将无法像此代码中那样计算大圆距离,这需要两个经纬度位置。
    • 您的 C# 变体中存在错误代码,例如:public override string ToString(),仅出于一个目的覆盖这样的全局方法是非常糟糕的,最好添加另一种方法,然后覆盖标准方法,可用于应用程序的其他部分,不适用于 gis 精确...
    • 这里是 Jan 的 GeoLocaiton 类的 PHP 端口的更新链接:github.com/anthonymartin/GeoLocation.php
    【解决方案4】:

    我改编了一个我发现的 PHP 脚本来执行此操作。您可以使用它来查找围绕某个点(例如 20 公里)的框的角。我的具体示例是针对 Google Maps API:

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

    【讨论】:

    • -1 OP正在寻找的是:给定一个参考点(纬度,经度)和距离,找到最小的盒子,使得所有点都
    • 好吧,碰巧,这正是我所需要的。所以谢谢你,即使它不能完全回答 this 问题:)
    • 嗯,我很高兴它可以帮助别人!
    【解决方案5】:

    在这里,我已将 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));
    }
    

    【讨论】:

    • 谢谢,这对我有用。不得不手动测试代码,不知道如何为此编写单元测试,但它会产生准确的结果,达到我需要的准确度
    • 这里的 haldSideinKm 是什么?不明白...我必须在这个参数中传递什么,地图中两点之间的半径还是什么?
    • @GeloVolro:这是你想要的边界框的一半长度。
    • 您不一定要编写自己的 MapPoint 类。 System.Device.Location 中有一个以纬度和经度为参数的 GeoCoordinate 类。
    • 这很好用。我非常感谢 C# 端口。
    【解决方案6】:

    很简单,只要去panoramio网站,然后从panoramio网站打开世界地图。然后去指定的位置,需要经纬度。

    然后你在地址栏中找到了纬度和经度,例如在这个地址中。

    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>
    

    【讨论】:

      【解决方案7】:

      我正在研究边界框问题,作为查找静态 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 部分,我遗漏了一些位置点。

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

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

      【讨论】:

        【解决方案8】:

        我写了一个 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()
          ];
        };
        

        【讨论】:

        • 这段代码根本不起作用。我的意思是,即使修复了 minLon = void 0;maxLon = MAX_LON; 这样的明显错误,它仍然不起作用。
        • @aroth,我刚刚测试过,没有问题。请记住,centerPoint 参数是一个由两个坐标组成的数组。例如,getBoundingBox([42.2, 34.5], 50)void 0 是“未定义”的 CoffeeScript 输出,不会影响代码的运行能力。
        • 此代码不起作用。 degLat.degToRad 不是函数
        • 代码在 Node 和 Chrome 中按原样工作,直到我将其放入我正在处理的项目中并开始收到“degToRad 不是函数”错误。从来没有发现原因,但Number.prototype. 对于这样的实用程序函数不是一个好主意,所以我将它们转换为普通的本地函数。同样重要的是要注意返回的框是 [LNG,LAT,LNG,LAT] 而不是 [LAT,LNG,LAT,LNG]。我在使用它时修改了返回函数以避免混淆。
        【解决方案9】:

        @Jan Philip Matuschek 出色解释的插图。(请投票给他的答案,而不是这个;我添加这个是因为我花了一点时间来理解原始答案)

        优化寻找最近邻居的边界框技术需要为距离 d 处的点 P 推导出最小和最大纬度、经度对。落在这些之外的所有点肯定与该点的距离大于 d。 这里要注意的一件事是交叉点纬度的计算,正如 Jan Philip Matuschek 解释中强调的那样。交点的纬度不在点 P 的纬度上,而是稍微偏离它。在确定距离 d 的点 P 的正确最小和最大边界经度时,这是一个经常被忽略但很重要的部分。这在验证中也很有用。

        P的(交点纬度,经度高)到(纬度,经度)的半正弦距离等于距离d。

        这里是 Python 要点https://gist.github.com/alexcpn/f95ae83a7ee0293a5225

        【讨论】:

          【解决方案10】:

          如果有人感兴趣,我在这里将 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));
          }
          ?>
          

          【讨论】:

          • 对于 PHP 7.4,我不得不稍微修改代码以使其正常工作:` # Semi-axes of WGS-84 geoidal reference $WGS84_a = floatval( 6378137.0 ); # 主要半轴 [m] $WGS84_b = floatval( 6356752.3 ); # 次半轴 [m] `
          【解决方案11】:

          由于我需要一个非常粗略的估计,所以为了在 elasticsearch 查询中过滤掉一些不需要的文档,我使用了以下公式:

          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 = 指定位置所需的公里数。对于您的情况 N=10

          不准确但方便。

          【讨论】:

          • 确实,不准确,但仍然有用且很容易实现。
          【解决方案12】:

          这是一个使用 javascript 的简单实现,它基于将纬度转换为 kms 其中1 degree latitude ~ 111.2 km

          我正在根据给定的纬度、经度和半径(以千米为单位)计算地图的边界。

          function getBoundsFromLatLng(lat, lng, radiusInKm){
               var lat_change = radiusInKm/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;
          }
          

          【讨论】:

            【解决方案13】:

            感谢@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 应该换成一个类,我用过只是在这里分享一下。

            【讨论】:

              【解决方案14】:

              以上所有答案都只是部分正确。特别是在像澳大利亚这样的地区,他们总是包括极点并计算一个非常大的矩形,即使是 10 公里。

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

              我发现罗彻斯特理工学院的 Drupal API 地球算法在极点周围和其他地方都表现得更好,并且更容易实施。

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

              使用上述算法中的earth_latitude_rangeearth_longitude_range计算边界矩形

              并使用谷歌地图记载的距离计算公式计算距离

              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阅读我的详细回答

              【讨论】:

                【解决方案15】:

                这是 Federico Ramponi 在围棋中的答案。注意:没有错误检查:(

                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),
                    }
                }
                

                【讨论】:

                  【解决方案16】:

                  这是用于获取基于纬度/经度和距离的边界框坐标的 javascript 代码。经过测试,工作正常。

                  Number.prototype.degreeToRadius = function () {
                      return this * (Math.PI / 180);
                  };
                  
                  Number.prototype.radiusToDegree = function () {
                      return (180 * this) / Math.PI;
                  };
                  
                  function getBoundingBox(fsLatitude, fsLongitude, fiDistanceInKM) {
                  
                      if (fiDistanceInKM == null || fiDistanceInKM == undefined || fiDistanceInKM == 0)
                          fiDistanceInKM = 1;
                      
                      var MIN_LAT, MAX_LAT, MIN_LON, MAX_LON, ldEarthRadius, ldDistanceInRadius, lsLatitudeInDegree, lsLongitudeInDegree,
                          lsLatitudeInRadius, lsLongitudeInRadius, lsMinLatitude, lsMaxLatitude, lsMinLongitude, lsMaxLongitude, deltaLon;
                      
                      // coordinate limits
                      MIN_LAT = (-90).degreeToRadius();
                      MAX_LAT = (90).degreeToRadius();
                      MIN_LON = (-180).degreeToRadius();
                      MAX_LON = (180).degreeToRadius();
                  
                      // Earth's radius (km)
                      ldEarthRadius = 6378.1;
                  
                      // angular distance in radians on a great circle
                      ldDistanceInRadius = fiDistanceInKM / ldEarthRadius;
                  
                      // center point coordinates (deg)
                      lsLatitudeInDegree = fsLatitude;
                      lsLongitudeInDegree = fsLongitude;
                  
                      // center point coordinates (rad)
                      lsLatitudeInRadius = lsLatitudeInDegree.degreeToRadius();
                      lsLongitudeInRadius = lsLongitudeInDegree.degreeToRadius();
                  
                      // minimum and maximum latitudes for given distance
                      lsMinLatitude = lsLatitudeInRadius - ldDistanceInRadius;
                      lsMaxLatitude = lsLatitudeInRadius + ldDistanceInRadius;
                  
                      // minimum and maximum longitudes for given distance
                      lsMinLongitude = void 0;
                      lsMaxLongitude = void 0;
                  
                      // define deltaLon to help determine min and max longitudes
                      deltaLon = Math.asin(Math.sin(ldDistanceInRadius) / Math.cos(lsLatitudeInRadius));
                  
                      if (lsMinLatitude > MIN_LAT && lsMaxLatitude < MAX_LAT) {
                          lsMinLongitude = lsLongitudeInRadius - deltaLon;
                          lsMaxLongitude = lsLongitudeInRadius + deltaLon;
                          if (lsMinLongitude < MIN_LON) {
                              lsMinLongitude = lsMinLongitude + 2 * Math.PI;
                          }
                          if (lsMaxLongitude > MAX_LON) {
                              lsMaxLongitude = lsMaxLongitude - 2 * Math.PI;
                          }
                      }
                  
                      // a pole is within the given distance
                      else {
                          lsMinLatitude = Math.max(lsMinLatitude, MIN_LAT);
                          lsMaxLatitude = Math.min(lsMaxLatitude, MAX_LAT);
                          lsMinLongitude = MIN_LON;
                          lsMaxLongitude = MAX_LON;
                      }
                  
                      return [
                          lsMinLatitude.radiusToDegree(),
                          lsMinLongitude.radiusToDegree(),
                          lsMaxLatitude.radiusToDegree(),
                          lsMaxLongitude.radiusToDegree()
                      ];
                  };
                  

                  使用下面的getBoundingBox函数来绘制一个边界框。

                  var lsRectangleLatLong = getBoundingBox(parseFloat(latitude), parseFloat(longitude), lsDistance);
                              if (lsRectangleLatLong != null && lsRectangleLatLong != undefined) {
                                  latLngArr.push({ lat: lsRectangleLatLong[0], lng: lsRectangleLatLong[1] });
                                  latLngArr.push({ lat: lsRectangleLatLong[0], lng: lsRectangleLatLong[3] });
                                  latLngArr.push({ lat: lsRectangleLatLong[2], lng: lsRectangleLatLong[3] });
                                  latLngArr.push({ lat: lsRectangleLatLong[2], lng: lsRectangleLatLong[1] });
                              }
                  

                  【讨论】:

                  • 像魅力一样工作......谢谢。
                  猜你喜欢
                  • 2016-12-16
                  • 2020-05-04
                  • 2013-10-27
                  • 2016-11-21
                  • 1970-01-01
                  • 1970-01-01
                  • 2011-03-17
                  • 2011-07-29
                  • 1970-01-01
                  相关资源
                  最近更新 更多