grib2转arl过程中遇到问题，求指点！

良辰

1. #--------------------------------------------------------
   
2. # Author: Yaqiang Wang
   
3. # Date: 2016-2-4
   
4. # Purpose: Convert GRIB data to ARL data
   
5. # Note: Sample
   
6. #-----------------------------------------------------------
   
7. #---- Set data folder
   
8. datadir = 'D:/Temp'
   
9. #---- Set output data file
   
10. outfn = os.path.join(datadir, 'arl/fnl_20150218_06_00.grib2.arl')
    
11. #if os.path.exists(outfn):
    
12. # os.remove(outfn)
    
13. #---- Read a GRIB data file
    
14. infn = os.path.join(datadir, 'grib/fnl_20150218_06_00.grib2')
    
15. print infn
    
16. inf = addfile(infn)
    
17. print 'GRIB data file has been opened...'
    
18. #---- Set output ARL data file
    
19. arlf = addfile(outfn, 'c', dtype='arl')
    
20. #---- Set variable and level list
    
21. gvar2d = ['Pressure_surface','Temperature_surface','u-component_of_wind_height_above_ground',\
    
22. 'v-component_of_wind_height_above_ground']
    
23. gvar3d = ['Geopotential_isobaric','Temperature_isobaric','Pressure_Vertical_velocity_isobaric',\
    
24. 'u-component_of_wind_isobaric','v-component_of_wind_isobaric','Specific_humidity_isobaric']
    
25. avar2d = ['PRSS','T02M','U10M','V10M']
    
26. avar3d = ['HGTS','TEMP','WWND','UWND','VWND','SPHU']
    
27. gv = inf['Geopotential_isobaric']
    
28. nx = gv.dimlen(gv.ndim - 1)
    
29. ny = gv.dimlen(gv.ndim - 2)
    
30. levels = gv.dimvalue(gv.ndim - 3)[::-1]
    
31. nz = len(levels)
    
32. arlf.setlevels(levels)
    
33. arlf.set2dvar(avar2d)
    
34. for l in levels:
    
35. arlf.set3dvar(avar3d)
    
36. 
    
37. #---- Write ARL data file
    
38. arlf.setx(gv.dimvalue(gv.ndim - 1))
    
39. arlf.sety(gv.dimvalue(gv.ndim - 2))
    
40. tNum = inf.timenum()
    
41. fhour = 0
    
42. for t in range(0, tNum):
    
43. print 'Time index: ' + str(t)
    
44. atime = inf.gettime(t)
    
45. print atime.strftime('%Y-%m-%d %H:00')
    
46. dhead = arlf.getdatahead(inf.proj, 'RSMC', 2, fhour)
    
47. #Pre-write index record without checksum - will be over-write latter
    
48. arlf.writeindexrec(atime, dhead)
    
49. #Checksum list
    
50. ksumlist = []
    
51. # Write 2d variables
    
52. ksums = []
    
53. for avname,gvname in zip(avar2d, gvar2d):
    
54. print avname + ' ' + gvname
    
55. if avname == 'U10M' or avname == 'V10M':
    
56. gdata = inf[gvname][t,0,:,:]
    
57. else:
    
58. gdata = inf[gvname][t,:,:]
    
59. if avname == 'PRSS':
    
60. gdata = gdata * 0.01
    
61. ksum = arlf.writedatarec(atime, 0, avname, fhour, 99, gdata)
    
62. ksums.append(ksum)
    
63. ksumlist.append(ksums)
    
64. # Write 3d variables
    
65. for lidx in range(0, nz):
    
66. ksums = []
    
67. llidx = nz - lidx - 1
    
68. print lidx
    
69. print llidx
    
70. for avname,gvname in zip(avar3d, gvar3d):
    
71. print avname + ' ' + gvname
    
72. gdata = inf[gvname][t,llidx,:,:]
    
73. if avname == 'WWND':
    
74. gdata = gdata * 0.01
    
75. elif avname == 'SPHU':
    
76. gdata = gdata * 1000.
    
77. ksum = arlf.writedatarec(atime, lidx + 1, avname, fhour, 99, gdata)
    
