[ncl-talk] Plotting Lat/Lon Points on a map
Mary Haley
haley at ucar.edu
Tue May 30 13:08:42 MDT 2017
My guess is that the issue is with these two lines:
water = a->VAR_209_3_34_P0_L102_GLL0(:,:)
hail = water *0.0393701;
When you read "water", it has all the necessary metadata (that is,
coordinate arrays) for plotting.
However, when you do a calculation using "water"
hail = water *0.0393701;
this causes metadata to be stripped off, and "hail" has no lat/lon arrays
attached to it. This means there's no way for NCL to plot it correctly over
a map.
To preserve metadata, use this little "trick":
hail = water ; this makes an exact copy of "water", including all
metadata
hail = hail *0.0393701 ; the metadata won't get stripped from "hail"
Generally, and this is just good programming practice, you want to update
the metadata when you do a calculation like this.
hail at long_name = "Hail Size"
hail at units = "Inches"
--Mary
On Sun, May 28, 2017 at 10:13 AM, Tim Melino <melino33 at gmail.com> wrote:
> Thanks Barry and Dennis for your help on this!
>
> Dennis, I added your my other question was why is the point in NCL
> plotted south of the North Carolina border with NCL? The actual location is
> in Virginia, I noticed this was occurring in the image that you sent me
> too.
>
> [image: Inline image 1]
>
> On Sun, May 28, 2017 at 12:00 PM, Barry Lynn <barry.h.lynn at gmail.com>
> wrote:
>
>> Hi Dennis:
>>
>> Thank you for helping. Sorry I was unavailable to (try to) do more.
>>
>> Barry
>>
>> On Sun, May 28, 2017 at 6:25 PM, Dennis Shea <shea at ucar.edu> wrote:
>>
>>> See attached. Please look carefully at the output produced by
>>> printVarSummary and printMinMax
>>>
>>> Note:
>>>
>>> [1] GRIB is a moving target. 6.3.0 did not have the look-up tabl
>>> builtin while NCL 6.4.0 does.
>>>
>>> [2] the grib file has 'no report' coded as as -1 ... the script sets
>>> those values to _FillValue
>>>
>>> [3] The grid is *rectilinear* :
>>> In netCDF speak, a one dimensional variable where the variable name is
>>> the same as the variable dimension name is a coordinate variable.
>>> lat_0(lat_0) and lon_0(lon_0)
>>>
>>> Good Luck
>>>
>>>
>>>
>>> On Sun, May 28, 2017 at 10:21 AM, Barry Lynn <barry.h.lynn at gmail.com>
>>> wrote:
>>>
>>>> Tim:
>>>>
>>>> If you can wait a day or two, I suggest you write to:
>>>>
>>>> Carrie.Langston at noaa.gov
>>>>
>>>> and ask Him/Her what is the format of the gridded data and if they have
>>>> an NCL program to read it.
>>>>
>>>> It's looks like very useful data once you "crack the code."
>>>>
>>>> Barry
>>>>
>>>> On Sun, May 28, 2017 at 5:09 PM, Tim Melino <melino33 at gmail.com> wrote:
>>>>
>>>>> Barry,
>>>>>
>>>>> The source of the data is here:
>>>>>
>>>>> http://www.nssl.noaa.gov/projects/mrms/operational/tables.php
>>>>>
>>>>>
>>>>>
>>>>> Tim
>>>>>
>>>>>
>>>>>
>>>>> On Sun, May 28, 2017 at 10:03 AM, Barry Lynn <barry.h.lynn at gmail.com>
>>>>> wrote:
>>>>>
>>>>>> Tim:
>>>>>>
>>>>>> If you write out the latitude and longitude information is is
>>>>>> rectilinear?
>>>>>>
>>>>>> When trying to sort this out, you should first check if the lat/lon
>>>>>> can be plotted by a single coordinate or two. If two, then there is a
>>>>>> chance that your grid is not rectilinear, and you will have to account for
>>>>>> the distortion in the grid.
>>>>>>
>>>>>> If there is distortion, then you need to account for it (e.g.,
>>>>>> Lambert, curvilinear.
>>>>>>
>>>>>> I can only offer advice at the moment...
>>>>>>
>>>>>> What is the data source for the data?
>>>>>>
>>>>>> Barry
>>>>>>
>>>>>> On Sun, May 28, 2017 at 4:52 PM, Tim Melino <melino33 at gmail.com>
>>>>>> wrote:
>>>>>>
>>>>>>> Hey Barry,
>>>>>>>
>>>>>>> The only information about the projection was: *grid_type :
>>>>>>> Latitude/Longitude. *I attached the file!
>>>>>>>
>>>>>>> Tim
>>>>>>>
>>>>>>> On Sun, May 28, 2017 at 9:00 AM, Barry Lynn <barry.h.lynn at gmail.com>
>>>>>>> wrote:
>>>>>>>
>>>>>>>> Tim:
>>>>>>>>
>>>>>>>> Please remind me what the name of the file is.
>>>>>>>>
>>>>>>>> Also, was there any information at the top of the printing about
>>>>>>>> the projection?
>>>>>>>>
>>>>>>>> Yes, the markers would (could) look off.
>>>>>>>>
>>>>>>>> On Sun, May 28, 2017 at 3:32 PM, Tim Melino <melino33 at gmail.com>
>>>>>>>> wrote:
>>>>>>>>
>>>>>>>>> Hi Barry,
>>>>>>>>>
>>>>>>>>> Here is the information when printing (a). Would the map
>>>>>>>>> projection also cause a marker to be printed in the wrong location?
>>>>>>>>>
>>>>>>>>> Variable: a
>>>>>>>>> Type: file
>>>>>>>>> filename: test
>>>>>>>>> path: ./test.grib2
>>>>>>>>> file global attributes:
>>>>>>>>> dimensions:
>>>>>>>>> lat_0 = 30
>>>>>>>>> lon_0 = 40
>>>>>>>>> variables:
>>>>>>>>> float VAR_209_3_34_P0_L102_GLL0 ( lat_0, lon_0 )
>>>>>>>>> center : US NOAA Office of Oceanic and Atmospheric
>>>>>>>>> Research
>>>>>>>>> production_status : Research products
>>>>>>>>> long_name : unknown variable name
>>>>>>>>> units : unknown
>>>>>>>>> _FillValue : 1e+20
>>>>>>>>> grid_type : Latitude/longitude
>>>>>>>>> parameter_discipline_and_category : 209, 3
>>>>>>>>> parameter_template_discipline_category_number : (
>>>>>>>>> 0, 209, 3, 34 )
>>>>>>>>> level_type : Specific altitude above mean sea level (m)
>>>>>>>>> level : 500
>>>>>>>>> forecast_time : 0
>>>>>>>>> forecast_time_units : minutes
>>>>>>>>> initial_time : 07/09/2016 (04:46)
>>>>>>>>>
>>>>>>>>> float lat_0 ( lat_0 )
>>>>>>>>> long_name : latitude
>>>>>>>>> grid_type : Latitude/Longitude
>>>>>>>>> units : degrees_north
>>>>>>>>> Dj : 0.01
>>>>>>>>> Di : 0.009999676
>>>>>>>>> Lo2 : 280.735
>>>>>>>>> La2 : 36.695
>>>>>>>>> Lo1 : 280.345
>>>>>>>>> La1 : 36.405
>>>>>>>>>
>>>>>>>>> float lon_0 ( lon_0 )
>>>>>>>>> long_name : longitude
>>>>>>>>> grid_type : Latitude/Longitude
>>>>>>>>> units : degrees_east
>>>>>>>>> Dj : 0.01
>>>>>>>>> Di : 0.009999676
>>>>>>>>> Lo2 : 280.735
>>>>>>>>> La2 : 36.695
>>>>>>>>> Lo1 : 280.345
>>>>>>>>> La1 : 36.405
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> On Sun, May 28, 2017 at 12:48 AM, Barry Lynn <
>>>>>>>>> barry.h.lynn at gmail.com> wrote:
>>>>>>>>>
>>>>>>>>>> Hi:
>>>>>>>>>>
>>>>>>>>>> When the data is skewed, it is often because the projection
>>>>>>>>>> you're reading is not the one your plotting.
>>>>>>>>>>
>>>>>>>>>> Do you know the projection of the data you're using?
>>>>>>>>>>
>>>>>>>>>> If you print(a) where a is the file you specify from addfile, you
>>>>>>>>>> can see the variables (and hopefully projection) indicated within.
>>>>>>>>>>
>>>>>>>>>> Barry
>>>>>>>>>>
>>>>>>>>>> On Sat, May 27, 2017 at 4:20 PM, Tim Melino <melino33 at gmail.com>
>>>>>>>>>> wrote:
>>>>>>>>>>
>>>>>>>>>>> Hi,
>>>>>>>>>>>
>>>>>>>>>>> I am attempting to plot a subset of the MRMS data from NCEP.
>>>>>>>>>>> The latitude and longitude I am plotting is :
>>>>>>>>>>>
>>>>>>>>>>> LAT: 36.5459089
>>>>>>>>>>> LON: -79.4585847
>>>>>>>>>>>
>>>>>>>>>>> I attached an image of the output along with the grib file I am
>>>>>>>>>>> using for input. I am noticing a couple things first the data appears to
>>>>>>>>>>> be skewed, so I am not sure if it is an issue with the map projection.
>>>>>>>>>>> Second the point which should be plotted in Virginia is actually plotted in
>>>>>>>>>>> North Carolina, the only thing that I have been able to do if offset the
>>>>>>>>>>> marker to compensate for the issue. Does anyone have any ideas on how to
>>>>>>>>>>> fix these issues? Any ideas would be greatly appreciated!
>>>>>>>>>>>
>>>>>>>>>>> Thanks!
>>>>>>>>>>> Tim
>>>>>>>>>>>
>>>>>>>>>>> Code:
>>>>>>>>>>>
>>>>>>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>>>>>>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>>>>>>>>>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>>>>>>>>>>>
>>>>>>>>>>> begin
>>>>>>>>>>>
>>>>>>>>>>> lati=LAT
>>>>>>>>>>> loni=LON
>>>>>>>>>>> date=DATE
>>>>>>>>>>> address=ADDRESS
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ; Convert the Lat/Lon to a Float
>>>>>>>>>>> latc= stringtofloat(lati)
>>>>>>>>>>> lonc = stringtofloat(loni)
>>>>>>>>>>>
>>>>>>>>>>> ; Set the Map Bounds
>>>>>>>>>>> minlat = latc - .12
>>>>>>>>>>> maxlat = latc + .12
>>>>>>>>>>> minlon = lonc - .15
>>>>>>>>>>> maxlon = lonc + .15
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ;========================
>>>>>>>>>>> ; get list of all files and open as "one big file"
>>>>>>>>>>> ;========================
>>>>>>>>>>> a = addfile( "./" + lati + "_" + loni + ".grib2","r")
>>>>>>>>>>>
>>>>>>>>>>> type = "png"
>>>>>>>>>>> wks = gsn_open_wks(type,"./" + lati + "" + loni + ".png")
>>>>>>>>>>>
>>>>>>>>>>> ; Set some Basic Plot options
>>>>>>>>>>> res = True
>>>>>>>>>>> res at gsnDraw = False ; do not draw the
>>>>>>>>>>> plot
>>>>>>>>>>> res at gsnFrame = False ; do not advance
>>>>>>>>>>> the frame
>>>>>>>>>>> pltres = True
>>>>>>>>>>> pltres at NoTitles = True
>>>>>>>>>>> pltres at CommonTitle = True
>>>>>>>>>>> pltres at PlotTitle = ""
>>>>>>>>>>> pltres at PanelPlot = True
>>>>>>>>>>> pltres at FramePlot = False
>>>>>>>>>>>
>>>>>>>>>>> ; Set the Colormap
>>>>>>>>>>> gsn_define_colormap(wks,"hail")
>>>>>>>>>>> lat = a->lat_0
>>>>>>>>>>> lon = a->lon_0
>>>>>>>>>>> water = a->VAR_209_3_34_P0_L102_GLL0(:,:)
>>>>>>>>>>> hail = water *0.0393701;
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> ; Plotting options for Hail
>>>>>>>>>>> res at gsnMaximize = True
>>>>>>>>>>> res at mpGeophysicalLineThicknessF = 4.0
>>>>>>>>>>> ;res at mpProjection = "LambertConformal"
>>>>>>>>>>> res at mpGridLineThicknessF = 0.5
>>>>>>>>>>> res at mpLimbLineThicknessF = 4.0
>>>>>>>>>>> res at mpNationalLineThicknessF = 4.0
>>>>>>>>>>> res at mpUSStateLineThicknessF = 4.0
>>>>>>>>>>> res at mpProvincialLineThicknessF = 4.0
>>>>>>>>>>> res at mpOutlineBoundarySets = "AllBoundaries"
>>>>>>>>>>> res at mpDataBaseVersion = "Ncarg4_1"
>>>>>>>>>>> res at mpDataSetName = "Earth..4"
>>>>>>>>>>> res at mpFillOn = False ; fill continents
>>>>>>>>>>> res at mpPerimOn = False ; no box around map
>>>>>>>>>>> res at tiMainOn = True
>>>>>>>>>>> res at tiMainFontColor = "Black"
>>>>>>>>>>> res at tiMainString =""
>>>>>>>>>>> res at tiXAxisString = ""
>>>>>>>>>>> res at tiYAxisString = ""
>>>>>>>>>>> res at gsnLeftString = "Date: " + date ; add the
>>>>>>>>>>> gsn titles
>>>>>>>>>>> res at gsnRightString = ""
>>>>>>>>>>> res at gsnStringFontHeightF = 0.016
>>>>>>>>>>> res at cnInfoLabelOrthogonalPosF = 0.07 ; offset second label
>>>>>>>>>>> information
>>>>>>>>>>> res at gsnContourLineThicknessesScale = 0.001
>>>>>>>>>>> res at cnFillOn = True
>>>>>>>>>>> res at lbTitleOn = True ; remove
>>>>>>>>>>> field name from label bar
>>>>>>>>>>>
>>>>>>>>>>> ; Set the Map Bounds
>>>>>>>>>>> res at sfYArray = lat
>>>>>>>>>>> res at sfXArray = lon
>>>>>>>>>>> res at mpLimitMode = "LatLon"
>>>>>>>>>>> res at mpMinLatF = minlat
>>>>>>>>>>> res at mpMaxLatF = maxlat
>>>>>>>>>>> res at mpMinLonF = minlon
>>>>>>>>>>> res at mpMaxLonF = maxlon
>>>>>>>>>>> res at tmXBOn = False
>>>>>>>>>>> res at tmYLOn = False
>>>>>>>>>>> res at tmXTOn = False
>>>>>>>>>>> res at tmYROn = False
>>>>>>>>>>> res at gsnAddCyclic = False
>>>>>>>>>>> res at tfDoNDCOverlay = False ; do
>>>>>>>>>>> not transform data
>>>>>>>>>>> res at cnLevelSelectionMode = "ExplicitLevels" ; set
>>>>>>>>>>> explicit contour levels
>>>>>>>>>>> res at cnLevels = (/ 0.5,0.75, 1.00, 1.25, 1.50, 1.75,
>>>>>>>>>>> 2.00, 2.50, 3.00, 3.50, 4.00, 4.50, 5.00/) ; set levels
>>>>>>>>>>> res at tiMainFontHeightF = 0.018 ; font
>>>>>>>>>>> height of title
>>>>>>>>>>> res at lbAutoManage = False ; we
>>>>>>>>>>> control label bar
>>>>>>>>>>> res at pmLabelBarDisplayMode = "Always" ; turns
>>>>>>>>>>> on label bar
>>>>>>>>>>> res at lbOrientation = "Horizontal" ; ncl
>>>>>>>>>>> default is vertical
>>>>>>>>>>> res at pmLabelBarSide = "Bottom" ; default
>>>>>>>>>>> is right
>>>>>>>>>>> res at lbLabelStride = 1 ; skip
>>>>>>>>>>> every other label
>>>>>>>>>>> res at pmLabelBarWidthF = 0.7 ; default
>>>>>>>>>>> is shorter
>>>>>>>>>>> res at pmLabelBarHeightF = 0.1 ; default
>>>>>>>>>>> is taller
>>>>>>>>>>> res at lbLabelFontHeightF = .012 ; default
>>>>>>>>>>> is HUGE
>>>>>>>>>>> res at lbPerimOn = False ; default
>>>>>>>>>>> has box
>>>>>>>>>>> res at lbTitleString = "Hail Size (Inches)"
>>>>>>>>>>> res at lbTitleFontHeightF = 0.012
>>>>>>>>>>> plot = gsn_csm_contour_map(wks,hail,res)
>>>>>>>>>>>
>>>>>>>>>>> ; Setup the point
>>>>>>>>>>> mstring1 = "z"
>>>>>>>>>>> fontnum1 = 35
>>>>>>>>>>> xoffset1 = .4
>>>>>>>>>>> yoffset1 = .7
>>>>>>>>>>> ratio1 = 0.0
>>>>>>>>>>> size1 = 1.0
>>>>>>>>>>> angle1 = 0.0
>>>>>>>>>>>
>>>>>>>>>>> Tm_index = NhlNewMarker(wks, mstring1, fontnum1, xoffset1,
>>>>>>>>>>> yoffset1, ratio1, size1, angle1)
>>>>>>>>>>>
>>>>>>>>>>> ; Plot the sites and label them
>>>>>>>>>>> gres = True
>>>>>>>>>>> gres at gsMarkerIndex = Tm_index ; 16
>>>>>>>>>>> gres at gsMarkerSizeF = 0.015 ; .009
>>>>>>>>>>> gres at gsMarkerColor = "black"
>>>>>>>>>>> gres at txFont ="helvetica-bold"
>>>>>>>>>>> gres at txFontHeightF = "25"
>>>>>>>>>>>
>>>>>>>>>>> site = (/address/)
>>>>>>>>>>> lats = (/latc/)
>>>>>>>>>>> lons = (/lonc/)
>>>>>>>>>>>
>>>>>>>>>>> dot = gsn_add_polymarker(wks,plot,lons,lats,gres)
>>>>>>>>>>> tlats =lats + .03
>>>>>>>>>>> tlons =lons
>>>>>>>>>>> text = gsn_add_text(wks,plot,site,tlons,tlats ,gres)
>>>>>>>>>>>
>>>>>>>>>>> ; Draw the graphics
>>>>>>>>>>> draw(plot)
>>>>>>>>>>> frame(wks)
>>>>>>>>>>>
>>>>>>>>>>> end
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>> _______________________________________________
>>>>>>>>>>> ncl-talk mailing list
>>>>>>>>>>> ncl-talk at ucar.edu
>>>>>>>>>>> List instructions, subscriber options, unsubscribe:
>>>>>>>>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>>>>>>>>
>>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> --
>>>>>>>>>> Barry H. Lynn, Ph.D
>>>>>>>>>> Senior Lecturer,
>>>>>>>>>> The Institute of the Earth Science,
>>>>>>>>>> The Hebrew University of Jerusalem,
>>>>>>>>>> Givat Ram, Jerusalem 91904, Israel
>>>>>>>>>> Tel: 972 547 231 170
>>>>>>>>>> Fax: (972)-25662581
>>>>>>>>>>
>>>>>>>>>> C.E.O, Weather It Is, LTD
>>>>>>>>>> Weather and Climate Focus
>>>>>>>>>> http://weather-it-is.com
>>>>>>>>>> Jerusalem, Israel
>>>>>>>>>> Local: 02 930 9525
>>>>>>>>>> Cell: 054 7 231 170
>>>>>>>>>> Int-IS: x972 2 930 9525
>>>>>>>>>> US 914 432 3108 <(914)%20432-3108>
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> --
>>>>>>>> Barry H. Lynn, Ph.D
>>>>>>>> Senior Lecturer,
>>>>>>>> The Institute of the Earth Science,
>>>>>>>> The Hebrew University of Jerusalem,
>>>>>>>> Givat Ram, Jerusalem 91904, Israel
>>>>>>>> Tel: 972 547 231 170
>>>>>>>> Fax: (972)-25662581
>>>>>>>>
>>>>>>>> C.E.O, Weather It Is, LTD
>>>>>>>> Weather and Climate Focus
>>>>>>>> http://weather-it-is.com
>>>>>>>> Jerusalem, Israel
>>>>>>>> Local: 02 930 9525
>>>>>>>> Cell: 054 7 231 170
>>>>>>>> Int-IS: x972 2 930 9525
>>>>>>>> US 914 432 3108 <(914)%20432-3108>
>>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> Barry H. Lynn, Ph.D
>>>>>> Senior Lecturer,
>>>>>> The Institute of the Earth Science,
>>>>>> The Hebrew University of Jerusalem,
>>>>>> Givat Ram, Jerusalem 91904, Israel
>>>>>> Tel: 972 547 231 170
>>>>>> Fax: (972)-25662581
>>>>>>
>>>>>> C.E.O, Weather It Is, LTD
>>>>>> Weather and Climate Focus
>>>>>> http://weather-it-is.com
>>>>>> Jerusalem, Israel
>>>>>> Local: 02 930 9525
>>>>>> Cell: 054 7 231 170
>>>>>> Int-IS: x972 2 930 9525
>>>>>> US 914 432 3108 <(914)%20432-3108>
>>>>>>
>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> Barry H. Lynn, Ph.D
>>>> Senior Lecturer,
>>>> The Institute of the Earth Science,
>>>> The Hebrew University of Jerusalem,
>>>> Givat Ram, Jerusalem 91904, Israel
>>>> Tel: 972 547 231 170
>>>> Fax: (972)-25662581
>>>>
>>>> C.E.O, Weather It Is, LTD
>>>> Weather and Climate Focus
>>>> http://weather-it-is.com
>>>> Jerusalem, Israel
>>>> Local: 02 930 9525
>>>> Cell: 054 7 231 170
>>>> Int-IS: x972 2 930 9525
>>>> US 914 432 3108 <(914)%20432-3108>
>>>>
>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk at ucar.edu
>>>> List instructions, subscriber options, unsubscribe:
>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>>
>>>>
>>>
>>
>>
>> --
>> Barry H. Lynn, Ph.D
>> Senior Lecturer,
>> The Institute of the Earth Science,
>> The Hebrew University of Jerusalem,
>> Givat Ram, Jerusalem 91904, Israel
>> Tel: 972 547 231 170
>> Fax: (972)-25662581
>>
>> C.E.O, Weather It Is, LTD
>> Weather and Climate Focus
>> http://weather-it-is.com
>> Jerusalem, Israel
>> Local: 02 930 9525
>> Cell: 054 7 231 170
>> Int-IS: x972 2 930 9525
>> US 914 432 3108 <(914)%20432-3108>
>>
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20170530/ebc5f7f5/attachment.html
-------------- next part --------------
A non-text attachment was scrubbed...
Name: sample_point.png
Type: image/png
Size: 95010 bytes
Desc: not available
Url : http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20170530/ebc5f7f5/attachment.png
More information about the ncl-talk
mailing list