需求

  1. 已知两点经纬度,求两点间的距离
  2. 已知两点经纬度,求两点间的方位角
  3. 已知一点经纬度为中心、r为半径,求出该经纬度区间
  4. 已知一点的经纬度、距离,方位角,求另一个点的经纬度

公式

  1. 已知两点经纬度,求两点间的距离:
    AB=R*arccos[sin(wA)sin(wB)+cos(wA)cos(wB)*cos(jA-jB)]
  2. 已知两点经纬度,求两点间的方位角
    r*arccos[cos(y1)cos(y2)cos(x1-x2)+sin(y1)*sin(y2)]
    这里得到的是弧度,我们需要将弧度转换成我们需要的方位角,这里需要考虑象限问题
    纬度是横轴,经度时纵轴
    B点在第一象限,Azimuth=A;
    B点在第四象限,Azimuth=360+A;
    B点在第二三象限,Azimuth=180-A
  3. 已知一点经纬度为中心、r为半径,求出该经纬度区间
  4. 已知一点的经纬度、距离,方位角,求另一个点的经纬度

实现代码

import java.util.HashMap;
import java.util.Map;

public class DistanceUtil {     

    /*
    * 大地坐标系资料WGS-84 长半径a=6378137 短半径b=6356752.3142 扁率f=1/298.2572236
    */
    /** 长半径a=6378137 */
    private static double EARTH_RADIUS = 6378137;
    /** 短半径b=6356752.3142 */
    private static double b = 6356752.3142;
    /** 扁率f=1/298.2572236 */
    private static double f = 1 / 298.2572236;
    
    private static double rad(double d) {
           return d * Math.PI / 180.0;
    }
    
    /**
     * 弧度换成度
     *
     * @param x
     *弧度
     * @return 度
     */
    private static double deg(double x) {
         return x * 180 / Math.PI;
    }
    
    /**
     * 根据两个位置的经纬度,来计算两地的距离(单位为KM)
     * 参数为String类型
     * @param lat1 用户经度
     * @param lon1 用户纬度
     * @param lat2 商家经度
     * @param lon2 商家纬度
     * @return
     */
    public static double getDistance(int lat1, int lon1, int lat2, int lon2) {
        double tmplat1 = lat1/1000000.0;
        double tmplon1 = lon1/1000000.0;
        double tmplat2 = lat2/1000000.0;
        double tmplon2 = lon2/1000000.0;;
        
        double radLat1 = rad(tmplat1);
        double radLat2 = rad(tmplat2);
        double difference = radLat1 - radLat2;
        double mdifference = rad(tmplon1) - rad(tmplon2);
        double distance = 2 * Math.asin(Math.sqrt(Math.pow(Math.sin(difference / 2), 2)
        + Math.cos(radLat1) * Math.cos(radLat2)
        * Math.pow(Math.sin(mdifference / 2), 2)));
        distance = distance * EARTH_RADIUS;
        distance = Math.round(distance * 10000) / 10000;
        
        return distance;
    }
    
    /**
     * 获取当前用户一定距离以内的经纬度值,距离单位是米
     * 单位米 return minLat
     * 最小经度 minLon
     * 最小纬度 maxLat
     * 最大经度 maxLon
     * 最大纬度 minLat
     */
    public static Map<String,Integer> getAround(int lat, int lon, double raidus) {
        Map<String,Integer> map = new HashMap();
        
        double latitude = lat/1000000.0;// 传值给经度
        double longitude = lon/1000000.0;// 传值给纬度
        
        Double degree = (24901 * 1609) / 360.0; // 获取每度
        double raidusMile = raidus;
        
        Double mpdLon = Double.parseDouble((degree * Math.cos(latitude * (Math.PI / 180))+"").replace("-", ""));
        Double dpmLon = 1 / mpdLon;
        Double radiusLon = dpmLon * raidusMile;
        //获取最小经度
        double minLat = latitude - radiusLon;
        // 获取最大经度
        double maxLat = latitude + radiusLon;
        
        Double dpmLat = 1 / degree;
        Double radiusLat = dpmLat * raidusMile;
        // 获取最小纬度
        double minLon = longitude - radiusLat;
        // 获取最大纬度
        double maxLon = longitude + radiusLat;
        
        map.put("minLat", new Double(minLat*1000000).intValue());
        map.put("maxLat", new Double(maxLat*1000000).intValue());
        map.put("minLon", new Double(minLon*1000000).intValue());
        map.put("maxLon", new Double(maxLon*1000000).intValue());
        
        return map;
    }
    