78. ksums.append(ksum)
    
79. ksumlist.append(ksums)
    
80. #Re-write index record with checksum
    
81. arlf.writeindexrec(atime, dhead, ksumlist)
    
82. fhour += 1
    
83. arlf.close()
    
84. print 'Finished!'

复制代码

王老师，我在进行grib2转arl过程中，使用您的脚本，出现如下错误，求指点！！！
>>> run script...
D:/Temp\grib/fnl_20150218_06_00.grib2
GRIB data file has been opened...
Traceback (most recent call last):
  File "<iostream>", line 28, in <module>
AttributeError: DimVariable instance has no attribute 'ndim'

MeteoInfo

你先要看看你的数据里的变量名是什么，不一定和示例脚本中的一样

长春市气象台

MeteoInfo 发表于 2016-8-1 10:59
你先要看看你的数据里的变量名是什么，不一定和示例脚本中的一样

王老师，昨天晚上我调试的时候看过文件名了！是一致的！附件的txt是用meteoinfo导出的文件信息！

MeteoInfo

长春市气象台 发表于 2016-8-1 11:54
王老师，昨天晚上我调试的时候看过文件名了！是一致的！附件的txt是用meteoinfo导出的文件信息！

你是“良辰”？数据里有“Geopotential_isobaric”变量吗？

良辰

MeteoInfo 发表于 2016-8-1 12:03
你是“良辰”？数据里有“Geopotential_isobaric”变量吗？

是我，王老师！我在单位值班，有的单位电脑！账号没切换！对不起，王老师，我看的不仔细，Geopotential_height_isobaric是这个名字！给王老师带来了无用的工作量！深表歉意！谢谢王老师！叩谢！！！

良辰

MeteoInfo 发表于 2016-8-1 12:03
你是“良辰”？数据里有“Geopotential_isobaric”变量吗？

1. 
   
2. <P> #--------------------------------------------------------
   
3. # Author: Yaqiang Wang
   
4. # Date: 2016-2-4
   
5. # Purpose: Convert GRIB data to ARL data
   
6. # Note: Sample
   
7. #-----------------------------------------------------------
   
8. #---- Set data folder
   
9. datadir = 'D:/Temp'
   
10. #---- Set output data file
    
11. outfn = os.path.join(datadir, 'arl/fnl_20150218_06_00.grib2.arl')
    
12. #if os.path.exists(outfn):
    
13. #    os.remove(outfn)
    
14. #---- Read a GRIB data file
    
15. infn = os.path.join(datadir, 'grib/fnl_20150218_06_00.grib2')
    
16. print infn
    
17. inf = addfile(infn)
    
18. print 'GRIB data file has been opened...'
    
19. #---- Set output ARL data file
    
20. arlf = addfile(outfn, 'c', dtype='arl')
    
21. #---- Set variable and level list
    
22. gvar2d = ['Pressure_surface','Temperature_surface','u-component_of_wind_height_above_ground',\
    
23.     'v-component_of_wind_height_above_ground']
    
24. gvar3d = ['Geopotential_height_isobaric','Temperature_isobaric',\
    
25.     'u-component_of_wind_isobaric','v-component_of_wind_isobaric',\
    
26.     'Vertical_velocity_pressure_isobaric','Relative_humidity_isobaric']
    
27. avar2d = ['PRSS','T02M','U10M','V10M']
    
28. avar3d = ['HGTS','TEMP','UWND','VWND','WWND','RELH']
    
29. gv = inf['Geopotential_height_isobaric']
    
30. nx = gv.dimlen(gv.ndim - 1)
    
31. ny = gv.dimlen(gv.ndim - 2)
    
32. levels = gv.dimvalue(gv.ndim - 3)[::-1]
    
33. nz = len(levels)
    
34. arlf.setlevels(levels)
    
35. arlf.set2dvar(avar2d)
    
36. for l in levels:
    
37.     arlf.set3dvar(avar3d)</P>
    
38. <P>#---- Write ARL data file
    
39. arlf.setx(gv.dimvalue(gv.ndim - 1))
    
