格点插值为站点数据批量处理

lovechang1314

MeteoInfo 发表于 2015-7-19 21:45
脚本里有对时间的循环

怪我没看细心。。谢谢

wsjqxjy

谢谢，很好很实用！

大侠MVP

您好：我这个上来就报错了。。1.2.8版本

Traceback (most recent call last):
  File "<iostream>", line 2, in <module>
NameError: name 'readtable' is not defined

求解惑

MeteoInfo

大侠MVP 发表于 2015-9-8 11:49
您好：我这个上来就报错了。。1.2.8版本

Traceback (most recent call last):

在MeteoInfoLab里运行，不是在MeteoInfo里。

大侠MVP

MeteoInfo 发表于 2015-9-8 12:51
在MeteoInfoLab里运行，不是在MeteoInfo里。

依然报错

>>> run script...
Longitude
Latitude
T2
Traceback (most recent call last):
  File "<iostream>", line 16, in <module>
  File "D:\MeteoInfo\pylib\mipylib\dimarray.py", line 407, in tostation
    gdata = self.asgriddata()
  File "D:\MeteoInfo\pylib\mipylib\dimarray.py", line 30, in asgriddata
    xdata = self.dims[1].getDimValue()
IndexError: index out of range: 1
>>>

MeteoInfo

大侠MVP 发表于 2015-9-8 15:29
依然报错

>>> run script...

你把脚本程序都贴出来看看

大侠MVP

0000000000

大侠MVP

MeteoInfo 发表于 2015-9-8 16:40
你把脚本程序都贴出来看看

折腾一整天了，麻烦大神帮忙看下，另外meteoinfo里一个时次一个时次手点是可以的，脚本却不行

脚本
#Open station file as table data
stdata = readtable('E:/STATION/2015050101.csv', delimiter=',', format='%s%f%f')
x = stdata['Longitude']
y = stdata['Latitude']
#Add data file
f = addfile('E:/STATION/wrfout_d02.ctl')
#Get data variable
psv = f['T2']
#Get time dimension length
tn = psv.dimlen(0)
#Loop
for i in range(0, tn):
    #Get dimension array
    ps = psv[i,(0,197),(0,227)]
    #Interpolate to stations
    ps_st = ps.tostation(x, y)
    #Add column to table data
    colname = 'T2' + str(i)
    stdata.addcol(colname, '%d', ps_st)
#Save table data to a file
fn = 'E:/STATION/test_st.csv'
stdata.savefile(fn)
print 'Finish...'

