GPS坐标转换为XYZ坐标

#include <cstdio>
#include <cmath>
#include <algorithm>
#include <iostream>
#include <fstream>
#include <iomanip>
#include <sstream>
#include <string>
// double wgs84a;
// double wgs84b;
// double dsq;
// double eccsq;
// double clat;
// double slat;
// double clon;
// double slon;
// double esq;
// double  x,y,z;
//  double latitude,  longitude, altitude;
// using namespace std;
// double gpsPoint2DescartesPointX(double latitude, double longitude, double altitude);
// double gpsPoint2DescartesPointY(double latitude, double longitude,double altitude);
// // double latitude = 30.337705166;
// //  double longitude = 121.238201404;
// //  double altitude = 7.321;


// double gpsPoint2DescartesPointX(double latitude, double longitude, double  altitude)
// {
// 	//wgs84 WGS84 Earth Constants
// 	double wgs84a = 6378.137;
// 	double wgs84f = 1.0 / 298.257223563;
// 	double wgs84b = wgs84a * (1.0 - wgs84f);
 
// 	//earthcon
// 	double f = 1 - wgs84b / wgs84a;
// 	double eccsq = 1 - (wgs84b* wgs84b) / (wgs84a * wgs84a);
// 	double ecc = sqrt(eccsq);
// 	double esq = ecc * ecc;
 
// 	//llhxyz
// 	double dtr = M_PI / 180.0;
// 	//qDebug() << dtr << gpsPoint.latitude << endl;
// 	double clat = cos(dtr * latitude);
// 	double slat = sin(dtr * latitude);
// 	double clon = cos(dtr * longitude);
// 	double slon = sin(dtr * longitude);
// 	//qDebug() << clat << slon << endl;
 
// 	//radcur compute the radii at the geodetic latitude lat (in degrees)
// 	double  dsq = 1.0 - eccsq * slat *slat;
// 	double  d = sqrt(dsq);
// 	//qDebug() << d;
// 	double rn = wgs84a / d;
// 	double rm = rn * (1.0 - eccsq) / dsq;
 
// 	double rho = rn * clat;
// 	double zz = (1.0 - eccsq) * rn *slat;
// 	double rsq = rho * rho + zz*zz;
// 	double r = sqrt(rsq);
 
// 	double x = (rn + altitude) * clat * clon;
// // 	double y = (rn + altitude) * clat * slon;
// // 	double z = ((1 - esq)*rn + altitude) * slat;
// return x;
// }


// double gpsPoint2DescartesPointY(double latitude, double longitude, double altitude)
// {
// 	//wgs84 WGS84 Earth Constants
// 	double wgs84a = 6378.137;
// 	double wgs84f = 1.0 / 298.257223563;
// 	double wgs84b = wgs84a * (1.0 - wgs84f);
 
// 	//earthcon
// 	double f = 1 - wgs84b / wgs84a;
// 	double eccsq = 1 - (wgs84b* wgs84b) / (wgs84a * wgs84a);
// 	double ecc = sqrt(eccsq);
// 	double esq = ecc * ecc;
 
// 	//llhxyz
// 	double dtr = M_PI / 180.0;
// 	//qDebug() << dtr << gpsPoint.latitude << endl;
// 	double clat = cos(dtr * latitude);
// 	double slat = sin(dtr * latitude);
// 	double clon = cos(dtr * longitude);
// 	double slon = sin(dtr * longitude);
// 	//qDebug() << clat << slon << endl;
 
// 	//radcur compute the radii at the geodetic latitude lat (in degrees)
// 	double  dsq = 1.0 - eccsq * slat *slat;
// 	double  d = sqrt(dsq);
// 	//qDebug() << d;
// 	double rn = wgs84a / d;
// 	double rm = rn * (1.0 - eccsq) / dsq;
 
// 	double rho = rn * clat;
// 	double zz = (1.0 - eccsq) * rn *slat;
// 	double rsq = rho * rho + zz*zz;
// 	double r = sqrt(rsq);
 
// //double x = (rn + altitude) * clat * clon;
// double y = (rn + altitude) * clat * slon;
// // 	double z = ((1 - esq)*rn + altitude) * slat;
// return  y;
// }



// int  main()
// {
//     double a = 30.337705166;
//     double b= 121.238201404;
//     double c = 7.321;

//     x = gpsPoint2DescartesPointX( a,   b,  c);
//     y = gpsPoint2DescartesPointY(a,  b,  c);
//     // double a = 30.337705166;
//     // double b= 121.238201404;
//     // double c = 7.321;
//     printf ("X是%f",x);
//     printf ("Y是%f",y);




//     return 0;
// }