Re: Trying to get correct map over rotated lat-lon data

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Fri Nov 30 2012 - 13:36:38 MST

Yes. That is my understanding.

Maybe, after you can modify Example 5 at:
http://www.ncl.ucar.edu/Applications/eof.shtml

To complete the example?

On 11/30/12 12:42 PM, David Brown wrote:
> Hi David,
> The issue with sfXArray and sfYArray is that, by default, the gsn_csm
> routines set the resource gsnAddCyclic to True. This appends a copy of
> the first column of the data to the end of the array resulting in the X
> dimension of the data array increasing by one element. I have to note
> here that there are 2 ways to supply 2d coordinate arrays to the gsn
> routines. You can 1) attach the special attributes lon2d and lat2d to
> the data array, or 2) you can set sfXArray and sfYArray. In the first
> case, the gsn routines set sfXArray and sfYArray for you using the
> values supplied as attributes, and also make sure to add the extra
> column to both the coordinate arrays. Unfortunately they do not adjust
> the sfXArrary and sfYArray values if you set them directly.
> The purpose of the additional column is to eliminate the gap that can
> appear in plots of global grids at the cyclic longitude point.
>
> In your case, since you do not have a global grid, you can fix the
> problem by setting gsnAddCyclic to False. You could instead change your
> code to assign the attributes lon2d and lat2d to the psfc variable
> rather than assigning to the resource values sfXArray and sfYArray. That
> seems to work okay with either setting of gsnAddCyclic. Be sure to turn
> off tfDoNDCOverlay whichever way you do it. This does result in a better
> fit of the landmask to the map outlines. Maybe the pole_lat and
> stand_lon values are not exactly correct.
> -dave
>
>
>
>
> On Nov 30, 2012, at 8:12 AM, David Jones wrote:
>
>> Dear Dave,
>>
>> thank you very much for your response. I really appreciate your help.
>>
>> Thanks for the file you sent. I made some slight changes (i.e. exact
>> corners, etc.) (attached) and was able to get it to work reasonably
>> well, but it's not really completely representing the borders in
>> places. For example, if I plot the variable LANDMASK instead of PSFC,
>> there are gaps near the edges of, for example, the Arabian Peninsula.
>> (Ignore the strange looking landmask in the polar regions, Greenland
>> etc., as that is an artifact of the reanalysis data i'm using). I
>> think that the method you've provided could well be the best that it
>> can work using the non-wrf routines in NCL.
>>
>> I've uploaded the file I'm using on the ftp (it's called
>> wrfout_d01_2009-12-01_00:00:00) so you can see this yourself if you
>> like. I've reduced it to just one time step to make the file smaller.
>>
>> When I try to use the 2D arrays XLAT and XLONG to define where the map
>> borders should be drawn I get an error message like
>> "warning:ScalarFieldSetValues: 2d coordinate array sfXArray has an
>> incorrect dimension size: defaulting sfXArray". Note that those 2D
>> arrays give the geographical lat and lon for each grid point in the
>> domain. My hope was that there would be some way to do this using the
>> WRF routines, as WRF output is not fully NetCDF compliant - which
>> could cause problems for some of the routines maybe. I thought that
>> maybe this could be done through attributes like LatLonOverlay,
>> sfXArray and sfYArray.
>>
>> FYI, ref_lat and ref_lon define the geographical latitude/longitude of
>> a specific grid point in the domain. By default it is the middle of
>> the domain, however in my case I have shifted it northwards to the
>> grid point defined by (ref_x,ref_y) as listed in the comments at the
>> bottom of my original script (re-attached).
>>
>> Many thanks for your assistance,
>>
>> Dave.
>>
>> ------------------------------------------------------------------------
>> *From:* David Brown <dbrown@ucar.edu <mailto:dbrown@ucar.edu>>
>> *To:* David Jones <jonesd647@yahoo.ca <mailto:jonesd647@yahoo.ca>>
>> *Cc:* ncl-talk Talk <ncl-talk@ucar.edu <mailto:ncl-talk@ucar.edu>>;
>> wrfhelp@ucar.edu <mailto:wrfhelp@ucar.edu>
>> *Sent:* Friday, November 30, 2012 5:56:38 AM
>> *Subject:* Re: [ncl-talk] Trying to get correct map over rotated
>> lat-lon data
>>
>> Hi David,
>>
>> I think that rotated lat/lon projections in WRF is a confusing
>> subject, especially in the context of NCL, because I suspect that the
>> NCL WRF routines may not always map the attributes defined in the WRF
>> file correctly into the resource settings required for NCL (in the
>> rotated lat/lon case only). I will add the caveat that I do not fully
>> understand the WRF attributes defined in your file either. What I do
>> know is based on looking through this web page:
>> http://www.mmm.ucar.edu/wrf/users/docs/user_guide_V3/users_guide_chap3.htm#_Description_of_GEOGRID.TBL,
>> which may or may not be the most current description of the WRF
>> projections and corresponding namelist parameters.
>>
>> Since I do not have your data file I cannot run your script for
>> myself, but I have attempted to approximate the map projection of your
>> data based on the WRF parameters you list in the script, the
>> documentation listed above, and the appearance of the surface pressure
>> field plot you attached. I am attaching my test script along with the
>> output. I think it is fairly close to the projection of your data, but
>> of course it would be better if I had the locations of the lower left
>> and upper right corners of your data. If you use gsn_csm_contour_map
>> instead of wrf_map_overlays you should be able to plot the correct map
>> projection for your data. One thing the documentation states is that
>> CEN_LON and CEN_LAT do not apply for rotated lat/lon projections, but
>> I have to wonder why they are set at all in that case. I am also not
>> sure what the role of ret_lat and ref_lon are. I do note that ref_lat
>> - pole_lat is almost equal to CEN_LAT and - STAND_LON is not too far
>> from CEN_LON, but I don't know if that has any significance.
>>
>> It would be helpful if you would upload your data file. It uses a
>> projection and a domain we have not often encountered. If you need
>> instructions on how to do that, see
>> http://www.ncl.ucar.edu/report_bug.shtml. Thanks,
>> -dave
>>
>>
>> <test.png>
>>
>>
>> On Nov 28, 2012, at 11:50 AM, David Jones wrote:
>>
>>> Dear NCL talk help,
>>>
>>> I would like to kindly ask for some advice regarding the plotting of
>>> rotated lat-lon output. I am currently trying to plot the output from
>>> a rotated lat-lon grid using NCL...
>>>
>>> When I make the contour plot of my variable (psfc), all of the
>>> contoured data appears where I expect and want it to be, the land-sea
>>> boundaries can even be detected as the value of the variable changes.
>>> So everything seems fine in that respect.
>>>
>>> However, when I try to plot the national/geographical borders using a
>>> map overlay, the borders are all in the wrong places. Please see the
>>> attached image, in which it is clearly visible.
>>>
>>> The variable I am plotting (PSFC) is made up of (time, west_east,
>>> south_north), with west_east and south_north corresponding to a point
>>> on the computational grid (i.e. they're defined by a lat and lon from
>>> the computational grid).
>>>
>>> I also have 2D arrays (XLAT and XLONG) giving the geographical
>>> latitude and longitude for every point on the computational grid. I
>>> tried to use these to plot the geographical boundaries in the correct
>>> location, but sadly to no avail.
>>>
>>> I have looked through the NCL talk archives and on the NCL website,
>>> but unfortunately nothing I've tried has worked so far.
>>>
>>> I'll include my script below and also attach a copy of the output
>>> from an ncl_filedump on my data. You'll notice in the script below
>>> that lots of the map options are commented out; that's because
>>> several tries showed them to make no difference whatsoever to the
>>> plot I received, whether they were included or commented out (I did
>>> try plotting with them uncommented though, in several different
>>> combinations). The only one of those options that did make a
>>> difference was setting mpres@mpCenterLatF = a@CEN_LAT. But that only
>>> changed it from one incorrect plot to another incorrect one, albeit
>>> the new looks slightly better (attached).
>>>
>>> If you had any guidance as to where I might be going wrong I would be
>>> most appreciative.
>>>
>>> Many thanks for any advice,
>>> Dave.
>>>
>>> *----------------------SCRIPT-----------------------------*
>>>
>>> ; Basic script for plotting model output
>>>
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>>> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
>>>
>>> begin
>>>
>>> a = addfile("wrfout_d01_2009-12-01_00:00:00.nc","r")
>>>
>>> wks = gsn_open_wks("x11","PSFC")
>>>
>>> ;---------------------------------------------------------------
>>> ; Set basic resources
>>> res = True
>>> res@MainTitle = "PSFC" ; Give a main title
>>> pltres = True ; Plotting resources
>>> mpres = True ; Map resources
>>>
>>> ;---------------------------------------------------------------
>>> ; Get variables
>>> ;---1. WRF------------------
>>>
>>> psfc = a->PSFC(0,:,:)
>>> psfc_size = dimsizes(psfc(:,:))
>>> print("dimsizes(psfc) = " + psfc_size)
>>>
>>> ;---2. Lat/Lon arrays-------
>>>
>>> lat2d = a->XLAT(0,:,:)
>>> printVarSummary(lat2d)
>>> lon2d = a->XLONG(0,:,:)
>>> printVarSummary(lon2d)
>>>
>>> ;--------------------------------------------------------------
>>> ; Plotting options
>>> opts = res ; Add basic resources
>>> opts@FieldTitle = "PSFC"
>>> opts@cnFillOn = True ; Shaded plot
>>> opts@gsnSpreadColorEnd = -8 ;
>>>
>>> contour_psfc = wrf_contour(a,wks,psfc,opts)
>>>
>>> ; Set some map options
>>>
>>> ; mpres@tfDoNDCOverlay= True
>>>
>>> ; pltres@sfXArray = lon2d
>>> ; pltres@sfYArray = lat2d
>>> ; pltres@LatLonOverlay = True
>>>
>>> ; mpres@mpCenterLonF = a@CEN_LON ;= -16.84654
>>> mpres@mpCenterLatF = a@CEN_LAT ;= 53.4331
>>>
>>> ; mpres@mpProjection = "CylindricalEquidistant"
>>> ; mpres@mpLimitMode = "Corners"
>>> ; mpres@mpLeftCornerLatF = lat2d(0,0) ; bottom left
>>> ; mpres@mpLeftCornerLonF = lon2d(0,0) ; bottom left
>>> ; mpres@mpRightCornerLatF= lat2d(268,568) ; top right
>>> ; mpres@mpRightCornerLonF= lon2d(268,568) ; top right
>>>
>>> mpres@mpGeophysicalLineColor = "Black"
>>> mpres@mpGridLineColor = "Black"
>>> mpres@mpLimbLineColor = "Black"
>>> mpres@mpNationalLineColor = "Black"
>>> mpres@mpUSStateLineColor = "Black"
>>> mpres@mpPerimLineColor = "Black"
>>> mpres@mpGeophysicalLineThicknessF = 1.0
>>> mpres@mpGridLineThicknessF = 1.0
>>> mpres@mpLimbLineThicknessF = 1.0
>>> mpres@mpNationalLineThicknessF = 1.0
>>> mpres@mpUSStateLineThicknessF = 1.0
>>>
>>> ; Namelist map projection information:
>>> ; map_proj : "lat-lon" (cylindrical equidistant)
>>> ; ref_lat : 88.5 (geographical lat of grid point (ref_x,ref_y)
>>> ; ref_lon : -5.0 (geographical lon of grid point (ref_x,ref_y)
>>> ; ref_x : 285 (total x grid points = 570)
>>> ; ref_y : 265 (total y grid points = 270)
>>> ; pole_lat : 35.0 (lat of geog N hem pole on comp grid)
>>> ; pole_lon : 180.0(lon of geog N hem pole on comp grid)
>>> ; stand_lon: 17.5 (rotation about geographic poles)
>>>
>>> ; MAKE PLOTS
>>>
>>> plot = wrf_map_overlays(a,wks,(/contour_psfc/),pltres,mpres)
>>>
>>> delete(res)
>>>
>>> end
>>>
>>> <ncl_filedump.txt><psfc_rotlatlon.ncl><PSFC.000001.png>_______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>>
>> <test_rlatlon.ncl><psfc_rotlatlon.ncl>
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Nov 30 13:36:48 2012

This archive was generated by hypermail 2.1.8 : Fri Dec 07 2012 - 13:30:06 MST