站点降水数据绘制等值线图新脚本

MeteoInfo · 发表于 2013-7-12 12:09:15

登录后查看更多精彩内容~

您需要登录才可以下载或查看，没有帐号？立即注册

x

之前为了简化，在MIApp类中做了很多常用的方法，比如打开各种数据，但这样做使得脚本程序不是很清晰，且用户难以对类库的结构有深入认识，所以在写脚本的时候更好的方式是只把MIApp类作为一个数据显示的窗体，数据打开、处理都由类库里相应的类的实现。基于这样方式重新写了用micaps第一类数据绘制降水等值线图的脚本。

# This Python file uses the following encoding: GB2312
#---- 引入类库
import clr
clr.AddReferenceByPartialName("System")
clr.AddReferenceByPartialName("System.Windows.Forms")
clr.AddReferenceByPartialName("System.Drawing")
from System import *
from System.Windows.Forms import *
from System.Drawing import *
clr.AddReference("MeteoInfoC.dll")
from MeteoInfoC import *
from MeteoInfoC.Data.MapData import *
from MeteoInfoC.Data.MeteoData import *
from MeteoInfoC.Layout import *
from MeteoInfoC.Legend import *
from MeteoInfoC.Projections import *
from MeteoInfoC.Shape import *
#---- 设置路径变量
BaseDir = "C:\\Program Files (x86)\\MeteoInfo\\"
MapDir = BaseDir + "Map\\"
LegendDir = BaseDir + "Legend\\"
DataDir = BaseDir + "Sample\\MICAPS\\"
#---- 打开图层
bou2Layer = MapDataManage.OpenLayer(MapDir + "bou2_4p.shp")
lb = bou2Layer.LegendScheme.LegendBreaks[0]
lb.Color = Color.Yellow
lb.OutlineColor = Color.Gray
lb.OutlineSize = 1
lb.DrawFill = False
bou1Layer = MapDataManage.OpenLayer(MapDir + "bou1_4l.shp")
lb = bou1Layer.LegendScheme.LegendBreaks[0]
lb.Color = Color.Blue
chinaLayer = MapDataManage.OpenLayer(MapDir + "china.shp")
chinaLayer.Visible = False
res1Layer = MapDataManage.OpenLayer(MapDir + "res1_4m.shp")
lb = res1Layer.LegendScheme.LegendBreaks[0]
lb.Color = Color.Red
lb.OutlineColor = Color.Black
lb.Size = 5
res1Layer.LabelSet.FieldName = "NAME";
res1Layer.LabelSet.LabelFont = Font("楷体", 10);
res1Layer.LabelSet.YOffset = 15;
res1Layer.AddLabels();
#---- 创建MIApp类的对象
myApp = MIApp()
mapLayout = myApp.MapLayout
layoutMap = mapLayout.ActiveLayoutMap
mapFrame = layoutMap.MapFrame
mapView = mapFrame.MapView
#---- 添加图层
mapFrame.AddLayer(bou2Layer)
mapFrame.AddLayer(bou1Layer)
mapFrame.AddLayer(chinaLayer)
mapFrame.AddLayer(res1Layer)
#---- Lambert投影
projInfo = ProjectionInfo("+proj=lcc+lat_1=25+lat_2=47+lon_0=105")
mapView.ProjectLayers(projInfo)
#---- 按照经纬度范围缩放地图
mapView.ZoomToExtentLonLat(78,130,15,53)
#---- 设置屏蔽图层（只绘制中国境内图形）
mapView.MaskOut.SetMaskLayer = True
mapView.MaskOut.MaskLayer = chinaLayer.LayerName
#---- 设置ActiveLayoutMap（图层显示）
layoutMap.DrawGridLine = False
layoutMap.DrawNeatLine = False
layoutMap.DrawGridLabel = False
layoutMap.DrawGridTickLine = False
layoutMap.Left = 10
layoutMap.Top = 10
layoutMap.Width = 620
layoutMap.Height = 450
#---- 设置南海脚图（添加Map Frame）
aMapFrame = MapFrame()
aMapFrame.IsFireMapViewUpdate = True
aLayoutMap = LayoutMap(aMapFrame)
aLayoutMap.DrawGridLabel = False
aLayoutMap.DrawGridTickLine = False
aLayoutMap.Left = 40
aLayoutMap.Top = 350
aLayoutMap.Width = 85
aLayoutMap.Height = 109
mapLayout.AddElement(aLayoutMap)
aMapFrame.AddLayer(bou1Layer)
aProjInfo = ProjectionInfo("+proj=lcc+lat_1=25+lat_2=47+lon_0=105")
aMapFrame.MapView.ProjectLayers(aProjInfo)
aMapFrame.MapView.ZoomToExtentLonLat(106.5,122.5,1,23)
#---- 设置图形标题
aText = mapLayout.AddText("Temp",220,50,"黑体",12)
#---- 设置图例
aLegend = mapLayout.AddLegend(550,220)
aLegend.LegendStyle = LegendStyles.Normal
aLegend.DrawNeatLine = False
aLegend.Title = "降水量（毫米）"
aLegend.BackColor = Color.White
aLegend.Font = Font("Arial", 8)
#---- 设置起始结束时间
sTime = DateTime.Parse("2010-10-14 14:00")
eTime = DateTime.Parse("2010-10-14 20:00")
#---- 站点数据插值为格点数据的设置
interpSet = InterpolationSetting(60,140,-20,60,160,160,"IDW_Radius",2,1)
#---- 设置图例
ls = LegendScheme(ShapeTypes.Polygon)
ls.ImportFromXMLFile(LegendDir + "rain1.lgs")
#---- 设置MeteoDataInfo
mid = MeteoDataInfo()
#---- 循环
while sTime <= eTime:
inFile = sTime.ToString("yyMMddHH") + ".000"
print inFile
#---- 打开MICAPS数据文件
mid.OpenMICAPSData(DataDir + inFile)
#---- 获取6小时降水量站点数据
stData = mid.GetStationData("Precipitation6h")
#---- 将站点数据差值为格点数据
gData = DrawMeteoData.InterpolateData(stData, interpSet)
#---- 生成6小时降水量图层
rainLayer = DrawMeteoData.CreateShadedLayer(gData, ls, "Precipitation", "Rain")
rainLayer.IsMaskout = True
#---- 添加图层
mapFrame.AddLayer(rainLayer)
print "Display finished"
#---- 调整图层顺序（以避免压盖）
mapFrame.MoveLayerToBottom(rainLayer)
#---- 设置标题名称
title = " 全国降水量实况图\n6小时降水量 (" + sTime.ToString("yyyy-MM-dd HH:00") + ")"
aText.SetLabelText(title)
#---- 设置图例名称
aLegend.Title = "降水量（毫米）"
#---- 绘制图形
mapLayout.PaintGraphics()
#---- 输出图形为文件
outFile = "Prec_" + sTime.ToString("yyyyMMddHH") + ".png"
print outFile
#mapLayout.ExportToPicture("E:\\Temp\\" + outFile)
#---- 删除6小时降水量图层
if sTime < eTime:
mapFrame.RemoveLayer(rainLayer)
#---- 时间加6小时
sTime = sTime.AddHours(6)
#---- 显示程序窗体（只是为了看效果，自动运行时不需要）
#Application.Run(myApp)
myApp.Show()

