怎么把119度带的坐标转换为120度带...

3度带和6度带坐标转换|3度带和6度带高斯-克吕格投影正反转换程序下载 绿色版_ - pc6下载站坐标转换经纬度
坐标转换经纬度
范文一:84的经纬度坐标转换84的经纬度坐标转换摘要:关于这个坐标系的转化网上有很多文章探讨了各种转换的方法。通过自己的学习,我自己做了一下总结,同时给出了其中要遇到的部分术语和数据,方便以后查阅使用。主要介绍的是:3参数(七参数)转换法,三参数坐标纠正法一:3参数(七参数)转换法从本质上来说,转换的步骤应该大致遵循这样的规则:首先,将84的经纬度坐标(B84,L84,H84)转换为以地心为中心点的大地坐标(X84,Y84,Z84);然后根据七参数法(或3参数法)将其转换为54下的地心坐标(X54,Y54,Z54);然后根据54下的椭球参数,将第二步得到的地心坐标转换为大地坐标(B54,L54,H54);最后根据工程需要以及各种投影(如高斯克吕格)规则进行投影得到对应的投影坐标。只有在第二步的时候涉及到七个参数的计算,其他的步骤都有现成的公式可供计算,稍后我会将各种论文贴上来。如果这里涉及到您的利益还请跟我联系,我将马上删除下载链接,我本意只是用于学习使用。其实如果在公司或者做项目的时候,当对这起个参数要求的很急的时候,我们可以从政府部门或者通过坐标转换软件求出这七个参数或者三个参数,这个可以大大提高效率,节省时间。这些坐标转换软件有:坐标转换大师(这个不错),coorconvert.exe(一般),COORD.exe(这个不错)。一旦求出了七个参数,可以进行坐标转换的软件除了上述这些小软件可以进行转换外,一些比较有名的GIS开发软件或者开发平台都提供了利用七个参数转换整个数据的功能或者提供了转换单个点的功能,这些在ARC GIS,superMap,mapGis中都有。二:三参数坐标纠正法这个方法是这次我在实践中得出来的。因为求出七个参数太过麻烦,所以选用了本方法。本方法的使用范围为:大比例尺地形图比较适用,如县范围等。具体方法:1.从测区取出适量的坐标控制点,坐标控制点是些这样的点,他们拥有84下的经纬度坐标,同时也拥有54下的投影坐标;2.取出后利用将经纬度坐标在esupermap平台中编写程序将其转成84下的高斯克吕格投影坐标(可以看成是一种虚假的投影);3.由2步中得到的投影坐标和原54下的投影坐标相比较得到一个差值p1(x1,y1,z1),并将其保存起来;4.重复第二步一直到把所有的点都计算完,计算完后将差值进行汇总并得到一个平均值p(x,y,z).通过此方法得到的三个参数经过测试和验证,他的精度在厘米或者亚米级的进度,这个对于一般的定位来说已经足够了。他正宗的参数法法的精度还要高。总结:当然,如果要提高精度的话,最好还是用七参数法,他的定位精度基本上都在厘米或者毫米级。54大地参数:参考椭球体:Krasovsky_1940长半轴:6378245短半轴:8扁率:298.384大地参数: 参考椭球体:WGS 84
长半轴:6378137
扁率:298.25722480大地参数: 参考椭球体:IAG 75
长半轴:6378140
短半轴:2 扁率:298.257000原文地址:
范文二:UTM和经纬度坐标的转换轻松实现坐标转换(UTM和经纬度坐标的转换) 19:15来源:developerWorks 中国很多流行的应用程序可以提供基于位置的服务,但是计算机如何识别真实世界中的位置呢?很多方法都涉及到地理坐标系统,并且在实际应用中存在不同的此类系统。在本文中,应用程序架构师 Sami Salkosuo 演示了使用 Java(TM) 代码在两种流行的系统之间转换位置数据:人们较为熟悉的经纬度系统和统一横轴墨卡托投影(Universal Transverse Mercator)系统。位置服务 ---- 包括基于 GPS 的导航系统和地图站点(如 Google Maps 和 Yahoo! Maps)---- 现在深受客户欢迎。很多企业已经利用了某些位置感知服务,而更多的用户将加入到这个行列中来,因为他们已认识到该服务带来的优势和潜能。在 2006 年,Garter 就曾表示,“位置感知服务在未来两到五年内将成为主流”,并且已经有 “越来越多的组织部署了位置感知移动业务应用程序。”(请参阅 参考资料,获得该报告的链接)。 当企业决定实现某种位置感知应用程序时,编写此类应用程序的任务最终都落在开发人员的身上。构建位置感知服务涉及多种任务,或大或小,其中一项任务(相对较小)可能要将一种系统坐标转换为另一种系统坐标。本文将演示执行此类转换的代码,从而帮助您节省大量的工作。两种不同的坐标系统在详细研究本文代码之前,首先需要讨论即将处理的代码所属的坐标系统:较为熟悉的经纬度系统和统一横轴墨卡托投影系统(Universal Transverse Mercator,UTM)。我们还要提到以 UTM 为基础的军事格网参考系 (MGRS)。经纬度系统经纬度系统可能是最为人熟知的地理坐标设计方法。它使用两个数值表示位置。纬度 表示从地球中心到地球表面东西方向线之间的角度。经度 指从地球中心到地球表面南北方向线之间的角度。经纬度可以表示为十进制角度(DD),或表示为度、分、和秒(DMS);后者的格式可表示为诸如 49°30'00地球以赤道(0° 纬线)为界,分为南半球和北半球,又以 0° 经线(从南极到北极的假想线,通过英国的格林威治市)为界分为东西半球。北半球的纬度从 0 度到 90 度,而南半球的纬度从 0 度到 -90 度。东半球的经度范围从 0 度到 180 度,西半球的经度范围为 0 度到 -180 度。举例说明,坐标 61.44,25.40(使用 DD 单位)或 61°26'24''N,25°23'60''E(使用 DMS 单位)位于芬兰南部。坐标 -47.04, -73.48(使用 DD 单位)或 47°02'24''S,73°28'48''W(使用 DMS 单位)位于智利南部。统一横轴墨卡托投影UTM 坐标系统使用基于网格的方法表示坐标。UTM 系统将地球分为 60 个区,每个区基于横轴墨卡托投影。绘图法中的地图投影方法可以在平面中表示一个两维的曲面,例如一个标准地图。UTM 经度区范围为 1 到 60;其中 58 个区的东西跨度为 6°(稍后详细讨论另外两个区)。经度区涵盖了地球中纬度范围从 80°S 到 84°N 之间的所有区域。一共有 20 个 UTM 纬度区,每个区的南北跨度为 8°;使用字母 C 到 X 标识(其中没有字母 I 和O)。A、B、Y、Z 区不在系统范围以内;它们覆盖了南极和北极区。两个非标准的经度区:32V 区被扩展为覆盖整个挪威的南部,而 31V 区被缩小,所以只覆盖了一片汪洋大海。UTM 坐标的表示格式为:经度区纬度区以东以北,其中以东 表示从经度区的中心子午线的投影距离,而以北表示距离赤道的投影距离。这个两个值的单位均为米。举例来说,使用 UTM 表示经/纬度坐标 61.44,25.40的结果就是 35 V 2844;而经/纬度坐标 -47.04,-73.48 的表示结果为 18 G 9269。军事格网参考系MGRS 是北约(NATO)军事组织使用的标准坐标系统。MGRS 以 UTM 为基础并进一步将每个区划分为 100km × 100 km 的小方块。这些方块使用两个相连的字母标识:第一个字母表示经度区的东西位置,而第二个字母表示南北位置。例如,UTM 点 35 V 2844 等价于 MGRS 点 35VMJ。该 MGRS 点精度为米,使用 15个字符表示,其中最后 10 个字符表示指定网格中的以东和以北的值。可以使用 15 个字符表示 MGRS 值(如前例),也可表示为 13、11、9 或 7 个字符;因此,所表示的值的精度分别为 1 米、10 米、100 米、1,000米或 10,000 米。本文并未对 MGRS 进行详细说明,但是本文的下载代码包含了经纬度坐标和 MGRS 坐标之间的转换。坐标转换确定地球上某个位置的经度和纬度坐标的最低需求是,你至少能够看到星星和太阳,并具备一个六分仪和能够显示 GMT 时间的时钟 T。根据空中某个物体与地平线之间的角度可以确定纬度,然后根据地球旋转计算出经度。本文并未详细讨论这些细节(想要了解更多请参阅 参考资料),相反,我们假设您已经具有DD、DMS 或 UTM 格式的坐标。在十进制角度和度/分/秒格式之间进行转换DD 和 DMS 坐标格式之间的转换非常简单。下面给出了 DD 到 DMS 的转换公式:DD: dd.ffDMS: dd mm ssdd=ddmm.gg=60*ffss=60*gg这里的 gg 代表计算的小数部分。负纬度表示位于南半球(S)的位置而负经度表示西半球(W)的位置。例如,假设您具有一个 DD 格式的坐标 61.44,25.40。按照下面的公式将其转换:lat dd=61lat mm.gg=60*0.44=26.4lat ss=60*0.4=24以及: lon dd=25lon mm.gg=60*0.40=24.0lon ss=60*0.0=0因此,转换为 DMS 格式的坐标变成了 61°26'24''N 25°24'00''E。将 DMS 转换为 DD 格式的公式如下所示:DD: dd.ffDMS: dd mm ssdd.ff=dd + mm/60 + ss/3600注意,南半球(S)的位置为负纬度,西半球(W)位置为负经度。现在将 DMS 格式坐标 47°02'24''S 和 73°28'48''W 转换为 DD 格式的坐标:lat dd.ff= - (47 + 2/60 + 24/3600 )=-47.04lon dd.ff= - (73 + 28/60 + 48/3600)=-73.48转换后的 DD 格式的坐标为 -47.04 和 -73.48。在经纬度和 UTM 坐标之间进行转换十进制坐标可通过一个六分仪和一个记时计确定,与此不同的是,必须通过计算才能确定 UTM 坐标。虽然这些计算无非是最基本的三角形和代数计算,但是所使用的公式非常复杂。如果您阅读了 “The Universal Grids:Universal Transverse Mercator (UTM) and Universal Polar Stereographic (UPS)”(参阅 参考资料 获得链接),就知道它有多复杂了。本文没有给出 UTM 转换公式,但是可从下面一节中给出的源代码中窥探一二,更多信息请参阅 参考资料 提供的链接。使用 Java 代码转换坐标本节介绍了执行坐标转换(十进制角度和 UTM)的库类的源代码。该 Java 类名为com.ibm.util.CoordinateConversion;其思想是构建一个提供转换方法的类。该类包含实际执行转换的内部类;如果需要的话,可以从 CoordinateConversion 类中重构内部类,从而创建一个库包或向现有包添加类。该类执行的转换精度低于 1 米。CoordinateConversion 的源代码包含大约 750 行代码,因此本文没有全部显示。以下小节描述了有关方法,本文的 下载 小节中附带了完整的源代码。CoordinateConversionCoordinateConversion 是主类,它被实例化为在需要是执行坐标转换。清单 1 展示了相关的公共方法,以及CoordinateConversion 类中包含的私有内部类:清单 1. CoordinateConversionpublic class CoordinateConversion{public CoordinateConversion(){}public double[] utm2LatLon(String UTM){UTM2LatLon c = new UTM2LatLon(); return c.convertUTMToLatLong(UTM); }public String latLon2UTM(double latitude, double longitude){LatLon2UTM c = new LatLon2UTM(); return c.convertLatLonToUTM(latitude, longitude);}//..implementation omittedprivate class LatLon2UTM{public String convertLatLonToUTM(double latitude, double longitude){//..implementation omitted}//..implementation omitted}private class LatLon2MGRUTM extends LatLon2UTM{public String convertLatLonToMGRUTM(double latitude, double longitude){//..implementation omitted}//..implementation omitted}private class MGRUTM2LatLon extends UTM2LatLon{public double[] convertMGRUTMToLatLong(String mgrutm){//..implementation omitted}//..implementation omitted}private class UTM2LatLon{public double[] convertUTMToLatLong(String UTM){//..implementation omitted}//..implementation omitted}private class Digraphs{//used to get digraphs when doing conversion between//lat/long and MGRS//..implementation omitted}private class LatZones{//include methods to determine latitude zones//..implementation omitted}下一节将进一步探讨经纬度与 UTM 之间的转换。将经纬度转换为 UTM将经纬度坐标转换为 UTM 坐标需要使用 String latLon2UTM(double latitude, double longitude) 方法。该方法的实现创建了一个新的内部类 LatLon2UTM c = new LatLon2UTM(); 实例,并将 UTM 坐标返回为由 15 个字符组成的字符串(即精度为 1 米)。LatLon2UTM 方法的实现如清单 2 所示:清单 2. public String convertLatLonToUTM(double latitude, double longitude)public String convertLatLonToUTM(double latitude, double longitude){validate(latitude, longitude); String UTM =setVariables(latitude, longitude);String longZone = getLongZone(longitude); LatZones latZones = new LatZones(); String latZone = latZones.getLatZone(latitude);double _easting = getEasting(); double _northing = getNorthing(latitude);UTM = longZone +return UTM; }该方法执行转换的方法为:调用各种方法获得经纬度区,然后计算以东和以北值,等等。使用 validate() 方法对输入进行验证;如果 (latitude < -90.0 || latitude > 90.0 || longitude < -180.0 || longitude >= 180.0)子句为真,将抛出一个 IllegalArgumentException。清单 3 中的 setVariables() 方法设置计算转换所需的各种变量(请查看 “The Universal Grids” 获取更多信息;可从 参考资料 获取链接):清单 3. protected void setVariables(double latitude, double longitude)protected void setVariables(double latitude, double longitude){latitude = degreeToRadian(latitude); rho = equatorialRadius * (1 - e * e) / POW(1 - POW(e * SIN(latitude), 2), 3 / 2.0);nu = equatorialRadius / POW(1 - POW(e * SIN(latitude), 2), (1 / 2.0));double var1; if (longitude < 0.0){var1 = ((int) ((180 + longitude) / 6.0)) + 1; }else{var1 = ((int) (longitude / 6)) + 31; }double var2 = (6 * var1) - 183; double var3 = longitude - var2; p = var3 * 3600 / 10000;S = A0 * latitude - B0 * SIN(2 * latitude) + C0 * SIN(4 * latitude) - D0* SIN(6 * latitude) + E0 * SIN(8 * latitude);K1 = S * k0; K2 = nu * SIN(latitude) * COS(latitude) * POW(sin1, 2) * k0 * () / 2; K3 = ((POW(sin1, 4) * nu * SIN(latitude) * Math.pow(COS(latitude), 3)) / 24)* (5 - POW(TAN(latitude), 2) + 9 * e1sq * POW(COS(latitude), 2) + 4* POW(e1sq, 2) * POW(COS(latitude), 4))* k0* (00000L);K4 = nu * COS(latitude) * sin1 * k0 * 10000;K5 = POW(sin1 * COS(latitude), 3) * (nu / 6)* (1 - POW(TAN(latitude), 2) + e1sq * POW(COS(latitude), 2)) * k0* 0L;A6 = (POW(p * sin1, 6) * nu * SIN(latitude) * POW(COS(latitude), 5) / 720)* (61 - 58 * POW(TAN(latitude), 2) + POW(TAN(latitude), 4) + 270* e1sq * POW(COS(latitude), 2) - 330 * e1sq* POW(SIN(latitude), 2)) * k0 * (1E+24);}清单 4 中的 getLongZone() 方法和 LatZones 类(可从 源代码 获得)用来获得经纬度区。经度区通过longitude 参数计算而来,而纬度区很难使用 LatZones 类中的数组进行编码。清单 4. protected String getLongZone(double longitude)protected String getLongZone(double longitude){double longZone = 0; if (longitude < 0.0){longZone = ((180.0 + longitude) / 6) + 1; }else{longZone = (longitude / 6) + 31; }String val = String.valueOf((int) longZone); if (val.length() == 1){val =getNorthing() 方法(清单 5)和 getEasting() 方法(清单 6)计算以北和以东的值。两种方法都使用 清单 3中的 setVariables() 方法设置的变量。清单 5. protected double getNorthing(double latitude)protected double getNorthing(double latitude){double northing = K1 + K2 * p * p + K3 * POW(p, 4); if (latitude < 0.0){northing =
+ } }清单 6. protected double getEasting()protected double getEasting(){return 500000 + (K4 * p + K5 * POW(p, 3)); }清单 7 包含了一些示例输出,包括一些经纬度坐标和对应的 UTM 坐标:清单 7. Latitude/longitude-to-UTM 测试值( 0.0 )将 UTM 坐标转换为经纬度坐标UTM 坐标到经纬度坐标的转换要比相反的转换过程容易一些。同样,“The Universal Grids”(请参阅 参考资料)提供了转换公式。清单 8 展示了 convertUTMToLatLong() 方法的代码。该方法返回一个双数组,其中第一个元素(数组索引 [0])表示纬度,而第二个元素(数组索引 [1])表示经度。由于 UTM 字符串参数的精度为 1 米,因此经纬度坐标具有与之相同的精度。清单 8. public double[] convertUTMToLatLong(String UTM)public double[] convertUTMToLatLong(String UTM){double[] latlon = { 0.0, 0.0 }; String[] utm = UTM.split(if (hemisphere.equals(latitude = 180 * (phi1 - fact1 * (fact2 + fact3 + fact4)) / Math.PI;if (zone > 0){zoneCM = 6 * zone - 183.0; }else{zoneCM = 3.0; }longitude = zoneCM - _a3; if (hemisphere.equals(latlon[0] = latlon[1] =}convertUTMToLatLong() 方法将传入的 UTM 字符串(格式为 34 G 2631)***,并使用getHemisphere() 方法确定字符串表示的位置所在的半球。这种确定非常简单:纬度区A、C、D、E、F、G、H、J、K、L 和 M 位于南半球,而其余区位于北半球。清单 9 所示的 setVariables() 方法将设置计算所需的变量,然后立即计算纬度值。经度值则通过经度区计算。清单 9. protected void setVariables()protected void setVariables(){arc = northing / k0; mu = arc/ (a * (1 - POW(e, 2) / 4.0 - 3 * POW(e, 4) / 64.0 - 5 * POW(e, 6) / 256.0));ei = (1 - POW((1 - e * e), (1 / 2.0)))/ (1 + POW((1 - e * e), (1 / 2.0)));ca = 3 * ei / 2 - 27 * POW(ei, 3) / 32.0;cb = 21 * POW(ei, 2) / 16 - 55 * POW(ei, 4) / 32; cc = 151 * POW(ei, 3) / 96; cd = 1097 * POW(ei, 4) / 512; phi1 = mu + ca * SIN(2 * mu) + cb * SIN(4 * mu) + cc * SIN(6 * mu) + cd* SIN(8 * mu);n0 = a / POW((1 - POW((e * SIN(phi1)), 2)), (1 / 2.0));r0 = a * (1 - e * e) / POW((1 - POW((e * SIN(phi1)), 2)), (3 / 2.0)); fact1 = n0 * TAN(phi1) / r0;_a1 = 500000 - dd0 = _a1 / (n0 * k0);fact2 = dd0 * dd0 / 2;t0 = POW(TAN(phi1), 2); Q0 = e1sq * POW(COS(phi1), 2); fact3 = (5 + 3 * t0 + 10 * Q0 - 4 * Q0 * Q0 - 9 * e1sq) * POW(dd0, 4) / 24;fact4 = (61 + 90 * t0 + 298 * Q0 + 45 * t0 * t0 - 252 * e1sq - 3 * Q0* Q0)* POW(dd0, 6) / 720;lof1 = _a1 / (n0 * k0); lof2 = (1 + 2 * t0 + Q0) * POW(dd0, 3) / 6.0; lof3 = (5 - 2 * Q0 + 28 * t0 - 3 * POW(Q0, 2) + 8 * e1sq + 24 * POW(t0, 2))* POW(dd0, 5) / 120; _a2 = (lof1 - lof2 + lof3) / COS(phi1); _a3 = _a2 * 180 / Math.PI;}setVariables() 使用以东和以北值设置所需的变量。这些都是类变量并且在 convertUTMToLatLong(StringUTM) 方法中进行设置(参见 清单 8)。其他方法源代码 还提供了其他公共和私有方法以及类。例如,提供了可对经纬度和 MGRS 进行坐标转换的方法和类,以及执行度和弧度之间转换的辅助方法,还提供了各种数学操作函数(例如 POW、SIN、COS 和 TAN)。 结束语本文简单介绍了有关世界坐标系统的一些知识,并提供了执行坐标转换的 Java 类。尽管没有详细介绍所有的坐标转换公式,您可以从 参考资料 小节了解详细内容。一般来说,日常的开发工作并不需要了解这些理论---- 只有极少数情况下需要,正如我最近遇到的坐标转换任务一样。我需要在经纬度、UTM 和 MGRS 之间进行坐标转换,因此我做了些基础研究并使用 Java 类实现了转换。开发工作花费了我好几个小时的时间,我希望本文能帮助您在执行其他任务时节省时间,并且CoordinateConversion 类能为您提供帮助。关于作者Sami Salkosuo 从 1999 年起一直在 IBM 工作。他是 Sun 认证的 Java 程序员,IBM 认证的 WebSphereMessage Broker 解决方案开发人员,还是 IBM 认证的 WebSphere MQ 的解决方案设计师。阅读详情:
范文三:大地坐标转经纬度和经纬度换大地坐标地质工作中常要对进行大地坐标转经纬度和经纬度换大地坐标,以下步骤请大家熟记:一、大地座标→经纬度(地理坐标)1、在文本文件中输入大地坐标数据,格式为Y空格X。如下,原始的大地坐标由一个8位的Y和一个7位的X组成,“新建文本文档.txt -记事本”显示如下:3500350075007500这组坐标数据中的Y的前两位为31,是分带号,一般使用的分带有三分带,六分带,这里的坐标是三分带的,记下Y前的这两位数,在原始数据中去除掉,
现在数据变为:Y--6位,X--7位。“新建文本文档.txt -记事本”显示如下:
3500350075007500保存这个TXT的文本文件。2、打开MAPGIS,启动坐标投影变形程序如果是MAPGIS6.7版,请选择“实用服务→投影变换系统→用户文件投影转换”→点击打开文件,打开刚才的大地坐标的文本文件。“指定数据起始位置”中出现刚才的的文本文档,显示如下:3500350075007500在设置用户文件选项中,一般选:按行读取数据,X→Y顺序,生成点。最后点击确定。3、设置输入数据的格式,点击用户投影参数,并完成设置。坐标系类型----大地坐标系投影类型----5:高斯克吕格投影比例尺分母----1椭球面高程----0投影面高程----0投影带类型----3度带或6度带投影带序号----31X,Y的平移均设0这里我们的大地座标为3度带的第31带,注意填好,坐标单位为米接着为:设置输出的格式,我们要求输出的是经纬度,点结果转换参数,完成设置。4、输入投影参数坐标系类型----地理坐标系我们输出的经纬度的单位应该是DDDMMMSS。SS注意点写到文件,保存就大功告成了,注意:保存的文件要写上.TXT的后缀最后,在文本文件中计算出的结果如下:x= 560000 y= 4503500
yp=x= 565000 y=4503500
yp=x=565000 y=4507500
yp=X=568500 y=4507500
yp=xp为经度,7就是123度42分34。357秒,YP为纬度,5就是40度39分50。255秒(纬度没有最多90,所以没有三位数)二、经纬度→大地座标同样,输入文本文件格式如下,1234234
4039501234607
4039481234608
4041581234837
403157这里面的数据前面的为经度,格式为DDDMMSS,后面的为纬度,格式为DDMMSS接下来的转换过程和大地坐标转换一样,只要将刚才的用户转换参数和结果转换参数交换即可,要注意分带号的确定,如果你不知道分带号,就应该先计算分带号,算法是
经度/3得到的整数为三度带的分带号经度/6得到的整数为六度带的分带号计算所得的结果格式如下x=1234234
yp=。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。其中的xp为地图上的y坐标,记得在前面加上带号,其中的yp为地图上的x坐标多练习几次,很快熟悉哦!阅读详情:
范文四:3、Mapgis大地坐标转经纬度和经纬度换大地坐标地质工作中常要对进行大地坐标转经纬度和经纬度换大地坐标,以下步骤请大家熟记:一、大地坐标→经纬度(地理坐标)1、在文本文件中输入大地坐标数据,格式为Y空格X。如下,原始的大地坐标由一个8位的Y和一个7位的X组成,“新建文本文档.txt -记事本”显示如下:3500350075007500这组坐标数据中的Y的前两位为31,是分带号,一般使用的分带有三分带,六分带,这里的坐标是三分带的,记下Y前的这两位数,在原始数据中去除掉, 现在数据变为:Y--6位,X--7位。“新建文本文档.txt -记事本”显示如下:
3500350075007500保存这个TXT的文本文件。2、打开MAPGIS,启动坐标投影变形程序如果是MAPGIS6.7版,请选择“实用服务→投影变换系统→用户文件投影转换”,点击打开文件,打开刚才的大地坐标的文本文件。“指定数据起始位置”中出现刚才的的文本文档,显示如下:在设置用户文件选项中,一般选:按行读取数据,X→Y顺序,生成点。 接着设置输入数据的格式,点击用户投影参数,并完成设置:坐标系类型----大地坐标系投影类型----5:高斯克吕格投影比例尺分母----1椭球面高程----0投影面高程----0投影带类型----3度带或6度带投影带序号----31X,Y的平移均设0这里我们的大地坐标为3度带的第31带,注意填好,坐标单位为米。填好以后设置输出的格式,我们要求输出的是经纬度,点击“结果投影参数”,完成以下设置。坐标系类型----地理坐标系我们输出的经纬度的单位应该是DDDMMMSS.SS(度.分.秒),注意点击投影变换、写到文件,保存就大功告成了,注意:保存的文件要写上.TXT的后缀
最后,在文本文件中计算出的结果如下:x= 560000 y= 4503500
yp=x= 565000 y=4503500
yp=x=565000 y=4507500
yp=X=568500 y=4507500
yp=xp为经度,就是93度42分34.357秒,yp为纬度,就是40度39分50.255秒(纬度最大为90,所以没有三位数)二、经纬度→大地座标同样,输入文本文件格式如下,1234234
4039501234607
4039481234608
4041581234837
403157这里面的数据前面的为经度,格式为DDDMMMSS.SS,后面的为纬度,格式为DDDMMMSS.SS,接下来的转换过程和大地坐标转换一样,只要将刚才的用户转换参数和结果转换参数交换即可,要注意分带号的确定,如果你不知道分带号,就应该先计算分带号,算法是经度/3得到的整数为三度带的分带号(经度+3)/6得到的整数为六度带的分带号计算所得的结果格式如下x=1234234
yp=。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。。其中的xp为地图上的y坐标,记得在前面加上带号,其中的yp为地图上的x坐标阅读详情:
范文五:MAPGIS经纬度坐标转换为大地坐标MAPGIS经纬度坐标转换为大地坐标1 打开MAPGIS【实用服务】-【投影变换】2 【投影变换】-【P输入单点投影坐标】后,出现如下界面:打开【P输入单点投影坐标】以西藏某点为例(15度带)X: 492791,Y: 3707650为例(X转换为经度,Y为纬度):设置原始投影参数和结果投影参数:原始投影参数:结果投影参数:修改后点击【投影点】,投影后如图所示:MAPGIS把经纬度坐标转换为大地坐标投影变换下的【投影转换】菜单下【输入单点投影转换】。设置“原始投影参数”和“结果投影参数”,并将已知点输进去“投影点”,影转换模块,投影转换菜单下,输入单点投影变换功能。设置当前投影:地理坐标系,单位可以是度,分,秒或ddmmss格式。根据数据决定。如数据是98.78度,那么你的单位就是度,依次类推。设置目的投影:投影平面直角坐标系,高斯投影,比例尺分母是1,单位是米,根据你的经度范围输入中央经度。其他不用设置,点击投影点按钮,在右边就计算出该点的大地坐标。阅读详情:
范文六:北京54坐标和经纬度坐标转换算法(C++)北京54坐标和经纬度坐标转换算法c++代码//坐标正算lbxy(double l, double b, double *x, double *y, int l0){double sa,sb,sep,sn,sy2,st,sm,sx,double xx,yy,hd,//判断值的范围if (l > 360 || l
360 || b{*x =*y =}l = l-l0;sa = 6378245;sb = ;sep= 0.;hd = b*PI;hb = hd/180.0;st = tan(hb);sn=pow(sa,(double)2)/sqrt(pow(sa,(double)2)*pow(cos(hb),(double)2)+pow(sb,(double)2)*pow(sin(hb),(double)2));sy2=sep*pow(cos(hb),(double)2);sd = cos(hb)*l*PI;sm = sd/180.0;sx = *b-(32005.78*sin(hb)+133.924*pow(sin(hb),(double)3)+0.697*pow(sin(hb),(double)5))*cos(hb);xx = sx+sn*st*(0.5*pow(sm,(double)2)+1.0/24.0*(5.0-pow(st,(double)2)+9.0*sy2)*pow(sm,(double)4));yy = sn*(sm+1.0/6.0*(1.0-pow(st,(double)2)+sy2)*pow(sm,(double)3)+1.0/120.0*(5.0-18.0*pow(st,(double)2)+pow(st,(double)4))*pow(sm,(double)5));*x =*y = yy+500000;}//坐标反算xylb(double l0, double x, double y, double *l, double *b){double bf,vf,nf,ynf,tf,yf2,double sa,sb,se2,sep2,double w1,w2,w,w3,w4;double pi = 3.1415926;x = x/;y = y - ;bf = 9.*x-0.*pow(x,2.0)-0.*pow(x,3.0)-0.*pow(x,4.0)+0.*pow(x,5.0)-0.*pow(x,6.0);hbf = bf * pi/ 180.0;sa = ;sb = ;se2 = 0.;sep2 = 0.;w1 = sin(hbf);w2 = 1.0 - se2 * pow(w1,(double)2);w = sqrt(w2);mf = sa*(1.0-se2)/pow(w,(double)3);w3 = cos(hbf);w4 = pow(sa,(double)2)*pow(w3,(double)2) + pow(sb,(double)2)*pow(w1,(double)2);nf = pow(sa,(double)2) / sqrt(w4);ynf = y/vf = nf/tf = tan(hbf);yf2 = sep2 * pow(w3, (double)2);*b = bf - 1.0/2.0 * vf * tf * (pow(ynf,(double)2)-1.0/12.0*(5.0+3.0*pow(tf,(double)2)+yf2-9.0*yf2*pow(tf,(double)2))*pow(ynf,(double)4))*180.0/*l = 1.0/w3*ynf*(1.0-1.0/6.0*(1.0+2.0*pow(tf,(double)2)+yf2)*pow(ynf,(double)2)+1.0/120.0*(5.0+28.0*pow(tf,(double)2)+24.0*pow(tf,(double)2)+6.0*yf2+8.0*yf2*pow(tf,(double)2))*pow(ynf,(double)4))*180.0/*l = l0 + *l;}delphi:unit Unit10;interfaceusesWindows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,StdCtrls, Buttons, ExtCtypeTForm10 = class(TForm)RadioGroup1: TRadioGGroupBox1: TGroupBLabel1: TLLabel2: TLLabel3: TLLabel4: TLLabel5: TLLabel6: TLLabel7: TLLabel8: TLLabel9: TLLabel10: TLLabel11: TLLabel12: TLEdit1: TEEdit2: TEEdit3: TEEdit4: TEEdit5: TEEdit6: TEEdit7: TEEdit8: TEEdit9: TEEdit10: TEButton1: TBButton2: TBButton3: TBprocedure Button1Click(Sender: TObject);procedure FormShow(Sender: TObject);procedure RadioGroup1Click(Sender: TObject);function fc(data1:doucs:integer):procedure calsub1();procedure calsub2();procedure calsub3();procedure calsub4();procedure calsub5(xzid:integer);procedure Button2Click(Sender: TObject);procedure Button3Click(Sender: TObject);private{ Private declarations }public{ Public declarations }varForm10: TForm10;///////////////全局变量var x,y,bb,w:var x1,b,p,p1,p2,t,f,n,p3,p4,p5,m:var x2,yyy,p6:var k,l1,l2,l3,b1,b2,b3,xs,ks,ys,kk:var jh,xzb,kzb,yzb,l1zb,l2zb,l3zb,b1zb,b2zb,b3zb://////////////////////implementation{$R *.DFM}function TForm10.fc(data1:cs:integer):var i:fhz:begin//计算方次(乘方)fhz:=1.0;for i:=1 to cs dofhz:=fhz*data1;fc:=procedure TForm10.calsub1();begin//sub1x1:=1*b-*sin(2*b*bb);x1:=x1+16.8281*sin(4*b*bb);x1:=x1-0.022*sin(6*b*bb);x1:=x1+0.00003*sin(8*b*bb);procedure TForm10.calsub2();begin//sub2p:=1/180*3.1415926;p1:=n*sin(b*bb)*cos(b*bb)*p*p/2;p2:=(5-t*t+9*f+4*sqr(f))*n*sqr(p)*sqr(p)*sin(b*bb)*(cos(b*bb)*sqr(cos(b*bb)))/24;p3:=(61-58*fc(t,6))*n*fc(p,6)*sin(b*bb)*fc(cos(b*bb),5)/720;p4:=n*cos(b*bb)*p+(1-sqr(t)+f)*n*fc(cos(b*bb),3)*fc(p,3)/6;p5:=(5-18*sqr(t)+fc(t,4)+14*f-58*sqr(t)*f)*n*fc(cos(b*bb),5)*fc(p,5)/120;procedure TForm10.calsub3();begin//sub3w:=sqrt(1-0.*sqr(sin(b*bb)));m:=./fc(w,3);procedure TForm10.calsub4();begin//sub4n:=6378245/w;f:=0.*sqr(cos(b*bb));t:=sin(b*bb)/cos(b*bb);procedure TForm10.calsub5(xzid:integer);var l,lb,zz:begin//sub5bb:=3.;if xzid=3 then //由经纬度计算三、六带坐标beginb:=b1+b2/60+b3/3600;l:=l1+l2/60+l3/3600;if ll:=l+3elsel:=l-3;calsub1;calsub3;calsub4;calsub2;x:=x1+p1+p2+p3;y:=p4+p5+500000;xs:=int(x*100+0.5)/100;ks:=int((l1+3)/6+0.5);ys:=int(y*100+0.5)/100;kk:=int(l1/3+0.5);if (xzid=1) or (xzid=2) then //由三、六带坐标计算 经纬度beginxs:=x;ks:=k;ys:=y;y:=y-500000;b:=x/1;calsub1;x2:=x-x1;b1:=b;while abs(x2)>=0.001 dobegincalsub3;b:=b+x2*(180/3.1415926)/m;calsub1; //调用calsub1后,返回x1的值x2:=x-x1;b1:=b;calsub4; //调用calsub4后,返回n,f,t的值p1:=t*sqr(y)/(m*n*2);p2:=(5+3*sqr(t)+f-9*f*sqr(t))*t*fc(y,4)/(24*m*fc(n,3));yyy:=fc(y,6)/(720*m*fc(n,5));p3:=(61+90*sqr(t)+45*fc(t,4))*t*p4:=y/(n*cos(b*bb));p5:=(1+2*sqr(t)+f)*fc(y,3)/(6*fc(n,3)*cos(b*bb));p6:=(5+28*sqr(t)+2*4*fc(t,4)+6*f+8*f*sqr(t))*fc(y,5)/(120*fc(n,5)*cos(b*bb));b:=b1-(p1-p2+p3)*180/3.1415926;l:=(p4-p5+p6)*180/3.1415926;zz:=l;if l>0 thenl:=l-3elsel:=l+3;calsub1;calsub3;calsub4;calsub2;x:=int((x1+p1+p2+p3)*100+0.5)/100;y:=int((p4+p5+0+0.5)/100;if (xzid=1) or (xzid=2) thenbeginl1:=int(l);b1:=int(b);lb:=(l-l1)*60;bb:=(b-b1)*60;l2:=int(lb);b2:=int(bb);l3:=(lb-l2)*60;b3:=(bb-b2)*60;l3:=int(l3*100+0.5)/100;b3:=int(b3*100+0.5)/100;if (xzid=1) and (zzl1:=l1+(k-1)*3;if (xzid=1) and (zz>=0) thenl1:=l1+(k+1)*3;if (xzid=2) and (zzl1:=l1+(k-1)*6;if (xzid=2) and (zz>=0) thenl1:=l1+(k)*6;if xzid=1 thenkk:=int((l1+3)/6+0.5);if xzid=2 thenkk:=int(l1/3+0.5);procedure TForm10.Button1Click(Sender: TObject);//var a1,a2,a3:begin//计算case radiogroup1.ItemIndex of0: begin////////////////////////////////////// inputxzb:=trim(edit1.Text);kzb:=copy(trim(edit2.Text),1,2);yzb:=copy(trim(edit2.Text),3,length(trim(edit2.Text))-2);x:=strtofloat(xzb);y:=strtofloat(yzb);k:=strtofloat(kzb);//////////////////////////////////////////calculatecalsub5(1); //三度带时,调用calsub5过程////////////////////////////////////////// outputedit3.Text:=floattostr(x);edit4.Text:=trim(floattostr(kk))+floattostr(y);edit5.Text:=floattostr(l1); //经度的度edit6.Text:=floattostr(l2); //经度的分edit7.Text:=floattostr(l3); //经度的秒edit8.Text:=floattostr(b1); //纬度的度edit9.Text:=floattostr(b2); //纬度的分edit10.Text:=floattostr(b3); //纬度的秒1: beginxzb:=trim(edit3.Text);kzb:=copy(trim(edit4.Text),1,2);yzb:=copy(trim(edit4.Text),3,length(trim(edit4.Text))-2);x:=strtofloat(xzb);y:=strtofloat(yzb);k:=strtofloat(kzb);//////////////////////////////////////////calculatecalsub5(2); //六度带时,调用calsub5////////////////////////////////////////// outputedit1.Text:=floattostr(x);edit2.Text:=trim(floattostr(kk))+floattostr(y);edit5.Text:=floattostr(l1); //经度的度edit6.Text:=floattostr(l2); //经度的分edit7.Text:=floattostr(l3); //经度的秒edit8.Text:=floattostr(b1); //纬度的度edit9.Text:=floattostr(b2); //纬度的分edit10.Text:=floattostr(b3); //纬度的秒2: beginl1zb:=trim(edit5.Text);l2zb:=trim(edit6.Text);l3zb:=trim(edit7.Text);b1zb:=trim(edit8.Text);b2zb:=trim(edit9.Text);b3zb:=trim(edit10.Text);l1:=strtoint(l1zb);l2:=strtoint(l2zb);l3:=strtoint(l3zb);b1:=strtoint(b1zb);b2:=strtoint(b2zb);b3:=strtoint(b3zb);//////////////////////////////////////////calculutecalsub5(3);////////////////////////////////////////// outputedit1.Text:=floattostr(x);edit2.Text:=trim(floattostr(kk))+floattostr(y);edit3.Text:=floattostr(xs);edit4.Text:=trim(floattostr(ks))+floattostr(ys);procedure TForm10.FormShow(Sender: TObject);begin//初始化radiogroup1.ItemIndex:=0;edit1.Text:='';edit2.Text:='';edit3.Text:='';edit4.Text:='';edit5.Text:='';edit6.Text:='';edit7.Text:='';edit8.Text:='';edit9.Text:='';edit10.Text:='';edit1.Enabled:=edit2.Enabled:=edit3.Enabled:=edit4.Enabled:=edit5.Enabled:=edit6.Enabled:=edit7.Enabled:=edit8.Enabled:=edit9.Enabled:=edit10.Enabled:=edit1.SetFprocedure TForm10.RadioGroup1Click(Sender: TObject);begin//选择项case radiogroup1.ItemIndex of0: beginedit1.Text:='';edit2.Text:='';edit3.Text:='';edit4.Text:='';edit5.Text:='';edit6.Text:='';edit7.Text:='';edit8.Text:='';edit9.Text:='';edit10.Text:='';edit1.Enabled:=edit2.Enabled:=edit3.Enabled:=edit4.Enabled:=edit5.Enabled:=edit6.Enabled:=edit7.Enabled:=edit8.Enabled:=edit9.Enabled:=edit10.Enabled:=edit1.SetF1: beginedit1.Text:='';edit2.Text:='';edit3.Text:='';edit4.Text:='';edit5.Text:='';edit6.Text:='';edit7.Text:='';edit8.Text:='';edit9.Text:='';edit10.Text:='';edit3.Enabled:=edit4.Enabled:=edit1.Enabled:=edit2.Enabled:=edit5.Enabled:=edit6.Enabled:=edit7.Enabled:=edit8.Enabled:=edit9.Enabled:=edit10.Enabled:=edit3.SetF2: beginedit1.Text:='';edit2.Text:='';edit3.Text:='';edit4.Text:='';edit5.Text:='';edit6.Text:='';edit7.Text:='';edit8.Text:='';edit9.Text:='';edit10.Text:='';edit1.Enabled:=edit2.Enabled:=edit3.Enabled:=edit4.Enabled:=edit5.Enabled:=edit6.Enabled:=edit7.Enabled:=edit8.Enabled:=edit9.Enabled:=edit10.Enabled:=edit5.SetFprocedure TForm10.Button2Click(Sender: TObject);beginprocedure Tform10.Button3Click(Sender: TObject);beginedit1.Text:='';edit2.Text:='';edit3.Text:='';edit4.Text:='';edit5.Text:='';edit6.Text:='';edit7.Text:='';edit8.Text:='';edit9.Text:='';edit10.Text:='';case radiogroup1.ItemIndex of0: edit1.SetF1: edit3.SetF2: edit5.SetFend阅读详情:
范文七:北京54坐标和经纬度坐标转换算法(C++)//坐标正算lbxy(double l, double b, double *x, double *y, int l0){double sa,sb,sep,sn,sy2,st,sm,sx,double xx,yy,hd,//判断值的范围if (l > 360 || l
360 || b{*x =*y =}l
= l-l0;sa = 6378245;sb = ;sep= 0.;hd = b*PI;hb = hd/180.0;st = tan(hb);sn=pow(sa,(double)2)/sqrt(pow(sa,(double)2)*pow(cos(hb),(double)2)+pow(sb,(double)2)*pow(sin(hb),(double)2));sy2=sep*pow(cos(hb),(double)2);sd = cos(hb)*l*PI;sm = sd/180.0;sx = *b-(32005.78*sin(hb)+133.924*pow(sin(hb),(double)3)+0.697*pow(sin(hb),(double)5))*cos(hb);xx = sx+sn*st*(0.5*pow(sm,(double)2)+1.0/24.0*(5.0-pow(st,(double)2)+9.0*sy2)*pow(sm,(double)4));yy = sn*(sm+1.0/6.0*(1.0-pow(st,(double)2)+sy2)*pow(sm,(double)3)+1.0/120.0*(5.0-18.0*pow(st,(double)2)+pow(st,(double)4))*pow(sm,(double)5));*x =*y = yy+500000;}//坐标反算xylb(double l0, double x, double y, double *l, double *b){double bf,vf,nf,ynf,tf,yf2,double sa,sb,se2,sep2,double w1,w2,w,w3,w4;double pi = 3.1415926;x = x/;y = y - ;bf = 9.*x-0.*pow(x,2.0)-0.*pow(x,3.0)-0.*pow(x,4.0)+0.*pow(x,5.0)-0.*pow(x,6.0);hbf = bf * pi/ 180.0;sa = ;sb = ;se2 = 0.;sep2 = 0.;w1 = sin(hbf);w2 = 1.0 - se2 * pow(w1,(double)2);w = sqrt(w2);mf = sa*(1.0-se2)/pow(w,(double)3);w3 = cos(hbf);w4 = pow(sa,(double)2)*pow(w3,(double)2) + pow(sb,(double)2)*pow(w1,(double)2);nf = pow(sa,(double)2) / sqrt(w4);ynf = y/vf = nf/tf = tan(hbf);yf2 = sep2 * pow(w3, (double)2);*b = bf - 1.0/2.0 * vf * tf * (pow(ynf,(double)2)-1.0/12.0*(5.0+3.0*pow(tf,(double)2)+yf2-9.0*yf2*pow(tf,(double)2))*pow(ynf,(double)4))*180.0/*l = 1.0/w3*ynf*(1.0-1.0/6.0*(1.0+2.0*pow(tf,(double)2)+yf2)*pow(ynf,(double)2)+1.0/120.0*(5.0+28.0*pow(tf,(double)2)+24.0*pow(tf,(double)2)+6.0*yf2+8.0*yf2*pow(tf,(double)2))*pow(ynf,(double)4))*180.0/*l = l0 + *l;}阅读详情:
范文八:经纬度坐标转换成屏幕坐标经纬度坐标转换成屏幕坐标地理坐标定义规则:X轴(代表经度)向右递增,Y轴(纬度)向上递增,就好比小学学过的平面坐标。向左、向下的规则。屏幕坐标定义规则:X轴向右递增,Y轴向下递增。可以看出,地理坐标和屏幕坐标的区别仅仅只是在于Y轴递增方向是相反的(这就是不同)。这里强调一点的就是为了保证精度,地理坐标的度*3600换算成秒,所有的取值用double来计算,最后的结果再转换成int。1 已知道屏幕的高(y)和宽(h),地理坐标区域的范围(maxLon,minLon,maxLat,minLat),这里我们知道了这些已知的参数。 2 我们可以算出每像素所代表的经度和纬度(有人称这个为比例因子)。公式:scaleX = ((maxLon-minLon)*3600)/h ----------X轴上每像素代表的经度秒数;公式:scaleY = ((maxLat-minLat)*3600)/y ----------Y轴上每像素代表的纬度秒数;这两个比例因子就是两个坐标系之间的关系。3 很简单的一步了,那就是算出该地理坐标区域中的任何一点(lon,lat)在屏幕上的坐标了。公式:screenX = lon*3600/scaleX;---------屏幕坐标X轴坐标 公式:screenY = lat*3600/scaleY; ---------屏幕坐标Y轴坐标还有最后一步,那就是我们要把该地理区域占满占个屏幕该怎么办呢? 4 接着我们需要该地理区域占满占个屏幕该怎么办呢公式:minX = minLon*3600/scaleX;区域左边置最左端公式:minY = minLat*3600/scaleY; 区域上面置最上端5 当地地理范围区域占满整个屏幕时,我们需要用到第三步计算出来的 screenX和screenY两个参数,该区域中的任何一点的公式如下:
公式:X = screenX - minX = (lon - minLon)*3600/scaleX;
公式:Y = screenMaxLat - screenLat = (maxLat - lat)*3600/scaleY; 6 总结:经纬度转屏幕坐标的最终公式如下:公式:X = (lon - minLon)*3600/scaleX;公式:Y = (maxLat - lat)*3600/scaleY;接着我们由上面的公式可以推出屏幕坐标转经纬度坐标公式如下:
公式:lon = X * scaleX/3600 + minLon;公式:lat = maxLat - y* scaleY/3600阅读详情:
范文九:像平面坐标到经纬度坐标的转换地理空间信息GEOSPATIALINFORMATION坐标转换是几乎每个遥感、GIS、航测等软件都具备的一个基本功能[1],在功能代码实现过程中会遇到大规模的矩阵运算问题,例如一幅图像为甚至更多个像元,一般在用一些开发软件进行编码的时候显得计算量巨大,从而导致计算运行速度慢的问题,也是大型商用软件最忌讳的问题。IDL(InteractiveDataLanguage)交互式数据语言是由Kodak公司的全资子公司RSI开发并投向市场的、面向矩阵的且用于数据可视化研究与应用开发的第四代计算机语言。其行业涉及到航天、军事、医学、地球科学、天文学、商用开发等各个领域。IDL的长期成功来源于它的灵活多样的操作模式。它可以快速地分析数据和实现可视化。它还拥有强大的程序环境以及成熟的跨平台终端应用商业软件产品,例如高级遥感影像处理系统(ENVI)就是由IDL编写的。大型的数组矩阵运算是IDL的强势,在IDL中可以为任何IDL数据类型创建数组。紧凑的数组语法能保证数组运行时不用循环的操作。另外,数组操作的优点充分体现在运行速度上,将数组作为整体进行数组操作,其速度远远大于对数组元素的循环操作。一收稿日期:项目来源:国家863计划资助项目();国防科技工业民用专项科研技术研究资助项目(KJSX0601);64地理空间信息GEOSPATIALINFORMATION个成为高效的IDL程序员的诀窍就是灵活有效地利用数组[2]。选择IDL做为开发工具去实现坐标转换这一涉及大量矩阵数组的运算不仅能大大加快运行速度,而且代码通俗易读。1算法分析设计1.1算法分析目前,大多数的坐标转换都是用四参数或者七参数方法实现的,传统的方法适用于大多数情况。而本文所要做的像平面坐标到经纬度坐标的转换不需要用这种传统的方法,而是用一种简单的数学思维去实现此功能一种方法。[3-4]在进行算法设计之前,应该考虑如下几个问题:1)地球是不规则的球体,两极扁,赤道鼓。不能简单地做为圆球体来计算;[5]2)每纬度在地理位置上的距离是一致的(约111km),但是每经度之间的距离是不同的(越往两极,东西方向上的距离就越小),不能直接转换为直角坐标系;3)给出的经纬度坐标发生变化后,需要动态的缩放所显示地图的大小。1.2算法流程设计算法结构图如图1所示,其中输入参数有lat0,lon0,m,n,col,row,divisor:lat0代表已知点的纬度;lon0代表已知点的经度;m代表已知点在像平面坐标系中的行坐标;n代表代表已知点在像平面坐标系中的列坐标;col代表图像的像元列数;row代表图像的像元行数;divisor代表图像的分辨率。返回的是一个col*row的结构体数组,结构体中每个字段都是col*row的数组。图1算法基本流程2算法的详细设计1)创建假定的经纬度坐标系统。分别创建一个代表纬度的col*row数组和一个代表经度的row*col数Feb.,2010Vol.8,No.1组,在这里所建立的这两个数组构成了一个假定的像平面坐标系。而这个假定的坐标系统就是我们以后计算的核心数据。建立的两个数组如图2所示,假定的是col=100,row=100,即假定图像为100*100的。这两个数组看起来虽然简单,但却蕴含了北半球的经纬度的走向规律,以后作为运算的一部分。创建经纬度数组的代码如下:(此代码只能实现行列数相同的情况,如果图像行列数不相同,可都定义为行和列数较大的那个数来建立方阵,然后根据需要的col列row行截取出来经纬度数组,此代码略。)图2建立的索引数组(左为纬度数组,右为经度数组)IDL>arrlat=indgen(col);建立一个有col个元素的col*1的索引数组IDL>arrlat=reform(arrlat,1,col);变为1*col的索引列数组IDL>arrlat=rebin(arrlat,col,row,/sample);扩展为col*row的索引数组IDL>arrlon=transpose(arrlat);转置上面的数组,经度数组创建完毕IDL>arrlat=rotate(arrlat,2);旋转数组,纬度数组创建完毕由此两个数组确定了一个假定的经纬度坐标系,如图3所示。图3假定的经纬度坐标系2)计算所有像元点到横轴和纵轴方向上的地理距离。IDL>Dx=(arrlat-m)*divisorFeb.,2010Vol.8,No.1地理空间信息GEOSPATIALINFORMATION可以明显地看出,当矩阵运算不是很大的时候,数组整体运算同For循环运算比较优势不是很大,但是当矩阵运算越来越大的时候,数组整体运算的优势性就明显表现出来了,速度要比For循环操作快了1066地理空间信息GEOSPATIALINFORMATIONFeb.,2010Vol.8,No.1表2精度的测试转换后的经纬度坐标(35°22′07.23″N,119°34′22.21″E)(35°22′46.39″N,119°33′56.48″E)(35°22′24.21″N,119°33′17.27″E)(39°22′34.63″N,119°32′59.56″E)多倍甚至更多,且随着运算量的增加而递增。3.2精度的验证特征点位GoogleEarth的经纬度坐标P1P2P3P4(35°22′06.81″N,119°34′21.76″E)(35°22′45.58″N,119°33′55.27″E)(35°22′23.50″N,119°33′16.84″E)(35°22′33.95″N,119°32′59.01″E)图4上为经过坐标转换的航空摄影正射图像下为有GoogleEarth的QuickBrid影像如图4所示,准备两幅有相同区域的遥感影像,一幅是经过坐标转换的航空摄影正射影像(中国科学院遥感应用研究所于2008年12月在烟台航空摄影所得),另外一幅有正确的经纬度坐标的QuickBrid影像。从每幅影像中均匀地挑出几个特征明显的点进行比较,比对结果如表2所示。阅读详情:
范文十:经纬度坐标与高斯坐标的转换代码/* 功能说明: 将绝对高斯坐标(y,x)转换成绝对的地理坐标(wd,jd)。
输入参数: 高斯坐标的横坐标,以米为单位//
输入参数: 高斯坐标的纵坐标,以米为单位// short
输入参数: 带号,表示上述高斯坐标是哪个带的// double *L;
输出参数: 指向经度坐标的指针,其中经度坐标以秒为单位// double *B;
输出参数: 指向纬度坐标的指针,其中纬度坐标以秒为单位void GaussToGeo(double y, double x, short DH, double *L, double *B, double LP){double l0;
tf = tg(Bf0),注意要将Bf转换成以弧度为单位
n = y * sqrt( 1 + etf ** 2) / c, 其中etf = e'**2 * cos(Bf0) ** 2double t_l0;
l0,经差,以度为单位double t_B0;
B0,纬度,以度为单位double Bf0;
etf,其中etf = e'**2 * cos(Bf0) ** 2double X_3 ;double PI=3.79;double b_e2=0.7;double b_c=78271;X_3 = x /
// 以兆米(1000000)为单位// 对于克拉索夫斯基椭球,计算Bf0Bf0 = 27. + 9. * X_3 - 0. * pow(X_3,2)
- 0. * pow(X_3,3) + 0. * pow(X_3,4)
+ 0. * pow(X_3,5) - 0. * pow(X_3,6) ;
tf = tan(Bf0*PI/180);
tf = tg(Bf),注意这里将Bf转换成以弧度为单位etf = b_e2 * pow(cos(Bf0*PI/180),2);
etf = e'**2 * cos(Bf) ** 2nf = y * sqrt( 1 + etf ) / b_c;
n = y * sqrt( 1 + etf ** 2) / c// 计算纬度,注意这里计算出来的结果是以度为单位的t_B0 = Bf0 - (1.0+etf) * tf / PI * (90.0 * pow(nf,2)- 7.5 * (5.0 + 3 * pow(tf,2) + etf - 9 * etf * pow(tf,2)) * pow(nf,4)
+ 0.25 * (61 + 90 * pow(tf,2) + 45 * pow(tf,4)) * pow(nf,6))
// 计算经差,注意这里计算出来的结果是以度为单位的t_l0 = (180 * nf - 30 * ( 1 + 2 * pow(tf,2) + etf ) * pow(nf,3)
+ 1.5 * (5 + 28 * pow(tf,2) + 24 * pow(tf,4)) * pow(nf,5))
/ ( PI * cos(Bf0*PI/180) ) ;l0 = (t_l0 * 3600.0);
将经差转成秒if (LP == -1000){*L = (double)((DH * 6 - 3) * 3600.0 + l0);
// 根据带号计算出以秒为单位的绝对经度,返回指针}else{*L = LP * 3600.0 + l0;
// 根据带号计算出以秒为单位的绝对经度,返回指针
}//----------------------------------*B = (double)(t_B0 * 3600.0) ;
将纬差转成秒,并返回指针}/* 功能说明: (1)将地理坐标(wd,jd)转换成绝对的高斯坐标(y,x)(2)本函数支持基于六度带(或三度带)、克拉索夫斯基椭球进行转换
*//* 适用范围: 本函数适用于将地球东半球中北半球(即东经0度到东经180度,北纬0度至90度)范围内所有地理坐标到高斯坐标的转换
*//* 使用说明: 调用本函数后返回的结果应在满足精度的条件下进行四舍五入
输入参数: 地理坐标的经度,以秒为单位//
输入参数: 地理坐标的纬度,以秒为单位// short
输入参数: 三度带或六度带的带号/*
六度带(三度带)的带号是这样得到的:从东经0度到东经180度自西向东按每6度(3度)顺序编号(编号从1开始),这个顺序编号就称为六度带(三度带)的带号。因此,六度带的带号的范围是1-30,三度带的带号的范围是1-60。如果一个点在图号为TH的图幅中,那麽该点所处的六度带的带号就可以这样得到:将该图号的第3、4位组成的字符串先转换成数字,再减去30。例如某点在图幅中,该点所在的带号就是49-30,即19。如果调用本函数去进行一般的从地理坐标到基于六度带高斯坐标的变换(非邻带转换),则参数DH的选取按前一段的方法去确定。如果调用本函数去进行基于六度带邻带转换,则参数DH的选取先按上述方法去确定,然后看是往前一个带还是后一个带进行邻带转换再确定是加1还是减1。
*/void GeoToGauss(double jd, double wd, short DH, short DH_width, double *y, double *x, double LP){
t=tgBdouble L;
中央经线的经度double l0;
经差double jd_hd,wd_
将jd、wd转换成以弧度为单位double et2;
et2 = (e' ** 2) * (cosB ** 2)double N;
N = C / sqrt(1 + et2)double X;
克拉索夫斯基椭球中子午弧长
m = cosB * PI/180 * l0double tsin,
sinB,cosBdouble PI=3.79;double b_e2=0.7;double b_c=78271;jd_hd = jd / 3600.0 * PI / 180.0 ;
// 将以秒为单位的经度转换成弧度
wd_hd = wd / 3600.0 * PI / 180.0 ;
// 将以秒为单位的纬度转换成弧度// 如果不设中央经线(缺省参数: -1000),则计算中央经线,// 否则,使用传入的中央经线,不再使用带号和带宽参数//L = (DH - 0.5) * DH_
// 计算中央经线的经度if (LP == -1000){L = (DH - 0.5) * DH_
// 计算中央经线的经度}else{L = LP ;}l0 = jd / 3600.0 - L
// 计算经差tsin = sin(wd_hd);
// 计算sinBtcos = cos(wd_hd);
// 计算cosB// 计算克拉索夫斯基椭球中子午弧长XX = 1 / 3600.0 * wd - ( * tsin + 133.9238 * pow(tsin,3)
+ 0.6976 * pow(tsin,5) + 0.0039 * pow(tsin,7) ) *et2 = b_e2 * pow(tcos,2) ;
et2 = (e' ** 2) * (cosB ** 2)
= b_c / sqrt( 1 + et2 ) ;
N = C / sqrt(1 + et2)t
= tan(wd_hd);
= PI/180 * l0 *
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 );
}阅读详情:

参考资料

 

随机推荐