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/WRFUserARW.ncl" begin ; ; The WRF ARW input file. ; This needs to have a ".nc" appended, so just do it. dir = "/media/schaefer_hd/" base_folder = "BASE/" fire_folder = "FIRE/" sensible_folder = "SENSIBLE/" latent_folder = "LATENT/" files = systemfunc("ls " + dir + base_folder + "wrfout_d03_2012-08-30*") nfiles = dimsizes(files) inner_dom = addfile(dir + fire_folder + "wrfout_d03_2012-08-29_19:00:00.nc","r") outer_dom = addfile(dir + fire_folder + "wrfout_d03_2012-08-29_19:00:00.nc","r") inner_lat = wrf_user_getvar(inner_dom, "lat", 0) inner_lon = wrf_user_getvar(inner_dom, "lon", 0) inner_nlat = dimsizes(inner_lat) inner_nlon = dimsizes(inner_lon) outer_lat = wrf_user_getvar(outer_dom, "lat", 0) outer_lon = wrf_user_getvar(outer_dom, "lon", 0) outer_nlat = dimsizes(outer_lat) outer_nlon = dimsizes(outer_lon) fire_cells_lat = wrf_user_getvar(inner_dom, "FXLAT", 0) fire_cells_lon = wrf_user_getvar(inner_dom, "FXLONG", 0) points = (/ (/ 66,100,0 /), (/ 66,140,1 /), (/ 66,180,2 /),\ (/ 116,100,3 /), (/ 116,140,4 /), (/ 116,180,5 /),\ (/ 166,100,6 /), (/ 166,140,7 /), (/ 166,180,8 /) /) pt_strs = (/ "P1", "P2", "P3", "P4", "P5", "P6", "P7", "P8", "P9" /) ; type = "x11" ; type = "pdf" ; type = "ps" ; type = "ncgm" type = "png" type@wkHeight = 2500 type@wkWidth = 2500 wks = gsn_open_wks(type,"points_fire_perim") setvalues NhlGetWorkspaceObjectId() "wsMaximumSize" : 200000000 end setvalues gsn_define_colormap(wks,"WhBlGrYeRe") ; gsn_define_colormap(wks,"radar") ; Set some Basic Plot options res = True res@tfDoNDCOverlay = True res@gsnAddCyclic = False res@gsnSpreadColors = True res@gsnMaximize = True res@tiMainFont = 25 res@tiMainString = "" res@tiMainOffsetYF = 0.02 res@tiMainFontHeightF = 0.0225 res@tiMainFontThicknessF = 40 res@gsnDraw = False res@gsnFrame = False pltres = True projres = res ;wrf_map_resources(a_2,res) projres@mpProjection = "LambertConformal" projres@mpRightCornerLatF = outer_lat(0,0) projres@mpRightCornerLonF = outer_lon(0,0) projres@mpLeftCornerLatF = outer_lat(outer_nlon-1,outer_nlat-1) projres@mpLeftCornerLonF = outer_lon(outer_nlon-1,outer_nlat-1) projres@mpLambertParallel1F = outer_dom@TRUELAT1 projres@mpLambertParallel2F = outer_dom@TRUELAT2 projres@mpLambertMeridianF = outer_dom@STAND_LON projres@mpCenterLonF = outer_dom@CEN_LON projres@mpCenterLatF = outer_dom@CEN_LAT projres@mpLimitMode = "Corners" projres@pmTickMarkDisplayMode = "Always" projres@mpDataBaseVersion = "Ncarg4_1" projres@mpDataSetName = "Earth..4" projres@mpOutlineBoundarySets = "AllBoundaries" ; all boundaries projres@mpNationalLineThicknessF = 1.0 projres@mpGeophysicalLineThicknessF = 1.0 projres@mpUSStateLineThicknessF = 2.0 projres@mpFillOn = False ;projres@mpDefaultFillColor = "white" ;projres@mpFillColor = "white" projres@tiMainString = "Points for Time-Height Plots" ;projres@lgAutoManage = False ;projres@pmLegendSide = "Bottom" ;projres@pm ;projres@lgLegendOn = True ;projres@lgLineColors = (/ "Black","Red","Blue" /) ;projres@lgOrientation = "Horizontal" ;projres@lgPerimOn = False ;projres@lgLabelFontHeightF = .08 ;projres@vpWidthF = 0.15 ;projres@vpHeightF = 0.1 ;projres@lgMonoDashIndex = True ;projres@pmLegendDisplayMode = "Always" ;projres@lgLabelString = (/ "FIRE","LATENT","SENSIBLE" /) do nf = 0 ,nfiles-1 outer_1 = addfile(files(nf) + ".nc" ,"r") ;;; These three plots have the fire with them file_2 = str_sub_str(files(nf),"BASE","FIRE") outer_2 = addfile(file_2 + ".nc" ,"r") file_3 = str_sub_str(files(nf),"BASE","LATENT") outer_3 = addfile(file_3 + ".nc" ,"r") file_4 = str_sub_str(files(nf),"BASE","SENSIBLE") outer_4 = addfile(file_4 + ".nc" ,"r") ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ; What times and how many time steps are in the data set? times = wrf_user_list_times(outer_1) ; get times in the file ntimes = dimsizes(times) ; number of times in the file angle = 0.0 pii = 3.14159 aspect_ratio = .7 ; This is the big loop over all of the time periods to process. do it =0,0;ntimes-1 time = it ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; z = wrf_user_getvar(outer_1, "z",time) ; grid point height z@lat2d = outer_lat z@lon2d = outer_lon perim = wrf_user_getvar(outer_2, "FIRE_AREA",time) perim@lat2d = fire_cells_lat perim@lon2d = fire_cells_lon perim_2 = wrf_user_getvar(outer_3, "FIRE_AREA",time) perim_2@lat2d = fire_cells_lat perim_2@lon2d = fire_cells_lon perim_3 = wrf_user_getvar(outer_4, "FIRE_AREA",time) perim_3@lat2d = fire_cells_lat perim_3@lon2d = fire_cells_lon conres = True conres@gsnDraw = False conres@gsnFrame = False conres@cnLinesOn = True conres@cnLevelSelectionMode = "ManualLevels" conres@cnMinLevelValF = .9 conres@cnMaxLevelValF = 1 conres@cnLevelSpacingF = .1 conres@cnLineThicknessF = 3 conres@cnInfoLabelOn = False conres@gsnRightString = "" conres@gsnLeftString = "" conres_2 = True conres_2@gsnDraw = False conres_2@gsnFrame = False conres_2@cnLinesOn = True conres_2@cnLevelSelectionMode = "ManualLevels" conres_2@cnMinLevelValF = .9 conres_2@cnMaxLevelValF = 1 conres_2@cnLevelSpacingF = .1 conres_2@cnLineThicknessF = 3 conres_2@cnInfoLabelOn = False conres_2@gsnRightString = "" conres_2@gsnLeftString = "" cont = gsn_csm_contour(wks,perim,conres) projres@tmXTOn = False projres@tmYROn = False delete(projres@gsnSpreadColors) delete(projres@gsnAddCyclic) ;projres@tiMainString = conres@tiMainString map = gsn_csm_map(wks,projres) overlay(map,cont) conres_2@cnLineColor = "Red" cont_2 = gsn_csm_contour(wks,perim_2,conres_2) overlay(map,cont_2) conres_2@cnLineColor = "Blue" cont_3 = gsn_csm_contour(wks,perim_3,conres_2) overlay(map,cont_3) ;------------Legend---------------------; lgres = True lgres@lgAutoManage = False lgres@lgLegendOn = True lgres@lgLineColors = (/ "Black","Red","Blue" /) lgres@lgOrientation = "Horizontal" lgres@lgPerimOn = False lgres@lgLabelFontHeightF = .02 lgres@vpWidthF = 0.5 lgres@vpHeightF = 0.1 lgres@lgMonoDashIndex = True lgres@lgLabelPosition = "Center" lgres@lgLineThicknessF = 3 ;lgres@lgLabelOffsetF = 0.1 ;lgres@lgLabelJust = "BottomCenter" lg_id = gsn_create_legend(wks,3,(/ "FIRE","LATENT","SENSIBLE" /), lgres) amres = True amres@amParallelPosF = 0 amres@amOrthogonalPosF = 0.45 anno_id = gsn_add_annotation(map,lg_id,amres) ;---------Plot Crawford -----------------; ; crawford_lat =42.72 ; crawford_lon =-102.97 mkres = True mkres@gsMarkerIndex = 16 mkres@gsMarkerSizeF = 0.015 mkres@gsMarkerColor = "Black" ; mkid_1 = gsn_add_polymarker(wks,plot,crawford_lon,crawford_lat,mkres) ; mkid_2 = gsn_add_polymarker(wks,plot(1),crawford_lon,crawford_lat,mkres) ; mkid_3 = gsn_add_polymarker(wks,plot(2),crawford_lon,crawford_lat,mkres) ; mkid_4 = gsn_add_polymarker(wks,plot(3),crawford_lon,crawford_lat,mkres) pt_mkrs = new(9,graphic) txtres = True txtres@txFontHeightF = 0.03 do pt = 0,8 point = points(pt,:) ll = wrf_user_ij_to_ll(inner_dom,point(0)+12,point(1),True) ll_2 = wrf_user_ij_to_ll(inner_dom,point(0),point(1),True) text = gsn_add_text(wks,map,pt_strs(pt),ll(0),ll(1),txtres) pt_mkrs(pt) = gsn_add_polymarker(wks,map,ll_2(0),ll_2(1),mkres) end do maximize_output(wks,True) ;draw(map) ;frame(wks) end do end do end