# This Python file uses the following encoding: GB2312 
#---- 引入类库 
import clr 
clr.AddReferenceByPartialName("System") 
clr.AddReferenceByPartialName("System.Windows.Forms") 
clr.AddReferenceByPartialName("System.Drawing") 
from System import * 
from System.Windows.Forms import * 
from System.Drawing import * 
clr.AddReference("MeteoInfoC.dll") 
from MeteoInfoC import * 
from MeteoInfoC.Data.MapData import * 
from MeteoInfoC.Data.MeteoData import * 
from MeteoInfoC.Layout import * 
from MeteoInfoC.Legend import * 
from MeteoInfoC.Projections import * 
from MeteoInfoC.Shape import * 
 
#---- 设置路径变量 
BaseDir = "C:\\Program Files (x86)\\MeteoInfo\\" 
MapDir = BaseDir + "Map\\" 
LegendDir = BaseDir + "Legend\\" 
DataDir = BaseDir + "Sample\\MICAPS\\" 
 
#---- 打开图层 
bou2Layer = MapDataManage.OpenLayer(MapDir + "bou2_4p.shp") 
lb = bou2Layer.LegendScheme.LegendBreaks[0] 
lb.Color = Color.Yellow 
lb.OutlineColor = Color.Gray 
lb.OutlineSize = 1 
lb.DrawFill = False 
bou1Layer = MapDataManage.OpenLayer(MapDir + "bou1_4l.shp") 
lb = bou1Layer.LegendScheme.LegendBreaks[0] 
lb.Color = Color.Blue 
chinaLayer = MapDataManage.OpenLayer(MapDir + "china.shp") 
chinaLayer.Visible = False 
res1Layer = MapDataManage.OpenLayer(MapDir + "res1_4m.shp") 
lb = res1Layer.LegendScheme.LegendBreaks[0] 
lb.Color = Color.Red 
lb.OutlineColor = Color.Black 
lb.Size = 5 
res1Layer.LabelSet.FieldName = "NAME"; 
res1Layer.LabelSet.LabelFont = Font("楷体", 10); 
res1Layer.LabelSet.YOffset = 15; 
res1Layer.AddLabels(); 
 
