最全的坐标系转换--经验证均可放心使用

  • 北京54/国家80坐标(x,y)转经纬度(大地坐标(B,L))
  • 大地坐标转换北京54/国家80坐标
  • 大地坐标(B,L)转火星坐标(GCj02)转百度坐标(BD-09)
  • 方法一:
  • 大地转火星
  • 火星转百度
  • 方法二:
  • 大地转火星
  • 火星转百度
  • 火星转大地
  • 百度转火星
  • 例子:


北京54/国家80坐标(x,y)转经纬度(大地坐标(B,L))

// 高斯反算:由54/80坐标反算成经纬度
    //, double *longitude, double *latitude)
    // X实质是测量坐标的Y,除带号,整数共6位,Y是y+500000的结果
   	//Y实质是测量坐标的x
    public static double[] GaussToBL(double X, double Y) {
        int ProjNo; int ZoneWide; 带宽
        double[] output = new double[2];
        double longitude1,latitude1, longitude0, X0,Y0, xval,yval;//latitude0,
        double e1,e2,f,a, ee, NN, T,C, M, D,R,u,fai, iPI;
        iPI = 0.0174532925199433; 3.1415926535898/180.0;
        //a = 6378245.0; f = 1.0/298.3; //54年北京坐标系参数
        a=6378140.0; f=1/298.257; //80年西安坐标系参数
        ZoneWide = 6; 6度带宽
        ProjNo = (int)(X/1000000L) ; //查找带号
        longitude0 = (ProjNo-1) * ZoneWide + ZoneWide / 2;
        //longitude0 = longitude0 * iPI ; //中央经线
        longitude0 = 117 * iPI ; //中央经线


        X0 = ProjNo*1000000L+500000L;
        Y0 = 0;
        xval = X-X0; yval = Y-Y0; //带内大地坐标
        e2 = 2*f-f*f;
        e1 = (1.0-Math.sqrt(1-e2))/(1.0+Math.sqrt(1-e2));
        ee = e2/(1-e2);
        M = yval;
        u = M/(a*(1-e2/4-3*e2*e2/64-5*e2*e2*e2/256));
        fai = u+(3*e1/2-27*e1*e1*e1/32)*Math.sin(2*u)+(21*e1*e1/16-55*e1*e1*e1*e1/32)*Math.sin(
                4*u)
                +(151*e1*e1*e1/96)*Math.sin(6*u)+(1097*e1*e1*e1*e1/512)*Math.sin(8*u);
        C = ee*Math.cos(fai)*Math.cos(fai);
        T = Math.tan(fai)*Math.tan(fai);
        NN = a/Math.sqrt(1.0-e2*Math.sin(fai)*Math.sin(fai));
        R = a*(1-e2)/Math.sqrt((1-e2*Math.sin(fai)*Math.sin(fai))*(1-e2*Math.sin(fai)*Math.sin(fai))*(1-e2*Math.sin
                (fai)*Math.sin(fai)));
        D = xval/NN;
        //计算经度(Longitude) 纬度(Latitude)
        longitude1 = longitude0+(D-(1+2*T+C)*D*D*D/6+(5-2*C+28*T-3*C*C+8*ee+24*T*T)*D
                *D*D*D*D/120)/Math.cos(fai);
        latitude1 = fai -(NN*Math.tan(fai)/R)*(D*D/2-(5+3*T+10*C-4*C*C-9*ee)*D*D*D*D/24
                +(61+90*T+298*C+45*T*T-256*ee-3*C*C)*D*D*D*D*D*D/720);
        //转换为度 DD
        output[0] = longitude1 / iPI;
        output[1] = latitude1 / iPI;
        return output;
        //*longitude = longitude1 / iPI;
        //*latitude = latitude1 / iPI;
    }

大地坐标转换北京54/国家80坐标