    /**
     * 已知两个点的坐标,求两个的方位角
     * @param lat1 起始点
     * @param lon1
     * @param lat2 结束点
     * @param lon2
     * @return
     */
    public static double getAngle(int lat1, int lon1, int lat2, int lon2) {
        double tmplat1 = lat1/1000000.0;
        double tmplon1 = lon1/1000000.0;
        double tmplat2 = lat2/1000000.0;
        double tmplon2 = lon2/1000000.0;
        double dRotateAngle = Math.atan2(Math.abs(tmplat1 - tmplat2), Math.abs(tmplon1 - tmplon2));
        
        
        /*
        纬度是横轴,经度时纵轴
         */
        if (tmplat2 >= tmplat1){
            if (tmplon2 >=tmplon1){
                  //第一象限
            } else {
                  // 第四象限
                  dRotateAngle =2 * Math.PI - dRotateAngle;
            }
        } else{
            if (tmplon2 >= tmplon1) {
                 //第二象限
                  dRotateAngle = Math.PI - dRotateAngle;
            } else {
                 //第三象限
                 dRotateAngle = Math.PI + dRotateAngle;
            }
        }
        
        dRotateAngle = dRotateAngle * 180 / Math.PI;
        return dRotateAngle;
    }
    
    /**
     * 求两经纬度方向角
     *
     *@param lat1
     *第一点的纬度
     * @param lon1
     *第一点的经度
     * @param lat2
     *第二点的纬度
      @param lon2
     *第二点的经度
     * @return 方位角,角度(单位:°)
     * */
    public static double getAzimuth(int lat1,int lon1,int lat2,int lon2) {
        double tmplat1 = lat1/1000000.0;
        double tmplon1 = lon1/1000000.0;
        double tmplat2 = lat2/1000000.0;
        double tmplon2 = lon2/1000000.0;
        tmplat1 = rad(tmplat1);
        tmplat2 = rad(tmplat2);
        tmplon1 = rad(tmplon1);
        tmplon2 = rad(tmplon2);
        double azimuth = Math.sin(tmplat1) * Math.sin(tmplat2) + Math.cos(tmplat1)
        * Math.cos(tmplat2) * Math.cos(tmplon2 - tmplon1);
        azimuth = Math.sqrt(1 - azimuth * azimuth);
        azimuth = Math.cos(tmplat2) * Math.sin(tmplon2 - tmplon1) / azimuth;
        azimuth = Math.asin(azimuth) * 180 / Math.PI;
        
        if (Double.isNaN(azimuth)) {
            if (lon1 < lon2) {
                azimuth = 90.0;
            } else {
                azimuth = 270.0;
            }
        }
        
        /*
        纬度是横轴,经度时纵轴
        B点在第一象限,Azimuth=A;
        B在第四象限,Azimuth=360+A;
        B在第二三象限,Azimuth=180-A
         */
        
        if (tmplat2 >= tmplat1){
            if (tmplon2 >=tmplon1){
                //第一象限
            } else {
                // 第四象限
                azimuth =360 + azimuth;
            }
        } else{
            if (tmplon2 >= tmplon1) {
                 //第二象限
                azimuth = 180 - azimuth;
            } else {
                //第三象限
                 azimuth = 180 - azimuth;
            }
        }
        
        return azimuth;
    }
    
