OMI Ozone 数据转为GrADS二进制数据

MeteoInfo

卫星臭氧观测数据OMI Ozone Level 3 Data可以在此网站上下载：http://ozoneaq.gsfc.nasa.gov/OMIOzone.md。数据的分辨率是1度*1度或0.25度*0.25度，每天一个数据文件，文件是文本格式，这里以1度*1度的数据为例。为了更方便地处理该数据，可以将多个时次的文本数据文件合并为一个GrADS二进制文件，脚本程序如下（需要MeteoInfo最新文件，见置顶帖子）：

1. #--------------------------------------------------------        
   
2. # Author: Yaqiang Wang                                          
   
3. # Date: 2012-5-30                                                
   
4. # Purpose: Convert OMI ozon data to GrADS binary data  
   
5. # Note: Sample                                                   
   
6. #-----------------------------------------------------------     
   
7. import clr
   
8. from System import *
   
9. clr.AddReference("MeteoInfoC.dll")
   
10. from MeteoInfoC import *
    
11. from MeteoInfoC.Data import *
    
12. from MeteoInfoC.Data.MeteoData import *
    
13. 
    
14. import os.path
    
15. 
    
16. #---- Set data folders
    
17. baseDir = 'F:\\Temp\\satellite\\'
    
18. dataDir = baseDir + 'OMI_ozone_201201\\'
    
19. 
    
20. #---- Set output data file
    
21. outCtlFile = baseDir + 'omi_ozone_201201.ctl'
    
22. outDataFile = baseDir + 'omi_ozone_201201.dat'
    
23. 
    
24. #---- Set start/end time
    
25. sTime = DateTime(2012, 1, 1)
    
26. eTime = DateTime(2012, 1, 31)
    
27. tNum = (eTime - sTime).Days + 1
    
28. print 'Start time: ' + sTime.ToString("yyyy-MM-dd")
    
29. print 'End time: ' + eTime.ToString("yyyy-MM-dd")
    
30. print 'Time Number: ' + str(tNum)
    
31. 
    
32. #---- Define and set output GrADS data info
    
33. gDataInfo = GrADSDataInfo()
    
34. gDataInfo.FileName = outCtlFile
    
35. gDataInfo.DSET = outDataFile
    
36. gDataInfo.TITLE = 'OMI Ozone data'
    
37. gDataInfo.DTYPE = 'GRIDDED'
    
38. gDataInfo.XDEF.Type = 'LINEAR'
    
39. gDataInfo.XDEF.XNum = 360
    
40. gDataInfo.XDEF.XMin = -179.5
    
41. gDataInfo.XDEF.XDelt = 1.0
    
42. gDataInfo.YDEF.Type = 'LINEAR'
    
43. gDataInfo.YDEF.YNum = 180
    
44. gDataInfo.YDEF.YMin = -89.5
    
45. gDataInfo.YDEF.YDelt = 1.0
    
46. gDataInfo.ZDEF.Type = 'LINEAR'
    
47. gDataInfo.ZDEF.ZNum = 1
    
48. gDataInfo.ZDEF.SLevel = 1
    
49. gDataInfo.ZDEF.ZDelt = 1
    
50. gDataInfo.TDEF.Type = 'LINEAR'
    
51. gDataInfo.TDEF.TNum = tNum
    
52. gDataInfo.TDEF.STime = sTime
    
53. gDataInfo.TDEF.TDelt = '1dy'
    
54. aVar = VAR()
    
55. aVar.VName = 'ozone'
    
56. aVar.LevelNum = 0
    
57. aVar.Units = '99'
    
58. aVar.Description = 'OMI Ozone'
    
59. gDataInfo.VARDEF.AddVar(aVar)
    
60. 
    
61. #---- Write GrADS data file
    
62. gDataInfo.CreateDataFile(outDataFile)
    
63. aTime = sTime
    
64. for t in range(0, tNum):
    
65.         print aTime.ToString("yyyy-MM-dd")
    
66.         inFilePath = dataDir + 'L3_ozone_omi_' + aTime.ToString("yyyyMMdd") + '.txt'
    
67.         if os.path.isfile(inFilePath):
    
68.                 print inFilePath
    
69. #---- Read grid data form OMI Ozone data file
    
70.                 gData = Array.CreateInstance(Double, 180, 360)
    
71.                 dataFile = open(inFilePath)
    
72.                 lineNum = 0
    
73.                 dataStr = ''
    
74.                 yIdx = 0
    
75.                 xIdx = 0
    
76.                 for aLine in dataFile.readlines():
    
77.                         lineNum += 1
    
78.                         if lineNum <= 3:
    
79.                                 continue
    
80.                                
    
81.                         aLine = aLine[1:]
    
82.                         dataStr = dataStr + aLine
    
83.                         for i in range(0, 25):
    
84.                                 if xIdx == 360:
    
85.                                         xIdx = 0
    
86.                                         yIdx += 1
    
87.                                         break
    
88.                                
    
89.                                 dStr = aLine[i * 3:i * 3 + 3]
    
90.                                 if dStr == '  0':
    
91.                                         gData[yIdx, xIdx] = gDataInfo.UNDEF
    
92.                                 else:
    
93.                                         gData[yIdx, xIdx] = float(dStr) / 10
    
94.                                 xIdx += 1
    
95.                        
    
96.                 dataFile.close()
    
97.                 gDataInfo.WriteGridData(gData)
    
98.         else:
    
99.                 gDataInfo.WriteGridData_Null()
    
100.                 print 'Null'
     
101.        
     
102.         aTime = aTime.AddDays(1)
     
103. 
     
104. gDataInfo.CloseDataFile()
     
105. 
     
106. #---- Write GrADS control file
     
107. gDataInfo.WriteGrADSCTLFile()
     
108. 
     
109. print 'Finished...'

复制代码

转换后的GrADS二进制数据可以用MeteoInfo绘图，或用脚本进行月平均之类的统计计算，也可以插值到站点和地面观测进行对比分析，图示如下：

mofangbao

木有用过  支持一下

michael_grace

王老师，太有用了，谢谢！

♂雨已~○

没用过这个数据，不过看起来挺漂亮的，牛逼~~~~王老师~~~

qw5622390

强烈支持！！！

2292163582

谢谢老师

2292163582

谢谢老师