// 高斯正算:由经纬度反算成高斯投影坐标
 	// ProjNo:带号
 	// 最后计算的x是测量中的y,并且格式为带号+(y值+500000);
 	// 最后计算的y是测量中的x
    public void GaussToBLToGauss(double longitude, double latitude) {
        int ProjNo=20; int ZoneWide; 带宽
        double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval;
        double a,f, e2,ee, NN, T,C,A, M, iPI;
        iPI = 0.0174532925199433; 3.1415926535898/180.0;
        ZoneWide = 6; 6度带宽
        //a=6378245.0; f=1.0/298.3; //54年北京坐标系参数
        a=6378140.0; f=1/298.257; //80年西安坐标系参数
        ProjNo = (int)(longitude / ZoneWide) ;
        longitude0 = ProjNo * ZoneWide + ZoneWide / 2;
        longitude0 = longitude0 * iPI ;
        latitude0 = 0;
        System.out.println(latitude0);
        longitude1 = longitude * iPI ; //经度转换为弧度
        latitude1 = latitude * iPI ; //纬度转换为弧度
        e2=2*f-f*f;
        ee=e2*(1.0-e2);
        NN=a/Math.sqrt(1.0-e2*Math.sin(latitude1)*Math.sin(latitude1));
        T=Math.tan(latitude1)*Math.tan(latitude1);
        C=ee*Math.cos(latitude1)*Math.cos(latitude1);
        A=(longitude1-longitude0)*Math.cos(latitude1);
        M=a*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude1-(3*e2/8+3*e2*e2/32+45*e2*e2
                *e2/1024)*Math.sin(2*latitude1)
                +(15*e2*e2/256+45*e2*e2*e2/1024)*Math.sin(4*latitude1)-(35*e2*e2*e2/3072)*Math.sin(6*latitude1));
        xval = NN*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120);
        yval = M+NN*Math.tan(latitude1)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24
                +(61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720);
        X0 = 1000000L*(ProjNo+1)+500000L;
        Y0 = 0;
        xval = xval+X0; yval = yval+Y0;
        //*X = xval;
        //*Y = yval;
        System.out.println("x:"+xval);
        System.out.println("y:"+yval);

    }

大地坐标(B,L)转火星坐标(GCj02)转百度坐标(BD-09)

方法一:

大地转火星
火星转百度
//方法一
    //坐标转换
    /**
     * WGS84转GCj02
     * GCj02转BD-09
     * 火星坐标系 (GCJ-02) 与百度坐标系 (BD-09) 的转换
     * @param lng
     * @param lat
     * @returns {*[]}
     */
    public static void wgs84togcj02tobd09(double lng,double lat) {
        double xPI = 3.14159265358979324 * 3000.0 / 180.0;
        double PI = 3.1415926535897932384626;
        double a = 6378245.0;
        double ee = 0.00669342162296594323;
        // WGS84转GCj02
        double dlat = transformlat(lng - 105.0, lat - 35.0);
        double dlng = transformlng(lng - 105.0, lat - 35.0);
        double radlat = lat / 180.0 * PI;
        double magic = Math.sin(radlat);
        magic = 1 - ee * magic * magic;
        double sqrtmagic = Math.sqrt(magic);
        dlat = (dlat * 180.0) / ((a * (1 - ee)) / (magic * sqrtmagic) * PI);
        dlng = (dlng * 180.0) / (a / sqrtmagic * Math.cos(radlat) * PI);
        double mglat = lat + dlat;
        double mglng = lng + dlng;
        System.out.println("火星坐标: lng:"+mglng+",lat:"+mglat);
        // 火星坐标系 (GCJ-02) 与百度坐标系 (BD-09) 的转换
        double z = Math.sqrt(mglng * mglng + mglat * mglat) + 0.00002 * Math.sin(mglat * xPI);
        double theta = Math.atan2(mglat, mglng) + 0.000003 * Math.cos(mglng * xPI);
        double bdlng = z * Math.cos(theta) + 0.0065;
        double bdlat = z * Math.sin(theta) + 0.006;
        // return [bdlng, bdlat]
        //return {lng: bdlng, lat: bdlat}
        System.out.println("百度坐标:  lng:"+bdlng+",lat:"+bdlat);
    }

    public static Double transformlat(Double lng,Double lat){
        double PI = 3.1415926535897932384626;
        double ret = -100.0 + 2.0 * lng + 3.0 * lat + 0.2 * lat * lat + 0.1 * lng * lat + 0.2 * Math.sqrt(Math.abs(lng));
        ret += (20.0 * Math.sin(6.0 * lng * PI) + 20.0 * Math.sin(2.0 * lng * PI)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(lat * PI) + 40.0 * Math.sin(lat / 3.0 * PI)) * 2.0 / 3.0;
        ret += (160.0 * Math.sin(lat / 12.0 * PI) + 320 * Math.sin(lat * PI / 30.0)) * 2.0 / 3.0;
        return ret;
    }

    public static Double transformlng(Double lng,Double lat){
        double PI = 3.1415926535897932384626;
        double ret = 300.0 + lng + 2.0 * lat + 0.1 * lng * lng + 0.1 * lng * lat + 0.1 * Math.sqrt(Math.abs(lng));
        ret += (20.0 * Math.sin(6.0 * lng * PI) + 20.0 * Math.sin(2.0 * lng * PI)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(lng * PI) + 40.0 * Math.sin(lng / 3.0 * PI)) * 2.0 / 3.0;
        ret += (150.0 * Math.sin(lng / 12.0 * PI) + 300.0 * Math.sin(lng / 30.0 * PI)) * 2.0 / 3.0;
        return ret;
    }

