c语言 坐标系转换 axistoradius,GIS坐标系、投影和转换

整理资料信息如下：

坐标转换工具：几个常用坐标系，在投影坐标和高斯坐标之间相互转换。

空间直角坐标系转换：使用七参数转换不同的坐标系。

GIS坐标系和ArcGIS坐标系介绍：GIS地理坐标系(大地坐标系)、投影坐标系(主要是高斯克吕格投影、测量坐标系)、地理坐标和投影坐标之间转换。同一个坐标系的经纬度坐标和投影坐标之间可以直接转换，不需要其它参数。

高斯克吕格投影：高斯克吕格投影坐标系中X轴和Y轴的定义，投影坐标值的构成方式、3度分带和6度分带的划分方式。

Yaw(偏航角)、Pitch(倾斜角)、Roll(旋转角)：三维坐标系中这个三个角度的定义。

Java使用七参数做坐标变换：

package gis.util;

public class XYTransformor {

private double dx = 0 ;

private double dy = 0 ;

private double dz = 0 ;

private double wx = 0 ;

private double wy = 0 ;

private double wz = 0 ;

private double m = 0 ;

private double targetX = 0 ;

private double targetY = 0 ;

private double targetZ = 0 ;

/**

* 设置坐标转换的七参数

* @param dx ：x方向平移量

* @param dy ：y方向平移量

* @param dz ：z方向平移量

* @param wx ：x方向旋转量

* @param wy ：y方向旋转量

* @param wz ：z方向旋转量

* @param m ：比例因子

*/

public XYTransformor(double dx,double dy,double dz,double wx,double wy,double wz,double m){

this.dx = dx ;

this.dy = dy ;

this.dz = dz ;

this.wx = wx ;

this.wy = wy ;

this.wz = wz ;

this.m = m ;

}

//........set and get method

public boolean transform(double x,double y ,double z){

targetX = dx + (1+m)*x + wz*y - wy*z ;

targetY = dy + (1+m)*y - wz*x + wx*z ;

targetZ = dz + (1+m)*z + wy*x - wx*y ;

return true ;

}

}

坐标投影、距离计算、三维坐标角度计算、经纬度与弧度相互转换帮助类：

package gis.util;

public class GISUtil {

/**

* 默认地球长半轴长度

*/

public static double EARTH_RADIUS = 6378140.0;

/**

* 弧度转度参数

*/

public static double RADIAN_CONSTANT = 180.0/Math.PI ;

/**

* 度转弧度参数

*/

public static double ANGLE_CONSTANT = Math.PI/180.0 ;

/**

* 计算地球上两点之间的实际距离，单位米

* @param longitude：第一个点经度，单位度

* @param latitude：第一个点纬度，单位度

* @param longitude2：第二个点经度，单位度

* @param latitude2：第二个点纬度，单位度

*/

public static double getGeographicDistance(double longitude,double latitude,double longitude2,double latitude2,double semimajorAxis){

double tempLongitude = angleToRadian(longitude);

double tempLatitude = angleToRadian(latitude);

double tempLongitude2 = angleToRadian(longitude2);

double tempLatitude2 = angleToRadian(latitude2);

double angleValue = Math.sin(tempLatitude)*Math.sin(tempLatitude2)+Math.cos(tempLatitude)*Math.cos(tempLatitude2)*Math.cos(Math.abs(tempLongitude-tempLongitude2)) ;

return semimajorAxis * Math.acos(angleValue) ;

}

/**

* angle to radian

* @param angle

* @return the length ,unit is meter

*/

public static double angleToRadian(double angle){

return angle * ANGLE_CONSTANT ;

}

/**

* 获取偏航角

* @param x：单位米

* @param y：单位米

* @param x2：单位米

* @param y2：单位米

* @return

*/

public static double getYawAngle(double x, double y, double x2, double y2){

return Math.atan2(y2-y,x2-x) * RADIAN_CONSTANT ;

}

/**

* 获取倾斜角

* @param x：第一个点经度，单位度

* @param y：第一个点纬度，单位度

* @param z：第一个点高程，单位米

* @param x2：第二个点经度，单位度

* @param y2：第二个点纬度，单位度

* @param z2：第二个点高程，单位米

* @return

*/

public static double getPitchAngle(double x, double y, double z, double x2, double y2, double z2,double semimajorAxis){

double length = getGeographicDistance(x,y,x2,y2,semimajorAxis);

return Math.atan2(-length,z2-z)* RADIAN_CONSTANT ;

//return Math.atan2(z2-z,length)* RADIAN_CONSTANT ;

}

/**

* @param longitude：经度

* @param latitude：纬度

* @param zoneWide：经度分带类型

* @param semimajorAxis：地球长半轴长度，单位米

* @param inverseFlattening：扁心率的倒数

*/

public static double[] gaussKrugerProject(double longitude, double latitude,int zoneWide,double semimajorAxis,double inverseFlattening){

//默认为6的分带

if(zoneWide != 3&&zoneWide != 6){

zoneWide = 6 ;

}

计算所在位置分带的中央经线

double primeMeridian = 0 ;

int projNo= (int)((longitude-1.5) / zoneWide) ;

if(zoneWide==3){

primeMeridian = projNo * zoneWide + zoneWide/2.0 + 1.5;

}else{

primeMeridian = projNo * zoneWide + zoneWide/2.0 ;

}

//度转弧度

primeMeridian = angleToRadian(primeMeridian);

longitude = angleToRadian(longitude);

latitude = angleToRadian(latitude);

//第一扁心率

double e2=(2-1/inverseFlattening)/inverseFlattening ;

//第二扁心率

double ee=e2*(1.0-e2);

double N=semimajorAxis/Math.sqrt(1.0-e2*Math.sin(latitude)*Math.sin(latitude));

double T = Math.pow(Math.tan(latitude),2);

double C = ee*Math.pow(Math.cos(latitude),2);

double A=(longitude-primeMeridian)*Math.cos(latitude);

double M=semimajorAxis*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude-(3*e2/8+3*e2*e2/32+45*e2*e2

*e2/1024)*Math.sin(2*latitude)+(15*e2*e2/256+45*e2*e2*e2/1024)*Math.sin(4*latitude)-(35*e2*e2*e2/3072)*Math.sin(6*latitude));

double x = N*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120);

double y = M+N*Math.tan(latitude)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24 + (61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720);

//正北方向为x轴，正东方向为y轴，为了避免x轴出现负值，将坐标原点向西偏移500千米，并添加前缀分带编号

x += 1000000L*(projNo+1)+500000L;

return new double[]{y,x} ;

}

}