    /**
     * 计算另一点经纬度
     *
     * @param lon
     *经度
     * @param lat
     *维度
     * @param brng
     *方位角
     * @param dist
     *距离(米)
     */
    public static Map<String,Integer> computerThatLonLat(int lat,int lon, double brng, double dist) {
    
        double tmplon = lon / 1000000.0;
        double tmplat = lat / 1000000.0;
        
        double alpha1 = rad(brng);
        double sinAlpha1 = Math.sin(alpha1);
        double cosAlpha1 = Math.cos(alpha1);
        
        double tanU1 = (1 - f) * Math.tan(rad(tmplat));
        double cosU1 = 1 / Math.sqrt((1 + tanU1 * tanU1));
        double sinU1 = tanU1 * cosU1;
        double sigma1 = Math.atan2(tanU1, cosAlpha1);
        double sinAlpha = cosU1 * sinAlpha1;
        double cosSqAlpha = 1 - sinAlpha * sinAlpha;
        double uSq = cosSqAlpha * (EARTH_RADIUS * EARTH_RADIUS - 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 cos2SigmaM=0;
        double sinSigma=0;
        double cosSigma=0;
        double sigma = dist / (b * A), sigmaP = 2 * Math.PI;
        while (Math.abs(sigma - sigmaP) > 1e-12) {
        cos2SigmaM = Math.cos(2 * sigma1 + sigma);
        sinSigma = Math.sin(sigma);
        cosSigma = Math.cos(sigma);
        double deltaSigma = B * sinSigma * (cos2SigmaM + B / 4 * (cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)
        - B / 6 * cos2SigmaM * (-3 + 4 * sinSigma * sinSigma) * (-3 + 4 * cos2SigmaM * cos2SigmaM)));
        sigmaP = sigma;
        sigma = dist / (b * A) + deltaSigma;
        }
        
        double tmp = sinU1 * sinSigma - cosU1 * cosSigma * cosAlpha1;
        double lat2 = Math.atan2(sinU1 * cosSigma + cosU1 * sinSigma * cosAlpha1,
        (1 - f) * Math.sqrt(sinAlpha * sinAlpha + tmp * tmp));
        double lambda = Math.atan2(sinSigma * sinAlpha1, cosU1 * cosSigma - sinU1 * sinSigma * cosAlpha1);
        double C = f / 16 * cosSqAlpha * (4 + f * (4 - 3 * cosSqAlpha));
        double L = lambda - (1 - C) * f * sinAlpha
        * (sigma + C * sinSigma * (cos2SigmaM + C * cosSigma * (-1 + 2 * cos2SigmaM * cos2SigmaM)));
        
        double revAz = Math.atan2(sinAlpha, -tmp); // final bearing
        
        Map<String,Integer> map = new HashMap<String, Integer>();
        map.put("lon",+new Double((tmplon+deg(L))*1000000).intValue());
        map.put("lat",+new Double(deg(lat2)*1000000).intValue());
        return map;
    }
    
    /**
     * 判断是否漂移
     * @return
     */
    /**
     * 判断是否漂移,如果漂移则返回true,否则返回false
     * @param lat1 单位是10^-6
     * @param lon1 单位是10^-6
     * @param raidus 半径单位是米
     * @param lat2
     * @param lon2
     * @return
     */
    public static boolean isDrift(int lat1,int lon1,int raidus,int lat2,int lon2){
    
        //中国纬度范围:4°~53°
        if(lat1<4000000 || lat1 >53000000 ){
        return true;
        }
        //中国经度范围:73°~135°
        if(lon1<73000000 || lon1 >135000000){
        return true;
        }
        
        Map<String,Integer> map = getAround(lat1,lon1,raidus);
        
        if(lat2>map.get("maxLat") || lat2<map.get("minLat")){
            return true;
        }
        if(lon2>map.get("maxLon") || lon2<map.get("minLon")){
            return true;
        }

        return false;
    }

}