根据经纬度判断点是否位于多边形内以及到多边形的距离

最近，在业务中遇到一个场景，需要根据经纬度，判断一个点是否位于多边形区域内，如果不在多边形内部，计算点到多边形的最短距离，即到多边形各个边的距离的最短值。

经过研究和参考网上的文章，最终实现了该功能。

一、定义多边形和点

// 定义多边形，沿顺时针或逆时针方向，依次添加点的经纬度
List<double[]> polygon = new ArrayList<>();
polygon.add(new double[]{121.553793, 29.808325});
polygon.add(new double[]{121.553793, 29.808329});
polygon.add(new double[]{121.553799, 29.808329});
polygon.add(new double[]{121.553799, 29.808325});

// 定义点
double[] point = new double[]{121.553792, 29.808325};

二、判断点是否在多边形内部

点在多边形边上，认为点不在多边形内部

/**
    * 判断一个点是否在多边形内部；在多边形边上也不算是在多边形内部
    * @param polygon 经纬度点组成的多边形
    * @param point 经纬度点
    * @return 点是否在多边形内
    */
public static boolean isInside(List<double[]> polygon, double[] point) {
    boolean inside = false;
    for (int i = 0, j = polygon.size() - 1; i < polygon.size(); j = i++) {
        if (((polygon.get(i)[1] > point[1]) != (polygon.get(j)[1] > point[1])) &&
                (point[0] < (polygon.get(j)[0] - polygon.get(i)[0]) * (point[1] - polygon.get(i)[1]) /
                        (polygon.get(j)[1] - polygon.get(i)[1]) + polygon.get(i)[0])) {
            inside = !inside;
        }
    }
    return inside;
}

三、计算外部点到多边形最短距离

计算外部点到多边形的距离，分为以下几个步骤：

• 点到多边形各个边的距离

        多边形各个边由多边形的点顺时针或者逆时针取相邻两个点，组成线段

• 取点到边的距离的最小值

/**
    * 计算点到多边形的最短距离
    * @param polygon 多边形各个点的经纬度
    * @param point 点
    * @return 点到多边形最短距离
    */
public static double shortestDistance(List<double[]> polygon, double[] point) {
    double shortestDistance = Double.POSITIVE_INFINITY;
    for (int i = 0, j = polygon.size() - 1; i < polygon.size(); j = i++) {
        double[] p1 = polygon.get(i);
        double[] p2 = polygon.get(j);
        double d = calcP2L(point, p1, p2);
        if (d < shortestDistance) {
            shortestDistance = d;
        }
    }
    return shortestDistance;
}

/**
    * 计算点到线段的距离
    * @param point 点的经度和纬度
    * @param point1 线段其中一个端点的经度和纬度
    * @param point2 线段另一个端点的经度和纬度
    * @return 点到线段的最短距离
    */
public static double calcP2L(double[] point, double[] point1, double[] point2) {
    double x = point[0], y = point[1];
    double x1 = point1[0], y1 = point1[1];
    double x2 = point2[0], y2 = point2[1];
    double d1 = (x2 - x1) * (x - x1) + (y2 - y1) * (y - y1);
    if (d1 < 0) {
        return calcP2P(point, point1);
    }
    double d2 = (x2 - x1) * (x2 - x1) + (y2 - y1) * (y2 - y1);
    if (d1 >= d2) {
        return calcP2P(point, point2);
    }
    double r = d1 / d2;
    double px = x1 + (x2 - x1) * r;
    double py = y1 + (y2 - y1) * r;
    double[] ppoint = new double[]{px, py};
    return calcP2P(point, ppoint);
}

/**
    * 计算两点间的球面距离
    * @param apoint a点经纬度
    * @param bpoint b点经纬度
    * @return 点到点的距离
    */
public static double calcP2P(double[] apoint, double[] bpoint) {
    double alat = apoint[0];
    double alon = apoint[1];
    double blat = bpoint[0];
    double blon = bpoint[1];

    double R = 6378.137 * 1000;// 地球半径(米)
    // 判断点是在北半球还是南半球，本程序中输入的数据若为负则表示在南边球
    double distance = 0.0D;
    double _alat = alat * (Math.PI / 180);//弧度
    double _alon = alon * (Math.PI / 180);
    double _blat = blat * (Math.PI / 180);
    double _blon = blon * (Math.PI / 180);
    double c = Math.sin(_alat) * Math.sin(_blat) + Math.cos(_alat) * Math.cos(_blat) * Math.cos(_alon - _blon);// Java中三角函数角度以弧度制表示
    if (c > 1) {
        c = 1;
    }

    distance = Math.acos(c) * R; // 弧长公式：弧长 = 弧度 * 半径
    if (distance <= 0.01) { // GPS误差
        distance = 0.0D;
    }
    return distance * 1.9088; // 经过测试，distance和真实距离差 1.9088 倍，在这里扩大 1.9088 倍进行调整
}

在主方法shortestDistance中，依次获取多边形的边p1、p2，调用calcP2L计算点到线段的距离；在calcP2L方法中，再调用calcP2P计算点到点的距离，获取点到点的距离，进而获取点到线段的距离，最终获取点到多边形的距离。

细节理论，后续补充。。。。。。