格点数据CTL
dset E:\STATION\wrfout_d02.dat
options  byteswapped
undef 1.e30
title  OUTPUT FROM WRF V3.5.1 MODEL
pdef  198 228 lcc  39.023  116.560   99.500  114.500  39.02700  39.02700  115.57000   3000.000   3000.000
xdef  548 linear  112.89408   0.01351351
ydef  478 linear   35.75130   0.01351351
zdef   31 linear 1 1  
tdef  193 linear 00Z30APR2015      60MN      
VARS    9
XLAT           1  0  LATITUDE, SOUTH IS NEGATIVE (degree_north)
XLONG          1  0  LONGITUDE, WEST IS NEGATIVE (degree_east)
T2             1  0  TEMP at 2 M (K)
PSFC           1  0  SFC PRESSURE (Pa)
U10            1  0  U at 10 M (m s-1)
V10            1  0  V at 10 M (m s-1)
HGT            1  0  Terrain Height (m)
td2            1  0  Dewpoint Temperature at 2m (C)
rh2            1  0  Relative Humidity at 2m (%)
ENDVARS
@ global String comment TITLE =  OUTPUT FROM WRF V3.5.1 MODEL
@ global String comment START_DATE = 2015-04-30_00:00:00
@ global String comment SIMULATION_START_DATE = 2015-04-30_00:00:00
@ global String comment WEST-EAST_GRID_DIMENSION =   199
@ global String comment SOUTH-NORTH_GRID_DIMENSION =   229
@ global String comment BOTTOM-TOP_GRID_DIMENSION =    32
@ global String comment DX =      3000.00
@ global String comment DY =      3000.00
@ global String comment STOCH_FORCE_OPT =     0
@ global String comment GRIDTYPE = C
@ global String comment DIFF_OPT =     1
@ global String comment KM_OPT =     4
@ global String comment DAMP_OPT =     0
@ global String comment DAMPCOEF =         0.20
@ global String comment KHDIF =         0.00
@ global String comment KVDIF =         0.00
@ global String comment MP_PHYSICS =     0
@ global String comment RA_LW_PHYSICS =     1
@ global String comment RA_SW_PHYSICS =     1
@ global String comment SF_SFCLAY_PHYSICS =     1
@ global String comment SF_SURFACE_PHYSICS =     2
@ global String comment BL_PBL_PHYSICS =     1
@ global String comment CU_PHYSICS =     0
@ global String comment SURFACE_INPUT_SOURCE =     1
@ global String comment SST_UPDATE =     0
@ global String comment GRID_FDDA =     0
@ global String comment GFDDA_INTERVAL_M =     0
@ global String comment GFDDA_END_H =     0
@ global String comment GRID_SFDDA =     0
@ global String comment SGFDDA_INTERVAL_M =     0
@ global String comment SGFDDA_END_H =     0
@ global String comment HYPSOMETRIC_OPT =     2
@ global String comment SF_URBAN_PHYSICS =     0
@ global String comment SHCU_PHYSICS =     0
@ global String comment MFSHCONV =     0
@ global String comment FEEDBACK =     1
@ global String comment SMOOTH_OPTION =     0
@ global String comment SWRAD_SCAT =         1.00
@ global String comment W_DAMPING =     0
@ global String comment DT =        18.00
@ global String comment RADT =        30.00
@ global String comment BLDT =         0.00
@ global String comment CUDT =         5.00
@ global String comment SWINT_OPT =     0
@ global String comment MOIST_ADV_OPT =     1
@ global String comment SCALAR_ADV_OPT =     1
@ global String comment TKE_ADV_OPT =     1
@ global String comment DIFF_6TH_OPT =     0
@ global String comment DIFF_6TH_FACTOR =         0.12
@ global String comment OBS_NUDGE_OPT =     0
@ global String comment BUCKET_MM =        -1.00
@ global String comment BUCKET_J =        -1.00
@ global String comment PREC_ACC_DT =         0.00
@ global String comment SF_OCEAN_PHYSICS =     0
@ global String comment ISFTCFLX =     0
@ global String comment ISHALLOW =     0
@ global String comment DFI_OPT =     0
@ global String comment WEST-EAST_PATCH_START_UNSTAG =     1
@ global String comment WEST-EAST_PATCH_END_UNSTAG =   198
@ global String comment WEST-EAST_PATCH_START_STAG =     1
@ global String comment WEST-EAST_PATCH_END_STAG =   199
@ global String comment SOUTH-NORTH_PATCH_START_UNSTAG =     1
@ global String comment SOUTH-NORTH_PATCH_END_UNSTAG =   228
@ global String comment SOUTH-NORTH_PATCH_START_STAG =     1
@ global String comment SOUTH-NORTH_PATCH_END_STAG =   229
@ global String comment BOTTOM-TOP_PATCH_START_UNSTAG =     1
@ global String comment BOTTOM-TOP_PATCH_END_UNSTAG =    31
@ global String comment BOTTOM-TOP_PATCH_START_STAG =     1
@ global String comment BOTTOM-TOP_PATCH_END_STAG =    32
@ global String comment GRID_ID =     2
@ global String comment PARENT_ID =     1
@ global String comment I_PARENT_START =    27
@ global String comment J_PARENT_START =    26
@ global String comment PARENT_GRID_RATIO =     3
@ global String comment CEN_LAT =        39.02
@ global String comment CEN_LON =       116.56
@ global String comment TRUELAT1 =        39.03
@ global String comment TRUELAT2 =        39.03
@ global String comment MOAD_CEN_LAT =        39.03
@ global String comment STAND_LON =       115.57
@ global String comment POLE_LAT =        90.00
@ global String comment POLE_LON =         0.00
@ global String comment GMT =         0.00
@ global String comment JULYR =  2015
@ global String comment JULDAY =   120
@ global String comment MAP_PROJ =     1
@ global String comment MAP_PROJ_CHAR = Lambert Conformal
@ global String comment MMINLU = USGS
@ global String comment NUM_LAND_CAT =    24
@ global String comment ISWATER =    16
@ global String comment ISLAKE =    -1
@ global String comment ISICE =    24
@ global String comment ISURBAN =     1
@ global String comment ISOILWATER =    14

MeteoInfo

大侠MVP 发表于 2015-9-8 20:21
折腾一整天了，麻烦大神帮忙看下，另外meteoinfo里一个时次一个时次手点是可以的，脚本却不行

脚本

两个问题：
1、ps = psv[i,(0,197),(0,227)]  语句的问题，如果要取整个区域用此语句： ps = psv[i,:,:]
2、数据是Lambert投影，因此在插值前需要将站点经纬度数组投影到Lambert投影。
    lx, ly = project(x, y, toproj=f.proj)
    ps_st = ps.tostation(lx, ly)

大侠MVP

MeteoInfo 发表于 2015-9-8 20:43
两个问题：
1、ps = psv  语句的问题，如果要取整个区域用此语句： ps = psv
2、数据是Lambert投影，因 ...

老兄，我真不知道该说啥了。。太厉害了

这个Lambert投影在meteoinfo里插值也要先转换站点经纬度的吧？
我手动在meteoinfo里转换的时候数据都不对应该与此有关。

非常感谢您！