<div dir="ltr"><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">Hi<br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">I want to plot surface station data on top of a WRF contour map. I want the polymarkers to be circled with a black border - there's no polymarker that can do this at the moment. So, I'm trying to plot gsmarker=4 type hollow markers on top of solid coloured polymarkers (gsmarker=16) but unable to achieve the desired result.<br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">Here's my script:<br><span style="font-family:comic sans ms,sans-serif"><i><font size="1">load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"<br>load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"<br>load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"<br><br>;----------------------------------------------------------------------<br>; Procedure to attach colored markers to an NCL map, given 1D<br>; data, and lat/lon arrays, a set of levels, and a color map.<br>;----------------------------------------------------------------------<br>procedure add_obs_markers_to_map(wks,plot,levels,colormap,data_1d,lat_1d,lon_1d)<br>local mkres, nlevels, colors, nlevels, n, ii<br>begin<br><br>;---Set resources for the markers<br>  mkres = True<br>  mkres@gsMarkerIndex = 16<br>;<br>; Based on the levels we had for the contour plot and the colormap used,<br>; get an array of colors that span the colormap. This uses the same<br>; algorithm that the contour used, so the colors will be identical.<br>;<br>  nlevels = dimsizes(levels)<br>  colors = span_color_rgba(colormap,nlevels+1)   ; Need one more color than number of levels<br><br>;<br>; Loop through each level, gather all the data values in the given level<br>; range, and add colored markers to the existing map.<br>;<br>  do n=0,nlevels<br><br>;---These "if" statements are how color contours are handled in NCL.<br>    if(n.eq.0) then<br>      ii := ind(data_1d.lt.levels(n))<br>    else if(n.eq.nlevels) then<br>      ii := ind(data_1d.ge.levels(n-1))<br>    else<br>      ii := ind(data_1d.ge.levels(n-1).and.data_1d.lt.levels(n))<br>    end if<br>    end if<br>    if(any(ismissing(ii))) then<br>      continue<br>    end if<br><br>;---Add the markers<br>    mkres@gsMarkerColor = colors(n,:)    ; colors is an N x 4 array<br>    plot@$unique_string("markers")$ = gsn_add_polymarker(wks,plot,lon_1d(ii),lat_1d(ii),mkres)<br>  end do<br>end<br><br>;----------------------------------------------------------------------<br>; Main code<br>;----------------------------------------------------------------------<br>begin<br><br>;---Open WRF output file and read data<br><br>  DATADir = "/data2/tabish/control-run-so4-ECMWF/"<br>FILES = systemfunc (" ls -1 " + DATADir + "wrfout_d03_2014-10* ")<br>numFILES = dimsizes(FILES)<br><br>s=1<br>a = addfile(FILES(0)+".nc","r"); THIS IS JUST FOR INITIALIZATION OF PM25 ARRAY, IT WILL BE SUBSTRACTED FROM THE SUM IN THE END<br>pm25init =  a->PM2_5_DRY(15,0,:,:)<br>pm25sum =  pm25init<br>times = wrf_user_getvar(a,"times",-1)  ; get all times in the file<br>ntimes = dimsizes(times)<br>do it = 16,ntimes-1,1             ; TIME LOOP<br>    s=s+1<br>    print("Working on time: " + times(it) )<br>    pm25sum = pm25sum + a->PM2_5_DRY(it,0,:,:)<br>end do ; END OF TIME LOOP<br><br>do ifil = 1,numFILES-2             ; FILE LOOP<br>  a = addfile(FILES(ifil)+".nc","r")   ; Open the next file<br>  times = wrf_user_getvar(a,"times",-1)  ; get all times in the file<br>  ntimes = dimsizes(times)<br>  do it = 0,ntimes-1,1             ; TIME LOOP<br>    s=s+1<br>    print("Working on time: " + times(it) )<br>    pm25sum = pm25sum + a->PM2_5_DRY(it,0,:,:)<br>  end do ; END OF TIME LOOP<br>end do<br><br>a = addfile(FILES(numFILES-1)+".nc","r")   ; Open the next file<br>  times = wrf_user_getvar(a,"times",-1)  ; get all times in the file<br>  ntimes = dimsizes(times)<br>  do it = 0,16,1             ; TIME LOOP<br>    s=s+1<br>    print("Working on time: " + times(it) )<br>    pm25sum = pm25sum + a->PM2_5_DRY(it,0,:,:)<br>  end do ; END OF TIME LOOP<br><br>pm25_avg= pm25sum/s<br>;pm25_avg@description = "NO2 conc in ppbv"<br><br>  lat2d = wrf_user_getvar(a,"XLAT",it)   ; latitude/longitude<br>  lon2d = wrf_user_getvar(a,"XLONG",it)<br><br>  pm25_avg@lat2d = lat2d    ; for plotting<br>  pm25_avg@lon2d = lon2d<br><br>;---Will use this for contours and filled markers<br>  colormap = "BlAqGrYeOrReVi200"<br>  levels   = ispan(0,150,15)<br>;  levels   = (/0,10,20,30,40,50,60,70,80,90,100,125,150,175,200,225,250/) * 0.01<br><br>  wks = gsn_open_wks("x11","wrf_obs_d03_12-31octmean_pm25")<br><br>;---Common resources shared by both contour plot and marker plot<br>  res                        = True<br>  res@gsnMaximize            = False<br>  res@gsnLeftString          = ""<br>  res@gsnRightString         = ""<br><br>  res@mpDataBaseVersion      = "MediumRes"<br> ; res@mpFillOn               = False<br> ; res@mpMinLatF              = min(pm25_avg@lat2d)<br> ; res@mpMaxLatF              = max(pm25_avg@lat2d)<br> ; res@mpMinLonF              = min(pm25_avg@lon2d)<br> ; res@mpMaxLonF              = max(pm25_avg@lon2d)<br> ; res@pmTickMarkDisplayMode  = "Always"      ; better map tickmarks<br><br> res@mpLimitMode       = "Corners"            ; choose range of map<br> res@mpLeftCornerLatF  = 36.05806<br> res@mpLeftCornerLonF  = 113.3356<br> res@mpRightCornerLatF = 42.4514<br> res@mpRightCornerLonF = 120.2583<br> res@mpProjection         = "LambertConformal"<br> res@mpLambertParallel1F = 20<br> res@mpLambertParallel2F = 50<br> res@mpLambertMeridianF  = 110<br><br>;---Resources for filled contour plot<br>  cnres                      = res<br>  cnres@cnFillOn             = True<br>  cnres@cnLevelSelectionMode = "ExplicitLevels"<br>  cnres@cnLevels             = levels<br>  cnres@cnFillPalette        = colormap<br>  cnres@cnLinesOn            = False<br>  cnres@cnLineLabelsOn       = False<br>;  cnres@lbLabelBarOn         = False   ; will add in panel later<br>  cnres@pmTitleZone          = 2       ; move title down<br>  cnres@gsnAddCyclic         = False<br><br>  cnres@mpLimitMode       = "Corners"            ; choose range of map<br>  cnres@mpLeftCornerLatF  = 36.05806<br>  cnres@mpLeftCornerLonF  = 113.3356<br>  cnres@mpRightCornerLatF = 42.4514<br>  cnres@mpRightCornerLonF = 120.2583<br>  cnres@mpProjection         = "LambertConformal"<br>  cnres@mpLambertParallel1F = 20<br>  cnres@mpLambertParallel2F = 50<br>  cnres@mpLambertMeridianF  = 110<br>  cnres@mpFillOn                    =  True<br>  cnres@mpOutlineDrawOrder          = "PostDraw"<br>  cnres@mpFillDrawOrder             = "PreDraw"<br> ; cnres@mpOutlineBoundarySets       = "National"<br>  cnres@mpOutlineBoundarySets       = "Allboundaries"<br>  cnresmpNationalLineColor =    "Black"<br>  cnresmpGeophysicalLineColor = "Black"<br>  cnres@mpUSStateLineDashPattern    = 2<br>  cnres@mpOutlineOn           = True<br>  cnres@mpDataBaseVersion        = "MediumRes"<br>  cnres@mpDataSetName            = "Earth..4"      ; U.S. counties<br>  cnres@mpGridAndLimbOn = True<br>  cnres@mpGridLineColor = "Black"<br>  cnres@mpFillColors = (/"transparent","transparent","gray","transparent"/)<br>  cnres@pmTickMarkDisplayMode = "Always"<br>  cnres@tiMainString = "WRF model data"<br>  cnres@tiDeltaF = 3.3<br>  cnres@tfDoNDCOverlay = True<br><br>  cnres@gsnDraw  = False                          ; don't draw<br>  cnres@gsnFrame = False                          ; don't advance frame<br>print(res)<br>;---Create smooth contour plot<br>  plot_wrf = gsn_csm_contour_map(wks,pm25_avg,cnres)<br><br>;----------------------------------------------------------------------<br>; Recreate similar plot, but by using filled markers drawn over a map.<br>;----------------------------------------------------------------------<br>  fname = "surfacedata_12-31_Octmean_improved.txt"<br>  lines = asciiread(fname,-1,"string")<br>  lon = tofloat(str_get_field(lines(:),2," "))<br>  lat = tofloat(str_get_field(lines(:),3," "))<br>  pm25 = tofloat(str_get_field(lines(:),4," "))<br><br>;---Create a map plot for adding markers to.<br>  res@gsnDraw         = False<br>  res@gsnFrame        = False<br>  res@tfPolyDrawOrder = "Draw"    ; Necessary to make sure filled markers drawn *under* map outlines<br>  res@tiMainString    = "'Observational' data"<br><br>  res@mpFillOn                    =  True<br>  res@mpOutlineDrawOrder          = "PostDraw"<br>  res@mpFillDrawOrder             = "PreDraw"<br> ; res@mpOutlineBoundarySets       = "National"<br>  res@mpOutlineBoundarySets       = "Allboundaries"<br>  res@mpUSStateLineDashPattern    = 2<br>  res@mpOutlineOn           = True<br>  res@mpDataBaseVersion        = "MediumRes"<br>  res@mpDataSetName            = "Earth..4"      ; U.S. counties<br> ; res@mpGridAndLimbOn = True<br>  res@mpFillColors = (/"transparent","transparent","gray","transparent"/)<br>  res@pmTickMarkDisplayMode = "Always"<br>  res@tiMainString = "Surface data"<br><br>  ;plot_obs = gsn_csm_map(wks,res)<br><br>;---Add colored markers to the map<br>  add_obs_markers_to_map(wks,plot_wrf,levels,colormap,pm25,lat,lon)<br><br>;*****************************************<br>; add black circles over obs data<br>;*****************************************<br>  polyres               = True          ; poly marker mods desired<br>  polyres@gsMarkerIndex = 4            ; choose circle as polymarker<br>  polyres@gsMarkerSizeF = 10.0           ; select size to avoid streaking<br>  polyres@gsMarkerColor = (/"black"/)   ; choose color<br>  gsn_polymarker(wks,plot_wrf,lat,lon,polyres)  ; draw polymarkers<br><br><br><br>;---This draws the map with the filled markers<br>  draw(plot_wrf)<br>  frame(wks)<br><br>end</font></i></span><br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif"><br><br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">Please let me know if that's possible.<br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">Thanks<br><br></div><div class="gmail_default" style="font-family:trebuchet ms,sans-serif">Tabish<br><br clear="all"></div><div><div class="gmail_signature"><div dir="ltr"><div><div><font size="1"><span style="font-family:tahoma,sans-serif">Tabish U Ansari<br></span></font></div><font size="1"><span style="font-family:tahoma,sans-serif">PhD student, Lancaster Environment Center<br></span></font></div><font size="1"><span style="font-family:tahoma,sans-serif">Lancaster Univeristy<br> <span>Bailrigg</span>, <span>Lancaster</span>, <br><span>LA1 4YW</span>, <span>United Kingdom</span></span></font><br></div></div></div>
</div>