Java空间大地坐标系坐标---高斯投影坐标系坐标(BLH-XY)正反算

/**算法是自己ç”究的，现在公开大家随便看，但是如果这份资料对ä½æœ‰ç”¨è¯·ä½ç»™ä¸€ä¸ªè¯„价，谢谢!*//***@Title:GeodeticForGauss*@Description:空间大地å��æ‡è½¬æ¢é«˜æ–¯æŠ•å½±å¹³é¢ç›´è§’叏æ‡*@param:@paramB纬度*@param:@paramLç»�度*@param:@paramtypeå‡åº¦åˆ†å¸¦*@param:@paramNumç›®æ‡å��æ‡ç³»*@return:void*/publicstaticvoidGeodeticForGauss(doubleB,doubleL,inttype,intNum){TuoQiuJiChuntuoQiuJiChun=null;switch(Num){case84://WGS-84å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(WGS-84å��æ‡ç³»,6378137,6356752.3142);break;case80://西安-80å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(西安-80å��æ‡ç³»,6378140,6356755.2881575287);break;case54://北京-54å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(北京-54å��æ‡ç³»,6378245,6356863.0187730473);break;}doublebeltNum;beltNum=Math.ceil((L-(type==3?1.5:0))/type);if(type==3&&beltNum*3==L-1.5){beltNum+=1;}doubleL0=L-(beltNum*type-(type==6?3:0));//中央经线的度数System.out.println(L0);doublerB,tB,m;rB=B*Math.PI/180;tB=Math.tan(rB);m=Math.cos(rB)*L0*Math.PI/180;doubleN=tuoQiuJiChun.getM_Long()/Math.sqrt(1-tuoQiuJiChun.getFirstE()*Math.sin(rB)*Math.sin(rB));doubleit2=tuoQiuJiChun.getSecondE()*Math.pow(Math.cos(rB),2);doublex=0.5*m*m+(double)1/24*(5-tB*tB+9*it2+4*it2*it2)*Math.pow(m,4)+(double)1/720*(61-58*tB*tB+Math.pow(tB,4))*Math.pow(m,6);doublem0=tuoQiuJiChun.getM_Long()*(1-tuoQiuJiChun.getFirstE());doublem2=(double)3/2*tuoQiuJiChun.getFirstE()*m0;doublem4=(double)5/4*tuoQiuJiChun.getFirstE()*m2;doublem6=(double)7/6*tuoQiuJiChun.getFirstE()*m4;doublem8=(double)9/8*tuoQiuJiChun.getFirstE()*m6;doublea0,a2,a4,a6,a8;a0=m0+(double)1/2*m2+(double)3/8*m4+(double)5/16*m6+(double)35/128*m8;a2=(double)1/2*m2+(double)1/2*m4+(double)15/32*m6+(double)7/16*m8;a4=(double)1/8*m4+(double)3/16*m6+(double)7/32*m8;a6=(double)1/32*m6+(double)1/16*m8;a8=(double)1/128*m8;//求孏又线弧长doubleX1=a0*rB-a2*Math.sin(2*rB)*0.5+(double)1/4*a4*Math.sin(4*rB)-(double)1/6*a6*Math.sin(6*rB)+(double)1/8*a8*Math.sin(8*rB);doubleX=X1+N*x*tB;doubleY=N*(m+(double)1/6*(1-tB*tB+it2)*Math.pow(m,3)+(double)1/720*(5-18*tB*tB+Math.pow(tB,4)+14*it2-58*tB*tB*it2)*Math.pow(m,5))+500000;}/***@Title:GaussForGeodetic*@Description:高斯投影平鏢直角叏æ‡è½¬æ�¢ç©ºé—´å¤§åœ°å��æ‡*@param:@paramX*@param:@paramY*@param:@paramtypeå‡åº¦åˆ†å¸¦*@param:@paramNumç›®æ‡å��æ‡ç³»*@param:@parambeltNum带å�·*@return:void*/publicstaticvoidGaussForGeodetic(doubleX,doubleY,inttype,intNum,doublebeltNum){doubleY1=Y-500000;TuoQiuJiChuntuoQiuJiChun=null;switch(Num){case84://WGS-84å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(WGS-84å��æ‡ç³»,6378137,6356752.3142);break;case80://西安-80å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(西安-80å��æ‡ç³»,6378140,6356755.2881575287);break;case54://北京-54å��æ‡ç³»tuoQiuJiChun=newTuoQiuJiChun(北京-54å��æ‡ç³»,6378245,6356863.0187730473);break;}if(Y11000000){intbeltnum=(int)Math.ceil(Y1/1000000)-1;Y1-=beltnum*1000000+500000;}doublem0=tuoQiuJiChun.getM_Long()*(1-tuoQiuJiChun.getFirstE());doublem2=(double)3/2*tuoQiuJiChun.getFirstE()*m0;doublem4=(double)5/4*tuoQiuJiChun.getFirstE()*m2;doublem6=(double)7/6*tuoQiuJiChun.getFirstE()*m4;doublem8=(double)9/8*tuoQiuJiChun.getFirstE()*m6;doublea0,a2,a4,a6,a8;a0=m0+1.0/2.0*m2+3.0/8.0*m4+5.0/16.0*m6+35.0/128.0*m8;a2=1.0/2.0*m2+1.0/2.0*m4+15.0/32.0*m6+7.0/16.0*m8;a4=1.0/8.0*m4+3.0/16.0*m6+7.0/32.0*m8;a6=1.0/32.0*m6+1.0/16.0*m8;a8=1.0/128.0*m8;doubleB1;doubleBf=X/a0;do{B1=Bf;Bf=(X+1.0/2.0*a2*Math.sin(2*B1)-1.0/4.0*a4*Math.sin(4*B1)+1.0/6.0*a6*Math.sin(6*B1))/a0;}while(Math.abs(B1-Bf)0.00000000001);doubletf=Math.tan(Bf);doubleit2=tuoQiuJiChun.getSecondE()*Math.cos(Bf);doubleNf=tuoQiuJiChun.getM_Long()/Math.sqrt(1-tuoQiuJiChun.getFirstE()*Math.sin(Bf)*Math.sin(Bf));doubleMf=Nf/(1+tuoQiuJiChun.getSecondE()*Math.cos(Bf)*Math.cos(Bf));doubleB1f=Bf-(1.0/2.0)*tf*Math.pow(Y1,2)/(Mf*Nf*Math.cos(Bf))+(1.0/24.0)*tf*(5+3*tf*tf+it2-9*it2*tf*tf)*Math.pow(Y1,4)/(Mf*Math.pow(Nf,3))-(1.0/720.0)*tf*(61+90*tf*tf+45*tf*tf*tf*tf)*Math.pow(Y1,6)/(Mf*Math.pow(Nf,5));doublel=Y1/(Nf*Math.cos(Bf))-(1.0/6.0)*(1+2*tf*tf+it2)*Math.pow(Y1,3)/(Math.pow(Nf,3)*Math.cos(Bf))+(1.0/120.0)*(5+28*tf*tf+24*tf*tf*tf*tf+6*it2+8*it2*tf*tf)*Math.pow(Y1,5)/(Math.pow(Nf,5)*Math.cos(Bf));l=l*180/Math.PI;doubleL=l+type*beltNum-((type==6)?3:0);doubleB=B1f*180/Math.PI;}