Re: interpolation to PBL height

All WRF questions should be sent to: wrfhelp@ucar.edu
This email is being forwarded as a cc

Given that it is Christmas Eve and tomorrow is Christmas Day,
you may not receive a reply for a few days.

Cheers
On 12/24/12 4:00 AM, BasitAli Khan wrote:
> Hi guys,
>
> I want to interpolate u, v, and w to the PBL height from my wrfoutput. Since PBL change/may change on every grid point, i therefore used wrf_interp_1d to interpolate every grid point value of u, v and w to each grid point value of PBL. I applied double do loop to accomplish this. Apparently it runs fine also however, this double do loop take almost one hour to process a single frame of my domain that comprises of 600 x 600 grid points. This makes wrf_interp_1d impractical for use.
>
> I would greatly appreciate if somebody could suggest a better and faster way to process these data. My ncl code is copied below.
>
> Best regards,
> Basit.
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>
> begin
>
> a = addfile("/project/k14/basit/wrf-sims/wrf_168/wrfout_d02_2006-05-31_06_60min.nc","r")
>
> ; We generate plots, but what kind do we prefer?
> type = "x11"
> wks = gsn_open_wks(type,"plt_HeightLevel")
> gsn_define_colormap(wks,"ViBlGrWhYeOrRe") ;WhiteYellowOrangeRed")
>
> ; Set some basic resources
> res = True
> res@MainTitle = "REAL-TIME WRF"
> res@Footer = False
>
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> times = wrf_user_getvar(a,"times",-1) ; get all times in the file
> ntimes = dimsizes(times) ; number of times in the file
> print(ntimes)
> lat2d = a->XLAT(0,:,:)
> lon2d = a->XLONG(0,:,:)
>
> pltres = True
> mpres = True
> mpres@gsnMaximize = True
> mpres@mpOutlineOn = True
> mpres@mpFillOn = False
> mpres@mpGridLineThicknessF = 0.5
> mpres@mpGridSpacingF = 3
> mpres@gsnDraw = False ; don't draw the plots
> mpres@gsnFrame = False ; don't advance the frame
> mpres@tmYROn = False
> mpres@tmXTOn = False
> mpres@mpGeophysicalLineColor = "Black"
> mpres@mpLimbLineColor = "Blue"
> mpres@mpGeophysicalLineThicknessF = 1.5
> mpres@mpDataBaseVersion = "HighRes" ; Use the high-res database
> mpres@mpDataResolution = "FinestResolution"
> mpres@pmTickMarkDisplayMode = "Always" ; Turn on map tickmarks.
> mpres@mpNationalLineColor = "Black"
> mpres@mpNationalLineThicknessF = 2.5
> mpres@mpGridLineColor = "Gray49" ; Grid line color
> mpres@mpGridLineDashPattern = 14.
> mpres@mpGridLineThicknessF = 0.8
> mpres@mpGridLatSpacingF = 1.
> mpres@mpGridLonSpacingF = 1.
> mpres@mpOutlineBoundarySets = "National"
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> do it = 10, ntimes-1, 1 ; TIME LOOP
>
> print("Working on time: " + times(it) )
> res@TimeLabel = times(it) ; Set Valid time to use on plots
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
> ; First get the variables we will need
>
> ua = wrf_user_getvar(a,"ua",it) ; u averaged to mass points
> va = wrf_user_getvar(a,"va",it) ; v averaged to mass points
> wa = wrf_user_getvar(a,"wa",it)
> z = wrf_user_getvar(a, "z",it) ; grid point height
> pbl= wrf_user_getvar(a,"PBLH",it)
>
> ua@lat2d = lat2d
> ua@lon2d = lon2d
> va@lat2d = lat2d
> va@lon2d = lon2d
> wa@lat2d = lat2d
> wa@lon2d = lon2d
> z@lat2d = lat2d
> z@lon2d = lon2d
> pbl@lat2d = lat2d
> pbl@lon2d = lon2d
>
> ; ==================================================================================
>
>
> dims_zZ = dimsizes(z)
> xn = dims_zZ(2)
> yn = dims_zZ(1)
>
> print(dims_zZ)
>
> ua_plane = new((/dims_zZ(1), dims_zZ(2)/), float)
> va_plane = new((/dims_zZ(1), dims_zZ(2)/), float)
> wa_plane = new((/dims_zZ(1), dims_zZ(2)/), float)
>
> pbl_ASL = z(0,:,:) + pbl ; pbl is above surface while z is above sea level (ASL).
> printVarSummary(pbl_with_ter)
>
> ua_plane@_FillValue = -999.99
> va_plane@_FillValue = -999.99
> wa_plane@_FillValue = -999.99
>
> do xi =0, xn-1, 1
> do yi = 0, yn-1, 1
> ua_plane(yi,xi) = wrf_interp_1d(ua(:,yi,xi), z(:,yi,xi), pbl_ASL(yi,xi))
> va_plane(yi,xi) = wrf_interp_1d(va(:,yi,xi), z(:,yi,xi), pbl_ASL(yi,xi))
> wa_plane(yi,xi) = wrf_interp_1d(wa(:,yi,xi), z(:,yi,xi), pbl_ASL(yi,xi))
> print(ua_plane(yi,xi)+"----"+z(0,yi,xi)+"---"+pbl_ASL(yi,xi))
> end do
> end do
>
> printVarSummary(ua_plane)
> wa_plane@description = "Vertical Velocity at PBL height"
> wa_plane@units = "~N~m s~S~-1~N~"
> ua_plane@units = "~N~m s~S~-1~N~"
> va_plane@units = "~N~m s~S~-1~N~"
>
> ; ==================================================================================
>
> ; Plotting options for Wa
> opts = res
> opts@cnFillOn = True
> opts@ContourParameters = (/-1., 1., 0.1/)
> cnWa = wrf_contour(a,wks,wa_plane,opts)
> delete(opts)
>
> ; -------------------------------------------
>
> ; Plotting options for Wind Vectors
>
> optsVc = res
> optsVc@FieldTitle = "Wind vector at PBL height"
> optsVc@vcRefAnnoOn = True ; to turn on/off vector/wbarb legend
> optsVc@gsnDraw = False ; Don't draw individual plot
> optsVc@gsnFrame = False ; Don't advance frame.
> optsVc@lbLabelStride = 2 ; every other color
> optsVc@gsnSpreadColors = True ; use full range of color map
> optsVc@vcRefMagnitudeF = 10.0 ; define vector ref mag
> optsVc@vcRefLengthF = 0.030 ; define length of vec ref
> optsVc@vcGlyphStyle = "CurlyVector" ;"CurlyVector" ; turn on curley vectors
> optsVc@vcMinDistanceF = 0.04 ; thin out vectors
> optsVc@vcVectorDrawOrder= "PostDraw"
>
> vector = wrf_vector(a,wks,ua_plane, va_plane, optsVc)
> delete(optsVc)
>
> ; MAKE PLOTS
> plot = wrf_map_overlays(a,wks,(/cnWa,vector/),pltres,mpres)
>
>
> ; end do ; END OF LEVEL LOOP
>
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
> delete(va_plane)
> delete(ua_plane)
> delete(wa_plane)
> delete(pbl_with_ter)
>
> end do ; END OF TIME LOOP
>
> end
>
>
>
>
>
>
> ________________________________
>
> This message and its contents including attachments are intended solely for the original recipient. If you are not the intended recipient or have received this message in error, please notify me immediately and delete this message from your computer system. Any unauthorized use or distribution is prohibited. Please consider the environment before printing this email.
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon Dec 24 08:42:30 2012