使用大圆公式计算纬度/经度边界框，弧度和度数不同

问：

我正在开发一个PHP类，以实现一些地理定位功能。为了帮助生成 SQL 查询以在数据库中查找业务，我编写了函数和 .这个想法是，您提供一个点的纬度和纬度，并为地球半径和要搜索的范围提供无单位值（您可以指定英里、公里或弗隆等）。此函数将返回一个包含最小和最大纬度和经度的边界框，您可以在高效查询中使用这些边界框从数据库中检索记录。从某种意义上说，查询具有简单的大于/小于比较，因此可以使用索引，而不是一些导致令人讨厌的表扫描的疯狂大圆距离计算。get_bounding_box_radsget_bounding_box_degrees

我的问题是边界框返回不同的经度（奇怪的是，到目前为止纬度似乎相同），用于以度为单位执行的计算与以弧度执行的计算。我认为这种差异出现在我从这篇维基百科文章中借来的圆距离计算中，其中说：

对于所有距离都准确的公式是具有相等长轴和短轴的椭球体的 Vincenty 公式的以下特例：

公式：

我的计算（参见下面的 PHP 脚本）返回不同的经度值，具体取决于我是使用弧度还是度数执行计算：

弧度min_lat 34.006625915236 max_lat 34.151356084764 min_lon-118.35241157423 max_lon-118.17519642577

度min_lat 34.006625915236 max_lat 34.151356084764 min_lon-118.35117374123 max_lon-118.17643425877

谁能告诉我差异的原因以及如何避免它？如果不是，度或弧度哪个更准确？我想知道这是由于 32 位还是 64 位舍入错误，或者我是否在某处犯了错误。

您应该能够简单地运行以下脚本并了解我的意思

<?php
/**
 * A class to encapsulate and standardize geography-related functions.
 *
 * @author sneakyimp
 */
class OGH_geography {

    /**
     * Radius of the earth in miles
     * @var float
     */
    const EARTH_RADIUS_MI = 3958.8;

    /**
     * outer limit of latitude in radians
     * @var float
     */
    const LAT_LIMIT_RAD = M_PI/2;

    /**
     * outer limit of latitude in degrees
     * @var float
     */
    const LAT_LIMIT_DEGREE = 90;

    /**
     * outer limit of longitude in radians
     * @var float
     */
    const LON_LIMIT_RAD = M_PI;

    /**
     * outer limit of longitude in degrees
     * @var float
     */
    const LON_LIMIT_DEGREE = 180;

    /**
     * Calculates great circle distance between two points on a sphere
     * @see https://en.wikipedia.org/wiki/Great-circle_distance
     * @param float $lat1 Point 1 latitude in rads
     * @param float $lon1 Point 1 longitude in rads
     * @param float $lat2 Point 2 latitude in rads
     * @param float $lon2 Point 2 longitude in rads
     * @param float $radius The unitless radius of the sphere (could be km or miles)
     * @return float unitless distance (depends on units in which radius is specified)
     */
    public static function great_circle_distance_rads($lat1, $lon1, $lat2, $lon2, $radius) : float {
        // broadly, distance = radius * central_angle

        // using the 'special case of the vincenty formula we get
        // central_angle = arctan(sqrt((cos(lat2) * sin(lon1)) ** 2 + ((cos(lat1) * sin(lat2)) - (sin(lat1) * cos(lat2) * cos(abs(lon2 - lon1)))) ** 2) / (sin(lat1) * sin(lat2)) + (cos(lat1) * cos(lat2) * cos(abs(lon2 - lon1))))
        $coslat1 = cos($lat1);
        $coslat2 = cos($lat2);
        $sinlat1 = sin($lat1);
        $sinlat2 = sin($lat2);
        $delta_lon = abs($lon2 - $lon1);
        $sin_delta_lon = sin($delta_lon);
        $cos_delta_lon = cos($delta_lon);

        $is1 = ($coslat2 * $sin_delta_lon) ** 2;
        $is2 = (($coslat1 * $sinlat2) - ($sinlat1 * $coslat2 * $cos_delta_lon)) ** 2;
        $sq_root = sqrt($is1 + $is2);
        $denom = ($sinlat1 * $sinlat2) + ($coslat1 * $coslat2 * $cos_delta_lon);
        //$central_angle_radians = atan($sq_root / $denom);
        $central_angle_radians = atan2($sq_root, $denom);

        return $central_angle_radians * $radius;
    }

