;---------------------------------------------------------------------- ; coneff_18.ncl ; ; Concepts illustrated: ; - Using functions for cleaner code ; - Setting contour line thicknesses and patterns ; - Using "setvalues" to set resource values ; - Using "getvalues" to retrieve resource values ; - Drawing partially transparent filled contours ;---------------------------------------------------------------------- ; ; These files are loaded by default in NCL V6.2.0 and newer ; 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" ;---------------------------------------------------------------------- ; Function that retrieves contour levels from a plot and changes the ; the requested contour levels to the requested patterns. ; ; This code is somewhat complicated because we are checking to make ; sure we don't override any line patterns that might have ; originally been set by the user. ; ; This code is very similar to the set_contour_line_thicknesses code. ;---------------------------------------------------------------------- undef("set_contour_line_patterns") procedure set_contour_line_patterns(plot,levels_to_change[*]:numeric,\ patterns[*]:integer) local ii, levels, nlevels, line_patterns, nchange, ncolor, n, \ mono_pattern, line_pattern, changed_a_level begin ;---Retrieve the original line patterns (or pattern) used for the plot. getvalues plot@contour "cnLevels" : levels "cnLineDashPatterns" : line_patterns "cnMonoLineDashPattern" : mono_pattern "cnLineDashPattern" : line_pattern end getvalues ;---------------------------------------------------------------------- ; If the original contour plot used a single dash pattern for ; all contour lines, then make sure we use that again for ; any contour lines that are not being changed. Otherwise, ; we assume the user set an array of dash patterns, and ; we'll use those. ;---------------------------------------------------------------------- if(mono_pattern) then line_patterns = line_pattern end if nlevels = dimsizes(levels) nchange = dimsizes(levels_to_change) npattern = dimsizes(patterns) if(npattern.ne.nchange) then print("set_contour_line_patterns: error: the contour line patterns must be an") print("array of the same length as the number of contour levels to change.") end if ;---Array to hold index values of contour levels that need to have a color applied changed_a_level = False do n=0,nchange-1 ii := ind(levels.eq.levels_to_change(n)) if(ismissing(ii(0))) then print("set_contour_line_patterns: warning, no contour level equal to " + levels_to_change(n)) else line_patterns(ii(0)) = patterns(n) changed_a_level = True end if end do if(changed_a_level) then setvalues plot@contour "cnMonoLineDashPattern" : False ; allows an array of line patterns to be set "cnLineDashPatterns" : line_patterns end setvalues end if end ;---------------------------------------------------------------------- ; Function that retrieves contour levels from a plot and changes the ; the requested contour levels to the requested thicknesses. ; ; This code is somewhat complicated because we are checking to make ; sure we don't override any line thicknesses that might have ; originally been set by the user. ; ; This code is very similar to the set_contour_line_patterns code. ;---------------------------------------------------------------------- undef("set_contour_line_thicknesses") procedure set_contour_line_thicknesses(plot,levels_to_change[*]:numeric,\ thicknesses[*]:numeric) local ii, levels, nlevels, line_thicknesses, nchange, nthick, \ mono_thickness, line_thickness,changed_a_level begin ;---Retrieve the original line thicknesses (or thickness) used for the plot. getvalues plot@contour "cnLevels" : levels "cnLineThicknesses" : line_thicknesses "cnMonoLineThickness" : mono_thickness "cnLineThicknessF" : line_thickness end getvalues ;---------------------------------------------------------------------- ; If the original contour plot used a single dash line thickness ; all contour lines, then make sure we use that again for ; any contour lines that are not being changed. Otherwise, ; we assume the user set an array of line thicknesses, and ; we'll use those. ;---------------------------------------------------------------------- ;---Check user set line thicknesses, and apply new thicknesses if needed if(mono_thickness) then line_thicknesses = line_thickness end if nlevels = dimsizes(levels) nchange = dimsizes(levels_to_change) nthick = dimsizes(thicknesses) if(nthick.ne.nchange) then print("set_contour_line_thicknesses: error: the contour line thicknesses must be an") print("array of the same length as the number of contour levels to change.") end if ;---Array to hold index values of contour levels that need to have a thickness applied changed_a_level = False do n=0,nchange-1 ii := ind(levels.eq.levels_to_change(n)) if(ismissing(ii(0))) then print("set_contour_line_thicknesses: warning, no contour level equal to " + levels_to_change(n)) else line_thicknesses(ii(0)) = thicknesses(n) changed_a_level = True end if end do if(changed_a_level) then setvalues plot@contour "cnMonoLineThickness" : False ; allows an array of line thicknesses to be set "cnLineThicknesses" : line_thicknesses end setvalues end if end undef("set_res1_list") function set_res1_list() begin res = True res@mpDataBaseVersion = "Ncarg4_1" ; choose more recent database res@mpDataSetName = "Earth..4" ; high resolution res@mpOutlineBoundarySets = "GeophysicalAndUSStates" ; National borders res@mpNationalLineColor = "Black" ; National borders color res@mpNationalLineThicknessF =16.0 res@mpGeophysicalLineThicknessF =12.0 res@mpUSStateLineThicknessF =12.0 ; res@gsnDraw = False ; don't draw ; res@gsnFrame = False ; don't advance frame res@cnFillOn = True ; turn on color res@cnLinesOn = False ; no contour lines ; res@gsnLeftString = "500 mb Wind Vectors" ; change left string res@gsnLeftString = "Cloud Top Temps" res@gsnRightString = "[~S~o~N~C]" ; assign right string res@mpFillOn = False ; no map fill res@cnLevelSelectionMode = "ExplicitLevels" ; set explicit contour levels res@cnLevelSelectionMode = "ExplicitLevels" ; use explicit levels res@cnLevels = ispan(-80,20,2) ; set the contour levels cmap = read_colormap_file("NCV_rainbow2") cmap := cmap(50:, :) cmap(0,:) = 1. ; Make first color white res@cnFillPalette = cmap ; set color map ; res@cnLevels = (/ -40.,-38.,-36.,-34.,-32.,-30.,-28.,\ ; -26.,-24.,-22.,-20.,-18.,-16.,-14.,\ ; -12.,-10.,-8.,-6.,-4.,-2.,0.,\ ; 2.,4.,6.,8.,10.,12.,14.,16.,18.,\ ;res@cnFillColors = (/ 10,12,14,15,17,19,22,25,31,35,39,43,48,53,\ ; 61,64,67,71,89,95,99,104,108,111,115,123,\ ; 130,135,143,148,150,155,159,163,165/) ; ; set the colors to be used ; res1 label bar res@cnInfoLabelOn = False ; turn off contour info label res@lbAutoManage = False ; control label bar res@pmLabelBarDisplayMode = "Always" ; turns on label bar res@lbOrientation = "Vertical" res@pmLabelBarSide = "Right" res@lbLabelAutoStride = True res@pmLabelBarWidthF = 0.20 res@pmLabelBarHeightF = 0.7 res@lbLabelFontHeightF = .025 res@lbPerimOn = False res@lbLabelFont = "Helvetica" ; label font res@lbTitleOn = True ; turn on title ; res@lbTitleString = "Heights [~S~o~N~C]" ; title string res@lbTitleString = "[~S~o~N~C]" ; title string res@lbTitlePosition = "Top" ; title position res@lbTitleFontHeightF= .020 ; make title smaller res@lbTitleDirection = "Across" ; title direction res@gsnAddCyclic = True ; data already has cyclic point return(res) end undef("set_res3_list") function set_res3_list() begin res = True res@gsnAddCyclic = False ; we will manually zoom in. res@gsnDraw = False ; do not draw the plot res@gsnFrame = False ; do not advance the frame res@gsnDraw = False ; do not draw the plot res@gsnFrame = False ; do not advance the frame res@cnConstFLabelPerimOn = False res@cnLevelSelectionMode = "ManualLevels" res@cnMinLevelValF = 4800 res@cnMaxLevelValF = 6200 res@cnLevelSpacingF = 100 res@cnLineLabelAngleF = 0.0 res@cnLineLabelFontHeightF = 0.025 res@cnLineLabelBackgroundColor = -1 res@gsnContourLineThicknessesScale = 8.0 res@tfDoNDCOverlay = True res@cnInfoLabelOn = False res@cnLineLabelDensityF = 1.5 ; increase the number of line labels/line res@cnLineLabelInterval = 1 ; labels for every line (default=2) res@gsnCenterString = "Increase Label Density + Every cn Line" return(res) end begin plots = new(6,graphic) in_file = "cld_top_temps_12utc.nc" in_file_18 = "cld_top_temps_18utc.nc" in_file_h_12 = "../heights_nomp_new_36h.nc" in_file_h_18 = "../heights_nomp_new_42h.nc" ;in_file = "cld_top_temps_c_7.nc" ;in_file = "cld_top_temps_c_7.nc" ;in_file = "heights_b.nc" ;in_file = "heights_cu1.nc" ;in_file = "heights_b_nomp.nc" ;in_file = "heights_7.nc" a = addfile(in_file,"r") a_18 = addfile(in_file_18,"r") a_h_12 = addfile(in_file_h_12,"r") a_h_18 = addfile(in_file_h_18,"r") xi = a->xi xi = xi yi = a->yi xo = a->xo x_orig = xo yo = a->yo fi = a->fi z_sat = a->z_nam xi_sat = a->xi_nam yi_sat = a->yi_nam xlat_sat = a->xlat_nam xlon_sat = a->xlon_nam ; fi_oth = a->fi_oth xo_min = xo xo_min = 180 xo = xo-xo_min xlat2d= a->xlat2d xlon2d= a->xlon2d z_gefs = a->z_gefs z_gefs_ave = a->z_gefs_ave z_sat_18 = a_18->z_nam fi_18 = a_18->fi z_gefs_18 = a_18->z_gefs_ave z_h_12 = a_h_12->z_gefs_ave z_h_18 = a_h_18->z_gefs_ave fi_h_12 = a_h_12->fi fi_h_18 = a_h_18->fi n_h_12 = a_h_12->z_nam n_h_18 = a_h_18->z_nam xlat_nam = a_h_18->xlat_nam xlon_nam = a_h_18->xlon_nam printMinMax(fi,False) ; print(xi) ; print("xi = " + xi) ; print("yi = " + yi) ; print("xo = " + xo) ; print("yo = " + yo) xoo = xo(180:250) yoo = yo(76:130:-1) ; print("xoo = " + xoo) ; print("yoo = " + yoo) outfile="cld_top_temps_maps" wtype = "png" ; wtype = "pdf" wtype@wkWidth = 7500 ; produces a better looking PNG wtype@wkHeight = 7500 wks = gsn_open_wks(wtype,outfile) gsn_define_colormap(wks,"BkBlAqGrYeOrReViWh200") res1 = set_res1_list() res1@tiMainString = "GEFS Mean (12 UTC)" res_slp = set_res3_list() res_slp@cnLowLabelsOn = True res_slp@cnHighLabelsOn = True ; res_slp@cnLowLabelFontHeightF = 0.15 res_slp@cnHighLabelFontHeightF = 0.05 res_slp@cnLowLabelBackgroundColor = -1 res_slp@cnHighLabelBackgroundColor = -1 res_slp@cnLowLabelFontColor = "White" ; res_slp@cnLowLabelFontColor = "Black" ;; Map Projection (min and max Lon and Lat) lat_beg = 27 lat_end = 47 lon_beg = -82 + 0 lon_end = -58 + 00 res1@mpMinLatF = lat_beg res1@mpMaxLatF = lat_end res1@mpMinLonF = lon_beg res1@mpMaxLonF = lon_end res1@mpLimitMode="LatLon" ;fi@long_sate = "height" ;fi@units = "(m)" ; yi@lat_sate = "observed latitudes" ; yi@units = "degrees_north" ; xi@long_sate = "observed longitudes" ; xi@units = "degrees_east" printVarSummary(xoo) printVarSummary(yoo) printVarSummary(fi) lon = xoo lat = yoo dims2d = dimsizes(z_gefs_ave) printVarSummary(xo) printVarSummary(yo) print(dims2d) dims_sat = dimsizes(z_sat) printVarSummary(z_sat) ; lat2d = conform_dims(dims2d,yo,0) ; lon2d = conform_dims(dims2d,xo,1) xlat2d_sat = conform_dims(dims_sat,xlat_sat,0) xlon2d_sat = conform_dims(dims_sat,xlon_sat,1) dims_nam = dimsizes(n_h_18) print("dims_nam = " + dims_nam) printVarSummary(n_h_12) ; lat2d = conform_dims(dims2d,yo,0) ; lon2d = conform_dims(dims2d,xo,1) printVarSummary(xlat_nam) printVarSummary(xlon_nam) ; xlat2d_nam = conform_dims(dims_nam,xlat_nam,0) ; xlon2d_nam = conform_dims(dims_nam,xlon_nam,1) ; printVarSummary(xlat2d_nam) ; printVarSummary(xlon2d_nam) copy_VarCoords(z_gefs,z_gefs_ave) copy_VarAtts(z_gefs,z_gefs_ave) dims2d = dimsizes(z_h_12) xo_o = where(x_orig.gt.180,x_orig-360,x_orig) ; lat2d = conform_dims(dims2d,yo,0) ; lon2d = conform_dims(dims2d,xo,1) lat2d_o = conform_dims(dims2d,yo,0) lon2d_o = conform_dims(dims2d,xo_o,1) h_gefs_12 = rcm2rgrid_Wrap(lat2d_o,lon2d_o,z_h_12,lat,lon,0) h_gefs_18 = rcm2rgrid_Wrap(lat2d_o,lon2d_o,z_h_18,lat,lon,0) z_fi_12 = rcm2rgrid_Wrap(xlat2d,xlon2d,fi_h_12,lat,lon,0) z_fi_18 = rcm2rgrid_Wrap(xlat2d,xlon2d,fi_h_18,lat,lon,0) printVarSummary(n_h_18) printVarSummary(xlat2d_sat) printVarSummary(xlon2d_sat) nam_h_18 = rcm2rgrid_Wrap(xlat_nam,xlon_nam,n_h_18,lat,lon,0) nam_h_12 = rcm2rgrid_Wrap(xlat_nam,xlon_nam,n_h_12,lat,lon,0) z_wrf_ave_o = rcm2rgrid_Wrap(xlat2d,xlon2d,fi,lat,lon,0) z_wrf_18 = rcm2rgrid_Wrap(xlat2d,xlon2d,fi_18,lat,lon,0) z_sat_o = rcm2rgrid_Wrap(xlat2d_sat,xlon2d_sat,z_sat,lat,lon,0) z_sat_18_o = rcm2rgrid_Wrap(xlat2d_sat,xlon2d_sat,z_sat_18,lat,lon,0) ;printVarSummary(z_wrf_ave_o) ;printVarSummary(xlat2d) ;copy_VarCoords(xlat2d,fi) ; z_wrf_ave_o!0 = "yoo" ; z_wrf_ave_o!1 = "xoo" ; z_wrf_ave_o&fi@units = "degrees_north" ; z_wrf_ave_o&fi@units = "degrees_east" printVarSummary(z_wrf_ave_o) ; z_wrf_ave_o@long_sate = "observations of some sort" ; z_wrf_ave_o@units = "not sure" ; z_wrf_ave_o!0 = "lat" ; name the dimension ; z_wrf_ave_o&lat = xi ; assign the coordinate variable ; z_wrf_ave_o!1 = "lon" ; name the dimension ; z_wrf_ave_o&lat = yi ; assign the coordinate variable ; lat_1@long_sate = "observed latitudes" ; lat_1@units = "degrees_north" ; lon_1@long_sate = "observed longitudes" ;printVarSummary(z_wrf_ave_o) printMinMax(z_wrf_ave_o,False) dims = dimsizes(z_wrf_ave_o) print(dims) do j= 0,dims(0)-1 do i= 0,dims(1)-1 check = ismissing(z_wrf_ave_o(j,i)) if (check.eq.False)then ;print("z, j,i = " + xoo(i) + " " + yoo(j)+" "+z_wrf_ave_o(j,i)) end if end do end do res1b = res_slp res1b@tiMainString = "GEFS Mean (18 UTC)" plot_1a = gsn_csm_contour_map(wks,z_gefs_ave,res1) res1b@gsnRightString = "" res1b@gsnLeftString = "" plot_1b = gsn_csm_contour(wks,h_gefs_12,res1b) set_contour_line_thicknesses(plot_1b,(/5500/),(/15/)) draw(plot1b) plots(2)= plot_1a overlay(plots(2),plot_1b) plot_4 = gsn_csm_contour_map(wks,z_gefs_18,res1) plot_1b = gsn_csm_contour(wks,h_gefs_18,res1b) set_contour_line_thicknesses(plot_1b,(/5500/),(/15/)) draw(plot1b) plots(3)= plot_4 overlay(plots(3),plot_1b) res2 = res1 ;res2@tiMainString = "WRF Mean (M-T)" res2@tiMainString = "WRF Mean (12 UTC)" res2@gsnAddCyclic = False res2b = res_slp res2b@tiMainString = "WRF Mean (18 UTC)" res2b@gsnRightString = "" res2b@gsnLeftString = "" plot_2 = gsn_csm_contour_map(wks,z_wrf_ave_o,res2) plots(4)= plot_2 plot_2b= gsn_csm_contour(wks,z_fi_12,res2b) overlay(plots(4),plot_2b) set_contour_line_thicknesses(plot_2b,(/5500/),(/15/)) draw(plot2b) plot_5 = gsn_csm_contour_map(wks,z_wrf_18,res2) plot_2b= gsn_csm_contour(wks,z_fi_18,res2b) set_contour_line_thicknesses(plot_2b,(/5500/),(/15/)) draw(plot2b) plots(5) = plot_5 overlay(plots(5),plot_2b) res3 = res1 ;res3@tiMainString = "WRF Mean (K-F)" ;res3@tiMainString = "WRF No Micro Mean" res3@tiMainString = "WRF No Micro Mean" res3@gsnAddCyclic = False ;plot_3 = gsn_csm_contour_map(wks,z_wrf_oth_o,res3) ;plots(3)= plot_3 res4 = res1 res4b = res_slp res4b@gsnRightString = "" res4b@gsnLeftString = "" res4@tiMainString = "Observed (12 UTC)" res4@gsnAddCyclic = False res4b = res_slp res4b@tiMainString = "Observed (18 UTC)" plot_4 = gsn_csm_contour_map(wks,z_sat_o,res4) plot_4b= gsn_csm_contour(wks,nam_h_12,res4b) set_contour_line_thicknesses(plot_4b,(/5500/),(/15/)) draw(plot4b) plots(0)= plot_4 overlay(plots(0),plot_4b) plot_8 = gsn_csm_contour_map(wks,z_sat_18_o,res4) plots(1)= plot_8 plot_4b= gsn_csm_contour(wks,nam_h_18,res4b) set_contour_line_thicknesses(plot_4b,(/5500/),(/15/)) draw(plot4b) overlay(plots(1),plot_4b) resP = True ; panel only resources resP@gsnMaximize = True ; maximize plots gsn_panel(wks,plots,(/3,2/),resP) end