#---- 创建MIApp类的对象 
myApp = MIApp() 
mapLayout = myApp.MapLayout 
layoutMap = mapLayout.ActiveLayoutMap 
mapFrame = layoutMap.MapFrame 
mapView = mapFrame.MapView 
 
#---- 添加图层 
mapFrame.AddLayer(bou2Layer) 
mapFrame.AddLayer(bou1Layer) 
mapFrame.AddLayer(chinaLayer) 
mapFrame.AddLayer(res1Layer) 
 
#---- Lambert投影 
projInfo = ProjectionInfo("+proj=lcc+lat_1=25+lat_2=47+lon_0=105") 
mapView.ProjectLayers(projInfo) 
 
#---- 按照经纬度范围缩放地图 
mapView.ZoomToExtentLonLat(78,130,15,53) 
 
#---- 设置屏蔽图层（只绘制中国境内图形） 
mapView.MaskOut.SetMaskLayer = True 
mapView.MaskOut.MaskLayer = chinaLayer.LayerName 
 
#---- 设置ActiveLayoutMap（图层显示） 
layoutMap.DrawGridLine = False 
layoutMap.DrawNeatLine = False 
layoutMap.DrawGridLabel = False 
layoutMap.DrawGridTickLine = False 
layoutMap.Left = 10 
layoutMap.Top = 10 
layoutMap.Width = 620 
layoutMap.Height = 450 
 
#---- 设置南海脚图（添加Map Frame） 
aMapFrame = MapFrame() 
aMapFrame.IsFireMapViewUpdate = True 
aLayoutMap = LayoutMap(aMapFrame) 
aLayoutMap.DrawGridLabel = False 
aLayoutMap.DrawGridTickLine = False 
aLayoutMap.Left = 40 
aLayoutMap.Top = 350 
aLayoutMap.Width = 85 
aLayoutMap.Height = 109 
mapLayout.AddElement(aLayoutMap) 
aMapFrame.AddLayer(bou1Layer) 
aProjInfo = ProjectionInfo("+proj=lcc+lat_1=25+lat_2=47+lon_0=105") 
aMapFrame.MapView.ProjectLayers(aProjInfo) 
aMapFrame.MapView.ZoomToExtentLonLat(106.5,122.5,1,23) 
 
#---- 设置图形标题 
aText = mapLayout.AddText("Temp",220,50,"黑体",12) 
 
#---- 设置图例 
aLegend = mapLayout.AddLegend(550,220) 
aLegend.LegendStyle = LegendStyles.Normal                
aLegend.DrawNeatLine = False 
aLegend.Title = "降水量（毫米）" 
aLegend.BackColor = Color.White 
aLegend.Font = Font("Arial", 8) 
 
#---- 设置起始结束时间 
sTime = DateTime.Parse("2010-10-14 14:00") 
eTime = DateTime.Parse("2010-10-14 20:00") 
 
#---- 站点数据插值为格点数据的设置 
interpSet = InterpolationSetting(60,140,-20,60,160,160,"IDW_Radius",2,1) 
 
#---- 设置图例 
ls = LegendScheme(ShapeTypes.Polygon) 
ls.ImportFromXMLFile(LegendDir + "rain1.lgs") 
 
#---- 设置MeteoDataInfo 
mid = MeteoDataInfo() 
 
