运用米勒投影原理Java实现：public class JWDUtil {private static DecimalFormat format = new DecimalFormat("###0.000");public static List<Double> MillierConvertion(double lon , double lat){ArrayList<Double> list = new ArrayList<>();double L = 6381372 * Math.PI * 2;//地球周长double W = L;// 平面展开后，x轴等于周长double H = L / 2;// y轴约等于周长一半double mill = 2.3;// 米勒投影中的一个常数，范围大约在正负2.3之间double x = lon * Math.PI / 180;// 将经度从度数转换为弧度double y = lat * Math.PI / 180;// 将纬度从度数转换为弧度y = 1.25 * Math.log(Math.tan(0.25 * Math.PI + 0.4 * y));// 米勒投影的转换// 弧度转为实际距离x = (W / 2) + (W / (2 * Math.PI)) * x;y = (H / 2) - (H / (2 * mill)) * y;Double v1 = (double) Math.round(x * 1000) / 1000000; //精确小数后三位，单位KMDouble v2 = (double) Math.round(y * 1000) / 1000000;String format1 = format.format(v1); //科学计数法下转换String format2 = format.format(v2);Double xDouble = Double.valueOf(format1);Double yDouble = Double.valueOf(format2);list.add(xDouble);list.add(yDouble);return list;}//XY坐标系转换成经纬度代码public static List<Double> XYTOJWD(double X,double Y){double L = 6381372 * Math.PI*2;double mill = 2.3;double JD = (X*1000-(L/2 ))*360/L; // 根据X轴计算经度double v = (L / 4 - Y*1000) * mill*2/ (L/2) ;double WD = (Math.atan(Math.pow(Math.E, (v / 1.25))) - (0.25 * Math.PI))/0.4 * 180 / Math.PI;//根据Y轴计算纬度ArrayList<Double> list = new ArrayList<>();list.add(JD);list.add(WD);return list;}}
Python 实现
import numpy as npdef MillierConvertion(lon ,lat):L=6381372*np.pi*2W=LH=L/2mill=2.3x=lon*np.pi/180y=lat*np.pi/180y=1.25*np.log(np.tan(0.25*np.pi+0.4*y))x=(W/2)+(W/(2*np.pi))*xy=(H/2)-(H/(2*mill))*ydou=[x/1000,y/1000]return doudef XYTOJWD(X,Y):L = 6381372 * np.pi * 2mill = 2.3JD = (X * 1000 - (L / 2)) * 360 / Lv = (L / 4 - Y * 1000) * mill * 2 / (L / 2)WD = (np.arctan(np.exp((v / 1.25))) - (0.25 * np.pi))/0.4 * 180 / np.pilist=[JD,WD]return listprint(MillierConvertion(120.7867915, 24.3729099))print(XYTOJWD(33500.41531086378,8133.060038815809))
运行结果
[33500.41531086378, 8133.060038815809]
[120.7867915, 24.372909900000003]
