Postgresql Postgis

Postgis：距离计算和排序优化

Posted by YI on 2022-11-20

Postgis：距离计算和排序优化

场景和数据

使用场景：使用Postgis寻找距离某个圆心一定距离内的所有企业，按照离圆心的距离排序。

地理数据说明：来自于高德地图地理编码api的经纬度，数据量40w，在Postgresql中使用geometry存储。注：高德地图采用GCJ20坐标系。

实现和优化思路

  1. 采用Postgis提供的函数st_dwithin圈选出范围

    Postgis St_Dwithin官方文档

    从文档中可以看出，我们可以用两种方式来描述距离，一是采用distance_of_srid，另一种是采用米。这两种方式对应了不同的计算逻辑。

    • 其中distance_of_srid的含义是srid定义的距离。srid是空间引用标识符，每个srid对应了特定椭圆体系的空间引用系统，可以理解为每个srid对应了一个对地球的建模方式，不同的空间引用系统会对地球进行不同的建模，比如srid=4326就对应了WGS 84空间坐标系统（1984世界大地坐标系），这个坐标系将地球描绘为长半轴6378137米，短半轴6356752米，扁率1/298.257223563的椭球，通过将经纬度投影到这个空间系统对应的平面坐标系中就可以计算距离了。不同的空间坐标系计算出的距离单位不同，比如WGS84计算出的距离是用度表示的。

    • 第二种是采用米的方式，使用这个距离要设定函数的参数use_spheroid，不论是设置为true还是fasle都会采用米来判断距离。使用了第二种方式的计算逻辑与前一种方式不同，这种方式不再将经纬度投影到平面坐标系中，而是采用球面计算的方式，直接在球体表面上计算距离，明显这种方式会增加计算的复杂度。use_spheroid设置为false时，会采用大圆距离（the great circle）的计算方式，也就是将地球视作一个规整的球体，而非椭圆，这样减少了计算的复杂度，也会比平面坐标系更加精准一些；如果将use_spheroid设置为true，则将地球视作椭球体，考虑长半轴、短半轴和扁率等参数，计算复杂度会增加，但也更加符合现实状况。然而不同的空间坐标系对地球的椭球体建模并不相同，所以计算出来的距离也会存在差异，在另一个计算距离的函数ST_DistanceSpheroid中，我们可以选择不同的空间坐标系来计算椭球表面的距离，如：

      ST_DistanceSpheroid(ST_Centroid(geom), ST_GeomFromText('POINT(-118 38)',4326), 'SPHEROID["WGS 84",6378137,298.257223563]')

      而在st_dwithin函数中却无法选择，只能使用默认的WGS84空间体系。

  2. 排序

    Postgis <->文档

    从文档中可以知道Postgis提供了支持索引的排序方式，在order by 中使用符号<->进行排序会使用gist索引进行加速，比使用st_distance排序速度更快。文档中提到这里采用的排序方式是knn最近邻，可以预见的是，当搜索范围变得足够大时，这个排序算法的时间消耗将会指数级增加，所以限制被排序数据的行数将会是优化排序速度的一个方向。

  3. 优化思路

    根据PG大佬digoal的研究，PG GiST距离排序操作符和过滤无法同时使用索引，所以我们的场景可能只有过滤的过程触发了索引，但个人实验时发现仅使用<->的执行计划中也没有显示走索引的情况，不知道是执行计划没有显示还是并没有触发。

    参考PG GiST距离排序操作符和过滤无法同时使用索引文章可以进行一些优化操作。

大圆距离

  1. 常规的大圆距离计算方式

    大圆距离是比较简单的对地球建模的方法，将地球视作一个规整的球体，并将地理经纬度直接作为球体的经纬度来计算距离。这个距离的算法如下：

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    private static double rad(double d) {
        return d * Math.PI / 180.0;
    }

    /**
    * 通过经纬度获取距离(单位：米)
    *
    * @return 距离(米)
    */
    public static double getDistance(MapPoint p1, MapPoint p2,Double radius) {
        double radLat1 = rad(p1.getLatitude());
        double radLat2 = rad(p2.getLatitude());
        double a = radLat1 - radLat2;
        double b = rad(p1.getLongitude()) - rad(p2.getLongitude());
        double s = 2 * Math.asin(Math.sqrt(Math.pow(Math.sin(a / 2), 2)
                                           + Math.cos(radLat1) * Math.cos(radLat2)
                                           * Math.pow(Math.sin(b / 2), 2)));
        s = s * radius;
        s = s * 1000;
        return s;
    }

    这里的s=s*radius步骤就是把角度转换成距离的一步，使用的radius一般是地球的平均半径（6371km），也有使用赤道半径（6378km）的情况，这通常会导致10米左右的差距。

高德地图API中的距离计算

  1. 高德地图的计算方式

    在postgis和大部分计算机语言的地理框架（比如python中geopy库的distance.great_circle和java的org.gavaghan.geodesy）中，大圆距离的计算都是采用地球平均半径的，但高德和腾讯地图（GCJ20坐标系）计算出来的距离跟采用赤道半径计算出来的大圆距离相同（个人测试，用上述算法和赤道半径计算距离与高德腾讯计算出来的距离在小数点6位前都能匹配上）。

    这将会导致一部分偏差，在追求高精度的情况下，这种偏差是需要注意的。

CATALOG
    FEATURED TAGS
    Postgresql Postgis
    FRIENDS