工具类-根据经纬度计算两点之间在地图上的距离和方向

---

前言

通过经纬度计算地图上两个点的距离和方向，有需要的自取。

工具类

import cn.hutool.core.util.ObjectUtil;

import java.math.BigDecimal;
import java.util.HashMap;
import java.util.Map;

import static java.lang.Math.*;
import static java.lang.Math.atan2;

/**
 * 坐标工具类
 *
 * @author qf
 */
public class CoordinateUtil {
    public static void main(String[] args) {
        double distance = getDistance(116.37047010, 39.83594606, 116.42088584, 39.87221456);
        System.out.println("两者之间的距离为："+distance);
        Address a = new Address();
        a.name = "a";
        a.lon = 116.37047010;
        a.lat = 39.83594606;
        Address b = new Address();
        b.name = "b";
        b.lon = 116.42088584;
        b.lat = 39.87221456;
        Map<String, Object> map = getDistanceAndDirection(a.lon, a.lat, b.lon, b.lat);
        System.out.println("两者之间的距离为："+map.get("distance"));
        System.out.println("b 在 a 的方向为："+map.get("direction"));
    }

    public static class Address{
        private String name;
        private double lon;
        private double lat;
    }

    //WGS84坐标系，长轴6378137.000m，短轴6356752.314，扁率1/298.257223563。
    public static final double f = 1 /  298.257223563;
    public static final double a= 6378137.0;
    public static final double b= 6356752.314245;

    /**
     * 调整给定的角度值 angle，使其保持在0到360度之间
     * @param angle
     * @return
     */
    public static double judge (double angle){
        if (angle < 0) {
            return angle + 360;
        } else {
            return angle;
        }
    }

