[ncl-talk] hovmoller diagram and reshape function
Ramchandra Karki
rammetro at hotmail.com
Sun Aug 7 11:16:17 MDT 2016
Dear NCL users,I have hourly wrf arw output and plotted the hovmoller diagramme along the cross section of two points for july as monthly average of each hour (average value of 0 hour, 1 hour ........23 hour of each day) precipitation for whole month ( 24 values ).ploted hovmoller for hourly time vs cross section grid of two points.i have used reshape function however i am confused whether its producing right output or not.(use of reshape function for my purpose). (array as 24, number of days, ncross_section_points and averaging for days furthermore, I would like to convert the UTC time to local time (5:45 hour ahead of UTC) in my plot but could not get idea to do so.I look forward for getting feedback from you.my script
; This script is used to plot ; get each hour average value for 24 hours; get data along two point section; plot hovmoller for time vs grid along section
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/wrf/WRF_contributed.ncl"load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
begin; The WRF ARW input file. ; This needs to have a ".nc" appended, so just do it. DATADir = "/geodata/data/private/fgvy024/Build_WRF/WRFV6/run/" FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d03_2015-07*") ; contain hourly time stamp numFILES = dimsizes(FILES) a = addfiles(FILES+".nc","r") type = "png" ; map type wks = gsn_open_wks(type,"july_hovmoller_hourly_ppt") ; Create a plot workstation rain_exp = wrf_user_getvar(a,"RAINNC",-1) rain_con = wrf_user_getvar(a,"RAINC",-1) rain = rain_exp + rain_con Remp = dimsizes(rain) nt = Remp(0) nnlat = Remp(1) nlon = Remp(2); get hourly rainfall Ree = rain(1:nt-1,:,:) - rain(0:nt-2,:,:) ; to get hourly precipitation ; get new variable with full time step to make it multiplable by 24 Re = new((/nt, nnlat, nlon/),float) Re(1:nt-1,:,:) = Ree(:,:, :) Re(0, :, :) = 0 ; first array will have zero value dims = dimsizes(Re) pA_lon = 86.2233 ; longitude of point A pA_lat = 27.88163333 ; lat of point A pB_lon = 86.541725 ; similar for point B pB_lat = 27.89950278 pointA = wrf_user_ll_to_ij(a, pA_lon, pA_lat, True) ; get array point pointB = wrf_user_ll_to_ij(a, pB_lon, pB_lat, True) xA = pointA(0) - 1 yA = pointA(1) - 1 xB = pointB(0) - 1 yB = pointB(1) - 1 opts = True plane = new(4,float) ; for section along two points plane = (/ yA,xA, yB,xB /)
; interpolate value along the section of A and B R_plane = wrf_user_intrp2d(Re,plane,0.,opts) Re_dims = dimsizes(R_plane) ntim = Re_dims(0) ncross = Re_dims(1) ndays = ntim/24 Re3d = reshape(R_plane,(/24,ndays, ncross/)) ; hours, days, crosssection copy_VarAtts(R_plane,Re3d) Re3d_hr_avg = dim_avg_n_Wrap(Re3d,1) ; hourly average (monthly, 0 hours average, 1 hour average, etc) res = True ; plot mods desired res at cnFillOn = True ; turn on color fill res at cnFillPalette = "BlWhRe" ; set color map res at tiMainString = "Hourly Precipitation" ; title res at cnLevelSelectionMode = "ManualLevels" ; manual contour levels res at cnMinLevelValF = 0. ; min level res at cnMaxLevelValF = 2 ; max level res at cnLevelSpacingF = 0.1 ; contour level spacing plot = gsn_csm_hov(wks, Re3d_hr_avg(:,:), res) end
RegardsRam
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160807/1d2a6523/attachment.html
More information about the ncl-talk
mailing list