    /**
     * Returns distance units per degree or latitude
     * latitude, unlike longitude has constant distance per unit, so doesn't depend on location
     * @param float $radius The unitless radius of the sphere
     * @return float distance/length units per degree
     */
    public static function distance_per_latitude_degree($radius) :float {
        // more subtle approach, use great circle distance calculation for points
        // of identical longitude, 1 deg of latitude apart
        $rads_per_deg = deg2rad(1);
        return self::great_circle_distance_rads(0, 0, $rads_per_deg, 0, $radius);
    }

    /**
     * Returns distance units per radian of latitude
     * latitude, unlike longitude has constant distance per unit, so doesn't depend on location
     * @param float $radius The unitless radius of the sphere
     * @return float distance/length units per radian
     */
    public static function distance_per_latitude_rad($radius) :float {
        // more subtle approach, use great circle distance calculation for points
        // of identical longitude, 1 rad of latitude apart
        return self::great_circle_distance_rads(0, 0, 1, 0, $radius);
    }

    /**
     * Returns distance units per radian of longitude
     * distance longitudinally varies with the latitude, so we need additional param of $lat
     * @param float $lat latitude, in radians, at which we want the longitudinal distance
     * @param float $radius The unitless radius of the sphere
     * @return float distance/length units per radian
     */
    public static function distance_per_longitude_rad($lat, $radius) :float {
        // use great circle distance calculation for points
        // of identical latitude, 1 rad of longitude apart
        return self::great_circle_distance_rads($lat, 0, $lat, 1, $radius);
    }

    /**
     * Returns distance units per degree of longitude
     * distance longitudinally varies with the latitude, so we need additional param of $lat
     * @param float $lat latitude, in degrees, at which we want the longitudinal distance
     * @param float $radius The unitless radius of the sphere
     * @return float distance/length units per degree
     */
    public static function distance_per_longitude_degree($lat, $radius) :float {
        // use great circle distance calculation for points
        // of identical latitude, 1 degree of longitude apart
        $rads_per_deg = deg2rad(1);
        $lat_rads = deg2rad($lat);
        return self::great_circle_distance_rads($lat_rads, 0, $lat_rads, $rads_per_deg, $radius);
    }

    /**
     * Given a point on the globe and a unitless range, will return upper and lower lat/lon in radians for
     * a roughly square bounding box, useful for constructing a query to find nearby points
     * @param float $radius unitless radius of globe
     * @param float $lat latitude of center in rads
     * @param float $lon longitude of center in rads
     * @param float $range unitless range distance from center
     * @return array associative array containing lat/lon floats in radians with keys min_lat, max_lat, min_lon, max_lon
     */
    public static function get_bounding_box_rads($lat, $lon, $radius, $range) {
        if (!is_numeric($lat) || $lat > self::LAT_LIMIT_RAD || $lat < -self::LAT_LIMIT_RAD) {
            throw new \Exception("$lat is not a valid latitude");
        }
        if (!is_numeric($lon) || $lon > self::LON_LIMIT_RAD || $lon < -self::LON_LIMIT_RAD) {
            throw new \Exception("$lon is not a valid longitude");
        }
        if (!is_numeric($radius)) {
            throw new \Exception("$radius is not a valid radius");
        }
        if (!is_numeric($range)) {
            throw new \Exception("$radius is not a valid range");
        }
        // sanity check...this prob shouldn't be used for more than 1/20 of the globe
        $max_range = $radius/20;
        if ($range < 0 || $range > $max_range) {
            throw new \Exception("$range is out of bounds. Must be between 0 and $max_range");
        }


        $lat_dist_per_rad = self::distance_per_latitude_rad($radius);
        $lon_dist_per_rad = self::distance_per_longitude_rad($lat, $radius);

        $lat_adjust = $range/$lat_dist_per_rad;
        $lon_adjust = $range/$lon_dist_per_rad;

        $min_lat = $lat - $lat_adjust;
        // watch out for min_lat out of bounds
        $min_lat = ($min_lat < -self::LAT_LIMIT_RAD) ? -self::LAT_LIMIT_RAD : $min_lat;
        $max_lat = $lat + $lat_adjust;
        // watch out for max_lat out of bounds
        $max_lat = ($max_lat > self::LAT_LIMIT_RAD) ? self::LAT_LIMIT_RAD : $max_lat;

        $min_lon = $lon - $lon_adjust;
        // watch out for min_lon out of bounds
        $min_lon = ($min_lon <= -self::LON_LIMIT_RAD) ? ($min_lon + (2 * M_PI)) : $min_lon;
        $max_lon = $lon + $lon_adjust;
        // watch out for max_lon out of bounds
        $max_lon = ($max_lon > self::LON_LIMIT_RAD) ? ($max_lon - (2 * M_PI)) : $max_lon;

        return [
            'min_lat' => $min_lat,
            'max_lat' => $max_lat,
            'min_lon' => $min_lon,
            'max_lon' => $max_lon
        ];
    }


