爱气象,爱气象家园! 

气象家园

 找回密码
 立即注册

QQ登录

只需一步,快速开始

新浪微博登陆

只需一步, 快速开始

搜索
查看: 28441|回复: 54

MeteoInfo脚本示例:MICAPS第4类数据转为NetCDF

[复制链接]

新浪微博达人勋

发表于 2012-12-11 16:17:52 | 显示全部楼层 |阅读模式

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

您需要 登录 才可以下载或查看,没有帐号?立即注册 新浪微博登陆

x
本帖最后由 MeteoInfo 于 2015-2-2 14:33 编辑

写NetCDF格式数据的基本步骤是:定义各种维(经度、纬度、高度、时间),定义变量(包括维变量和数据变量)以及变量的属性,定义全局属性,将各维变量的数据写入文件,将数据变量按照一定规则写入文件。

这里演示将多个MICAPS第4类数据文件写入一个NetCDF文件中,MICAPS第4类数据是简单的二维格点数据,通过读取不同高度、不同时间的数据文件形成4维数据,写入一个NetCDF文件。所需的MICAPS示例数据见此贴:MICAPS第四类格式数据转换的请教
http://bbs.06climate.com/forum.php?mod=viewthread&tid=11611&fromuid=106

需要更新MeteoInfo的最新文件(见置顶帖子)。

写出来的NC数据经过测试和原始MICAPS数据是吻合的。
Image00300.png
Image00299.png

