GIS坐标转换及其编程实现

最新推荐文章于 2024-06-29 11:41:21 发布
最新推荐文章于 2024-06-29 11:41:21 发布
阅读量6.3k 收藏 14
点赞数 1
分类专栏： GIS底层开发 GIS底层开发 文章标签： 编程 null c
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/zhouxuguang236/article/details/6700706
版权

最近，在做这个坐标转换的东西，涉及到大地测量学等很深奥的东西，在这里我就不讲解那些难懂的理论了，在此，我将会把代码贴出来和大家分享，其实要编写出这个代码，还真得把大地测量相关的知识弄熟，否则是无法理解代码的。好了废话少说。

源代码

/**
  * 空间大地直角坐标->大地坐标
  */
 public GeoPoint XYZ_BLH(int ellipse, Point3D point)
 {
  GeoPoint geoPoint = new GeoPoint();
  if(geoPoint == null || point == null)
  {
   System.out.println("对象为空！");
  }
  double a = 0,b = 0,e1 = 0,e2 = 0; //a为长半轴，b为短半轴，e1为第一偏心率，e2为第一偏心率的平方
  double deta = 0.0000001;
  switch(ellipse) //选择椭球体
  {
   case 0:  //WGS84椭球体
   {
    a = 6378137.0; //长半轴
    b = 6356752.3142;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006694384999588;
    break;
   }
   case 1:  //北京54椭球
   {
    a = 6378245.0; //长半轴
    b = 6356863.0187730473;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006693421622966;
    break;
   }
   case 2:  //西安80椭球
   {
    a = 6378140.0; //长半轴
    b = 6356755.2881575287;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006694384999588;
    break;
   }
   default:break;
  }
  //double W = Math.sqrt(1.0-e1*e1*Math.pow(Math.sin(point.getX()), 2));
  //double N = a/W;
  double X = point.getX();
  double Y = point.getY();
  double Z = point.getZ();
  double m = Math.sqrt(Math.pow(X, 2)+Math.pow(Y, 2));
  double L = geoPoint.getL(); //大地经度
  double B = geoPoint.getB(); //大地纬度
  double H = geoPoint.getH(); //大地高
  L = Math.atan(Y/X);
  if (L < 0)
  {
   L += Math.PI;
  }
  double e2_ = e2/(1-e2);
  double c = a*Math.sqrt(1+e2_);
  double ce2 = c*e2;
  double k = 1 + e2_;
  double front = Z/m;
  double temp = front;
  int count = 0; //迭代次数
  do
  {
   front = temp;
   m = Math.sqrt(Math.pow(X, 2)+Math.pow(Y, 2));
   temp = Z/m + ce2*front/(m*Math.sqrt(k+ Math.pow(front, 2)));
  }
  while(Math.abs(front-temp)<deta&&count<100000); //是否在误差范围之内
  B = Math.atan(temp);//求纬度
  if (B<0)
  {
   B += Math.PI;
  }
  double W = Math.sqrt(1-e1*e1*Math.sin(B)*Math.sin(B));
  double N = a/W;
  //N = (a*m - c*c)/(2*b*Z);
  System.out.println("N = " + N);
  H = m/Math.cos(B)-N;//求高
  //H = Z/Math.sin(B)-N*(1.0-e1*e1);
  //H = X/(Math.cos(B)*Math.cos(L))-N;
  //转换为角度值
  L = Math.toDegrees(L);
  B = Math.toDegrees(B);
  //H = Math.toDegrees(H);
  geoPoint.setL(L);
  geoPoint.setB(B);
  geoPoint.setH(H);
  return geoPoint;
 }
 
 /**
  * 大地坐标->空间大地直角坐标
  * @param ellipse 椭球体
  * @param geoPoint 转换前的坐标
  * @return 返回空间直角坐标
  */
 public Point3D BLH_XYZ(int ellipse, GeoPoint geoPoint)
 {
  Point3D point1 = new Point3D(1,1,1);
  if(geoPoint == null || point1 == null)
  {
   System.out.println("对象为空！");
  }
  double a = 0,b,e1 = 0,e2 = 0; //a为长半轴，b为短半轴，e1为第一偏心率，e2为第一偏心率的平方
  switch(ellipse) //选择椭球体
  {
   case 0:  //WGS84椭球体
   {
    a = 6378137.0; //长半轴
    b = 6356752.3142;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006694384999588;
    break;
   }
   case 1:  //北京54椭球
   {
    a = 6378245.0; //长半轴
    b = 6356863.0187730473;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006693421622966;
    break;
   }
   case 2:  //西安80椭球
   {
    a = 6378140.0; //长半轴
    b = 6356755.2881575287;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = 0.006694384999588;
    break;
   }
   default:break;
  }
  double B = geoPoint.getB()*Math.PI/180; //转换为弧度制
  double L = geoPoint.getL()*Math.PI/180; 
  double H = geoPoint.getH();
  //计算卯酉圈曲率半径
  double N = a/Math.sqrt(1.0-e1*e1*Math.pow(Math.sin(B), 2));
  //计算空间直角坐标
  double X = (N + H)*Math.cos(B)*Math.cos(L);
  double Y = (N + H)*Math.cos(B)*Math.sin(L);
  double Z = (N*(1-e1*e1)+ H)*Math.sin(B);
  point1.X = X;
  point1.Y = Y;
  point1.Z = Z;
  return point1;
 }
 
 /**
  * 高斯-克吕格坐标正算（从大地坐标到平面直角坐标）
  * @param ellipse 椭球体
  * @param zoneWide 带宽(3度或6度)
  * @param geoPoint 大地坐标
  * @param point  直角坐标
  */
 public void BL_xy(int ellipse,int zoneWide,GeoPoint geoPoint,Point point)
 {
  if(geoPoint == null || point == null)
  {
   System.out.println("对象为空！");
  }
  double a = 0,b,e1 = 0,e2 = 0; //a为长半轴，b为短半轴，e1为第一偏心率，e2为第一偏心率的平方
  switch(ellipse) //选择椭球体
  {
   case 0:  //WGS84椭球体
   {
    a = 6378137.0; //长半轴
    b = 6356752.3142;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006694384999588;
    break;
   }
   case 1:  //北京54椭球
   {
    a = 6378245.0; //长半轴
    b = 6356863.0187730473;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006693421622966;
    break;
   }
   case 2:  //西安80椭球
   {
    a = 6378140.0; //长半轴
    b = 6356755.2881575287;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006694384999588;
    break;
   }
   default:break;
  }
  
  double B = geoPoint.getB();
  double L = geoPoint.getL();
  double centerL = 0.0;//投影带中央经线的度数
  int n = 0;   //投影带的带号（6°带是1-60，3°带是1-120）
  double l = 0;  //该地所在经度与中央经度的差值
  if (6 == zoneWide) //如果是6度分带
  {
   n = (int)L/6;
   if (L%6 > 0)
   {
    n = n + 1;
   }
   centerL = 6*n-3;
   l = L - centerL;
  }
  if (3 == zoneWide)//如果是3度分带
  {
   n = (int)(L-1.5)/3;
   if ((L-1.5)%3>0)
   {
    n = n +1;
   }
   centerL = 3*n;
   l = L - centerL;
  }
  //转化为弧度
  B = B*Math.PI/180;
  L = L*Math.PI/180;
  l = l*Math.PI/180;
  double X;//从赤道起算的子午线弧长
  //计算子午线弧长的系数
  double A0 = 1.0+3.0/4.0*Math.pow(e1, 2)+45.0/64.0*Math.pow(e1, 4)
  +350.0/512.0*Math.pow(e1, 6)+11025.0/16384.0*Math.pow(e1, 8);
  double A2 = -3.0/4*Math.pow(e1, 2)+60.0/64*Math.pow(e1, 4)+
  525.0/512*Math.pow(e1, 6)+17640.0/16384.0*Math.pow(e1, 8)/2.0;
  double A4 = 15.0/64.0*Math.pow(e1, 4)+210.0/512.0*Math.pow(e1, 6)+
  8820.0/16384.0*Math.pow(e1, 8)/4.0;
  double A6 = -35.0/512.0*Math.pow(e1, 6)+2520.0/16384.0*Math.pow(e1, 8)/6.0;
  double A8 = 315.0/16384.0*Math.pow(e1, 8)/8.0;
  //计算子午线弧长X
  X = a*(1.0-Math.pow(e1, 2))*A0*B + A2*Math.sin(2*B) + A4*Math.sin(4*B)
   + A6*Math.sin(6*B) + A8*Math.sin(8*B);
  double t = Math.tan(B);
  double anke = e2*Math.cos(B);
  double N = a/Math.sqrt(1.0-Math.pow(e1, 2)*Math.pow(Math.sin(B), 2));
  //坐标计算
  double x = point.getX();
  double y = point.getY();
  x = X + 1.0/2.0*N*t*Math.pow(Math.cos(B), 2)*Math.pow(l, 2)+
  1.0/24.0*N*t*(5.0-Math.pow(t, 2)+9*Math.pow(anke, 2)+4*Math.pow(anke, 4))*Math.pow(Math.cos(B), 4)*Math.pow(l, 4)
  + 1.0/720*N*t*(61.0-58*Math.pow(t, 2)+Math.pow(t, 4)+270.0*Math.pow(anke, 2)-330.0*Math.pow(anke, 2)*Math.pow(t, 2))*Math.pow(Math.cos(B), 6);
  
  y = N*Math.cos(B)*l + 1.0/6*N*(1.0-Math.pow(t, 2)+Math.pow(anke, 2))*Math.pow(Math.cos(B), 3)*Math.pow(l, 3)+
  1.0/120*N*(5.0-18*t*t+Math.pow(t, 4)+14*Math.pow(anke, 2)-58*anke*anke*t*t)*
  Math.pow(Math.cos(B), 3)*Math.pow(l, 5);
  y += 500000.0;//加上500KM
  point.setX(x);
  point.setY(y);
 }
 
 /**
  * 高斯平面直角坐标转换为大地坐标
  * @param ellipse 椭球体
  * @param zoneWide 带宽
  * @param point  平面坐标
  * @param geoPoint 大地坐标
  */
 public void xy_BL(int ellipse,double centerL,Point point,GeoPoint geoPoint)
 {
  if(geoPoint == null || point == null)
  {
   System.out.println("对象为空！");
  }
  double a = 0,b,e1 = 0,e2 = 0; //a为长半轴，b为短半轴，e1为第一偏心率，e2为第二偏心率
  switch(ellipse) //选择椭球体
  {
   case 0:  //WGS84椭球体
   {
    a = 6378137.0; //长半轴
    b = 6356752.3142;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006694384999588;
    break;
   }
   case 1:  //北京54椭球
   {
    a = 6378245.0; //长半轴
    b = 6356863.0187730473;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006693421622966;
    break;
   }
   case 2:  //西安80椭球
   {
    a = 6378140.0; //长半轴
    b = 6356755.2881575287;//短半轴
    e1 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/a;
    e2 = Math.sqrt(Math.pow(a, 2)-Math.pow(b, 2))/b;
    //e2 = 0.006694384999588;
    break;
   }
   default:break;
  }
  double x = point.getX();
  double y = point.getY()-500000;
  double Bf = 0.0;//底点纬度
  double B0 = 1.0+3.0/4.0*e1*e1 + 45.0/64.0*Math.pow(e1, 4)+350.0/312.0*Math.pow(e1, 6)+
     11025.0/16384.0*Math.pow(e1, 8) + 43659.0/65536.0*Math.pow(e1, 10);
  Bf = x/(a*(1.0-Math.pow(e1, 2))*B0);
  double tf = Math.tan(Bf);
  double Mf = a*(1.0-e1*e1)/Math.sqrt(Math.pow(1.0-Math.pow(e1, 2), 3));
  double Nf = (a/Math.sqrt(1.0-Math.pow(e1, 2)))/Math.sqrt(1.0+Math.pow(e2*Math.cos(Bf),2));
  double anke = e2*Math.cos(Bf);
  double B = geoPoint.getB();
  double L = geoPoint.getL();
  //开始坐标计算
  B = Bf - tf/(2*Mf*Nf*Math.cos(Bf))*Math.pow(y, 2) + tf/(24*Mf*Math.pow(Nf, 3))*
  (5.0+3.0*tf*tf+anke*anke-9.0*anke*anke*tf*tf)*Math.pow(y, 4)-
  1.0/(720.0*Math.pow(Nf, 5)*Math.cos(Bf))*(61.0+90.0*tf*tf+45*Math.pow(tf,4))*Math.pow(y, 6);
  L = y/(Nf*Math.cos(Bf)) - (1.0+2*tf*tf+anke*anke)*Math.pow(y, 3)/(6.0*Math.pow(Nf, 3)*Math.cos(Bf))
  + (5.0+28*tf*tf+24*Math.pow(tf, 4)+6*anke*anke+8*anke*anke*tf*tf)*Math.pow(y, 5)/(120.0*Math.pow(Nf, 5)*Math.cos(Bf));
  
  B = B*180/Math.PI;
  L = L*180/Math.PI;
  L = L + centerL;//中央经线+经度差
  geoPoint.setB(B);
  geoPoint.setL(L);
 }

