MeteoInfoLab脚本示例：同时加载多个文件

MeteoInfo

MeteoInfoLab增加了一个addfiles函数用来同时加载多个数据文件，参数就是一个包含多个文件名的列表，在提取数组的时候可以用类似单个文件的方式方便处理。需要注意的是多个文件必须都有时间维，除了时间维，其它维的设置是相同的，每个文件所包含的变量中要有共同变量。addfiles函数比较适合于处理时间递增的多个数据文件。

示例脚本：

1. datadir = 'D:/Temp/grib'
   
2. st = datetime.datetime(2017,1,1,0)
   
3. et = datetime.datetime(2017,1,1,18)
   
4. fns = []
   
5. while st <= et:
   
6.     fn = os.path.join(datadir, 'fnl_' + st.strftime('%Y%m%d_%H') + \
   
7.         '_00.grib2')
   
8.     print fn
   
9.     fns.append(fn)
   
10.     st = st + datetime.timedelta(hours=6)
    
11. 
    
12. fs = addfiles(fns)
    
13. v = fs['v-component_of_wind_tropopause']
    
14. data = v[:,:,:]
    
15. ndata = mean(data, axis=0)
    
16. 
    
17. #Plot
    
18. axesm()
    
19. lworld = shaperead('D:/Temp/map/country1.shp')
    
20. geoshow(lworld, edgecolor='k')
    
21. layer = imshowm(ndata)
    
22. colorbar(layer)
    
23. title('Mutiple grid data files example\nMean V (' + st.strftime('%Y-%m-%d %H') + ' - ' + \
    
24.     et.strftime('%Y-%m-%d %H') + ')')

datadir = 'D:/Temp/grib'<br /> st = datetime.datetime(2017,1,1,0)<br /> et = datetime.datetime(2017,1,1,18)<br /> fns = []<br /> while st <= et:<br />     fn = os.path.join(datadir, 'fnl_' + st.strftime('%Y%m%d_%H') + \<br />         '_00.grib2')<br />     print fn<br />     fns.append(fn)<br />     st = st + datetime.timedelta(hours=6)<br /> <br /> fs = addfiles(fns)<br /> v = fs['v-component_of_wind_tropopause']<br /> data = v[:,:,:]<br /> ndata = mean(data, axis=0)<br /> <br /> #Plot<br /> axesm()<br /> lworld = shaperead('D:/Temp/map/country1.shp')<br /> geoshow(lworld, edgecolor='k')<br /> layer = imshowm(ndata)<br /> colorbar(layer)<br /> title('Mutiple grid data files example\nMean V (' + st.strftime('%Y-%m-%d %H') + ' - ' + \<br />     et.strftime('%Y-%m-%d %H') + ')')

良辰

王老师，要是gdas数据，文件名没有固定六小时的怎么办？可以list吗？直接写文件名

MeteoInfo

良辰 发表于 2017-1-17 19:27
王老师，要是gdas数据，文件名没有固定六小时的怎么办？可以list吗？直接写文件名

可以

长春市气象台

王老师，最新1.4版本加载多个gdas资料遇到问题，求指点！

1. #Set working directory
   
2. trajDir = 'D:/Temp/HYSPLIT'
   
3. meteoDir = r'D:\Temp\arl'
   
4. # Open meteorological data file
   
5. print 'Open meteorological data file...'
   
6. fns = []
   
7. meteofn1 = os.path.join(meteoDir, 'gdas1.jul09.w1')
   
8. meteofn2 = os.path.join(meteoDir, 'gdas1.jul09.w2')
   
9. meteofn3 = os.path.join(meteoDir, 'gdas1.jul09.w3')
   
10. meteofn4 = os.path.join(meteoDir, 'gdas1.jul09.w4')
    
11. meteofn5 = os.path.join(meteoDir, 'gdas1.jul09.w5')
    
12. fns.append(meteofn1)
    
13. fns.append(meteofn2)
    
14. fns.append(meteofn3)
    
15. fns.append(meteofn4)
    
16. fns.append(meteofn5)
    
17. meteof = addfiles(fns)</P>
    
18. #print 'Meteorological file: ' +meteof
    
19. # Read data
    
20. print 'Read data...'
    
21. latlim = [10,60]
    
22. lonlim = [60,140]
    
23. rh = meteof['RELH'][:,:,latlim,lonlim]

复制代码