方法二:

大地转火星
火星转百度
火星转大地
//方法二
    private final static double a = 6378245.0;
    private final static double pi = 3.1415926535897932384626;
    private final static double ee = 0.00669342162296594323;

    /**
     * Description:大地转火星 WGS-84 to GCJ-02 <BR>
     * @param latitude 纬度
     * @param longitude 经度
     * @return [纬度,经度]
     */
    public static double[] toGCJ02Point(double longitude,double latitude) {
        double[] dev = calDev(latitude, longitude);
        double retLat = latitude + dev[0];
        double retLon = longitude + dev[1];
        System.out.println("火星坐标:  lng:"+retLon+",lat:"+retLat);
        GCJ02_To_BD09(retLon,retLat);
        return new double[] { retLat, retLon };
    }
     /**
     * 火星转百度 GCJ-02 to BD-09
     */
    private static Double x_pi = 3.14159265358979324 * 3000.0 / 180.0;  //圆周率转换量
    public static void  GCJ02_To_BD09(Double GCJ02_Lon,Double GCJ02_Lat){
        double x = GCJ02_Lon, y = GCJ02_Lat;
        double z = Math.sqrt(x * x + y * y) + 0.00002 * Math.sin(y * x_pi);
        double theta = Math.atan2(y, x) + 0.000003 * Math.cos(x * x_pi);
        double BD_09_Lon = z * Math.cos(theta) + 0.0065;
        double BD_09_Lat = z * Math.sin(theta) + 0.006;
        System.out.println("百度坐标:  lng:"+BD_09_Lon+",lat:"+BD_09_Lat);
        //return [BD_09_Lon, BD_09_Lat];
    }
    /**
     * Description:火星转大地 GCJ-02 to WGS-84 <BR>
     * @param latitude 纬度
     * @param longitude 经度
     * @return [纬度,经度]
     */
    public static double[] toWGS84Point(double latitude, double longitude) {
        double[] dev = calDev(latitude, longitude);
        double retLat = latitude - dev[0];
        double retLon = longitude - dev[1];
        dev = calDev(retLat, retLon);
        retLat = latitude - dev[0];
        retLon = longitude - dev[1];
        return new double[] { retLat, retLon };
    }
	
    private static double[] calDev(double wgLat, double wgLon) {
        if (isOutOfChina(wgLat, wgLon)) {
            return new double[] { 0, 0 };
        }
        double dLat = calLat(wgLon - 105.0, wgLat - 35.0);
        double dLon = calLon(wgLon - 105.0, wgLat - 35.0);
        double radLat = wgLat / 180.0 * pi;
        double magic = Math.sin(radLat);
        magic = 1 - ee * magic * magic;
        double sqrtMagic = Math.sqrt(magic);
        dLat = (dLat * 180.0) / ((a * (1 - ee)) / (magic * sqrtMagic) * pi);
        dLon = (dLon * 180.0) / (a / sqrtMagic * Math.cos(radLat) * pi);
        return new double[] { dLat, dLon };
    }

    private static boolean isOutOfChina(double lat, double lon) {
        if (lon < 72.004 || lon > 137.8347)
            return true;
        if (lat < 0.8293 || lat > 55.8271)
            return true;
        return false;
    }

    private static double calLat(double x, double y) {
        double ret = -100.0 + 2.0 * x + 3.0 * y + 0.2 * y * y + 0.1 * x * y + 0.2 * Math.sqrt(Math.abs(x));
        ret += (20.0 * Math.sin(6.0 * x * pi) + 20.0 * Math.sin(2.0 * x * pi)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(y * pi) + 40.0 * Math.sin(y / 3.0 * pi)) * 2.0 / 3.0;
        ret += (160.0 * Math.sin(y / 12.0 * pi) + 320 * Math.sin(y * pi / 30.0)) * 2.0 / 3.0;
        return ret;
    }

    private static double calLon(double x, double y) {
        double ret = 300.0 + x + 2.0 * y + 0.1 * x * x + 0.1 * x * y + 0.1 * Math.sqrt(Math.abs(x));
        ret += (20.0 * Math.sin(6.0 * x * pi) + 20.0 * Math.sin(2.0 * x * pi)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(x * pi) + 40.0 * Math.sin(x / 3.0 * pi)) * 2.0 / 3.0;
        ret += (150.0 * Math.sin(x / 12.0 * pi) + 300.0 * Math.sin(x / 30.0 * pi)) * 2.0 / 3.0;
        return ret;
    }