#---- 循环 
while sTime <= eTime:    
  inFile = sTime.ToString("yyMMddHH") + ".000" 
  print inFile 
  #---- 打开MICAPS数据文件 
  mid.OpenMICAPSData(DataDir + inFile)   
  #---- 获取6小时降水量站点数据 
  stData = mid.GetStationData("Precipitation6h") 
  #---- 将站点数据差值为格点数据 
  gData = DrawMeteoData.InterpolateData(stData, interpSet) 
  #---- 生成6小时降水量图层 
  rainLayer = DrawMeteoData.CreateShadedLayer(gData, ls, "Precipitation", "Rain") 
  rainLayer.IsMaskout = True 
  #---- 添加图层 
  mapFrame.AddLayer(rainLayer) 
  print "Display finished"   
  #---- 调整图层顺序（以避免压盖） 
  mapFrame.MoveLayerToBottom(rainLayer)   
  #---- 设置标题名称 
  title = "      全国降水量实况图\n6小时降水量 (" + sTime.ToString("yyyy-MM-dd HH:00") + ")" 
  aText.SetLabelText(title)      
  #---- 设置图例名称 
  aLegend.Title = "降水量（毫米）" 
  #---- 绘制图形 
  mapLayout.PaintGraphics() 
  #---- 输出图形为文件 
  outFile = "Prec_" + sTime.ToString("yyyyMMddHH") + ".png" 
  print outFile 
  #mapLayout.ExportToPicture("E:\\Temp\\" + outFile)   
  #---- 删除6小时降水量图层 
  if sTime < eTime: 
          mapFrame.RemoveLayer(rainLayer) 
  #---- 时间加6小时 
  sTime = sTime.AddHours(6) 
 
#---- 显示程序窗体（只是为了看效果，自动运行时不需要） 
#Application.Run(myApp) 
myApp.Show()

rceclx · 发表于 2013-7-12 17:01:54

本帖最后由 rceclx 于 2013-7-12 17:09 编辑

登录/注册后可看大图

！王老师智慧啊！

Jasonluo · 发表于 2013-7-12 18:46:37

王老师，你这个降雨量空间分布图是基于什么插值方法呢？

MeteoInfo · 发表于 2013-7-12 18:48:34

Jasonluo 发表于 2013-7-12 18:46

登录/注册后可看大图

王老师，你这个降雨量空间分布图是基于什么插值方法呢？

代码里有，IDW

honview · 发表于 2013-7-25 20:27:53

好东西，厉害呀

wbriver · 发表于 2013-7-28 23:29:42

MeteoInfo 发表于 2013-7-12 18:48

登录/注册后可看大图

代码里有，IDW

王老师您好，之前用的是2013.3.2版本的dll，今天更新了下，发现报错：（1）“MeteoInfoC.Data.StationData”不包含“UNDEF”的定义，这个改成stationData.UNDEF→stationData.MissingValue？
（2）“MeteoInfoC.Data.MeteoData.MICAPS4DataInfo”不包含“GetGridData”的定义，这个不知道如何改呢？谢谢！

            MICAPS4DataInfo aDataInfo = new MICAPS4DataInfo();
            aDataInfo.ReadDataInfo(openFileDialog.FileName);
            GridData rainTemp = aDataInfo.GetGridData();//这句话如何改呢？
            GridData rain = rainTemp;

MeteoInfo · 发表于 2013-7-28 23:34:46

wbriver 发表于 2013-7-28 23:29

登录/注册后可看大图

王老师您好，之前用的是2013.3.2版本的dll，今天更新了下，发现报错：（1）“MeteoInfoC.Data.StationDat ...

1、你的理解是对的。
2、
MeteoDataInfo aDataInfo = new MeteoDataInfo();
aDataInfo.OpenMICAPSData(openFileDialog.FileName);
GridData rainTemp = aDataInfo.GetGridData("var")

wbriver · 发表于 2013-7-29 11:19:00

MeteoInfo 发表于 2013-7-28 23:34

登录/注册后可看大图

1、你的理解是对的。
2、
MeteoDataInfo aDataInfo = new MeteoDataInfo();

问题解决了，谢谢！

微微一笑 · 发表于 2013-7-31 13:35:14

怎样画南海角图呀？把那段程序复制到script中总是出错…我也新建了map frame呀？

MeteoInfo · 发表于 2013-7-31 13:49:01

微微一笑发表于 2013-7-31 13:35

登录/注册后可看大图

怎样画南海角图呀？把那段程序复制到script中总是出错…我也新建了map frame呀？

出错信息是什么？我不能凭空瞎猜呀！

		自动登录	找回密码
密码			立即注册