详细脚本代码如下:
  1. #--------------------------------------------------------
  2. # Author: Yaqiang Wang
  3. # Date: 2014-3-18
  4. # Purpose: Convert MICAPS 4 grid data to NetCDF data format
  5. # Note: Sample
  6. #-----------------------------------------------------------
  7. import clr
  8. clr.AddReferenceByPartialName("System.Windows.Forms")
  9. clr.AddReferenceByPartialName("System.Drawing")
  10. from System import *
  11. from System.Windows.Forms import *
  12. from System.Drawing import *
  13. clr.AddReference("MeteoInfoC.dll")
  14. from MeteoInfoC import *
  15. from MeteoInfoC.Data import *
  16. from MeteoInfoC.Data.MeteoData import *

  17. import os.path
  18. import datetime
  19. import time

  20. #---- Set data folder
  21. dataDir = "D:\\Temp\\test\\R\\"

  22. #---- Open a MICAPS 4 Data file to get the information
  23. inDataInfo = MeteoDataInfo()
  24. infile = dataDir + "1000\\14031708.000.000"
  25. inDataInfo.OpenMICAPSData(infile)
  26. undef = inDataInfo.MissingValue
  27. lons = inDataInfo.GetX()
  28. lats = inDataInfo.GetY()
  29. lonNum = lons.Length
  30. latNum = lats.Length
  31. ls = [1000, 925, 850, 700, 500, 300, 100]
  32. levelNum = len(ls)
  33. levels = Array.CreateInstance(Double, levelNum)
  34. i = 0
  35. for l in ls:
  36.         levels = l;
  37.         i += 1

  38. #---- Create netCDF data info object
  39. outfilepath = dataDir + "test.nc"
  40. outDataInfo = NetCDFDataInfo()
  41. outDataInfo.FileName = outfilepath
  42. outDataInfo.MissingValue = undef
  43. outDataInfo.unlimdimid = 3

  44. #---- Add dimensions: lon, lat, level, time
  45. lonDim = outDataInfo.AddDimension("lon", lonNum)
  46. latDim = outDataInfo.AddDimension("lat", latNum)
  47. levelDim = outDataInfo.AddDimension("level", levelNum)
  48. timeDim = outDataInfo.AddDimension("time", -1)

  49. #---- Add variables
  50. outDataInfo.AddVariable("lon", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([lonDim]))
  51. outDataInfo.AddVariable("lat", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([latDim]))
  52. outDataInfo.AddVariable("level", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([levelDim]))
  53. outDataInfo.AddVariable("time", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim]))
  54. outDataInfo.AddVariable("temadv", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim, levelDim, latDim, lonDim]))
  55. #---- Add variable attributes
  56. outDataInfo.AddVariableAttribute("lon", "units", "degrees_east")
  57. outDataInfo.AddVariableAttribute("lon", "long_name", "longitude")
  58. outDataInfo.AddVariableAttribute("lat", "units", "degrees_north")
  59. outDataInfo.AddVariableAttribute("lat", "long_name", "latitude")
  60. outDataInfo.AddVariableAttribute("level", "units", "level");
  61. outDataInfo.AddVariableAttribute("level", "long_name", "level");
  62. outDataInfo.AddVariableAttribute("level", "axis", "z");
  63. outDataInfo.AddVariableAttribute("time", "units", "hours since 1800-1-1 00:00:00")
  64. outDataInfo.AddVariableAttribute("time", "long_name", "time")
  65. outDataInfo.AddVariableAttribute("temadv", "units", "No")
  66. outDataInfo.AddVariableAttribute("temadv", "long_name", "TemperatureAdvection")
  67. outDataInfo.AddVariableAttribute("temadv", "missing_value", undef)

  68. #---- Add global attributes
  69. outDataInfo.AddGlobalAttribute("title", "Temperature advection field data")
  70. outDataInfo.AddGlobalAttribute("description", "script sample data")

  71. #---- Create netCDF file
  72. outDataInfo.CreateNCFile(outfilepath)

  73. #---- Write lon, lat, level, time data
  74. outDataInfo.WriteVar("lon", lons)
  75. outDataInfo.WriteVar("lat", lats)
  76. outDataInfo.WriteVar("level", levels)

  77. tnum = 6
  78. timeArray = Array.CreateInstance(Double, tnum)
  79. times = []
  80. st = datetime.datetime(1800, 1, 1, 0, 0, 0)
  81. t = datetime.datetime(2014, 3, 17, 8, 0, 0)        
  82. for i in range(0, tnum):        
  83.         times.append(t)
  84.         delta = t - st
  85.         hours = (delta.days * 86400 + delta.seconds) / 3600.0
  86.         t = t + datetime.timedelta(hours=3)
  87.         timeArray = hours
  88.         print timeArray
  89.         
  90. startArray = Array.CreateInstance(int, 1)
  91. countArray = Array.CreateInstance(int, 1)
  92. startArray[0] = 0
  93. countArray[0] = len(times)
  94. outDataInfo.WriteVara("time", startArray, countArray, timeArray)

  95. #---- Write data
  96. tidx = 0
  97. st = times[0]
  98. for t in times:
  99.         lidx = 0        
  100.         for l in ls:
  101.                 tt = tidx * 3
  102.                 tstr = '.%03d'%tt
  103.                 infile = dataDir + str(l) + '\\' + st.strftime("%y%m%d%H") + tstr + tstr
  104.                 #print infile
  105.                 if os.path.isfile(infile):
  106.                         print infile
  107.                         inDataInfo.OpenMICAPSData(infile)
  108.                         gdata = inDataInfo.GetGridData("var")
  109.                         dataArray = gdata.ToOneDimData()        
  110.                         startArray = Array.CreateInstance(int, 4)
  111.                         countArray = Array.CreateInstance(int, 4)
  112.                         startArray[0] = tidx
  113.                         startArray[1] = lidx
  114.                         startArray[2] = 0
  115.                         startArray[3] = 0
  116.                         countArray[0] = 1
  117.                         countArray[1] = 1
  118.                         countArray[2] = latDim.DimLength
  119.                         countArray[3] = lonDim.DimLength
  120.                         outDataInfo.WriteVara("temadv", startArray, countArray, dataArray)
  121.                 lidx += 1
  122.         tidx += 1

  123. #---- Close netCDF file
  124. outDataInfo.CloseNCFile()

  125. print 'Finished...!'




评分

参与人数 1金钱 +20 贡献 +2 收起 理由
njzqxt + 20 + 2 谢谢王老师对基层气象人的悉心指点!

查看全部评分

密码修改失败请联系微信:mofangbao

新浪微博达人勋

 成长值: 19710
发表于 2012-12-11 18:10:41 | 显示全部楼层
王老师好神速~
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-12-18 11:31:26 | 显示全部楼层
本帖最后由 njzqxt 于 2012-12-18 11:36 编辑

    王老师,谢谢了。
    我这试了多次也无法正确转出数据啊。
    数据放在K:\physic\ttadv\下的分别为1000、 925、 850、 700、 500、 400、 300、 200、 150、 100的文件夹中。
    这个脚本中,要改动哪些呢?我改了以下几处。