百度转火星
//定义一些常量
	var x_PI = 3.14159265358979324 * 3000.0 / 180.0;
	var PI = 3.1415926535897932384626;
	var a = 6378245.0;
	var ee = 0.00669342162296594323;
	
	/**
	 * 百度坐标系 (BD-09) 与 火星坐标系 (GCJ-02)的转换
	 * 即 百度 转 谷歌、高德
	 * @param bd_lon
	 * @param bd_lat
	 * @returns {*[]}
	 */
	function bd09togcj02(bd_lon, bd_lat) {
	    var x_pi = 3.14159265358979324 * 3000.0 / 180.0;
	    var x = bd_lon - 0.0065;
	    var y = bd_lat - 0.006;
	    var z = Math.sqrt(x * x + y * y) - 0.00002 * Math.sin(y * x_pi);
	    var theta = Math.atan2(y, x) - 0.000003 * Math.cos(x * x_pi);
	    var gg_lng = z * Math.cos(theta);
	    var gg_lat = z * Math.sin(theta);
	    return [gg_lng, gg_lat]
	}

例子:

public class Demo {

    public static void main(String[] args) {
        //80坐标转大地坐标:
        // X实质是测量坐标的Y,除带号,整数共6位,Y是y+500000的结果
        //Y实质是测量坐标的x
        double[] doubles = GaussToBL(489822.82000000036,4241542.691);
        System.out.println("大地坐标:"+doubles[0] +","+ doubles[1]);
        //大地转80(转换结果:2.048982282000024E7,4241542.69101974)
        GaussToBLToGauss(116.88364353980904,38.306634937035255);

        //大地转火星转百度
        //方法一
        wgs84togcj02tobd09(116.95562589104205,38.2883423673599);
        //方法二:
        toGCJ02Point(116.95562589104205,38.2883423673599);
    }

