NCL绘制降雨

382702235

请问一下，我在网上搜到的这个NCL的降雨程序，我想用它画我需要的时间段应该怎么修改，为啥我只能读取一个时间点，1是我的图，2是他的图，程序如下：
;**********************************************************
  load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
  load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
  load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
  load "$NCARG_ROOT/lib/ncarg/nclscripts/cnmap/cnmap.ncl"
;**********************************************************

begin

begTime = get_cpu_time()

diri_constant= "./"

;; Open WRF constant file and read variables from it.
a  = addfile(diri_constant + "wrfout_d01_2020-06-11_00_00_00","r")

lon2d = a->XLONG(0,:,:)               
lat2d = a->XLAT(0,:,:)                 
lat = lat2d(:,0)              ;; Xarray of latitudes
lon = lon2d(0,:)              ;; Xarray of Longitudes
times = wrf_user_list_times(a)

minLat = min(lat)
maxLat = max(lat)
minLon = min(lon)
maxLon = max(lon)

;calculate the total precipitation
rain1 = wrf_user_getvar(a,"RAINNC",-1)
rain2 = wrf_user_getvar(a,"RAINC",-1)
rain3 = wrf_user_getvar(a,"RAINSH",-1)
dim_sizes=dimsizes(rain1)
NT = dim_sizes(0)
NY = dim_sizes(1)
NX = dim_sizes(2)

deltaRain1 = rain1(NT-1,:,:) - rain1(0,:,:)
deltaRain2 = rain2(NT-1,:,:) - rain2(0,:,:)
deltaRain3 = rain3(NT-1,:,:) - rain3(0,:,:)
rain_tot   = deltaRain1 + deltaRain2 + deltaRain3

delete(rain1)
delete(rain2)

rain_tot!0 = "lat"
rain_tot!1 = "lon"
rain_tot&lat = lat
rain_tot&lon = lon
rain_tot@description = "Precipitation"
rain_tot@units = "mm"

; create a window
fig_name = "precipitaion_"+times(0)+"-"+times(NT-1)+"UTC.png"
wks = gsn_open_wks("png",fig_name)

;define colormap
cmap =  (/ (/255, 255, 255/), \
            (/0, 0,0/), \
            (/242, 255, 255/), \
            (/242, 242, 242/), \
            (/154, 192, 205/), \
            (/178, 223, 238/), \
            (/191, 239, 255/), \
            (/  0, 235, 235/), \
            (/  0, 163, 247/), \
            (/  0, 255,   0/), \
            (/  0, 199,   0/), \
            (/  0, 143,   0/), \
            (/  0,  63,   0/), \
            (/255, 255,   0/), \
            (/255, 143,   0/), \
            (/255,   0,   0/), \
            (/215,   0,   0/), \
            (/191,   0,   0/), \
            (/255,   0, 255/), \
            (/155,  87, 203/), \
            (/ 92,  52, 176/) /) /255.0


gsn_define_colormap(wks,cmap)
;gsn_define_colormap(wks,"Rainbow")

res = True
res@gsnDraw       = False
res@gsnFrame      = False
res@gsnAddCyclic  = False
res@mpFillOn      = False

res@mpMaxLatF = maxLat                      ; specify the plot domain
res@mpMinLatF = minLat                      ;
res@mpMaxLonF = maxLon                      ;
res@mpMinLonF = minLon  
res@tiMainString = "Precipitation: "+times(0)+" to "+times(NT-1)+"UTC"
res@tiMainFontHeightF = 0.013
res@gsnLeftStringFontHeightF  = 0.012
res@gsnRightStringFontHeightF = 0.012


;res@cnLevelSelectionMode = "ExplicitLevels"
res@cnLevelSelectionMode = "ManualLevels"
res@cnLevelSpacingF      = 10
res@cnMinLevelValF       = 0
res@cnMaxLevelValF       = 150

res@cnFillOn = True
res@cnLinesOn = False
res@lbLabelAutoStride = True
res@cnFillDrawOrder = "PreDraw"
res@cnInfoLabelOn = False

res@lbLabelBarOn          = True        ; 色标
res@pmLabelBarHeightF     = 0.15
res@pmLabelBarWidthF      = 0.6
res@pmLabelBarOrthogonalPosF = 0.07
res@lbLabelFontHeightF    = 0.010
res@cnInfoLabelOn         = False            ; 去掉图底端的标签信息

lon_value = fspan(minLon, maxLon, 5)
lon_label = lon_value + "~S~o~N~E"
res@tmXBMode = "Explicit"
res@tmXBValues = lon_value
res@tmXBLabels = lon_label
res@tmXBLabelFontHeightF = 0.01
lat_value = fspan(minLat, maxLat , 5)
lat_label = lat_value + "~S~o~N~N"
res@tmYLValues = lat_value
res@tmYLLabels = lat_label

;; plot
plot    = gsn_csm_contour_map(wks,rain_tot,res)

; add China map
cnres           = True
cnres@china     = True       ;draw china map or not
cnres@river     = False       ;draw changjiang&huanghe or not
cnres@province  = True      ;draw province boundary or not
cnres@nanhai    = False       ;draw nanhai or not
cnres@diqu      = True       ; draw diqujie or not

chinamap = add_china_map(wks,plot,cnres)

draw(plot)
frame(wks)

print("Time Required " + (get_cpu_time() - begTime) + " seconds")

end

张如也

检查rain_tot这个变量，从图上看他是定值0，可能是计算错误。

382702235

张如也 发表于 2024-7-18 23:10
检查rain_tot这个变量，从图上看他是定值0，可能是计算错误。

谢谢您，但是我查看时间变量发现它好像只读取到了下面设置的这一个文件，我的文件是分开的每一个小时一个文件，应该首先是这里的问题吧，但是我改用addfiles直接读多个文件后面运行又会出错，请问这里怎么修改
diri_constant= "./"
;; Open WRF constant file and read variables from it.
a  = addfile(diri_constant + "wrfout_d01_2020-06-11_00_00_00","r")

张如也

addfiles的报错信息看不到。
授人以鱼不如授人以渔，cd_calendar函数，你学习后，可以想办法做循环读取。

382702235

张如也 发表于 2024-7-25 08:50
addfiles的报错信息看不到。
授人以鱼不如授人以渔，cd_calendar函数，你学习后，可以想办法做循环读取。

好的好的，谢谢您

kongcw

你好，问题解决了吗

kongcw

你好，问题解决了吗