    /**
     * Given a point on the globe and a unitless range, will return upper and lower lat/lon in degrees for
     * a roughly square bounding box, useful for constructing a query to find nearby points
     * @param float $radius unitless radius of globe
     * @param float $lat latitude of center in degrees
     * @param float $lon longitude of center in degrees
     * @param float $range unitless range distance from center
     * @return array associative array containing lat/lon floats in degrees with keys min_lat, max_lat, min_lon, max_lon
     */
    public static function get_bounding_box_degrees($lat, $lon, $radius, $range) {
        if (!is_numeric($lat) || $lat > self::LAT_LIMIT_DEGREE || $lat < -self::LAT_LIMIT_DEGREE) {
            throw new \Exception("$lat is not a valid latitude");
        }
        if (!is_numeric($lon) || $lon > self::LON_LIMIT_DEGREE || $lon < -self::LON_LIMIT_DEGREE) {
            throw new \Exception("$lon is not a valid longitude");
        }
        if (!is_numeric($radius)) {
            throw new \Exception("$radius is not a valid radius");
        }
        if (!is_numeric($range)) {
            throw new \Exception("$radius is not a valid range");
        }
        // sanity check...this prob shouldn't be used for more than 1/20 of the globe
        $max_range = $radius/20;
        if ($range < 0 || $range > $max_range) {
            throw new \Exception("$range is out of bounds. Must be between 0 and $max_range");
        }


        $lat_dist_per_degree = self::distance_per_latitude_degree($radius);
        $lon_dist_per_degree = self::distance_per_longitude_degree($lat, $radius);


        $lat_adjust = $range/$lat_dist_per_degree;
        $lon_adjust = $range/$lon_dist_per_degree;


        $min_lat = $lat - $lat_adjust;
        // watch out for min_lat out of bounds
        $min_lat = ($min_lat < -self::LAT_LIMIT_DEGREE) ? -self::LAT_LIMIT_DEGREE : $min_lat;
        $max_lat = $lat + $lat_adjust;
        // watch out for max_lat out of bounds
        $max_lat = ($max_lat > self::LAT_LIMIT_DEGREE) ? self::LAT_LIMIT_DEGREE : $max_lat;

        $min_lon = $lon - $lon_adjust;
        // watch out for min_lon out of bounds
        $min_lon = ($min_lon <= -self::LON_LIMIT_DEGREE) ? ($min_lon + 360) : $min_lon;
        $max_lon = $lon + $lon_adjust;
        // watch out for max_lon out of bounds
        $max_lon = ($max_lon > self::LON_LIMIT_DEGREE) ? ($max_lon - 360) : $max_lon;

        return [
                'min_lat' => $min_lat,
                'max_lat' => $max_lat,
                'min_lon' => $min_lon,
                'max_lon' => $max_lon
        ];
    }
} // class OGH_geography


$lat = 34.078991;
$lon = -118.263804;

$lat_rad = deg2rad($lat);
$lon_rad = deg2rad($lon);

$box = OGH_geography::get_bounding_box_rads($lat_rad, $lon_rad, 3958.8, 5);
var_dump($box);
foreach($box as $key => $val) {
    $val = rad2deg($val);
    echo "$key $val\n";
}

$box = OGH_geography::get_bounding_box_degrees($lat, $lon, 3958.8, 5);
var_dump($box);
foreach($box as $key => $val) {
    echo "$key $val\n";
}