这些代码可能有什么不合理的地方，还有就是由空间直角坐标转换为大地坐标的时候，大地高计算有很大偏差，一致找不出来是什么原因

周旭光
关注
  • 1
    点赞
  • 踩
  • 14
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
专栏目录
基于Matlab实现坐标转换及地图投影（源码）.rar
07-24
本资源"基于Matlab实现坐标转换及地图投影（源码）.rar"提供了一种利用MATLAB编程语言来处理这些问题的方法。MATLAB是一种强大的数学计算软件，常用于科学计算、数据分析和算法开发。以下是对这两个主题的详细解释：...
c#大地测量学各坐标系坐标转换（窗体应用程序）
yitian_00的博客
04-25 2288
【代码】c#大地测量学各坐标系坐标转换（窗体应用程序）
GIS数据格式坐标转换（地球坐标WGS84、GCJ-02、火星坐标、百度坐标BD-09、国家大地坐标系CGCS2000）
漏刻有时数据可视化大屏（PHP&ECHARTS智能化开源软件系统）
02-12 1万+
地理信息系统 (GIS) 是一个创建、管理、分析和绘制所有类型数据的系统。GIS 将数据连接到地图，将位置数据（事物所在位置）与所有类型的描述性信息（事物在该位置的情况）集成到一起。这可以为适用于自然科学和几乎所有行业的制图和分析提供基础。GIS 帮助用户了解模式、关系和地理环境。其优势包括改善沟通、提高效率以及更好地管理和决策。解决使用百度webAPI坐标转换时的问题：X→GPS根据法律规定，不支持将任何类型的坐标转换为GPS坐标；
ArcGIS中将测绘数据投影坐标（平面坐标）转地理坐标（球面经纬度坐标）
最新发布
AIGIS的博客
06-29 3218
本文介绍的测绘数据是西安80的投影坐标，即平面坐标，我们要做的转为经纬度坐标，即地理坐标，其他坐标系转换类似。
gis中的坐标系
hfhan的博客
05-17 4527
地球的三级逼近 地球的自然表面不是平整的，需要想办法用数学公式描述地球表面，只能设想一个近似的数学面。 大地水准面是地球表面的第一级逼近。假设当海水处于完全静止的平衡状态时，从海平面延伸到所有大陆下部，而与地球重力方向处处正交的一个连续、闭合的曲面，这就是大地水准面。 地球椭球体是地球表面的第二级逼近。大地水准面可以近似成一个规则成椭球体，但并不是完全规则，其形状接近一个扁率极小的椭圆绕短轴旋转所形成的规则椭球体，这个椭球体称为地球椭球体。 由一系列参数可确定一个椭球体，地球椭球体的基本参数如下 大地基准
GIS地图坐标系转换方法，快快收藏
layliangbo的博客
09-08 2537
开发地图的时最烦的就是不同的地图工具以及数据来源坐标系不同造成偏移甚至报错。以下整理了部分的地图坐标转换，一个文件直接搞定所有情况。
ArcGIS坐标系转化
小贝壳
08-09 9565
ArcGIS坐标系转化
坐标转换工具，需在ArcGIS 10.1及以上版本内使用
11-21
在GIS（地理信息系统）领域，坐标转换是一项至关重要的任务，尤其在处理来自不同来源或采用不同坐标系统的地理数据时。本文将围绕“坐标转换工具”这一主题，详细讲解相关知识点，以及如何在ArcGIS 10.1及以上版本中...
GIS中坐标转换程序的设计与实现
04-13
VC++作为一种成熟的开发语言，具有强大的计算能力和图形界面处理能力，是实现GIS坐标转换程序的理想选择。利用VC++开发的坐标转换程序不仅能够快速实现单点坐标转换，还可以处理批量数据转换，通过文件化管理坐标...
坐标转换与VC编程.pdf
10-01
《坐标转换与VC编程》是一本深入探讨地理坐标系统转换及其在Visual C++（简称VC）环境下实现的专著。本书旨在为读者提供一个全面理解坐标转换原理和技术，并结合实践运用到VC编程中的指南。 第一章绪论介绍了坐标...
七参数转换及大地坐标系与空间坐标转换Qt实现
01-21
总的来说，"七参数转换及大地坐标系与空间坐标转换Qt实现"项目是一个结合了地理信息技术和软件工程实践的优秀案例，不仅涵盖了基础的GIS理论，还体现了Qt编程的技巧，对学习者提升理论知识和编程能力都大有裨益。...
GIS基础教程之坐标系
2301_77413856的博客
09-25 2238
在ArcGIS中对数据定义坐标系,可以使用Catalog的数据管理功能，也可以使用ArcToolbox的定义投影（Define Projection）工具。在Catalog中查看数据的坐标系属性并修改定义投影工具注意，上面这两种方法都是直接修改数据的坐标系属性，适用条件是：数据坐标系信息缺失或坐标系信息错误如何走好GIS数据处理的第一步显示坐标系会直接影响地图出图的形状。专题要素到标准地图的整合——以中国5A级景区分布图为例。
ArcGis系列-坐标系转换
占星安啦的博客
06-19 6649
本文将介绍ArcGis三个需要坐标转换的场景： - Arcgis Pro设置项目坐标 - GP分析指定输出结果坐标系 - 使用python发布图层或要素服务时指定坐标系
GIS中的坐标系
weixin_46720482的博客
04-22 2142
一、坐标系类型 地理坐标系和投影坐标系两大类； 地理坐标系又可以分为参心和地心大地或者空间直角坐标系，投影坐标系按照高斯克吕格3°和6°分带； (1)大地坐标系 是大地测量的基本坐标系，是全世界公用的最方便的坐标系统，对于研究地球形状、编制地图和大地测量的计算等方面都有很大作用。 空间中任意点的位置采用大地纬度B、大地经度L和大地高H来表示。 大地纬度是指空间中某一点P的法线与赤道面的夹角，赤道以北叫做北纬；赤道以南叫做南纬。 大地经度是指P点所在的子午面与参考椭球的起始子午面所构成的夹角..
常见GIS坐标系的坐标转换
邢慧鹏的博客
04-26 1285
几种常见的 GIS 坐标系的简单介绍及坐标数据的相互转换方法。
ArcGis坐标转换
kupePoem的专栏
05-30 1201
一般不同大地测量系统坐标系之间的转换会涉及保密参数的使用，此参数为【三参数】或者【七参数】，均为国家保密参数，需要到当地的测绘部门或者国土部门，以单位名义签保密协议进行购买，此参数在各个区域是不同地。坐标系转换可以是地理坐标系与投影坐标系之间的转换，可以是地理坐标系与地理坐标系之间的转换，也可以是投影坐标系与投影坐标系之间转换。（注:坐标系转换前记得先进行坐标系的定义哦）。需要注意的是，有带号的坐标系都是8位，前两位表示的是带号。地理坐标系经过投影后变成投影坐标系，投影坐标系因此由地理坐标系和投影组成。
arcgis里怎么转换坐标系（转载）
热门推荐
KK的笔记本
04-02 4万+
1、首先我们要查看数据当前坐标系统是什么，我们点击“开始”，然后点击“ArcGIS”目录下的“ArcCatalog”，从而打开电脑里的打开ArcCatalog。 2、打开了软件后我们找到数据所在位置。 3、然后右键点击数据集，打开“图层属性”界面。 4、选择“XY坐标系”选项卡，查看数据的当前坐标系统。GCS_WGS_1984，即世界级地理坐标系(World Geodetic System 1984)，是为GPS全球定位系统使用而建立的坐标系统，也是全球所有卫星所用的坐标系统，主要参数是地球经纬度。 5、接
GIS中坐标系的基本概念
专注GIS软件开发。
12-02 2412
GIS坐标系的分类 1、地理坐标系（Geographic Coordinate System, GCS） 2、投影坐标系（Projected Coordinate System, PCS） 平面极坐标系，这个在初中的时候都学过。GIS中的坐标系可以从这个概念入门学习。可以理解为“为确定两个值之间相对变化而产生的确定值“。这两个值在平面上。 在理解一下球坐标系，“为确定三个值之间相对变化而产生的确定值“。这三个值球面上。 现在再来看下GIS中的地理坐标系，就是在GIS学科中的球坐标系。以地球的质心为圆心，
Arcgis 导入经纬度坐标、导出经纬度坐标
知识搬运者
03-24 1万+
此内容不用管，仅仅记录：（如果要解决表”没有 Object-ID 字段“的提示--将XY数据图层另存为要素类）>>>> 可以右击数据表图层，点击【数据】|【导出数据】。以下图表为例：度=MID(B3,1,3)+MID(B3,5,2)/60+MID(B3,8,2)/3600。其中MID（数据框（B1），从第几位数为开始，提取多少位数），度分秒也各占一位数。右键样本点的图层，选择数据，导出数据，即可导出shp文件。选择添加数据，找到自己的经纬度坐标的表格，点击添加。其中，X是经度，Y是纬度。
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值