dataDir = "K:\\physic\\ttadv\\"                       #要转K:\physic\ttadv\  下的10个文件夹中的数据
infile = dataDir + "1000\\12121408.000"       #每个文件夹下第一个数据的名称都是12121408.000
infile = dataDir + str(l) + "\\" + t.strftime("%y%m%d%H") + ".000"      这个怎么改呢?
运行时提示错误  Error:NetCDF:Not a valid ID
确定后在脚本框的下面出现图2中的数据。转出来的数据才1kb。

1.JPG
2.JPG

ttadv.rar

432.15 KB, 下载次数: 13, 下载积分: 金钱 -5

密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-12-18 12:01:28 | 显示全部楼层

你把所有代码贴出来看看
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-12-18 13:10:12 | 显示全部楼层
#--------------------------------------------------------
# Author: Yaqiang Wang
# Date: 2012-12-11
# Purpose: Convert MICAPS 4 grid data to NetCDF data format
# Note: Sample
#-----------------------------------------------------------
import clr
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 import *
from MeteoInfoC.Data.MeteoData import *

import os.path
import datetime
import time

#---- Set data folder
dataDir = "K:\\physic\\ttadv\\"

#---- Open a MICAPS 4 Data file to get the information
inDataInfo = MeteoDataInfo()
infile = dataDir + "1000\\12121408.000"
inDataInfo.OpenMICAPSData(infile)
undef = inDataInfo.UNDEF
lons = inDataInfo.GetX()
lats = inDataInfo.GetY()
lonNum = lons.Length
latNum = lats.Length
ls = [1000, 925, 850, 700, 500, 400, 300, 200, 150, 100]
levelNum = len(ls)
levels = Array.CreateInstance(Double, levelNum)
i = 0
for l in ls:
        levels[i] = l;
        i += 1

#---- Create netCDF data info object
outfilepath = "D:\\" + "test.nc"
outDataInfo = NetCDFDataInfo()
outDataInfo.fileName = outfilepath
outDataInfo.UNDEF = undef
outDataInfo.unlimdimid = 3

#---- Add dimensions: lon, lat, level, time
lonDim = outDataInfo.AddDimension("lon", lonNum)
latDim = outDataInfo.AddDimension("lat", latNum)
levelDim = outDataInfo.AddDimension("level", levelNum)
timeDim = outDataInfo.AddDimension("time", -1)

#---- Add variables
outDataInfo.AddVariable("lon", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([lonDim]))
outDataInfo.AddVariable("lat", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([latDim]))
outDataInfo.AddVariable("level", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([levelDim]))
outDataInfo.AddVariable("time", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim]))
outDataInfo.AddVariable("temadv", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim, levelDim, latDim, lonDim]))
#---- Add variable attributes
outDataInfo.AddVariableAttribute("lon", "units", "degrees_east")
outDataInfo.AddVariableAttribute("lon", "long_name", "longitude")
outDataInfo.AddVariableAttribute("lat", "units", "degrees_north")
outDataInfo.AddVariableAttribute("lat", "long_name", "latitude")
outDataInfo.AddVariableAttribute("level", "units", "level");
outDataInfo.AddVariableAttribute("level", "long_name", "level");
outDataInfo.AddVariableAttribute("level", "axis", "z");
outDataInfo.AddVariableAttribute("time", "units", "hours since 1800-1-1 00:00:00")
outDataInfo.AddVariableAttribute("time", "long_name", "time")
outDataInfo.AddVariableAttribute("temadv", "units", "No")
outDataInfo.AddVariableAttribute("temadv", "long_name", "TemperatureAdvection")
outDataInfo.AddVariableAttribute("temadv", "missing_value", undef)

#---- Add global attributes
outDataInfo.AddGlobalAttribute("title", "Temperature advection field data")
outDataInfo.AddGlobalAttribute("description", "script sample data")

#---- Create netCDF file
outDataInfo.CreateNCFile(outfilepath)

#---- Write lon, lat, level, time data
outDataInfo.WriteVar("lon", lons)
outDataInfo.WriteVar("lat", lats)
outDataInfo.WriteVar("level", levels)