40. arlf.sety(gv.dimvalue(gv.ndim - 2))
    
41. tNum = inf.timenum()
    
42. fhour = 0
    
43. for t in range(0, tNum):
    
44.     print 'Time index: ' + str(t)
    
45.     atime = inf.gettime(t)
    
46.     print atime.strftime('%Y-%m-%d %H:00')
    
47.     dhead = arlf.getdatahead(inf.proj, 'RSMC', 2, fhour)
    
48.     #Pre-write index record without checksum - will be over-write latter
    
49.     arlf.writeindexrec(atime, dhead)
    
50.     #Checksum list
    
51.     ksumlist = []
    
52.     # Write 2d variables
    
53.     ksums = []
    
54.     for avname,gvname in zip(avar2d, gvar2d):
    
55.         print avname + ' ' + gvname
    
56.         if avname == 'U10M' or avname == 'V10M':
    
57.             gdata = inf[gvname][t,0,:,:]
    
58.         else:
    
59.             gdata = inf[gvname][t,:,:]
    
60.         if avname == 'PRSS':
    
61.             gdata = gdata * 0.01
    
62.         ksum = arlf.writedatarec(atime, 0, avname, fhour, 99, gdata)
    
63.         ksums.append(ksum)
    
64.     ksumlist.append(ksums)
    
65.     # Write 3d variables
    
66.     for lidx in range(0, nz):
    
67.         ksums = []
    
68.         llidx = nz - lidx - 1
    
69.         print lidx
    
70.         print llidx
    
71.         for avname,gvname in zip(avar3d, gvar3d):
    
72.             print avname + ' ' + gvname
    
73.             gdata = inf[gvname][t,llidx,:,:]   
    
74.             if avname == 'WWND':
    
75.                gdata = gdata * 0.01      
    
76.             ksum = arlf.writedatarec(atime, lidx + 1, avname, fhour, 99, gdata)
    
77.             ksums.append(ksum)
    
78.         ksumlist.append(ksums)
    
79.     #Re-write index record with checksum
    
80.     arlf.writeindexrec(atime, dhead, ksumlist)
    
81. fhour += 1
    
82. arlf.close()
    
83. print 'Finished!'</P>

复制代码

王老师，还要求助您！运行上面的脚本，出现下面的错误！我感觉是WWND Vertical_velocity_pressure_isobaric有问题，或是llidx = nz - lidx - 1货lidx有问题，或是缩进的问题！python是个小白，给老师添了好多麻烦！求指点~~~
另：如果我要加多个grib2数据怎么处理，就是t变量增加！
>>> run script...
D:/Temp\grib/fnl_20150218_06_00.grib2
GRIB data file has been opened...
Time index: 0
2015-02-18 06:00
PRSS Pressure_surface
T02M Temperature_surface
U10M u-component_of_wind_height_above_ground
V10M v-component_of_wind_height_above_ground
0
25
HGTS Geopotential_height_isobaric
TEMP Temperature_isobaric
UWND u-component_of_wind_isobaric
VWND v-component_of_wind_isobaric
WWND Vertical_velocity_pressure_isobaric
Index out of range!
Traceback (most recent call last):
  File "<iostream>", line 75, in <module>
TypeError: unsupported operand type(s) for *: 'NoneType' and 'float'

MeteoInfo

良辰 发表于 2016-8-1 15:25
王老师，还要求助您！运行上面的脚本，出现下面的错误！我感觉是WWND Vertical_velocity_pressure_isobar ...

你看看你的数据里Vertical_velocity_pressure_isobaric有几层，是不是比其它变量层少一些。

良辰

MeteoInfo 发表于 2016-8-1 15:40
你看看你的数据里Vertical_velocity_pressure_isobaric有几层，是不是比其它变量层少一些。

是的，其他要素26层！Vertical_velocity_pressure_isobaric有21层！

MeteoInfo

良辰 发表于 2016-8-1 15:46
是的，其他要素26层！Vertical_velocity_pressure_isobaric有21层！

那只能在结果文件中写21层了。

良辰

好的，我改一下输入层次！谢谢