[ncl-talk] Masking
Mary Haley
haley at ucar.edu
Thu Nov 9 13:28:25 MST 2017
Muhammad,
A long email follows. Please read carefully, visit the links I've provided,
and look at the two scripts I've provided.
As I said in a previous email, you have curvilinear data, which means your
data is represented by two-dimensional latitude / longitude arrays. You
cannot use coordinate arrays with curvilinear data, as you are trying to do
with this code:
lat1d = fspan(minlat,maxlat,nlat)
lon1d = fspan(minlon,maxlon,nlon)
lat1d at units = "degrees_north"
lon1d at units = "degrees_east"
tc2a = wrf_user_getvar(a, "T2",0)
tc2a at _FillValue = -9999
tc2a!0 = "lat"
tc2a!1 = "lon"
tc2a&lat = lat1d
tc2a&lon = lon1d
When you use coordinate arrays via "&", this tells NCL that you have
rectilinear data. This will NOT work for your data. To help understand the
difference between curvilinear and rectilinear data, please see our
"Plotting data on a map" examples page:
http://www.ncl.ucar.edu/Applications/plot_data_on_map.shtml
In particular, look at these two images, where the actual lat/lon grid is
drawn with lines:
http://www.ncl.ucar.edu/Applications/Images/dataonmap_grid_1_lg.png
http://www.ncl.ucar.edu/Applications/Images/dataonmap_grid_3_lg.png
The dataonmap_grid1.ncl script shows how to plot rectilinear data, and the
dataonma_grid_3.ncl is plotting curvilinear data.
You said that the code I sent you a couple of days ago was working, so I'm
not sure why you went back to trying to use coordinate arrays.
Here's the correct way, which is just copied-and-pasted from the email I
sent you before.
nt = 0
tc2a
= wrf_user_getvar(a, "T2",
nt)
tc2a at lat2d = wrf_user_getvar(a, "lat", nt) ; This is necessary for
shapefile masking, but do not
tc2a at lon2d = wrf_user_getvar(a,
"lon",
nt) ; write these two attributes to a new NetCDF file!
tc2a at _FillValue = -9999
opt = True
opt at debug = True
tc2a_mask = shapefile_mask_data(tc2a,"Shapefile/malay.shp",opt)
If you are wanting to apply the same mask to all other variables on the
file, then what you can do is have shapefile_mask_data return a mask array,
which is an array with 0s and 1s, where 1s represent data that falls inside
the area, and 0s represent data that falls outside the area. You can then
use this mask array to mask all other variables that are defined on the
same lat/lon grid.
Please see the attached wrf_gsn_5.ncl script, which I adapted from the one
on the web. It takes a shapefile that has the Mississippi River Basin
outline, and uses this to mask WRF data that covers the east part of the
United States. The WRF data doesn't cover the full area of the Mississippi
River Basin, but I wanted to use this example because I can provide you
with the data files so you can try it yourself.
The wrf_gsn_5.ncl script writes the mask array to a NetCDF file, and then
the wrf_mask.ncl script shows how you can read this mask array back in and
apply it to other variables on the WRF file. This is MUCH faster than
calling shapefile_mask_data every single time. In this case,
shapefile_mask_data takes over a minute to create the mask array. Once you
have this array, masking data takes a fraction of a second.
Try running wrf_gsn_5.ncl first, which will create several plots (to an X11
window). It will also write the mask array to "wrf_mask.nc". This script
will take over a minute to run.
Then, run wrf_mask.ncl, which masks the same variables, but this time it
uses the mask array read off wrf_mask.nc to mask the data. This script
takes a few seconds to run. It will produce a single PNG image, which I've
attached.
I've attached the four mrb.xxx files which are the files that represent the
shapefile.
To get the WRF output fie, which is larger:
wget ftp://ftp.cgd.ucar.edu//incoming/wrfout_d01_2005-12-14_13:00:00
or:
ftp ftp.cgd.ucar.edu
[login as anonymous]
cd incoming
get wrfout_d01_2005-12-14_13:00:00
quit
--Mary
On Thu, Nov 9, 2017 at 6:39 AM, Muhammad Omer Mughal <
m.mughal1 at postgrad.curtin.edu.au> wrote:
> Hi Marry
>
> Thank you again for the reply. I found the problem after I used CDO the
> time dimension stripped out from the data. In order to correct this problem
> I use the following
> lat=a->XLAT
> lon=a->XLONG
> minlat=min(lat)
> maxlat=max(lat)
> minlon=min(lon)
> maxlon=max(lon)
>
> nlat=129
> nlon=210
>
> lat1d = fspan(minlat,maxlat,nlat)
> lon1d = fspan(minlon,maxlon,nlon)
> lat1d at units = "degrees_north"
> lon1d at units = "degrees_east"
> tc2a = wrf_user_getvar(a, "T2",0)
> tc2a at _FillValue = -9999
>
> tc2a!0 = "lat"
> tc2a!1 = "lon"
> tc2a&lat = lat1d
> tc2a&lon = lon1d
> tc2a_mask = shapefile_mask_data(tc2a,"/data/muhdomer/NSCC/scratch/
> Shapefiles/SGP_adm0.shp",opt)
>
> Is there a way to get all variables masked using a shapefile_mask_data
> command ?
>
>
>
>
> Muhammad Omer Mughal
> MSc BSc Mechanical Engineering
> PhD Research Scholar
> Remote Sensing and Satellite Research Group
> Department of Imaging and Applied Physics
> Curtin University
>
> Curtin University
> Tel | +61 8 9266 7962
> Fax | +61 8 9266 2377
> Mobile | 0470 237 525
>
> Email | m.mughal1 at postgrad.curtin.edu.au<mailto:m.lynch at curtin.edu.au>
> Web | http://curtin.edu.au<http://curtin.edu.au/>
>
> Curtin University is a trademark of Curtin University of Technology.
> CRICOS Provider Code 00301J (WA), 02637B (NSW)
>
>
> ________________________________
> From: Mary Haley <haley at ucar.edu>
> Sent: Thursday, 9 November 2017 12:30:26 AM
> To: Muhammad Omer Mughal
> Cc: ncl-talk at ucar.edu
> Subject: Re: [ncl-talk] Masking
>
> Muhammad,
>
> Something is not right with data2.
>
> In order to correctly mask curvilinear data, your lat/lon arrays must be
> the same size as your data. If you look at the output from data2:
>
> Variable: data2
> Type: float
> Total Size: 108360 bytes
> 27090 values
> Number of Dimensions: 2
> Dimensions and sizes: [y | 129] x [x | 210]
> Coordinates:
> Number Of Attributes: 9
> lon2d : <ARRAY of 210 elements>
> lat2d : <ARRAY of 210 elements>
> stagger :
> description : TEMP at 2 M
>
> The data2 variable is 129 x 210, but the lat2d and lon2d are each 1D
> arrays of length 210. lat2d and lon2d should be 2D arrays that are the same
> size as data2.
>
> However, and this is important, attaching lat2d and lon2d attributes to
> the data variable is an NCL specific thing for plotting, regridding, and
> masking, and should NOT be used when writing new data to a file. This is
> needlessly expensive, because you likely already have the lat, lon arrays
> written as separate variables on the file.
>
> You can look at your original WRF output file to see how curvilinear data
> is generally represented on a NetCDF file. Here's an example:
>
> float T2 ( Time, south_north, west_east )
> FieldType : 104
> MemoryOrder : XY
> description : TEMP at 2 M
> units : K
> stagger :
> coordinates : XLONG XLAT XTIME
>
> Note the "coordinates" attribute attached T2. This is what tells you that
> the lat/lon arrays associated with T2 are called "XLONG" and "XLAT" on the
> file, which is what I used earlier to get the lat/lon for masking.
>
> The data2 variable longer has this coordinates attribute, so either NCO
> didn't do something correctly, or this variable was modified in some way.
>
> Please double-check that you are using NCO correctly, and don't write
> those lat2d / lon2d attributes. Also make sure that the coordinates
> attribute is included, and that you have the correct XLAT and XLONG
> variables written to the file. Once you do this, you can then mask your
> new data in the exact same way that you did for the original T2 variable.
>
> --Mary
>
>
>
>
>
> On Wed, Nov 8, 2017 at 12:35 AM, Muhammad Omer Mughal <
> m.mughal1 at postgrad.curtin.edu.au<mailto:m.mughal1 at postgrad.curtin.edu.au>>
> wrote:
> Hi Mary
>
> Thank you for the reply. The adjustment works perfectly to mask the data.
> I have another question related to this. When I use the wrfoutput file and
> use the commands lat2d I get this output
>
> Variable: data3
> Type: float
> Total Size: 108360 bytes
> 27090 values
> Number of Dimensions: 2
> Dimensions and sizes: [south_north | 129] x [west_east | 210]
> Coordinates:
> Number Of Attributes: 8
> lon2d : <ARRAY of 27090 elements>
> lat2d : <ARRAY of 27090 elements>
> coordinates : XLONG XLAT XTIME
> stagger :
> units : K
> description : TEMP at 2 M
> MemoryOrder : XY
> FieldType : 104
>
>
> When I use this to a sliced file from nco which has changed the attributes
> of the file I get the following output.
>
> Variable: data2
> Type: float
> Total Size: 108360 bytes
> 27090 values
> Number of Dimensions: 2
> Dimensions and sizes: [y | 129] x [x | 210]
> Coordinates:
> Number Of Attributes: 9
> lon2d : <ARRAY of 210 elements>
> lat2d : <ARRAY of 210 elements>
> stagger :
> description : TEMP at 2 M
> MemoryOrder : XY
> FieldType : 104
> coordinates : XLONG XLAT
> units : K
> time : 40320
>
> Also when I use the sliced file to mask the data I get the error
> (0) shapefile_mask_data:Error: not a valid rectilinear,curvilinear, or
> unstructured grid.
>
> I will appreciate your help
> Muhammad Omer Mughal
> MSc BSc Mechanical Engineering
> PhD Research Scholar
> Remote Sensing and Satellite Research Group
> Department of Imaging and Applied Physics
> Curtin University
>
> Curtin University
> Tel | +61 8 9266 7962
> Fax | +61 8 9266 2377
> Mobile | 0470 237 525
>
> Email | m.mughal1 at postgrad.curtin.edu.au<mailto:m.lynch at curtin.edu.au>
> Web | http://curtin.edu.au<http://curtin.edu.au/>
>
> Curtin University is a trademark of Curtin University of Technology.
> CRICOS Provider Code 00301J (WA), 02637B (NSW)
>
>
> ________________________________
> From: Mary Haley <haley at ucar.edu<mailto:haley at ucar.edu>>
> Sent: Wednesday, 8 November 2017 12:42:12 AM
> To: Muhammad Omer Mughal
> Cc: ncl-talk at ucar.edu<mailto:ncl-talk at ucar.edu>
> Subject: Re: [ncl-talk] Masking
>
> Muhammad,
>
> It looks like you are trying to use a rectilinear lat/lon grid for masking
> rather than the original lat/lon data associated with your "T2" variable.
>
> If your file is indeed a WRF output file, then you should use the
> XLAT/XLONG variables on the file, rather than trying to create your own
> lat/lon arrays.
>
> Your code would look something like this:
>
> nt = 0
> data
>
> = wrf_user_getvar(a, "T2",
>
> nt)
> data at lat2d = wrf_user_getvar(a, "lat", nt)
> data at lon2d = wrf_user_getvar(a,
>
> "lon",
>
> nt)
> data at _FillValue = -9999
>
> opt = True
> opt at debug = True
> data_mask = shapefile_mask_data(data,"Shapefile/malay.shp",opt)
>
>
> For an example that masks a WRF output file, see example
> shapefile_14_mask.ncl
>
> at:
>
> http://www.ncl.ucar.edu/Applications/shapefiles.shtml#ex14
>
> Note that this example is specifying "Ohio" as the specific area in the
> shapefile to mask against, but if you simply want to use the whole
> shapefile, then you do not need to set "shape_var" or "shape_names" as this
> script is doing.
>
> --Mary
>
> On Tue, Nov 7, 2017 at 6:48 AM, Muhammad Omer Mughal <
> m.mughal1 at postgrad.curtin.edu.au<mailto:m.mughal1 at postgrad.curtin.edu.au>>
> wrote:
> Hi
> I am using a WRF output file and I did the following changes "ONLY" to the
> code[ https://www.ncl.ucar.edu/Applications/Scripts/mask_9.ncl] . I want
> to mask the temperature in Malaysia using the shape file malay.shp. Can
> you kindly tell me what else do I need to change as currently I am not
> able to see the data masked within Malaysia. I will be grateful for your
> help
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
>
> begin
> MASK_INSIDE = True
>
> a = addfile("./UTC_0.nc","r")
> lat=a->XLAT
> lon=a->XLONG
>
> minlat=min(lat)
> maxlat=max(lat)
> minlon=min(lon)
> maxlon=max(lon)
>
> nlat=129
> nlon=210
>
> lat1d = fspan(minlat,maxlat,nlat)
> lon1d = fspan(minlon,maxlon,nlon)
>
> data = wrf_user_getvar(a, "T2", 0)
> data at _FillValue = -9999
>
>
>
> lat1d = fspan(minlat,maxlat,nlat)
> lon1d = fspan(minlon,maxlon,nlon)
> lat1d at units = "degrees_north"
> lon1d at units = "degrees_east"
>
> ;---Attach lat/lon coordinate array information.
> data!0 = "lat"
> data!1 = "lon"
> data&lat = lat1d
> data&lon = lon1d
>
>
>
> f = addfile("Shapefile/malay.shp", "r")
> mrb_lon = f->x
> mrb_lat = f->y
> nmrb = dimsizes(mrb_lon)
>
> min_mrb_lat = min(mrb_lat)
> max_mrb_lat = max(mrb_lat)
> min_mrb_lon = min(mrb_lon)
> max_mrb_lon = max(mrb_lon)
>
>
> www.ncl.ucar.edu<https://www.ncl.ucar.edu/Applications/Scripts/mask_9.ncl>
> www.ncl.ucar.edu<http://www.ncl.ucar.edu>
> ;***** ; mask_9.ncl ; ; Concepts illustrated: ; - Drawing the Mississippi
> River Basin using data from a shapefile ; - Masking a data array based on a
> geographical ...
>
>
>
> www.ncl.ucar.edu<https://www.ncl.ucar.edu/Applications/Scripts/mask_9.ncl>
> www.ncl.ucar.edu<http://www.ncl.ucar.edu>
> ;***** ; mask_9.ncl ; ; Concepts illustrated: ; - Drawing the Mississippi
> River Basin using data from a shapefile ; - Masking a data array based on a
> geographical ...
>
>
>
>
> Muhammad Omer Mughal
> MSc BSc Mechanical Engineering
> PhD Research Scholar
> Remote Sensing and Satellite Research Group
> Department of Imaging and Applied Physics
> Curtin University
>
> Curtin University
> Tel | +61 8 9266 7962
> Fax | +61 8 9266 2377
> Mobile | 0470 237 525
>
> Email | m.mughal1 at postgrad.curtin.edu.au<mailto:m.lynch at curtin.edu.au>
> Web | http://curtin.edu.au<http://curtin.edu.au/>
>
> Curtin University is a trademark of Curtin University of Technology.
> CRICOS Provider Code 00301J (WA), 02637B (NSW)
>
>
>
> _______________________________________________
> ncl-talk mailing list
> ncl-talk at ucar.edu<mailto: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/20171109/52c13b43/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: wrf_gsn_5.ncl
Type: application/octet-stream
Size: 4692 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0006.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: wrf_mask.png
Type: image/png
Size: 317853 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0001.png>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: wrf_mask.ncl
Type: application/octet-stream
Size: 3152 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0007.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: mrb.dbf
Type: application/octet-stream
Size: 71 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0008.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: mrb.prj
Type: application/octet-stream
Size: 167 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0009.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: mrb.shp
Type: application/octet-stream
Size: 1417196 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0010.obj>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: mrb.shx
Type: application/octet-stream
Size: 108 bytes
Desc: not available
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20171109/52c13b43/attachment-0011.obj>
More information about the ncl-talk
mailing list