圆形经纬度数据_在地理坐标系中创建圆形

最近在ES升级项目(ES5升级到ES6、ES7)中遇到了问题，新版的ES Java API竟然不支持圆形查询了，没办法，只能用正多边形近似代替圆形。当前的数据的坐标系为CGCS2000地理坐标系，所以还不能直接用正余弦计算。目前有两种办法生成圆：

• 先把经纬度中心点高斯正算到平面坐标，用得到的平面坐标作为中心点，根据半径值和各个方向的角度生成圆上的点，最后把圆上的点高斯反算到地理坐标系

• 不经过高斯正反算，直接根据经纬度中心点、方位角、距离进行计算圆上点

个人感觉第二种计算比较直接点，上代码

package com.inspur.geocoding.search.util;/** * Created by 行远必自 on 2020/6/2. */import org.jetbrains.annotations.Contract;import org.jetbrains.annotations.NotNull;import org.locationtech.jts.geom.Coordinate;import java.util.ArrayList;import java.util.List;/** * 项目名:   geocoding * 包名:     com.行远必自.util   * 类名:     SpheroidUtil * 描述:     [一句话描述该类的功能] * 创建人:    行远必自    * 创建时间:  2020/6/2 17:10   * 修改人:    行远必自   * 修改时间:  2020/6/2 17:10   * 修改备注:  [说明本次修改内容]   * 版本:     v1.0   **/public class SpheroidUtil {    @Contract(pure = true)    private static double DEG2RAD(double x){return ((x)* Math.PI/180);}    @Contract(pure = true)    private static double RAD2DEG(double r){ return (180.0 * (r) / Math.PI);}    @Contract(pure = true)    private static double POW2(double x) {return ((x)*(x));}    @NotNull    public static Coordinate CalcEndPoint(Coordinate startPoint, double distance, double azimuth){        // ellipsoid        double a = 6378137.0;        double b = 6356752.314140356;        double f = ( a - b ) / a;        if ( ( ( a < 0 ) && ( b < 0 ) ) ||                ( ( startPoint.x < -180.0 ) || ( startPoint.x > 180.0 ) || ( startPoint.y < -85.05115 ) || ( startPoint.y > 85.05115 ) ) )        {            // latitudes outside these bounds cause the calculations to become unstable and can return invalid results            return new Coordinate( 0, 0 );        }        double radians_lat = DEG2RAD( startPoint.y );        double radians_long = DEG2RAD( startPoint.x );        double b2 = POW2( b ); // spheroid_mu2        double omf = 1 - f;        double tan_u1 = omf * Math.tan( radians_lat );        double u1 =  Math.atan( tan_u1 );        double sigma, last_sigma, delta_sigma, two_sigma_m;        double sigma1, sin_alpha, alpha, cos_alphasq;        double u2, A, B;        double lat2, lambda, lambda2, C, omega;        int i = 0;        if ( azimuth < 0.0 )        {            azimuth = azimuth + Math.PI * 2.0;        }        if ( azimuth > ( Math.PI * 2.0 ) )        {            azimuth = azimuth - Math.PI * 2.0;        }        sigma1 = Math.atan2( tan_u1, Math.cos( azimuth ) );        sin_alpha = Math.cos( u1 ) * Math.sin( azimuth );        alpha = Math.asin( sin_alpha );        cos_alphasq = 1.0 - POW2( sin_alpha );        u2 = POW2( Math.cos( alpha ) ) * ( POW2( a ) - b2 ) / b2; // spheroid_mu2        A = 1.0 + ( u2 / 16384.0 ) * ( 4096.0 + u2 * ( -768.0 + u2 * ( 320.0 - 175.0 * u2 ) ) );        B = ( u2 / 1024.0 ) * ( 256.0 + u2 * ( -128.0 + u2 * ( 74.0 - 47.0 * u2 ) ) );        sigma = ( distance / ( b * A ) );        do        {            two_sigma_m = 2.0 * sigma1 + sigma;            delta_sigma = B * Math.sin( sigma ) * ( Math.cos( two_sigma_m ) + ( B / 4.0 ) * ( Math.cos( sigma ) * ( -1.0 + 2.0 * POW2( Math.cos( two_sigma_m ) ) - ( B / 6.0 ) * Math.cos( two_sigma_m ) * ( -3.0 + 4.0 * POW2( Math.sin( sigma ) ) ) * ( -3.0 + 4.0 * POW2( Math.cos( two_sigma_m ) ) ) ) ) );            last_sigma = sigma;            sigma = ( distance / ( b * A ) ) + delta_sigma;            i++;        }        while ( i < 999 && Math.abs( ( last_sigma - sigma ) / sigma ) > 1.0e-9 );        lat2 = Math.atan2( ( Math.sin( u1 ) * Math.cos( sigma ) + Math.cos( u1 ) * Math.sin( sigma ) *            Math.cos( azimuth ) ), ( omf * Math.sqrt( POW2( sin_alpha ) +            POW2( Math.sin( u1 ) * Math.sin( sigma ) - Math.cos( u1 ) * Math.cos( sigma ) *            Math.cos( azimuth ) ) ) ) );        lambda = Math.atan2( ( Math.sin( sigma ) * Math.sin( azimuth ) ), ( Math.cos( u1 ) * Math.cos( sigma ) -                Math.sin( u1 ) * Math.sin( sigma ) * Math.cos( azimuth ) ) );        C = ( f / 16.0 ) * cos_alphasq * ( 4.0 + f * ( 4.0 - 3.0 * cos_alphasq ) );        omega = lambda - ( 1.0 - C ) * f * sin_alpha * ( sigma + C * Math.sin( sigma ) *            ( Math.cos( two_sigma_m ) + C * Math.cos( sigma ) * ( -1.0 + 2.0 * POW2( Math.cos( two_sigma_m ) ) ) ) );        lambda2 = radians_long + omega;        return new Coordinate( RAD2DEG( lambda2 ), RAD2DEG( lat2 ) );    }    public static List CreateCircle(Coordinate center, double radius, int sides) {        List coordinates = new ArrayList<>();        for (int i = sides; i >= 0; i--) {            double azimuth = ((double) i / (double) sides) * Math.PI * 2.0;            coordinates.add(SpheroidUtil.CalcEndPoint(center, radius , azimuth));        }        return coordinates;    }}

调用下面的语句即可

SpheroidUtil.CreateCircle(center,radius,sides);

经过计算得到的正十二边形如图所示，其中外接圆是用openlayer生成的正六十四边形，可以看出用上述方法生成的正十二边形顶点计算都是正确的。