地图两点距离的简单计算

写在前面

在日常涉及地图的业务中，通常会涉及点到点的距离计算。
一般基于余弦定理公式，Haversine公式，勾股定理，Vincenty 公式 实现，此类方法可满足常见的距离计算需求。

核心代码

<?php
namespace App\Service\Common;


class MapDistanceUtil
{

    /**
     * 地球半径（米）
     */
    const EARTH_RADIUS = 6370993.0;

    /**
     * 计算两点间直线距离（球面余弦定理公式）
     * 适用于各种距离计算，精度较高
     *
     * @param float $lng1 起点经度
     * @param float $lat1 起点纬度
     * @param float $lng2 终点经度
     * @param float $lat2 终点纬度
     * @return float 距离（米）
     */
    public function getDistance(float $lng1, float $lat1, float $lng2, float $lat2): float
    {
        // 将角度转换为弧度
        $radLat1 = deg2rad($lat1);
        $radLat2 = deg2rad($lat2);
        $radLng1 = deg2rad($lng1);
        $radLng2 = deg2rad($lng2);

        // 纬度差
        $diffLat = $radLat1 - $radLat2;
        // 经度差
        $diffLng = $radLng1 - $radLng2;

        // 球面余弦定理公式
        $distance = 2 * asin(sqrt(pow(sin($diffLat / 2), 2) + cos($radLat1) * cos($radLat2) * pow(sin($diffLng / 2), 2)));
        $distance = $distance * self::EARTH_RADIUS;

        return round($distance, 2);
    }

    /**
     * 计算两点间直线距离（简化版勾股定理方法）
     * 适用于短距离计算，速度快但长距离有误差
     *
     * @param float $lng1 起点经度
     * @param float $lat1 起点纬度
     * @param float $lng2 终点经度
     * @param float $lat2 终点纬度
     * @return float 距离（米）
     */
    public function getShortDistance(float $lng1, float $lat1, float $lng2, float $lat2): float
    {
        // 将角度转换为弧度
        $radLat1 = deg2rad($lat1);
        $radLng1 = deg2rad($lng1);
        $radLat2 = deg2rad($lat2);
        $radLng2 = deg2rad($lng2);

        // 经纬度差值
        $diffLat = $radLat1 - $radLat2;
        $diffLng = $radLng1 - $radLng2;

        // 勾股定理计算直线距离
        $distance = sqrt(pow($diffLat * self::EARTH_RADIUS, 2) + pow($diffLng * self::EARTH_RADIUS * cos($radLat1), 2));

        return round($distance, 2);
    }

    /**
     * 高精度计算两点间距离（椭球体模型 Vincenty 公式）
     *
     * @param float $lng1 起点经度
     * @param float $lat1 起点纬度
     * @param float $lng2 终点经度
     * @param float $lat2 终点纬度
     * @return float 距离（米）
     */
    public function getPreciseDistance(float $lng1, float $lat1, float $lng2, float $lat2): float
    {
        $f = 1 / 298.257223563; // 扁率
        $b = (1 - $f) * self::EARTH_RADIUS; // 极半径

        // 将角度转换为弧度
        $lat1 = deg2rad($lat1);
        $lng1 = deg2rad($lng1);
        $lat2 = deg2rad($lat2);
        $lng2 = deg2rad($lng2);

        $L = $lng2 - $lng1; // 两点经度差
        $U1 = atan((1 - $f) * tan($lat1)); // 归化纬度
        $U2 = atan((1 - $f) * tan($lat2));
        $sinU1 = sin($U1);
        $cosU1 = cos($U1);
        $sinU2 = sin($U2);
        $cosU2 = cos($U2);

        $lambda = $L;
        $iterLimit = 100;
        do {
            $sinLambda = sin($lambda);
            $cosLambda = cos($lambda);
            $sinSigma = sqrt(($cosU2 * $sinLambda) * ($cosU2 * $sinLambda) +
                ($cosU1 * $sinU2 - $sinU1 * $cosU2 * $cosLambda) *
                ($cosU1 * $sinU2 - $sinU1 * $cosU2 * $cosLambda));

            if ($sinSigma == 0) return 0; // 重合点

            $cosSigma = $sinU1 * $sinU2 + $cosU1 * $cosU2 * $cosLambda;
            $sigma = atan2($sinSigma, $cosSigma);
            $sinAlpha = $cosU1 * $cosU2 * $sinLambda / $sinSigma;
            $cosSqAlpha = 1 - $sinAlpha * $sinAlpha;
            $cos2SigmaM = $cosSigma - 2 * $sinU1 * $sinU2 / $cosSqAlpha;

            if (is_nan($cos2SigmaM)) $cos2SigmaM = 0; // 赤道线情况

            $C = $f / 16 * $cosSqAlpha * (4 + $f * (4 - 3 * $cosSqAlpha));
            $lambdaPrev = $lambda;
            $lambda = $L + (1 - $C) * $f * $sinAlpha *
                ($sigma + $C * $sinSigma *
                    ($cos2SigmaM + $C * $cosSigma * (-1 + 2 * $cos2SigmaM * $cos2SigmaM)));
        } while (abs($lambda - $lambdaPrev) > 1e-12 && --$iterLimit > 0);

        // 迭代未收敛，回退到Haversine公式
        if ($iterLimit == 0) {
            return $this->getDistance($lat1, $lng1, $lat2, $lng2);
        }

        $uSq = $cosSqAlpha * (pow(self::EARTH_RADIUS, 2) - pow($b, 2)) / pow($b, 2);
        $A = 1 + $uSq / 16384 * (4096 + $uSq * (-768 + $uSq * (320 - 175 * $uSq)));
        $B = $uSq / 1024 * (256 + $uSq * (-128 + $uSq * (74 - 47 * $uSq)));
        $deltaSigma = $B * $sinSigma * ($cos2SigmaM + $B / 4 *
                ($cosSigma * (-1 + 2 * $cos2SigmaM * $cos2SigmaM) -
                    $B / 6 * $cos2SigmaM * (-3 + 4 * $sinSigma * $sinSigma) *
                    (-3 + 4 * $cos2SigmaM * $cos2SigmaM)));

        $distance = $b * $A * ($sigma - $deltaSigma);

        return round($distance, 2);
    }

    /**
     * 将距离转换为更友好的显示格式
     *
     * @param float $distance 距离（米）
     * @return string 格式化后的距离
     */
    public function formatDistance(float $distance): string
    {
        if ($distance < 1000) {
            return round($distance, 2) . 'm';
        } else {
            return round($distance / 1000, 2) . 'km';
        }
    }


}

本作品采用《CC 协议》，转载必须注明作者和本文链接