[ncl-talk] Fwd: Help with creating a script for precipitation stored in shapefiles

[Carlos J. Valle-Diaz]

Forgot to Cc ncl-talk...I apologize!

---------- Forwarded message ----------
From: Carlos J. Valle-Diaz <cj.vallediaz at gmail.com>
Date: Wed, Oct 21, 2015 at 8:27 PM
Subject: Re: [ncl-talk] Help with creating a script for precipitation
stored in shapefiles
To: Rick Brownrigg <brownrig at ucar.edu>

Thank you very much Rick,

I've been playing around with the script to understand the functions and
resources used to make the plot and learn NCL from it. I managed to focus
it on my area of interest (Puerto Rico). At the moment, I've been trying to
put a legend and text. I've seen an example at the website "polyg_8.ncl",
but I have a question.

I know that the gsn_polymarker_ndc(wks, x, y, res) function lets me put the
legend markers and the gsn_text_ndc(wks, text, x, y, res) function lets me
put the text. I'm using "gsres" as a resource for the polymarkers and
"txres" as a resource for the text to go along with it. My question lies in
the numeric values (x and y). In the example script I mentioned above, the
code is set for an array set defined. In my case I would have to define it
for the Globalvalues (precipitation) in my shapefile. In summary, how do I
define x and y to make a legend at an even interval of my GlobalValues and
to play with the location of both the markers and text? (See script
attached)

Also, I've been trying to include black circles in the loop that generates
the colored polymarkers in the plot to make it more visible, but haven't
been succesful. I've tried the following sequence, for example, and marks
an error.  What am I doing wrong?

   do i=0, dimsizes(f->Globvalue) - 1

      indx = doubletoint( (f->Globvalue(i)-globMin) / globRange *
numColors) + 1  ; skip first color in table
      if (indx.ge.numColors) then
          indx = numColors - 1
      end if

       gsres at gsMarkerIndex = 4
       gsres at gsMarkerThicknessF = 2
       gsres at gsMarkerColor = "Black"
       gsn_polymarker(wks, plot, f->Lon(i), f->Lat(i), gsres)
       gsres at gsMarkerIndex = 1
       gsres at gsMarkerSizeF = 0.03
       gsres at gsMarkerColor = cmap(indx, :)
       gsn_polymarker(wks, plot, f->Lon(i), f->Lat(i), gsres)

   end do

Thanks again,

Carlos
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20151021/0a2a5386/attachment.html 
-------------- next part --------------
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"

begin
  wks  = gsn_open_wks("X11","NWS-Precip-2011071215")

  f = addfile("nws_precip_2011071215.shp","r")

  minLat = 17.30
  maxLat = 19
  latRange = maxLat - minLat
  minLon = 292
  maxLon = 295
  lonRange = maxLon - minLon

  globMin = min(f->Globvalue)
  globRange = max(f->Globvalue) - min(f->Globvalue)

  mres                       = True
  mres at mpLimitMode           = "Corners"     ; corner method of zoom
  mres at mpLeftCornerLatF      = minLat - (latRange * .01)
  mres at mpLeftCornerLonF      = minLon - (lonRange * .01)
  mres at mpRightCornerLatF     = maxLat + (latRange * .01)
  mres at mpRightCornerLonF     = maxLon + (lonRange * .01)

  mres at gsnDraw               = True          ; don't draw yet
  mres at gsnFrame              = False         ; don't advance frame yet
  mres at gsnMaximize           = False
  mres at pmTickMarkDisplayMode = "Always"      ; turn on tickmarks
  mres at tiMainString          = "NWS Precipitation July 12 2011 at 15 UTC"

  plot = gsn_csm_map(wks,mres)

  cmap = read_colormap_file("prcp_1")
  numColors = dimsizes(cmap(:,0)) - 1        ; skip the first color in the table

  print(globMin + " " + globRange)

  txres = True
  txres at txFontHeightF = 0.015

  gsres = True

   do i=0, dimsizes(f->Globvalue) - 1

      indx = doubletoint( (f->Globvalue(i)-globMin) / globRange * numColors) + 1  ; skip first color in table
      if (indx.ge.numColors) then
          indx = numColors - 1
      end if

       gsres at gsMarkerIndex = 1
       gsres at gsMarkerSizeF = 0.03
       gsres at gsMarkerColor = cmap(indx, :)
       gsn_polymarker(wks, plot, f->Lon(i), f->Lat(i), gsres)
       gsn_polymarker_ndc(wks, x, y, gsres)
       gsn_text_ndc(wks, text, x, y, txres)

   end do

 draw(plot)
 frame(wks)

end