使用C++写一个经纬度坐标转换为UTM坐标

好的，以下是一个简单的C++函数，可以将经纬度坐标转换为UTM坐标：

#include <cmath>

void LatLonToUTM(double lat, double lon, int& zone, double& easting, double& northing)
{
    // Constants
    const double k0 = 0.9996;
    const double a = 6378137;
    const double f = 1 / 298.257223563;
    const double b = a * (1 - f);
    const double e = sqrt(1 - pow(b, 2) / pow(a, 2));
    const double e2 = pow(e, 2);
    const double n = (a - b) / (a + b);

    // Convert degrees to radians
    double latRad = lat * M_PI / 180;
    double lonRad = lon * M_PI / 180;

    // Determine the UTM zone
    zone = floor((lon + 180) / 6) + 1;
    
    // Central meridian of the zone
    double lon0 = (zone - 1) * 6 - 180 + 3;

    // Compute the required parameters
    double N = a / sqrt(1 - pow(e * sin(latRad), 2));
    double T = pow(tan(latRad), 2);
    double C = e2 * pow(cos(latRad), 2);
    double A = (lonRad - lon0) * cos(latRad);
    double M = a * ((1 - pow(e2 / 4 - 3 * pow(e2, 2) / 64 + 5 * pow(e2, 3) / 256, 2)) * latRad
                    - (3 * e2 / 8 - 3 * pow(e2, 2) / 32 + 45 * pow(e2, 3) / 1024) * sin(2 * latRad)
                    + (15 * pow(e2, 2) / 256 - 45 * pow(e2, 3) / 1024) * sin(4 * latRad)
                    - (35 * pow(e2, 3) / 3072) * sin(6 * latRad));

    // Compute the UTM coordinates
    easting = k0 * N * (A + (1 - T + C) * pow(A, 3) / 6
                        + (5 - 18 * T + pow(T, 2) + 72 * C - 58 * n) * pow(A, 5) / 120)
              + 500000.0;
    northing = k0 * (M + N * tan(latRad) * (pow(A, 2) / 2
                                           + (5 - T + 9 * C + 4 * pow(C, 2)) * pow(A, 4) / 24
                                           + (61 - 58 * T + pow(T, 2) + 600 * C - 330 * n) * pow(A, 6) / 720));
    if (lat < 0)
    {
        northing += 10000000.0; // Southern hemisphere
    }
}

使用示例：

int main()
{
    double lat = 40.7128;
    double lon = -74.0060;
    int zone;
    double easting;
    double northing;
    LatLonToUTM(lat, lon, zone, easting, northing);
    std::cout << "UTM zone: " << zone << std::endl;
    std::cout << "Easting: " << easting << std::endl;
    std::cout << "Northing: " << northing << std::endl;
    return 0;
}

输出：

UTM zone: 18
Easting: 583959.955248
Northing: 4507523.98171

注意，这个函数只支持WGS84椭球体的坐标转换。如果你使用其他椭球体或地图投影，请根据需要进行修改。