    /**
     * 获得两个点之间的距离(单位：km)和相对角度
     * @param lon1
     * @param lat1
     * @param lon2
     * @param lat2
     * @return 两个点之间的距离(单位：km)和相对角度
     */
    public static Map<String,Double> getDistanceAndBearing(double lon1, double lat1, double lon2, double lat2) {
        //将经纬度坐标转换为弧度
        lon1 = Math.toRadians(lon1);
        lat1 = Math.toRadians(lat1);
        lon2 = Math.toRadians(lon2);
        lat2 = Math.toRadians(lat2);

        double L = lon2 - lon1;
        // 增加一个微小值到L以避免分母(两点的纬度相同)为0的情况
        if (L == 0) {
            L += 1.7453292588953673E-8; // 防止lon1和lon2相等时出现除以零的情况
        }

        double tanU1 = (1 - f) * tan(lat1);
        double cosU1 = 1 / sqrt((1 + tanU1 * tanU1));
        double sinU1 = tanU1 * cosU1;

        double tanU2 = (1 - f) * tan(lat2);
        double cosU2 = 1 / sqrt((1 + tanU2 * tanU2));
        double sinU2 = tanU2 * cosU2;

        double lambda = L;
        double lambda_ = 0;
        double iterationLimit = 100;

        double cosSq_alpha = 0;

        double sin_delta = 0;
        double cos2_deltaM = 0;
        double cos_delta = 0;
        double delta = 0;
        double sinlambda = 0;
        double coslambda = 0;
        double s = 0;

        while (abs(lambda - lambda_) > 1e-12 && iterationLimit > 0) {
            iterationLimit = iterationLimit - 1;
            sinlambda = sin(lambda);
            coslambda = cos(lambda);
            double sinSq_delta = (cosU2 * sinlambda) * (cosU2 * sinlambda) + (cosU1 * sinU2 - sinU1 * cosU2 * coslambda) * (cosU1 * sinU2 - sinU1 * cosU2 * coslambda);
            sin_delta = sqrt(sinSq_delta);
            if (sin_delta == 0) {

            }
            cos_delta = sinU1 * sinU2 + cosU1 * cosU2 * coslambda;
            delta = atan2(sin_delta, cos_delta);
            double sin_alpha = cosU1 * cosU2 * sinlambda / sin_delta;
            cosSq_alpha = 1 - sin_alpha * sin_alpha;
            cos2_deltaM = cos_delta - 2 * sinU1 * sinU2 / cosSq_alpha;
            double C = f / 16 * cosSq_alpha * (4 + f * (4 - 3 * cosSq_alpha));
            lambda_ = lambda;
            lambda = L + (1 - C) * f * sin_alpha * (delta + C * sin_delta * (cos2_deltaM + C * cos_delta * (-1 + 2 * cos2_deltaM * cos2_deltaM)));
        }

        double uSq = cosSq_alpha * (a * a - b * b) / (b * b);
        double A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)));
        double B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)));
        double delta_delta = B * sin_delta * (cos2_deltaM + B / 4 * (cos_delta * (-1 + 2 * cos2_deltaM * cos2_deltaM) - B / 6 * cos2_deltaM * (-3 + 4 * sin_delta * sin_delta) * (-3 + 4 * cos2_deltaM * cos2_deltaM)));
        s = b * A * (delta - delta_delta);                                        //距离m
        double fwdAz = atan2(cosU2 * sinlambda, cosU1 * sinU2 - sinU1 * cosU2 * coslambda); //初始方位角
        double revAz = atan2(cosU1 * sinlambda, -sinU1 * cosU2 + cosU1 * sinU2 * coslambda); //最终方位角
        double fwdAz_ = judge(Math.toDegrees(fwdAz));  //将角度归一化至0 - 360°内
        double revAz_ = judge(Math.toDegrees(revAz));  //将角度归一化至0 - 360°内
        s = s / 1000;  //距离km

        HashMap map = new HashMap();
        //初始方位角
        map.put("initialAzimuth",fwdAz_);
        //最终方位角
        map.put("finalAzimuth",revAz_);
        //两点相距 距离km
        map.put("distance",s);

        return map;
    }

    /**
     * 获得两个点之间的距离(单位：m)
     * @param lon1
     * @param lat1
     * @param lon2
     * @param lat2
     * @return
     */
    public static double getDistance(double lon1, double lat1, double lon2, double lat2) {
        //将经纬度坐标转换为弧度
        lon1 = Math.toRadians(lon1);
        lat1 = Math.toRadians(lat1);
        lon2 = Math.toRadians(lon2);
        lat2 = Math.toRadians(lat2);

        // 检查经度是否相等（防止lon1和lon2相等时出现除以零的情况）
        double L = lon2 - lon1;
        // 增加一个微小值到L以避免分母为0的情况
        if (L == 0) {
            L += 1.7453292588953673E-8; // 防止lon1和lon2相等时出现除以零的情况
        }

        double tanU1 = (1 - f) * tan(lat1);
        double cosU1 = 1 / sqrt((1 + tanU1 * tanU1));
        double sinU1 = tanU1 * cosU1;

        double tanU2 = (1 - f) * tan(lat2);
        double cosU2 = 1 / sqrt((1 + tanU2 * tanU2));
        double sinU2 = tanU2 * cosU2;

        double lambda = L;
        double lambda_ = 0;
        double iterationLimit = 100;

        double cosSq_alpha = 0;

        double sin_delta = 0;
        double cos2_deltaM = 0;
        double cos_delta = 0;
        double delta = 0;
        double sinlambda = 0;
        double coslambda = 0;
        double s = 0;

        while (abs(lambda - lambda_) > 1e-12 && iterationLimit > 0) {
            iterationLimit = iterationLimit - 1;
            sinlambda = sin(lambda);
            coslambda = cos(lambda);
            double sinSq_delta = (cosU2 * sinlambda) * (cosU2 * sinlambda) + (cosU1 * sinU2 - sinU1 * cosU2 * coslambda) * (cosU1 * sinU2 - sinU1 * cosU2 * coslambda);
            sin_delta = sqrt(sinSq_delta);
            if (sin_delta == 0) {

            }
            cos_delta = sinU1 * sinU2 + cosU1 * cosU2 * coslambda;
            delta = atan2(sin_delta, cos_delta);
            double sin_alpha = cosU1 * cosU2 * sinlambda / sin_delta;
            cosSq_alpha = 1 - sin_alpha * sin_alpha;
            cos2_deltaM = cos_delta - 2 * sinU1 * sinU2 / cosSq_alpha;
            double C = f / 16 * cosSq_alpha * (4 + f * (4 - 3 * cosSq_alpha));
            lambda_ = lambda;
            lambda = L + (1 - C) * f * sin_alpha * (delta + C * sin_delta * (cos2_deltaM + C * cos_delta * (-1 + 2 * cos2_deltaM * cos2_deltaM)));
        }

        double uSq = cosSq_alpha * (a * a - b * b) / (b * b);
        double A = 1 + uSq / 16384 * (4096 + uSq * (-768 + uSq * (320 - 175 * uSq)));
        double B = uSq / 1024 * (256 + uSq * (-128 + uSq * (74 - 47 * uSq)));
        double delta_delta = B * sin_delta * (cos2_deltaM + B / 4 * (cos_delta * (-1 + 2 * cos2_deltaM * cos2_deltaM) - B / 6 * cos2_deltaM * (-3 + 4 * sin_delta * sin_delta) * (-3 + 4 * cos2_deltaM * cos2_deltaM)));
        s = b * A * (delta - delta_delta);                                        //距离m
        double fwdAz = atan2(cosU2 * sinlambda, cosU1 * sinU2 - sinU1 * cosU2 * coslambda); //初始方位角
        double revAz = atan2(cosU1 * sinlambda, -sinU1 * cosU2 + cosU1 * sinU2 * coslambda); //最终方位角
        double fwdAz_ = judge(Math.toDegrees(fwdAz));  //将角度归一化至0 - 360°内
        double revAz_ = judge(Math.toDegrees(revAz));  //将角度归一化至0 - 360°内
        return s;
    }

    /**
     * 获得两个点之间的距离(单位：m)和方向(2在1的相对方向)
     * @param lon1
     * @param lat1
     * @param lon2
     * @param lat2
     * @return
     */
    public static Map<String,Object> getDistanceAndDirection(double lon1, double lat1, double lon2, double lat2) {
        Map<String, Double> map = getDistanceAndBearing(lon1, lat1, lon2, lat2);
        Double finalAzimuth = map.get("finalAzimuth");
        BigDecimal bigDecimal = new BigDecimal(finalAzimuth).setScale(0, BigDecimal.ROUND_HALF_UP);
        String direction = getNameByDirection(bigDecimal.intValue());
        Map<String, Object> mapFinal = new HashMap<>();
        mapFinal.put("finalAzimuth",finalAzimuth);//最终方位角
        mapFinal.put("direction",direction);//方向
        mapFinal.put("distance",map.get("distance") * 1000);//距离
        return mapFinal;
    }

    /**
     * 相对角度转换为方向
     * @param direction
     * @return
     */
    private static String getNameByDirection(Integer direction) {
        if(direction == null){
            return null;
        }
        if((direction >= 338 && direction <= 360) || direction >= 0 && direction <= 22){
            return "北";
        } else if (direction >= 23 && direction <= 67) {
            return "东北";
        }else if (direction >= 68 && direction <= 112) {
            return "东";
        }else if (direction >= 113 && direction <= 157) {
            return "东南";
        } else if (direction >= 158 && direction <= 202) {
            return "南";
        }else if (direction >= 203 && direction <= 247) {
            return "西南";
        }else if (direction >= 248 && direction <= 292) {
            return "西";
        }else if (direction >= 293 && direction <= 337) {
            return "西北";
        }
        return "风向异常";
    }


}

---