timeArray = Array.CreateInstance(Double, 8)
times = []
st = datetime.datetime(1800, 1, 1, 0, 0, 0)
t = datetime.datetime(2012, 12, 7, 8, 0, 0)        
for i in range(0, 8):        
        times.append(t)
        delta = t - st
        hours = (delta.days * 86400 + delta.seconds) / 3600.0
        t = t + datetime.timedelta(hours=12)
        timeArray[i] = hours
        print timeArray[i]
        
startArray = Array.CreateInstance(int, 1)
countArray = Array.CreateInstance(int, 1)
startArray[0] = 0
countArray[0] = len(times)
outDataInfo.WriteVara("time", startArray, countArray, timeArray)

#---- Write data
tidx = 0
for t in times:
        lidx = 0
        for l in ls:
                infile = dataDir + str(l) + "\\" + t.strftime("%y%m%d%H") + ".000"
                #print infile
                if os.path.isfile(infile):
                        print infile
                        inDataInfo.OpenMICAPSData(infile)
                        gdata = inDataInfo.GetGridData("var")
                        dataArray = gdata.ToOneDimData()        
                        startArray = Array.CreateInstance(int, 4)
                        countArray = Array.CreateInstance(int, 4)
                        startArray[0] = tidx
                        startArray[1] = lidx
                        startArray[2] = 0
                        startArray[3] = 0
                        countArray[0] = 1
                        countArray[1] = 1
                        countArray[2] = latDim.DimLength
                        countArray[3] = lonDim.DimLength
                        outDataInfo.WriteVara("temadv", startArray, countArray, dataArray)
                lidx += 1
        tidx += 1

#---- Close netCDF file
outDataInfo.CloseNCFile()

print 'Finished...!'
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-12-18 15:38:39 | 显示全部楼层
njzqxt 发表于 2012-12-18 13:10
#--------------------------------------------------------
# Author: Yaqiang Wang
# Date: 2012-12-1 ...

时间设置错了,改为:
t = datetime.datetime(2012, 12, 12, 8, 0, 0)
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-12-18 15:58:01 | 显示全部楼层
本帖最后由 njzqxt 于 2012-12-18 15:59 编辑

生成 的test.nc还是1kb,还要改哪儿呢?
1.JPG
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-12-18 16:16:24 | 显示全部楼层
njzqxt 发表于 2012-12-18 15:58
生成 的test.nc还是1kb,还要改哪儿呢?

你把MeteoInfo关掉,重新打开再试试
密码修改失败请联系微信:mofangbao

新浪微博达人勋

发表于 2012-12-18 19:48:30 | 显示全部楼层
重新打开,运行,也不行,错误提示和7楼一样。
#--------------------------------------------------------
# Author: Yaqiang Wang
# Date: 2012-12-11
# Purpose: Convert MICAPS 4 grid data to NetCDF data format
# Note: Sample
#-----------------------------------------------------------
import clr
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 import *
from MeteoInfoC.Data.MeteoData import *

import os.path
import datetime
import time

#---- Set data folder
dataDir = "D:\\ttadv\\"

#---- Open a MICAPS 4 Data file to get the information
inDataInfo = MeteoDataInfo()
infile = dataDir + "1000\\12121208.000"
inDataInfo.OpenMICAPSData(infile)
undef = inDataInfo.UNDEF
lons = inDataInfo.GetX()
lats = inDataInfo.GetY()
lonNum = lons.Length
latNum = lats.Length
ls = [1000, 925, 850, 700, 500, 400, 300, 200, 150, 100]
levelNum = len(ls)
levels = Array.CreateInstance(Double, levelNum)
i = 0
for l in ls:
        levels[i] = l;
        i += 1

#---- Create netCDF data info object
outfilepath = dataDir + "test.nc"
outDataInfo = NetCDFDataInfo()
outDataInfo.fileName = outfilepath
outDataInfo.UNDEF = undef
outDataInfo.unlimdimid = 3

#---- Add dimensions: lon, lat, level, time
lonDim = outDataInfo.AddDimension("lon", lonNum)
latDim = outDataInfo.AddDimension("lat", latNum)
levelDim = outDataInfo.AddDimension("level", levelNum)
timeDim = outDataInfo.AddDimension("time", -1)

