在Java中计算距离给定经纬度坐标一定距离的边界框

42

给定一个坐标(纬度,经度),我正在尝试计算距离该坐标一定距离(例如50公里)的正方形边界框。因此,输入为纬度、经度和距离,输出为两个坐标;一个是西南角(左下角),另一个是东北角(右上角)。在此之前,我见过几篇试图在Python中解决这个问题的答案,但我特别想要一个Java实现。

只是为了明确,我打算仅在地球上使用算法,因此不需要适应可变半径。

它不必非常准确(+/-20%可以接受),并且它只会用于在小范围内计算边界框(不超过150公里)。因此,我愿意为高效的算法牺牲一些准确性。非常感谢任何帮助。

编辑:我应该更清楚一些,我真的是想要一个正方形,而不是圆形。我知道正方形的中心与其周长各点之间的距离不像圆形那样是常量值。我想表达的是如果你从中心到周长上的任意一点画一条垂直于周长一侧的线,那么这4条线的长度相等。


如果你有Python实现,我不明白为什么不能将其转换为Java。算法应该是相同的。 - Alexandru Luchian
我找到的最好的Python示例被作者标记为未经测试。考虑到我需要的是一个相对宽容于准确性的东西,我想另一个算法可能更合适。 - Bryce Thomas
这个问题只涉及到纯三角函数;我建议删除Java和算法标签。(如果可能的话。) - Kevin Bourrillion
你只需要一个公式,我相信有人会为您提供,并希望给您提供参考链接。 - Kevin Bourrillion
地球的球形特性需要考虑。如果输入位置是北极,您想要什么答案? - MarkJ
6个回答

59

2
我对为什么这些球面逼近被视为标准解决方案感到困惑。这里涉及了可怕数量的三角函数。此外,所有GPS设备都基于WGS84椭球体生成纬度/经度坐标--如果您想要完全准确,您根本无法使用球面逼近... 调整这些球面逼近中的方程以考虑椭球体将导致难以想象的巨大计算复杂度。看起来方法必须是在三维空间中执行计算和弧长积分。 - Steven Lu
1
@Mecki 抱歉我没有表达清楚,需要我反思一下才能想起来我投诉的原因。基本上,“问题”在于对于计算机而言,三角函数运算不是很快(与比较简单的操作如大于/小于相比)。这个问题要求进行基本的边界框检查。OP甚至用多个句子指定他不需要极高的精度,实际上他需要的是速度。因此,在每个测试点上使用大量的三角函数来获得是否在给定距离内的精确结果是过度杀伤力的。 - Steven Lu
1
OP 想要做的事情(就像我在上一个项目中所做的那样,在那个项目中我们使用纬度/经度进行了工作)是应用一种基本的快速近似方法:你需要两个常量,即纬度一度对应的公里数和赤道经度一度对应的公里数。然后只需使用这个来计算在多少度的范围内您的点应该通过检查。是的,这种近似法在靠近极地的地方会降解。虽然对于大多数应用程序来说,这并不重要。地理兴趣点通常离极地很远。(实际上我有点过于简化了,它还有点复杂) - Steven Lu
1
这是关于编程的内容,翻译成中文:这是在每个点进行AABB检查时,误报率微不足道(如果您的屏幕显示区域是矩形的话,误报率甚至比您预期的还要低,事实上OP也明确指出了这一点),与每个点进行3或9个三角函数调用而没有误报率之间的区别。棘手的部分是,当您接近极点时,您必须测试更大的经度范围以安全地覆盖您的框。 - Steven Lu
1
我想说的是,只要有人真正阅读了链接的文章,他们就应该能够理解如何以超准确的方式实现它,并考虑如何调整它并放松条件,以使代码运行更少的三角函数。有时您可以保证始终有足够的时间为每个查询点运行大量的三角函数。 - Steven Lu
显示剩余5条评论

16
double R = 6371;  // earth radius in km

double radius = 50; // km