    //方法一
    //坐标转换
    /**
     * WGS84转GCj02
     * GCj02转BD-09
     * 火星坐标系 (GCJ-02) 与百度坐标系 (BD-09) 的转换
     * @param lng
     * @param lat
     * @returns {*[]}
     */
    public static void wgs84togcj02tobd09(double lng,double lat) {
        double xPI = 3.14159265358979324 * 3000.0 / 180.0;
        double PI = 3.1415926535897932384626;
        double a = 6378245.0;
        double ee = 0.00669342162296594323;
        // WGS84转GCj02
        double dlat = transformlat(lng - 105.0, lat - 35.0);
        double dlng = transformlng(lng - 105.0, lat - 35.0);
        double radlat = lat / 180.0 * PI;
        double magic = Math.sin(radlat);
        magic = 1 - ee * magic * magic;
        double sqrtmagic = Math.sqrt(magic);
        dlat = (dlat * 180.0) / ((a * (1 - ee)) / (magic * sqrtmagic) * PI);
        dlng = (dlng * 180.0) / (a / sqrtmagic * Math.cos(radlat) * PI);
        double mglat = lat + dlat;
        double mglng = lng + dlng;
        System.out.println("火星坐标: lng:"+mglng+",lat:"+mglat);
        // 火星坐标系 (GCJ-02) 与百度坐标系 (BD-09) 的转换
        double z = Math.sqrt(mglng * mglng + mglat * mglat) + 0.00002 * Math.sin(mglat * xPI);
        double theta = Math.atan2(mglat, mglng) + 0.000003 * Math.cos(mglng * xPI);
        double bdlng = z * Math.cos(theta) + 0.0065;
        double bdlat = z * Math.sin(theta) + 0.006;
        // return [bdlng, bdlat]
        //return {lng: bdlng, lat: bdlat}
        System.out.println("百度坐标:  lng:"+bdlng+",lat:"+bdlat);
    }

    public static Double transformlat(Double lng,Double lat){
        double PI = 3.1415926535897932384626;
        double ret = -100.0 + 2.0 * lng + 3.0 * lat + 0.2 * lat * lat + 0.1 * lng * lat + 0.2 * Math.sqrt(Math.abs(lng));
        ret += (20.0 * Math.sin(6.0 * lng * PI) + 20.0 * Math.sin(2.0 * lng * PI)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(lat * PI) + 40.0 * Math.sin(lat / 3.0 * PI)) * 2.0 / 3.0;
        ret += (160.0 * Math.sin(lat / 12.0 * PI) + 320 * Math.sin(lat * PI / 30.0)) * 2.0 / 3.0;
        return ret;
    }

    public static Double transformlng(Double lng,Double lat){
        double PI = 3.1415926535897932384626;
        double ret = 300.0 + lng + 2.0 * lat + 0.1 * lng * lng + 0.1 * lng * lat + 0.1 * Math.sqrt(Math.abs(lng));
        ret += (20.0 * Math.sin(6.0 * lng * PI) + 20.0 * Math.sin(2.0 * lng * PI)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(lng * PI) + 40.0 * Math.sin(lng / 3.0 * PI)) * 2.0 / 3.0;
        ret += (150.0 * Math.sin(lng / 12.0 * PI) + 300.0 * Math.sin(lng / 30.0 * PI)) * 2.0 / 3.0;
        return ret;
    }


    //方法二
    private final static double a = 6378245.0;
    private final static double pi = 3.1415926535897932384626;
    private final static double ee = 0.00669342162296594323;

    /**
     * Description: WGS-84 to GCJ-02 <BR>
     * @param latitude 纬度
     * @param longitude 经度
     * @return [纬度,经度]
     */
    public static double[] toGCJ02Point(double longitude,double latitude) {
        double[] dev = calDev(latitude, longitude);
        double retLat = latitude + dev[0];
        double retLon = longitude + dev[1];
        System.out.println("火星坐标:  lng:"+retLon+",lat:"+retLat);
        GCJ02_To_BD09(retLon,retLat);
        return new double[] { retLat, retLon };
    }

    /**
     * GCJ-02 to BD-09
     */
    private static Double x_pi = 3.14159265358979324 * 3000.0 / 180.0;  //圆周率转换量
    public static void  GCJ02_To_BD09(Double GCJ02_Lon,Double GCJ02_Lat){
        double x = GCJ02_Lon, y = GCJ02_Lat;
        double z = Math.sqrt(x * x + y * y) + 0.00002 * Math.sin(y * x_pi);
        double theta = Math.atan2(y, x) + 0.000003 * Math.cos(x * x_pi);
        double BD_09_Lon = z * Math.cos(theta) + 0.0065;
        double BD_09_Lat = z * Math.sin(theta) + 0.006;
        System.out.println("百度坐标:  lng:"+BD_09_Lon+",lat:"+BD_09_Lat);
        //return [BD_09_Lon, BD_09_Lat];
    }