#---- Add variables
outDataInfo.AddVariable("lon", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([lonDim]))
outDataInfo.AddVariable("lat", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([latDim]))
outDataInfo.AddVariable("level", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([levelDim]))
outDataInfo.AddVariable("time", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim]))
outDataInfo.AddVariable("temadv", NetCDF4.NcType.NC_DOUBLE, Array[Dimension]([timeDim, levelDim, latDim, lonDim]))
#---- Add variable attributes
outDataInfo.AddVariableAttribute("lon", "units", "degrees_east")
outDataInfo.AddVariableAttribute("lon", "long_name", "longitude")
outDataInfo.AddVariableAttribute("lat", "units", "degrees_north")
outDataInfo.AddVariableAttribute("lat", "long_name", "latitude")
outDataInfo.AddVariableAttribute("level", "units", "level");
outDataInfo.AddVariableAttribute("level", "long_name", "level");
outDataInfo.AddVariableAttribute("level", "axis", "z");
outDataInfo.AddVariableAttribute("time", "units", "hours since 1800-1-1 00:00:00")
outDataInfo.AddVariableAttribute("time", "long_name", "time")
outDataInfo.AddVariableAttribute("temadv", "units", "No")
outDataInfo.AddVariableAttribute("temadv", "long_name", "TemperatureAdvection")
outDataInfo.AddVariableAttribute("temadv", "missing_value", undef)

#---- Add global attributes
outDataInfo.AddGlobalAttribute("title", "Temperature advection field data")
outDataInfo.AddGlobalAttribute("description", "script sample data")

#---- Create netCDF file
outDataInfo.CreateNCFile(outfilepath)

#---- Write lon, lat, level, time data
outDataInfo.WriteVar("lon", lons)
outDataInfo.WriteVar("lat", lats)
outDataInfo.WriteVar("level", levels)

timeArray = Array.CreateInstance(Double, 8)
times = []
st = datetime.datetime(1800, 1, 1, 0, 0, 0)
t = datetime.datetime(2012, 12, 12, 8, 0, 0)        
for i in range(0, 8):        
        times.append(t)
        delta = t - st
        hours = (delta.days * 86400 + delta.seconds) / 3600.0
        t = t + datetime.timedelta(hours=12)
        timeArray[i] = hours
        print timeArray[i]
        
startArray = Array.CreateInstance(int, 1)
countArray = Array.CreateInstance(int, 1)
startArray[0] = 0
countArray[0] = len(times)
outDataInfo.WriteVara("time", startArray, countArray, timeArray)

#---- Write data
tidx = 0
for t in times:
        lidx = 0
        for l in ls:
                infile = dataDir + str(l) + "\\" + t.strftime("%y%m%d%H") + ".000"
                #print infile
                if os.path.isfile(infile):
                        print infile
                        inDataInfo.OpenMICAPSData(infile)
                        gdata = inDataInfo.GetGridData("var")
                        dataArray = gdata.ToOneDimData()        
                        startArray = Array.CreateInstance(int, 4)
                        countArray = Array.CreateInstance(int, 4)
                        startArray[0] = tidx
                        startArray[1] = lidx
                        startArray[2] = 0
                        startArray[3] = 0
                        countArray[0] = 1
                        countArray[1] = 1
                        countArray[2] = latDim.DimLength
                        countArray[3] = lonDim.DimLength
                        outDataInfo.WriteVara("temadv", startArray, countArray, dataArray)
                lidx += 1
        tidx += 1

#---- Close netCDF file
outDataInfo.CloseNCFile()

print 'Finished...!'
密码修改失败请联系微信:mofangbao

新浪微博达人勋

 楼主| 发表于 2012-12-18 20:52:31 | 显示全部楼层
njzqxt 发表于 2012-12-18 19:48
重新打开,运行,也不行,错误提示和7楼一样。
#------------------------------------------------------ ...

我这里是可以的,脚本运行的截图如下,正常运行应该没有错误提示,而且会打印出所有文件名。

Image00952.png

转成的NetCDF文件: test.nc (962.17 KB, 下载次数: 8)
密码修改失败请联系微信:mofangbao
您需要登录后才可以回帖 登录 | 立即注册 新浪微博登陆

本版积分规则

Copyright ©2011-2014 bbs.06climate.com All Rights Reserved.  Powered by Discuz! (京ICP-10201084)

本站信息均由会员发表,不代表气象家园立场,禁止在本站发表与国家法律相抵触言论

快速回复 返回顶部 返回列表