double x1 = lon - Math.toDegrees(radius/R/Math.cos(Math.toRadians(lat)));

double x2 = lon + Math.toDegrees(radius/R/Math.cos(Math.toRadians(lat)));

double y1 = lat + Math.toDegrees(radius/R);

double y2 = lat - Math.toDegrees(radius/R);

虽然我也建议使用JTS。


9
JTS代表什么? - Gili
1
Java Topology Suite,它是一个用于几何相关操作的API。 - prettyvoid
如何从xyz瓦片计算半径的值? - Reyan Chougle

5
import com.vividsolutions.jts.geom.Envelope;

...
Envelope env = new Envelope(centerPoint.getCoordinate());
env.expandBy(distance_in_degrees); 
...

现在env包含您的信封。它实际上不是一个“正方形”(无论在球体表面意味着什么),但它应该可以胜任。
您应该注意,度数距离将取决于中心点的纬度。在赤道,1度纬度约为111公里,但在纽约只有约75公里。
真正酷的事情是,您可以将所有点投入com.vividsolutions.jts.index.strtree.STRtree,然后使用它快速计算那个信封内的点。

3
所有之前的答案都只是部分正确,特别是在像澳大利亚这样的地区,他们总是包括极点,并为10公里甚至计算一个非常大的矩形。
特别是Jan Philip Matuschek的算法(http://janmatuschek.de/LatitudeLongitudeBoundingCoordinates#UsingIndex)针对澳大利亚几乎每个点都包括了一个非常大的矩形从(-37,-90,-180,180)。这会影响到数据库中的大量用户,需要为几乎一半的国家的所有用户计算距离。
我发现Rochester Institute of Technology的Drupal API Earth Algorithm在极点和其他地方都更好,并且更容易实现。

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

使用上述算法中的earth_latitude_rangeearth_longitude_range来计算边界矩形。
以下是Java的实现。
    /**
 * Get bouding rectangle using Drupal Earth Algorithm
 * @see https://www.rit.edu/drupal/api/drupal/sites%21all%21modules%21location%21earth.inc/7.54
 * @param lat
 * @param lng
 * @param distance
 * @return
 */
default BoundingRectangle getBoundingRectangleDrupalEarthAlgo(double lat, double lng, int distance) {
    lng = Math.toRadians(lng);
    lat = Math.toRadians(lat);
    double radius = earth_radius(lat);
    List<Double> retLats = earth_latitude_range(lat, radius, distance);
    List<Double> retLngs = earth_longitude_range(lat, lng, radius, distance);
    return new BoundingRectangle(retLats.get(0), retLats.get(1), retLngs.get(0), retLngs.get(1));
}


/**
 * Calculate latitude range based on earths radius at a given point
 * @param latitude
 * @param longitude
 * @param distance
 * @return
 */
default List<Double> earth_latitude_range(double lat, double radius, double distance) {
      // Estimate the min and max latitudes within distance of a given location.

      double angle = distance / radius;
      double minlat = lat - angle;
      double maxlat = lat + angle;
      double rightangle = Math.PI / 2;
      // Wrapped around the south pole.
      if (minlat < -rightangle) {
        double overshoot = -minlat - rightangle;
        minlat = -rightangle + overshoot;
        if (minlat > maxlat) {
          maxlat = minlat;
        }
        minlat = -rightangle;
      }
      // Wrapped around the north pole.
      if (maxlat > rightangle) {
        double overshoot = maxlat - rightangle;
        maxlat = rightangle - overshoot;
        if (maxlat < minlat) {
          minlat = maxlat;
        }
        maxlat = rightangle;
      }
      List<Double> ret = new ArrayList<>();
      ret.add((minlat));
      ret.add((maxlat));
      return ret;
    }

/**
 * Calculate longitude range based on earths radius at a given point
 * @param lat
 * @param lng
 * @param earth_radius
 * @param distance
 * @return
 */
default List<Double> earth_longitude_range(double lat, double lng, double earth_radius, int distance) {
      // Estimate the min and max longitudes within distance of a given location.
      double radius = earth_radius * Math.cos(lat);

      double angle;
      if (radius > 0) {
        angle = Math.abs(distance / radius);
        angle = Math.min(angle, Math.PI);
      }
      else {
        angle = Math.PI;
      }
      double minlong = lng - angle;
      double maxlong = lng + angle;
      if (minlong < -Math.PI) {
        minlong = minlong + Math.PI * 2;
      }
      if (maxlong > Math.PI) {
        maxlong = maxlong - Math.PI * 2;
      }

      List<Double> ret = new ArrayList<>();
      ret.add((minlong));
      ret.add((maxlong));
      return ret;
    }

/**
 * Calculate earth radius at given latitude
 * @param latitude
 * @return
 */
default Double earth_radius(double latitude) {
      // Estimate the Earth's radius at a given latitude.
      // Default to an approximate average radius for the United States.
      double lat = Math.toRadians(latitude);

      double x = Math.cos(lat) / 6378137.0;
      double y = Math.sin(lat) / (6378137.0 * (1 - (1 / 298.257223563)));

      //Make sure earth's radius is in km , not meters
      return (1 / (Math.sqrt(x * x + y * y)))/1000;
    }

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

观察 earth_radius 代码。确保地球半径以公里为单位,而不是米:返回 (1 / (Math.sqrt(x * x + y * y)))/1000; - Sacky San

1
根据IronMan的回应:
/**
 * Calculate the lat and len of a square around a point.
 * @return latMin, latMax, lngMin, lngMax
 */
public static double[] calculateSquareRadius(double lat, double lng, double radius) {
    double R = 6371;  // earth radius in km
    double latMin = lat - Math.toDegrees(radius/R);
    double latMax = lat + Math.toDegrees(radius/R);
    double lngMin = lng - Math.toDegrees(radius/R/Math.cos(Math.toRadians(lat)));
    double lngMax = lng + Math.toDegrees(radius/R/Math.cos(Math.toRadians(lat)));

    return new double[] {latMin, latMax, lngMin, lngMax};
}

0

这里有一个简单的解决方案,我用它来生成边界框坐标,然后与GeoNames citieJSON API一起使用,以获取从GPS十进制坐标附近的大城市。

这是我的GitHub存储库中的Java方法:FusionTableModifyJava

我有一个十进制GPS位置,需要找到“附近”的最大城市/州。我需要一个相对准确的边界框,以传递给citiesJSON GeoNames网络服务,以获取该边界框中最大的城市。我传递位置和我感兴趣的“半径”(以公里为单位),然后它会返回所需的北、南、东、西十进制坐标,以传递给citiesJSON。

(我在研究过程中发现了这些资源:

计算纬度/经度点之间的距离、方位等更多内容。

经度 - 维基百科)

它不是非常准确,但对于我使用的目的足够了:

    // Compute bounding Box coordinates for use with Geonames API.
    class BoundingBox
    {
        public double north, south, east, west;
        public BoundingBox(String location, float km)
        {
             //System.out.println(location + " : "+ km);
            String[] parts = location.replaceAll("\\s","").split(","); //remove spaces and split on ,

            double lat = Double.parseDouble(parts[0]);
            double lng = Double.parseDouble(parts[1]);

            double adjust = .008983112; // 1km in degrees at equator.
            //adjust = 0.008983152770714983; // 1km in degrees at equator.

            //System.out.println("deg: "+(1.0/40075.017)*360.0);


            north = lat + ( km * adjust);
            south = lat - ( km * adjust);

            double lngRatio = 1/Math.cos(Math.toRadians(lat)); //ratio for lng size
            //System.out.println("lngRatio: "+lngRatio);

            east = lng + (km * adjust) * lngRatio;
            west = lng - (km * adjust) * lngRatio;
        }

    }

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接