求助：经纬度坐标和高斯坐标转换

[Hawker101]

在将经纬度坐标转换成高斯坐标时(用的是北京1954参考球和6度带投影法), 发现在跨带时Y坐标的增量不是平滑的，而是跳跃式的，
比如：
(32.30, 118) -> (3486844, 20595004)
(32.30, 119) -> (3488144, 20690022) DelatY = 95018
(32.30, 120) -> (3490312, 20785066) DelatY = 95044
(32.30, 121) -> (3488144, 21309977) DelatY = 524911 (???)

为何会出现这种情况，应该是连续的值？
希有大虾告知,小弟不胜感激!

[klbt]

计算误差？关注.

[Hawker101]

这是小弟用的高斯投影（6度带）的代码，各位看看还有问题.

double CFlyMapView::N(double B)
{
double BB, AX, BX, E, ES;
AX = 6378245.0;
BX = 6356863.0188;
E = sqrt(AX * AX - BX * BX) / AX;
ES = E * sin(B);
BB = AX / sqrt(1 - ES * ES);
return BB;
}

double CFlyMapView::S(double B)
{
double SS, S1, S2, S3;
S1 = 6367558.497 * B - 16036.4803 * sin(2 * B);
S2 = 16.8281 * sin(4 * B)- 0.02198 * sin(6 * B);
S3 = 0.000031 * sin(8 * B);
SS = S1 + S2 + S3;
return SS;

}
/***********************************************************************

JwdToGauss - convert the GPS coordinate to the gauss coordinate

Entry:
fa - the latitude of GPS (in)
leda - the longitude of the GPS (in)
xx - the x vector of the gauss coordinate (out)
yy - the y vector of the gauss coordinate (out)

Exit: None

/**********************************************************************/
void CFlyMapView::JwdToGauss(double fa, double leda, double &xx, double &yy)
{
double pi = 3.1415926;
double AX, BX, E1, B, L,N1, S1, AT, T, X, Y;
int numOfProjectStrap = (int)(leda/6) + 1;
leda = leda + 3 - numOfProjectStrap * 6.0;

AX = 6378245.0;
BX = 6356863.0188;
E1 = (AX * AX- BX * BX) / pow(BX, 2.0);
B = fa * pi / 180.0;
L = leda * pi / 180.0;
N1 = N(B);
S1 = S(B);
AT = E1 * pow(cos(B), 2.0);
T = tan(B);
X = S1 + L* L* N1 * sin(B) * cos(B)/2 + pow(L, 4.0)
* N1 * sin(B) * pow(cos(B), 3.0) *
(5- T * T + 9 * AT + 4 * AT * AT) / 24;
Y = L * N1 * cos(B) + pow(L, 3.0) * N1 * pow(cos(B), 3.0)
* (1- T * T + AT) / 6 + pow(L, 5.0) * N1 *
pow(cos(B), 5.0) * (5 - 18 * T * T + pow(T, 4.0)) / 120;
Y = Y + 500 * 1000.0 + (numOfProjectStrap) * 1000.0 * 1000.0;

xx = X;
yy = Y;
}

void CFlyMapView::GaussToJwd(double x, double y, double &B, double &L)
{
double e2, e12, tmp, bf, nf, Nf, tf,fa;
double yNf,yy;
int l0,n;
double fa0 =(0.1570460641219e-6); double b2 =(0.25184647783e-2);
double b4 =(0.36998873e-5);
double b6 =(0.74449e-8);
double b8 =(0.1828e-10);
double a =6378245;
double b =6356863.02;
double pi =3.1415926;

e2=(a-b)*(a+b)/a/a; //¼ÆËãÆ«ÂÊ
e12=(a-b)*(a+b)/b/b;
n=(int)(y/1000.0/1000.0); //¼ÆËã´øºÅ
//n=(int)(y/1000.0);
yy=y-(n*1000.0+500.0)*1000.0; //¼ÆËãyÖµ
l0=n*6-3; //¼ÆËãÖÐÑë¾­¶È

fa=fa0*x; //¼ÆËãbf
bf=fa+b2*sin(2.0*fa)+b4*sin(4.0*fa)+b6*sin(6.0*fa)+b8*sin(8.0*fa);

tf=tan(bf); //¼ÆËãtf,nf,Nf
nf=sqrt(e12)*cos(bf);
Nf=sqrt(1.0-e2*pow(sin(bf),2.0));
Nf=a/Nf;

yNf=yy/Nf; //¼ÆËãB,L(ÒÔ¶ÈΪµ¥Î»)
tmp=1.0-yNf*yNf*(5.0+3.0*tf*tf+nf*nf-9.0*tf*tf*nf*nf)/12.0;
tmp=tmp+pow(yNf,4.0)*(61.0+90.0*tf*tf+45.0*pow(tf,4.0))/360.0;
tmp=tmp*yNf*yNf*(1.0+nf*nf)*tf/2.0;
B=(bf-tmp)*180.0/pi;
tmp=1.0-yNf*yNf*(1.0+2.0*tf*tf+nf*nf)/6.0;
tmp=tmp+pow(yNf,4.0)*(5.0+28.0*tf*tf+24.0*pow(tf,4.0)+6.0*nf*nf+8.0*nf*nf*tf*tf)/120.0;
L=(yy/(Nf*cos(bf))*tmp)*180.0/pi;
L=L+l0;

}