Arcpy脚本——基于Arcgis实现shp/gdb等到火星/百度坐标系的转换

  在生产中,会遇到将我们既有的坐标系转换到高德或者百度坐标系下,从而可以使用百度高德的基础地图。亦或者将获得的百度高德数据转换到普通的84坐标系下。本文基于Arcpy脚本实现了对shp或者gdb数据的坐标系之间的相互转换。

1、从常用坐标系到百度高德坐标系的转换

思路如下:
1)读取原始数据的坐标系;
2)利用Arcpy将原始坐标转换成WGS84经纬度;
3)对于线和面,利用Arcpy对数据进行修复;
4)逐个取点,将WGS84经纬度转换成高德或者百度坐标系。

# -*- coding:UTF-8 -*-
import arcpy, os
import ConvertLag_Lat as CVLL
import sys
reload(sys)
sys.setdefaultencoding('utf8')

#输入工作空间
in_features = arcpy.GetParameterAsText(0)
#输出的空间坐标系
outCSStr=arcpy.GetParameterAsText(1)
#输出工作空间
output_folder = arcpy.GetParameterAsText(2)

outCS=arcpy.SpatialReference(4326)
arcpy.BatchProject_management(in_features, output_folder, outCS)
arcpy.env.workspace = output_folder
transType=0
if(outCSStr=='wgs84'):
  arcpy.AddMessage("transform succeed")
elif(outCSStr=='gaode'):
  arcpy.AddMessage("target coords gaode")
  transType='3'
elif(outCSStr=='baidu'):
  transType='6'
  arcpy.AddMessage("target coords baidu")
arcpy.AddMessage("start transforming......")
featureclasses = arcpy.ListFeatureClasses()
for fc in featureclasses:
  cur= arcpy.UpdateCursor(fc)
  SR= arcpy.Describe(fc).spatialReference
  des = arcpy.Describe(fc)
  arcpy.RepairGeometry_management(fc)
  fc_count=int(arcpy.GetCount_management(fc).getOutput(0))
  arcpy.AddMessage(fc_count)
  arcpy.SetProgressor("step", "transform "+fc,0, fc_count, 1)
  if des.shapeType.upper()=='POINT':
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetPoint(geom,transType))
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur 
  elif des.shapeType.upper()=='POLYLINE':   
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetLine(geom,SR,transType))  
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur 
  elif des.shapeType.upper()=='POLYGON':   
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetPolygon(geom,SR,transType))
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur
  arcpy.AddMessage(fc+" transform finished")

2、百度高德与WGS84互转

import arcpy, os
import ConvertLag_Lat as CVLL
import sys
reload(sys)
sys.setdefaultencoding('utf8')

#input path of folder where all shapefiles reside in
in_fold = arcpy.GetParameterAsText(0)
output_folder = arcpy.GetParameterAsText(1)
outCSStr= arcpy.GetParameterAsText(2)
outCS=arcpy.SpatialReference(4326)

arcpy.env.workspace = in_fold
#copy data
featureclasses = arcpy.ListFeatureClasses()
for fc in featureclasses:
  arcpy.Copy_management(fc,output_folder+"\\"+fc)
arcpy.env.workspace=output_folder

transType='0'
if(outCSStr=='mars2baidu'):
  transType='1'

elif(outCSStr=='baidu2mars'):
  transType='2'

elif(outCSStr=='mars2baidu'):
  transType='3'

elif(outCSStr=='mars2wgs84'):
  transType='4'

elif(outCSStr=='baidu2wgs'):
  transType='5'

elif(outCSStr=='wgs2baidu'):
  transType='6'


arcpy.AddMessage("transfor type is "+str(transType))
featureclasses = arcpy.ListFeatureClasses()
for fc in featureclasses:
  cur= arcpy.UpdateCursor(fc)
  SR= arcpy.SpatialReference(4326)
  des = arcpy.Describe(fc)
  arcpy.RepairGeometry_management(fc)
  fc_count=int(arcpy.GetCount_management(fc).getOutput(0))
  arcpy.SetProgressor("step", "transform "+fc,0, fc_count, 1)
  if des.shapeType.upper()=='POINT':
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetPoint(geom,transType))
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur 
  elif des.shapeType.upper()=='POLYLINE':   
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetLine(geom,SR,transType))  
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur 
  elif des.shapeType.upper()=='POLYGON':   
    for r in cur:
      geom = r.getValue("SHAPE")
      r.setValue("SHAPE",CVLL.offsetPolygon(geom,SR,transType))
      cur.updateRow(r)
      arcpy.SetProgressorPosition()
    del r,cur
  arcpy.AddMessage(fc+" transform finished")

不足

这里有几个方面的不足
1)我使用的是git上给出的坐标转换方式,经测试转换后有米级的误差;
2)这里用的是取点在重组的方式,因此对于弧,多部件要素的支持不行,只能支持对简单的点线面进行转换。
这里给出源码,希望有需要的人可以继续优化完善。

上一篇:[R] 添加误差棒的分组折线图:geom_path: Each group consists of only one observation. Do you need to adjust the...


下一篇:R语言使用混合模型进行聚类