有关经纬度转化高斯坐标的问题

我想写一个c语言程序把经纬度转化为高斯坐标（经度12005.28270, 纬度3604.15877 ddmm.mmmm，3度带，第四十带）
代码如下

#include <stdio.h>
#include <math.h>
// 3°带，第40带，青岛，之后要改
#define MID_LNG 120 
#define EARTH_A 6378137
#define EARTH_B 6356752.3142
#define EARTH_C 6399593.6359
#define EARTH_E2 0.00669437999013   // 目前该程序有问题，极有可能是公式错误， 需要修改 
#define PI 3.14159265358979
double lng = 12005.28270, lat = 3604.15877, x, y;

double to_sss(double angle)//将经纬度转化为秒为单位
{
    double dd = (int)(angle / 100.0);
    double mm = (angle - dd * 100.0) / 100.0;
    return (60 * (dd * 60 + mm));
}
void GeoToXY(double jd, double wd, double *y, double *x)
{
 double t;     //  t=tgB
 double L = MID_LNG ;     //  中央经线的经度
 double l0;    //  经差
 double jd_hd,wd_hd;  //  将jd、wd转换成以弧度为单位
 double et2;    //  et2 = (e' ** 2) * (cosB ** 2)
 double N;     //  N = C / sqrt(1 + et2)
 double X;     //  克拉索夫斯基椭球中子午弧长
 double m;     //  m = cosB * PI/180 * l0 
 double tsin,tcos;   //  sinB,cosB


 double b_e2=0.0067385254147;
 double b_c=6399698.90178271;

 jd_hd = jd / 3600.0 * PI / 180.0 ;    // 将以秒为单位的经度转换成弧度
 wd_hd = wd / 3600.0 * PI / 180.0 ;    // 将以秒为单位的纬度转换成弧度
  
 
 l0 = jd / 3600.0 - L  ;       // 计算经差
 tsin = sin(wd_hd);        // 计算sinB
 tcos = cos(wd_hd);        // 计算cosB
             // 计算克拉索夫斯基椭球中子午弧长X
 X = 111134.8611 / 3600.0 * wd - (32005.7799 * tsin + 133.9238 * pow(tsin,3)  
      + 0.6976 * pow(tsin,5) + 0.0039 * pow(tsin,7) ) * tcos; 
 et2 = b_e2 * pow(tcos,2) ;      //  et2 = (e' ** 2) * (cosB ** 2)
 N  = b_c / sqrt( 1 + et2 ) ;      //  N = C / sqrt(1 + et2)
 t  = tan(wd_hd);         //  t=tgB
 m  = PI/180 * l0 * tcos;       //  m = cosB * PI/180 * l0 
 *x = X + N * t * ( 0.5 * pow(m,2)                                           
          + (5.0 - pow(t,2) + 9.0 * et2 + 4 * pow(et2,2)) * pow(m,4)/24.0      
          + (61.0 - 58.0 * pow(t,2) + pow(t,4)) * pow(m,6) / 720.0 ) ;
 *y = N * ( m + ( 1.0 - pow(t,2) + et2 ) * pow(m,3) / 6.0                    
                + ( 5.0 - 18.0 * pow(t,2) + pow(t,4) + 14.0 * et2              
                   - 58.0 * et2 * pow(t,2) ) * pow(m,5) / 120.0 );

}
int main()
{
    lng = to_sss(lng); lat = to_sss(lat);
    printf("%f, %f\n", lng, lat);
    GeoToXY(lng, lat, &y, &x);
    printf("%f, %f", x, y);
    return 0;
}

这个程序很大部分是粘贴自网络上搜索“经纬度转化为高斯坐标c语言实现”，然后略作修改。对照另一篇论文，我认为代码里的公式没有问题

现在的情况是 转化的x，y坐标为3985690.215217, 79.385255，但是我用一个软件叫gausscoor转化的只有一个跟我一样

恳请大神们赐教，已经卡住好几天了

double to_sss(double angle)//将经纬度转化为秒为单位
你这经纬度是什么单位的值呢？
double mm = (angle - dd * 100.0) / 100.0;
为何还要除以100呢？