DKRZ NCL example use markers to plot data on grid#

Example:

;-------------------------------------------------------------------------------------
; DKRZ NCL Example:  quality_points_per_grid_cell.ncl
;
; Description:   Draw two plots, first plot is a raster contour plot
;       and the second shows the data using filled circles
;       which are sized by a quality variable.
;
; First plot:    - draw contour plot
;
; Second plot:    - draw contour plot
;      - clear the contour plot
;      - determine marker colors
;      - determine marker size (depending on quality variable)
;      - draw markers (filled circles)
;      - use dummy data
; 30.11.18 kmf
;-------------------------------------------------------------------------------------
begin
  VALUES = True ;-- turn on/off value annotation of first plot
  GRID = False ;-- turn on/off the data grid lines of second plot

  ;-----------------------------------------------------------------------------------
  ;-- which country to be shown
  ;-----------------------------------------------------------------------------------
  mask_specs = (/"Germany"/)
  minlat = 47.0 ;-- minimum latitude of map
  maxlat = 55.0 ;-- maximum latitude of map
  minlon = 5.0 ;-- minimum longitude of map
  maxlon = 16.0 ;-- maximum longitude of map

  ;-----------------------------------------------------------------------------------
  ;-- generate dummy data and coordinates
  ;-- generate lat and lon arrays with named dimensions
  ;-----------------------------------------------------------------------------------
  nlat = 16
  nlon = 22
  lat = fspan(minlat, maxlat,nlat)
  lon = fspan(minlon, maxlon,nlon)
  lat@units = "degrees_north"
  lon@units = "degrees_east"

  ;-- generate dummy data with named dimensions
  tempmin = -2.0
  tempmax = 2.0
  tempint = 0.5
  temp = random_uniform(tempmin,tempmax,(/nlat,nlon/))
  temp!0 = "lat"
  temp!1 = "lon"
  temp&lat = lat
  temp&lon = lon
  temp1d = ndtooned(temp)
  ncells = dimsizes(temp1d)

  ;-- generate random dummy quality data
  minqual = 1
  maxqual = 4

  quality = floor(random_uniform(minqual,maxqual+0.5,(/nlat,nlon/)))
  quality!0 = "lat"
  quality!1 = "lon"
  quality&lat = lat
  quality&lon = lon
  quality1d = ndtooned(quality)

  ;-----------------------------------------------------------------------------------
  ;-- graphics
  ;-----------------------------------------------------------------------------------
  ;-- choose color map
  cmap = "BlueDarkRed18"

  ;-- open workstation
  wks = gsn_open_wks("png","plot_quality_per_cell")

  ;-----------------------------------------------------------------------------------
  ;-- contour plot resources
  ;-----------------------------------------------------------------------------------
  res = True
  res@gsnDraw = False
  res@gsnFrame = False
  res@gsnAddCyclic = False
  res@gsnCenterString = "Quality: min = "+min(quality)+" max = "+max(quality)
  res@gsnCenterStringFontHeightF = 0.008

  res@vpXF = 0.09 ;-- viewport x-position
  res@vpYF = 0.98 ;-- viewport y-position
  res@vpWidthF = 0.77 ;-- viewport width
  res@vpHeightF = 0.8 ;-- viewport height

  res@cnFillMode = "RasterFill" ;-- use raster fill for contours
  res@cnFillOn = True ;-- filled contour areas
  res@cnLinesOn = False
  res@cnInfoLabelOn = False
  res@cnFillPalette = cmap ;-- set color map
  res@cnLevelSelectionMode = "ManualLevels" ;-- set manual data levels
  res@cnMinLevelValF = tempmin
  res@cnMaxLevelValF = tempmax
  res@cnLevelSpacingF = tempint

  res@lbBoxMinorExtentF = 0.15 ;-- decrease height of labelbar

  res@mpOutlineBoundarySets = "National" ;-- draw national map outlines
  res@mpOceanFillColor = "gray90"
  res@mpLandFillColor = "gray75"
  res@mpInlandWaterFillColor = "gray75"
  res@mpDataBaseVersion = "MediumRes" ;-- alias to Ncarg4_1
  res@mpDataSetName = "Earth..4" ;-- higher map resolution
  res@mpMinLatF = minlat-1.0
  res@mpMaxLatF = maxlat+1.0
  res@mpMinLonF = minlon-1.0
  res@mpMaxLonF = maxlon+1.0

  ;-----------------------------------------------------------------------------------
  ;-- first, we have to draw a contour map in memory to get the level and color values
  ;-----------------------------------------------------------------------------------
  contour = gsn_csm_contour_map(wks,temp,res)

  ;-- for checking purposes: draw text and markers at data locations !! this takes long !!
  if(VALUES) then
     txres = True
     txres@gsFontColor = "black"
     txres@txFontHeightF = 0.01
     do j=0,nlat-1
        do i=0,nlon-1
           txid = unique_string("txid")
           contour@$txid$ = gsn_add_text(wks, contour,""+quality(j,i),lon(i),lat(j),txres)
        end do
     end do
  end if

  ;-- draw the first plot
  draw(contour)
  frame(wks)

  ;-----------------------------------------------------------------------------------
  ;-- create second plot
  ;-----------------------------------------------------------------------------------
  plot = gsn_csm_contour_map(wks,temp,res)

  ;-- retrieve used levels and colors
  getvalues plot@contour
    "cnFillColors" : colors ;-- get the colors used
    "cnLevels" : levels ;-- retrieve levels
  end getvalues

  ncolors = dimsizes(colors) ;-- number of used colors
  nlevels = dimsizes(levels) ;-- number of used levels

  print("ncells: "+ncells)
  print("ncolors: "+ncolors)
  print("nlevels: "+nlevels)

  ;-- clean the plot that we can see the markers instead of filled cells
  plot = setColorContourClear(plot,min(temp)-1.,max(temp)+1.) ;-- clear plot, but keep all the information

  ;-- draw data coordinate grid lines just to see the lon/lat lines at cell center
  if(GRID) then
     pres = True
     pres@gsnDraw = False
     pres@gsnFrame = False
     pres@gsnCoordsAsLines = True ;-- default: points
     pres@gsnCoordsAttach = True ;-- attach the coordinate grid lines
     gsn_coordinates(wks,plot,temp,pres) ;-- draw grid lines
  end if

  ;-- now, create the marker size for each cell; marker sizes for quality 1-4
  marker_sizes = fspan(0.01,0.04,4)
  ms_array = new(ncells,float) ;-- create array for marker sizes depending
  ;-- on quality1d

  do ll=minqual,maxqual
     indsize := ind(quality1d .eq. ll)
     if(any(.not. ismissing(indsize))) then
        ms_array(indsize) = marker_sizes(ll-1)
     end if
  end do

  ;-- marker resources
  plmres = True
  plmres@gsMarkerIndex = 16 ;-- filled circles

  ;-- now, create the color array for each cell from temp1d
  gscolors = new(ncells,integer)

  ;-- set color for data less than given minimum value
  vlow = ind(temp1d .lt. levels(0)) ;-- get the indices of values less levels(0)
  if(any(.not. ismissing(vlow))) then
     gscolors(vlow) = colors(0) ;-- choose color
  end if

  ;-- set colors for all cells in between tempmin and tempmax
  do i=1,nlevels-1
     vind := ind(temp1d .ge. levels(i-1) .and. temp1d .lt. levels(i)) ;-- get the indices of 'middle' values
     if(any(.not. ismissing(vind))) then
        gscolors(vind) = colors(i) ;-- choose the colors
     end if
  end do

  ;-- set color for data greater than given maximum
  vhgh := ind(temp1d .gt. levels(nlevels-1)) ;-- get indices of values greater levels(nl)
  if(any(.not. ismissing(vhgh))) then
     gscolors(vhgh) = colors(ncolors) ;-- choose color
  end if

  ;-- add the marker to the plot
  n=0
  do ii=0,nlat-1
     do jj=0,nlon-1
        k = jj+ii
        plmres@gsMarkerSizeF = ms_array(n) ;-- marker size
        plmres@gsMarkerColor = gscolors(n) ;-- marker color
        str = unique_string("plmark")
        plot@$str$ = gsn_add_polymarker(wks,plot,lon(jj),lat(ii),plmres) ;-- add marker to plot
        n = n + 1
     end do
  end do

  ;-----------------------------------------------------------------------------------
  ;-- add marker size legend
  ;-----------------------------------------------------------------------------------
  legres = True ;-- legend resources
  legres@gsMarkerIndex = 16 ;-- filled dots
  legres@gsMarkerColor = "gray50" ;-- legend marker color

  txres = True ;-- text resources
  txres@gsFontColor = "black"
  txres@txFontHeightF = 0.01
  txres@txFont = 30

  x = 0.93
  y = 0.47
  ik = 0

  do il=minqual,maxqual
     legres@gsMarkerSizeF = marker_sizes(ik)
     gsn_polymarker_ndc(wks,x,y,legres)
     gsn_text_ndc(wks,""+il,x+0.03,y,txres)
     y = y + 0.05
     ik = ik + 1
  end do

  txres@txFontHeightF = 0.012
  gsn_text_ndc(wks,"Quality",x,y,txres) ;-- legend title

  ;-- add title string
  getvalues plot
    "vpWidthF" : vpw
    "vpXF" : vpx
  end getvalues
  xpos = vpw/2 + vpx

  title1 = "Draw data values at lat/lon location as circles"
  title2 = "the size is defined by the quality variable"

  txres@txFont = 21
  txres@txFontHeightF = 0.018
  gsn_text_ndc(wks,title1,xpos,0.96,txres)
  gsn_text_ndc(wks,title2,xpos,0.935,txres)

  ;-----------------------------------------------------------------------------------
  ;-- draw the second plot and advance the frame
  ;-----------------------------------------------------------------------------------
  draw(plot)
  frame(wks)
end

Result:

../../../../../../_images/plot_quality_per_cell.000001_w400.png ../../../../../../_images/plot_quality_per_cell.000002_w400.png