    /**
     * Description:GCJ-02 to WGS-84 <BR>
     * @param latitude 纬度
     * @param longitude 经度
     * @return [纬度,经度]
     */
    public static double[] toWGS84Point(double latitude, double longitude) {
        double[] dev = calDev(latitude, longitude);
        double retLat = latitude - dev[0];
        double retLon = longitude - dev[1];
        dev = calDev(retLat, retLon);
        retLat = latitude - dev[0];
        retLon = longitude - dev[1];
        return new double[] { retLat, retLon };
    }

    private static double[] calDev(double wgLat, double wgLon) {
        if (isOutOfChina(wgLat, wgLon)) {
            return new double[] { 0, 0 };
        }
        double dLat = calLat(wgLon - 105.0, wgLat - 35.0);
        double dLon = calLon(wgLon - 105.0, wgLat - 35.0);
        double radLat = wgLat / 180.0 * pi;
        double magic = Math.sin(radLat);
        magic = 1 - ee * magic * magic;
        double sqrtMagic = Math.sqrt(magic);
        dLat = (dLat * 180.0) / ((a * (1 - ee)) / (magic * sqrtMagic) * pi);
        dLon = (dLon * 180.0) / (a / sqrtMagic * Math.cos(radLat) * pi);
        return new double[] { dLat, dLon };
    }

    private static boolean isOutOfChina(double lat, double lon) {
        if (lon < 72.004 || lon > 137.8347)
            return true;
        if (lat < 0.8293 || lat > 55.8271)
            return true;
        return false;
    }

    private static double calLat(double x, double y) {
        double ret = -100.0 + 2.0 * x + 3.0 * y + 0.2 * y * y + 0.1 * x * y + 0.2 * Math.sqrt(Math.abs(x));
        ret += (20.0 * Math.sin(6.0 * x * pi) + 20.0 * Math.sin(2.0 * x * pi)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(y * pi) + 40.0 * Math.sin(y / 3.0 * pi)) * 2.0 / 3.0;
        ret += (160.0 * Math.sin(y / 12.0 * pi) + 320 * Math.sin(y * pi / 30.0)) * 2.0 / 3.0;
        return ret;
    }

    private static double calLon(double x, double y) {
        double ret = 300.0 + x + 2.0 * y + 0.1 * x * x + 0.1 * x * y + 0.1 * Math.sqrt(Math.abs(x));
        ret += (20.0 * Math.sin(6.0 * x * pi) + 20.0 * Math.sin(2.0 * x * pi)) * 2.0 / 3.0;
        ret += (20.0 * Math.sin(x * pi) + 40.0 * Math.sin(x / 3.0 * pi)) * 2.0 / 3.0;
        ret += (150.0 * Math.sin(x / 12.0 * pi) + 300.0 * Math.sin(x / 30.0 * pi)) * 2.0 / 3.0;
        return ret;
    }