错误提示为：
Traceback (most recent call last):
  File "<iostream>", line 23, in <module>
  File "D:\Temp\MeteoInfo\pylib\mipylib\dataset\dimvariable.py", line 383, in __getitem__
    data = minum.concatenate([data, aa])
  File "D:\Temp\MeteoInfo\pylib\mipylib\numeric\minum.py", line 1273, in concatenate
    r = ArrayUtil.concatenate(ars, axis)
        at ucar.ma2.Index3D.set(Index3D.java:139)
        at org.meteoinfo.data.ArrayUtil.concatenate(ArrayUtil.java:1010)
        at sun.reflect.NativeMethodAccessorImpl.invoke0(Native Method)
        at sun.reflect.NativeMethodAccessorImpl.invoke(Unknown Source)
        at sun.reflect.DelegatingMethodAccessorImpl.invoke(Unknown Source)
        at java.lang.reflect.Method.invoke(Unknown Source)

java.lang.ArrayIndexOutOfBoundsException: java.lang.ArrayIndexOutOfBoundsException

MeteoInfo

长春市气象台 发表于 2017-1-18 08:05
王老师，最新1.4版本加载多个gdas资料遇到问题，求指点！

这个新功能测试得很少，有空了我会试试。

MeteoInfo

长春市气象台 发表于 2017-1-18 08:05
王老师，最新1.4版本加载多个gdas资料遇到问题，求指点！

对addfiles功能进行了更新，你在试试。可以在百度网盘上下载：http://pan.baidu.com/s/1kUENi4n

长春市气象台

MeteoInfo 发表于 2017-1-18 18:05
对addfiles功能进行了更新，你在试试。可以在百度网盘上下载：http://pan.baidu.com/s/1kUENi4n

我试一下，谢谢王老师！

良辰

王老师，多个文件加载成功，但是在进行涡度计算时报错，以前是可以的。另外DimArray也报错。代码如下。

1. # Set working directory
   
2. trajDir = 'D:/Temp/HYSPLIT'
   
3. meteoDir = r'D:\Temp\arl'
   
4. # Open meteorological data file
   
5. print 'Open meteorological data file...'
   
6. fns = []
   
7. meteofn1 = os.path.join(meteoDir, 'gdas1.jul09.w1')
   
8. meteofn2 = os.path.join(meteoDir, 'gdas1.jul09.w2')
   
9. meteofn3 = os.path.join(meteoDir, 'gdas1.jul09.w3')
   
10. meteofn4 = os.path.join(meteoDir, 'gdas1.jul09.w4')
    
11. meteofn5 = os.path.join(meteoDir, 'gdas1.jul09.w5')
    
12. fns.append(meteofn1)
    
13. fns.append(meteofn2)
    
14. fns.append(meteofn3)
    
15. fns.append(meteofn4)
    
16. fns.append(meteofn5)
    
17. meteof = addfiles(fns)
    
18. #print 'Meteorological file: ' + meteof
    
19. # Read data
    
20. print 'Read data...'
    
21. latlim = [10,60]
    
22. lonlim = [60,140]
    
23. rh = meteof['RELH'][:,:,latlim,lonlim]
    
24. nx = rh.dimlen(3)
    
25. ny = rh.dimlen(2)
    
26. nz = rh.dimlen(1)
    
27. nt = rh.dimlen(0)
    
28. lat = rh.dimvalue(2)
    
29. lev = rh.dimvalue(1)
    
30. t0 = meteof['TEMP'][:,:nz-1,latlim,lonlim]
    
31. uwnd = meteof['UWND'][:,:nz-1,latlim,lonlim]
    
32. vwnd = meteof['VWND'][:,:nz-1,latlim,lonlim]
    
33. #vort = hcurl(uwnd, vwnd)
    
34. #divg = hdivg(uwnd, vwnd)
    
35. prs = zeros([nt,nz,ny,nx])
    
36. prs = DimArray(prs, rh.dims)

复制代码

错误码为1：（涡度）
Traceback (most recent call last):
  File "<iostream>", line 33, in <module>
  File "D:\MeteoInfo\pylib\mipylib\numeric\minum.py", line 1564, in hcurl
    r = ArrayMath.hcurl(u.asarray(), v.asarray(), xdim.getDimValue(), ydim.getDimValue())
注释掉后错误码为2：（气压数组）
Traceback (most recent call last):
  File "<iostream>", line 36, in <module>
NameError: name 'DimArray' is not defined

MeteoInfo

良辰 发表于 2017-1-20 16:45
王老师，多个文件加载成功，但是在进行涡度计算时报错，以前是可以的。另外DimArray也报错。代码如下。
错 ...

简单测试了一下，vort = hcurl(uwnd, vwnd) 应该没有问题。

prs = DimArray(prs, rh.dims)
改为：
prs = dim_array(prs, rh.dims)

良辰

王老师，程序还是不通呀，能帮忙测试一下吗？