- 积分
- 7108
- 贡献
-
- 精华
- 在线时间
- 小时
- 注册时间
- 2016-7-4
- 最后登录
- 1970-1-1
|
登录后查看更多精彩内容~
您需要 登录 才可以下载或查看,没有帐号?立即注册
x
参考脚本:
http://bbs.06climate.com/forum.php?mod=viewthread&tid=28890
http://bbs.06climate.com/forum.php?mod=viewthread&tid=35437
资料是FNl的grib2,调了一晚上都没成功,跪求大神帮忙,红色是报错
oad "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
begin
;************************************************
; open file and read in variable
;***********************************************
;levels 1000 925 850 700 600 500 400 300
; 层数 0 1 2 3 4 5 6 7
;***********************************************
nlat=73
nlon=144
nyear=1
nlev=8
grb_file = addfile("/home/fuchao/datafile/FNL/fnl_20140809_00_00.grib2","r")
tt = grb_file->TMP_P0_L102_GLL0(:,:,:)
q = grb_file->SPFH_P0_2L108_GLL0(:,:)
u = grb_file->UGRD_P0_L100_GLL0(:,:,:)
v = grb_file->VGRD_P0_L100_GLL0(:,:,:)
slp = grb_file->PRES_P0_L1_GLL0(:,:)
qu1 = new((/8,73,144/),float)
qv1 = new((/8,73,144/),float)
do iz=0,7
do i=0,72
do j=0,143
qu1(iz,i,j) = q(i,j)*u(iz,i,j)
qv1(iz,i,j) = q(i,j)*v(iz,i,j)
end do
end do
end do
lev =(/1000,925,850,700,600,500,400,300/)
lat=fspan(-90,90,73)
lon=fspan(0,357.5,144)
lon!0 = "lon"
lon@units = "degrees_east"
lat!0 = "lat"
lat@units = "degrees_north"
lev!0 = "lev"
lev@units = "hPa"
qu1!0 = "lev"
qu1!1 = "lat"
qu1!2 = "lon"
qu1&lon = lon
qu1&lat = lat
qu1&lev = lev
qv1!0 = "lev"
qv1!1 = "lat"
qv1!2 = "lon"
qv1&lon = lon
qv1&lat = lat
qv1&lev = lev
;**************************************************************************************
qu11 = qu1(lat|:, lon|:, lev|:)
delete(qu1)
qv11 = qv1(lat|:, lon|:, lev|:)
delete(qv1)
p = (/ 1000.,925.,850.,700.,600.,500., \
400.,300.,250.,200.,150.,100., \
70.,50.,30.,20.,10. /)
linlog = 2
pbot = 1100.
ptop = 300.
linlog = 2
; 1为线性积分 2为对数积分
data11 = (vibeta(p,qu11,linlog,slp,1100,300))/9.8 ; 返回的数值是2维的(lat,lon)
data12 = (vibeta(p,qv11,linlog,slp,1100,300))/9.8
qu = new((/73,144/),float)
qv = new((/73,144/),float)
mm = new((/73,144/),float)
mm1 = sqrt(data11^2+data12^2)
do i=0,72
do j=0,143
mm(i,j) = mm1(i,j)
end do
end do
do i=0,72
do j=0,143
qu(i,j) = data11(i,j)
end do
end do
do i=0,72
do j=0,143
qv(i,j) = data12(i,j)
end do
end do
mm!0 = "lat"
mm!1 = "lon"
mm&lat = lat
mm&lon = lon
lon!0 = "lon"
lon@units = "degrees_east"
lat!0 = "lat"
lat@units = "degrees_north"
qu!0 = "lat"
qu!1 = "lon"
qu&lat = lat
qu&lon = lon
lon!0 = "lon"
lon@units = "degrees_east"
lat!0 = "lat"
lat@units = "degrees_north"
qv!0 = "lat"
qv!1 = "lon"
qv&lat = lat
qv&lon = lon
lon!0 = "lon"
lon@units = "degrees_east"
lat!0 = "lat"
lat@untis = "degrees_north"
wks =gsn_open_wks("pdf","fnl_vapour_flux")
gsn_define_colormap(wks,"BlAqGrYeOrReVi200")
vcres = True
vcres@gsnAddCyclic = False
vcres@gsnDraw = False ; don't draw yet
vcres@gsnFrame = False ; don't advance frame yet
vcres@mpDataSetName = "Earth..4"
vcres@mpDataBaseVersion = "MediumRes" ; or "Ncarg4_1"
vcres@mpAreaMaskingOn = True
vcres@mpMaskAreaSpecifiers = (/"China"/)
vcres@mpOutlineSpecifiers = (/"China","China:Provinces"/)
vcres@mpLandFillColor = "white"
vcres@mpInlandWaterFillColor = "white"
vcres@mpOceanFillColor = "white"
vcres@mpFillBoundarySets = "NoBoundaries"
vcres@mpOutlineBoundarySets = "NoBoundaries"
vcres@mpNationalLineColor = "black"
vcres@mpProvincialLineColor = "black"
vcres@mpGeophysicalLineColor = "black"
vcres@mpNationalLineThicknessF = 2
vcres@mpProvincialLineThicknessF = 1
vcres@mpProjection = "CylindricalEquidistant"
vcres@tiMainString =""
vcres@mpGeophysicalLineThicknessF = 0.35 ; double the thickness of geophysical boundaries
vcres@mpNationalLineThicknessF = 0.5 ; double the thickness of national boundaries
vcres@pmTickMarkDisplayMode = "Always"
vcres@mpMinLatF = 0 ;最小纬度
vcres@mpMaxLatF = 50 ;最大纬度
vcres@mpMinLonF = 60 ;最小经度
vcres@mpMaxLonF = 160 ;最大经度
vcres@gsnDraw = False ; don't draw yet
vcres@gsnFrame = False ; don't advance frame yet
vcres@lbLabelBarOn = True
vcres@lbLabelStride = 2 ; plot every other colar bar label
vcres@lbOrientation = "vertical" ; vertical label bars
vcres@vcRefMagnitudeF = 10.0 ; make vectors larger
vcres@vcRefLengthF = 0.018 ; ref vec length
vcres@vcGlyphStyle = "CurlyVector" ; turn on curly vectors
vcres@vcMinDistanceF = 0.017 ; thin out vectors
vcres@vcRefAnnoOn = False
vcres@vcMonoFillArrowFillColor = True
vcres@vcLineArrowThicknessF = 1.0
vcres@vcRefAnnoOn =False
vcres@vcVectorDrawOrder= "PostDraw"
vcres@gsnLeftString =""
vcres@gsnRightString =""
vcres@tiMainString =""
plot(0)=gsn_csm_vector_map(wks,qu,qv,vcres) ; create plot
txres = True
txres@txFontHeightF = 0.020 ; Set the font height
txres@txPerimOn = True
txres@txBackgroundFillColor =0
txres@txPerimColor = 0
;************************************************************************************************
res = True ; plot mods desired
res@gsnDraw = False ; don't draw yet
res@gsnFrame = False ; don't advance frame yet
res@cnFillOn = True ; turn on color
res@gsnSpreadColors = True ; use full colormap
res@gsnSpreadColorStart = 0
res@gsnSpreadColorEnd = 199
res@cnLinesOn = False ; turn off contour lines
res@cnLineLabelsOn = False ; tuen off line labels
res@cnInfoLabelOn =False
res@cnLevelSelectionMode="ExplicitLevels"
res@cnLevels=(/14,18,22,26,30,34,38,42,46,50,52/)
res@gsnLeftString =""
res@gsnRightString =""
res@tiMainString =""
res@lbLabelBarOn = False
res@lbLabelStride = 2 ; plot every other colar bar label
res@lbOrientation = "vertical" ; vertical label bars
mm=smth9(mm,0.5,0.25,False)
plot01 = gsn_csm_contour(wks,mm,res)
overlay(plot(0),plot01)
end
|
|