Re: NARR contour plots

Hi folks,

Adam asked Dave Brown or I to take a look at this.

First, you only want to set tfDoNDCOverlay to True if:

  - your data is already projected (as you said your data was)
  - you know the *exact* map projection for your data

If you don't have the exact map projection, the data may
look distorted.

Also, if you set res_at_tfDoNDCOverlay to True, then you do
*not* want to set the special lat2d/lon2d attributes, because
this is redundant. You don't need lat/lon data if the map
projection has been set up to exactly match your grid.

However, if you don't know the exact map projection, *or*
you want to put the data on a different map than what it
was originally on, then you can use lat2d/lon2d. In this case,
make sure res_at_tfDoNDCOverlay is not set at all, or set it
to False.

Assuming you have the first scenario, where you know the
exact map projection, then here's what you want to do:

Comment out any lat2d/lon2d attribute settings:

slp_at_lat2d = lat2d
slp_at_lon2d = lon2d

Remove lines like:

> slp!0="lat2d"
> slp!1="lon2d"
> slp&lat2d_at_units="degrees_north"
> slp&lon2d_at_units="degrees_east"

because these don't make sense with 2D lat/lon arrays.
In fact, this should have given you an error, because
you can't have 2D coordinate arrays.

If the above doesn't work, then try scenario #2. That is,
set res_at_tfDoNDCOverlay to False (or don't set
it at all) *and* set the lat2d/lon2d attributes for
each data set you plot.

slp_at_lat2d = lat2d
slp_at_lon2d = lon2d
etc.

Do *not* set anything else with regard to lat2d/lon2d. The above is
all you need.

If neither of these scenarios work, then it would be helpful to get
your data and script.

--Mary

On Mar 19, 2009, at 4:04 PM, David Small wrote:

> Thanks Adam.
>
> I think the problem is that NARR is on a projection, not a grid. The
> gsn_contour won't accept the lat and long coordinates. Will the
> contour
> program be able to contour projected data and not have it be
> distorted?
>
> I pasted the code below. It gives an error where I call
> gsn_contour_map().
>
> Thanks for all the help. I
>
> Thanks again.
>
> Dave
>
> ;*************************************************
> ; station_2.ncl
> ;*************************************************
> ;
> ;
> ;
> ;
> ;
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> load "/users/davidsmall/desktop/NCLwork/cc2g.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
> begin
> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>
>
> mon = 9
> yr = 1999
>
> f = addfile
> ("/users/davidsmall/desktop/Thesis_misc/WIND_EVENTS/prmsl."+yr
> +".nc", "r")
> fz = addfile
> ("/users/davidsmall/desktop/Thesis_misc/WIND_EVENTS/hgt."+yr+"0"+mon
> +".nc",
> "r")
>
> lat2d = f->lat
> lon2d = f->lon
> dimlc = dimsizes(lat2d)
> nlat = dimlc(0)
> mlon = dimlc(1)
>
> t = f->time
> ut = ut_calendar(t, 0)
> mo=ut(:,1)
> am=ind(mo.eq.mon)
>
>
> nlat2=277
> nlon2=349
> ;################################################################
> do i=172,176
>
> slp = short2flt(f->prmsl(am(i),:,:))
> time = f->time(am(i))
> utc = ut_calendar(time, 0)
>
> year = utc(:,0)
> month=utc(:,1)
>
> day = utc(:,2)
> hour = utc(:,3)
>
> slp = slp/100.0
> ; x@_FillValue=(max(x))
>
> slp!0="lat2d"
> slp!1="lon2d"
> slp&lat2d_at_units="degrees_north"
> slp&lon2d_at_units="degrees_east"
> slp_at_lat2d = lat2d
> slp_at_lon2d = lon2d
> ;##########################################################################
>
>
> ;************************************************
> ; create plot
> ;************************************************
> wks = gsn_open_wks ("eps", "slp_hgt."+year+"."+month+"."+day
> +"."+hour)
> ; open workstation
> gsn_define_colormap (wks,"BlAqGrYeOrRe") ; choose color map
>
> res = True ; plot mods desired
> for
> original grid
> res_at_cnFillOn = False ; color fill
> res_at_cnLinesOn = True ; no contour lines
> ; res_at_gsnSpreadColors = True ; use total
> colormap
> ; res_at_gsnSpreadColorStart = 4
> ; res_at_gsnSpreadColorEnd = -1
> res_at_mpGridAndLimbOn = True
> res_at_pmTickMarkDisplayMode = "Always" ; turn on tickmarks
> res_at_tmXTOn = False
> res_at_gsnAddCyclic = False ; regional data
>
> res_at_mpLimitMode = "Corners"
>
> ; choose range of map
> res_at_mpLeftCornerLatF = lat2d(150,30)
> res_at_mpLeftCornerLonF = lon2d(150,30)
> res_at_mpRightCornerLatF = lat2d(nlat-10,mlon-150)
> res_at_mpRightCornerLonF = lon2d(nlat-10,mlon-150)
>
> res_at_tfDoNDCOverlay = True ; direct plot
> res_at_mpProjection = "LambertConformal"
> res_at_mpLambertParallel1F = f-
> >Lambert_Conformal_at_standard_parallel(0)
> ; lat2d_at_mpLambertParallel1F
> res_at_mpLambertParallel2F = f-
> >Lambert_Conformal_at_standard_parallel(1)
> res_at_mpLambertMeridianF =
> f->Lambert_Conformal_at_longitude_of_central_meridian
> res_at_gsnMaximize = True
> res_at_pmTickMarkDisplayMode = "Never"
> res_at_mpFillOn= False
>
> ; res_at_lbLabelBarOn = False
> res_at_gsnLeftString = " "
> res_at_gsnRightString = " "
>
> res_at_cnLevelSelectionMode = "ManualLevels"
> res_at_cnInfoLabelOn = False
> res_at_cnMinLevelValF = 950 ; set min contour level
> res_at_cnMaxLevelValF = 1035 ; set max contour level
> res_at_cnLevelSpacingF = 4 ; set contour spacing
>
> res_at_mpDataBaseVersion = "MediumRes"
> res_at_mpGeophysicalLineThicknessF = 4.0
> ; res_at_cnLineLabelFont="times-roman"
> ; res_at_lbLabelBarOn = True
> res_at_gsnLeftString = " "
> res_at_gsnRightString = " "
>
>
> res_at_gsnDraw = False ; don't draw
> res_at_gsnFrame = False ; don't advance frame
>
>
>
> res_at_gsnLeftString = "" ; no titles
> res_at_gsnRightString = ""
> res_at_tiXAxisString = ""
> res_at_tiYAxisString = ""
> ; res_at_lbOrientation = "vertical"
> res_at_lbLabelAutoStride = True
> res_at_lbLabelFont= "times-roman"
>
>
> res_at_mpGridAndLimbOn = False
>
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;************************************************
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> res_at_gsnAddCyclic = False ; regional data
>
>
> if(day.lt.10)
> dd="0"+day
> end if
> if(day.gt.9)
> dd=day
> end if
> if(hour.lt.10)
> hh=" 0"+hour
> end if
> if(hour.gt.10)
> hh=" "+hour
> end if
>
> nme = year+"0"+month+dd+" "+hh+"Z"
>
> res_at_gsnLeftString = ""
> res_at_gsnRightString = nme
> res_at_gsnStringFont= "times-roman"
> res_at_vcRefAnnoOn = False
> res_at_gsnScalarContour = True
> res_at_cnFillOn = False
> res_at_cnLinesOn = True ; turn off contour lines
> res_at_gsnSpreadColors = True ; use full range of
> colormap
> res_at_gsnContourLineThicknessesScale = 3.
> res_at_cnLineLabelFontHeightF = 0.0075
> ; res_at_vcMinDistanceF = 0.017 ; thin out vectors
> ; plot = gsn_csm_vector_scalar_map(wks,uwnd,vwnd,slp,res)
> plot = gsn_csm_contour_map(wks,slp,res)
>
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ; Add thickness contours
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
> ;$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
> $$$$$$
>
> z = short2flt(fz->hgt((i),16,:,:))
> z0 = short2flt(fz->hgt((i),0,:,:))
>
> z = z-z0
> delete(z0)
>
> z!0="lat2d"
> z!1="lon2d"
> z&lat2d_at_units="degrees_north"
> z&lon2d_at_units="degrees_east"
> z_at_lat2d = lat2d
> z_at_lon2d = lon2d
>
> res_at_tfDoNDCOverlay = True ; direct plot
> res_at_cnLevelSelectionMode = "ManualLevels"
> res_at_cnInfoLabelOn = False
> res_at_cnMinLevelValF = 4800 ; set min contour level
> res_at_cnMaxLevelValF = 6000 ; set max contour level
> res_at_cnLevelSpacingF = 60 ; set contour spacing
>
> res_at_cnLineColor = "blue"
> res_at_cnMonoLineDashPattern = True
> res_at_cnLineDashPattern = 2
> res_at_cnLineLabelFontColor = "blue"
> plot2 = gsn_csm_contour(wks,z,res)
>
> ;************************************************
> overlay(plot,plot2)
>
> end do
>
> end
>
>
>
> On 3/19/09 5:42 PM, "Adam Phillips" <asphilli_at_cgd.ucar.edu> wrote:
>
>> Hi David,
>> Yes, you should be able to overlay contour plots, NARR (w/2D lats/
>> lons)
>> or otherwise. Can you send the group the code you are using?
>>
>> One thing to watch for: Make sure the 2nd plot going into overlay
>> is not
>> a "_map" plotting function. For instance, if the first plot was
>> created
>> with this:
>> plot = gsn_csm_contour_map(wks,data,res)
>>
>> then the second plot should be created with gsn_csm_contour:
>> plot2 = gsn_csm_contour(wks,data2,res)
>> overlay(plot,plot2)
>>
>> Adam
>>
>>
>> David Small wrote:
>>> Thank you for your help. When I try to superimpose two sets of
>>> contours
>>> from NARR data, thickness and sea level pressure, I get the
>>> following error:
>>>
>>> NhlAddOverlay: plot class mapPlotClass cannot be overlay plot member
>>>
>>> I can contour thickness or pressure separately (with or without
>>> winds).
>>> But I can’t lay the contours on top of each other. I’ve been told
>>> in
>>> the past that it is not possible to overlay contours using NARR in
>>> NCL.
>>> Other people have told me it is possible. Is it possible? I really
>>> need to know because otherwise I will have to switch back to GEMPAK
>>> (which I don’t want to do).
>>>
>>> Thanks again,
>>> Dave Small
>>>
>>> On 3/18/09 11:36 AM, "Mary Haley" <haley_at_ucar.edu> wrote:
>>>
>>> Hi Dave,
>>>
>>> See example 13 at:
>>>
>>> http://www.ncl.ucar.edu/Applications/panel.shtml#ex13
>>>
>>> This isn't NARR data, but it shows how to overlay contours and
>>> vectors on a map.
>>>
>>> For NARR data, since you have 2D lat/lon, this example would
>>> additionally need the special "lon2d" and "lat2d" attributes
>>> attached to
>>> the data, OR, you would need to know the exact map projection
>>> that
>>> your data is on, and set the special "tfDoNDCOverlay" resource
>>> to True. See the NARR examples for how to use this special
>>> resource:
>>>
>>> http://www.ncl.ucar.edu/Applications/narr.shtml
>>>
>>> Meanwhile, if you don't have the exact map projection, and you
>>> want
>>> to use the 2D lat/lon coords: assumed that the data in the
>>> panel 13
>>> example had 2D lat/lon coordinates called "lat" and "lon" on the
>>> file. Then, we would do something like this:
>>>
>>> a = addfile("uv300.nc","r")
>>> u = a->U(1,:,:)
>>> v = a->V(1,:,:)
>>> speed = sqrt(u^2+v^2)
>>> u_at_lat2d = a->lat
>>> u_at_lon2d = a->lon
>>> v_at_lat2d = u_at_lat2d
>>> v_at_lon2d = u_at_lon2d
>>> speed_at_lat2d = u_at_lat2d
>>> speed_at_lon2d = u_at_lon2d
>>> Everything else would be the same.
>>> For the vector plot, if you want wind barbs, set:
>>> vcres_at_vcGlpyhStyle = "WindBarb"
>>> --Mary
>>>
>>> On Mar 17, 2009, at 8:45 PM, David Small wrote:
>>>
>>> Hello everyone,
>>>
>>> The NARR dataset uses 2 dimensional latitude and
>>> longitude. Is
>>> it possible to lay several sets of contours over top of one
>>> another, possibly with vector winds? If so, can someone
>>> tell
>>> me how to do it?
>>>
>>> Thank you very much.
>>>
>>> Dave
>>>
>>> --
>>> David Small
>>> Graduate Research Assistant
>>> Department of Atmospheric and Oceanic Sciences
>>> McGill University
>>> Montreal, Quebec, Canada
>>>
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>>
>>>
>>>
>>> --
>>> David Small
>>> Graduate Research Assistant
>>> Department of Atmospheric and Oceanic Sciences
>>> McGill University
>>> Montreal, Quebec, Canada
>>>
>>>
>>> ------------------------------------------------------------------------
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
> --
> David Small
> Graduate Research Assistant
> Department of Atmospheric and Oceanic Sciences
> McGill University
> Montreal, Quebec, Canada
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Thu Mar 19 2009 - 23:15:26 MDT