    //高斯正反算
    // 高斯反算:由54/80坐标反算成经纬度
    //, double *longitude, double *latitude)
    public static double[] GaussToBL(double X, double Y) {
        int ProjNo; int ZoneWide; 带宽
        double[] output = new double[2];
        double longitude1,latitude1, longitude0, X0,Y0, xval,yval;//latitude0,
        double e1,e2,f,a, ee, NN, T,C, M, D,R,u,fai, iPI;
        iPI = 0.0174532925199433; 3.1415926535898/180.0;
        //a = 6378245.0; f = 1.0/298.3; //54年北京坐标系参数
        a=6378140.0; f=1/298.257; //80年西安坐标系参数
        ZoneWide = 6; 6度带宽
        ProjNo = (int)(X/1000000L) ; //查找带号
        longitude0 = (ProjNo-1) * ZoneWide + ZoneWide / 2;
        //longitude0 = longitude0 * iPI ; //中央经线
        longitude0 = 117 * iPI ; //中央经线

        X0 = ProjNo*1000000L+500000L;
        Y0 = 0;
        xval = X-X0; yval = Y-Y0; //带内大地坐标
        e2 = 2*f-f*f;
        e1 = (1.0-Math.sqrt(1-e2))/(1.0+Math.sqrt(1-e2));
        ee = e2/(1-e2);
        M = yval;
        u = M/(a*(1-e2/4-3*e2*e2/64-5*e2*e2*e2/256));
        fai = u+(3*e1/2-27*e1*e1*e1/32)*Math.sin(2*u)+(21*e1*e1/16-55*e1*e1*e1*e1/32)*Math.sin(
                4*u)
                +(151*e1*e1*e1/96)*Math.sin(6*u)+(1097*e1*e1*e1*e1/512)*Math.sin(8*u);
        C = ee*Math.cos(fai)*Math.cos(fai);
        T = Math.tan(fai)*Math.tan(fai);
        NN = a/Math.sqrt(1.0-e2*Math.sin(fai)*Math.sin(fai));
        R = a*(1-e2)/Math.sqrt((1-e2*Math.sin(fai)*Math.sin(fai))*(1-e2*Math.sin(fai)*Math.sin(fai))*(1-e2*Math.sin
                (fai)*Math.sin(fai)));
        D = xval/NN;
        //计算经度(Longitude) 纬度(Latitude)
        longitude1 = longitude0+(D-(1+2*T+C)*D*D*D/6+(5-2*C+28*T-3*C*C+8*ee+24*T*T)*D
                *D*D*D*D/120)/Math.cos(fai);
        latitude1 = fai -(NN*Math.tan(fai)/R)*(D*D/2-(5+3*T+10*C-4*C*C-9*ee)*D*D*D*D/24
                +(61+90*T+298*C+45*T*T-256*ee-3*C*C)*D*D*D*D*D*D/720);
        //转换为度 DD
        output[0] = longitude1 / iPI;
        output[1] = latitude1 / iPI;
        return output;
        //*longitude = longitude1 / iPI;
        //*latitude = latitude1 / iPI;
    }
    // 高斯正算:由经纬度反算成高斯投影坐标
    public static void GaussToBLToGauss(double longitude, double latitude){
        int ProjNo=20; int ZoneWide; 带宽
        double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval;
        double a,f, e2,ee, NN, T,C,A, M, iPI;
        iPI = 0.0174532925199433; 3.1415926535898/180.0;
        ZoneWide = 6; 6度带宽
        //a=6378245.0; f=1.0/298.3; //54年北京坐标系参数
        a=6378140.0; f=1/298.257; //80年西安坐标系参数
        ProjNo = (int)(longitude / ZoneWide) ;
        longitude0 = ProjNo * ZoneWide + ZoneWide / 2;
        longitude0 = longitude0 * iPI ;
        latitude0 = 0;
        System.out.println(latitude0);
        longitude1 = longitude * iPI ; //经度转换为弧度
        latitude1 = latitude * iPI ; //纬度转换为弧度
        e2=2*f-f*f;
        ee=e2*(1.0-e2);
        NN=a/Math.sqrt(1.0-e2*Math.sin(latitude1)*Math.sin(latitude1));
        T=Math.tan(latitude1)*Math.tan(latitude1);
        C=ee*Math.cos(latitude1)*Math.cos(latitude1);
        A=(longitude1-longitude0)*Math.cos(latitude1);
        M=a*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude1-(3*e2/8+3*e2*e2/32+45*e2*e2
                *e2/1024)*Math.sin(2*latitude1)
                +(15*e2*e2/256+45*e2*e2*e2/1024)*Math.sin(4*latitude1)-(35*e2*e2*e2/3072)*Math.sin(6*latitude1));
        xval = NN*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120);
        yval = M+NN*Math.tan(latitude1)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24
                +(61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720);
        X0 = 1000000L*(ProjNo+1)+500000L;
        Y0 = 0;
        xval = xval+X0; yval = yval+Y0;
        //*X = xval;
        //*Y = yval;
        System.out.println("x:"+xval);
        System.out.println("y:"+yval);
    }
}