Re: smooth boundaries of a masked area

Hi Sebastian,

to mask data along a coastline you could you shapefiles. See shapefile
example mask_12.ncl:

http://ncl.ucar.edu/Applications/Scripts/mask_12.ncl
http://ncl.ucar.edu/Applications/Images/mask_12_lg.png

To get the data within a choosen sub-region and masking everything else
with _FillValue
you can use the ind(..) function. The following user defined function
data_mask(..)
generates such a mask data array:

;--------------------------------------------------------------------
undef("mask_data")
;--------------------------------------------------------------------
function mask_data(data,minlat,maxlat,minlon,maxlon)
;--------------------------------------------------------------------
local dims,lat,lon,nlat,nlon,lat1d,lon1d,nlatlon,ind_latlon,dims
begin
   dims = dimsizes(data) ;--
dimension of data
   lat = data&lat ;--
retrieve latitudes
   lon = data&lon ;--
retrieve longitudes
   nlat = dims(0) ;--
dimension of lat
   nlon = dims(1) ;--
dimension of lon

;-- to use the ind function we need to convert the dimesnions lat and lon
;-- and the data into a 1D array
   lat1d = ndtooned(conform_dims((/nlat,nlon/),data&lat,0)) ;--
convert to 1D array
   lon1d = ndtooned(conform_dims((/nlat,nlon/),data&lon,1)) ;--
convert to 1D array
   nlatlon = dimsizes(lat1d) ;--
dimension of 1D array
   data_mask_1d = new(nlatlon,integer) ;--
assign new integer array

;-- retrieve the indices of the lat/lon's inside the sub-region
   ind_latlon =
ind(lat1d.ge.minlat.and.lat1d.le.maxlat.and.lon1d.ge.minlon.and.lon1d.le.maxlon)

   if(.not.isatt(data,"_FillValue")) then
     data@_FillValue = default_fillvalue(typeof(data)) ;-- make
sure "data" has a missing value
   end if

   data_1d = ndtooned(data) ;--
convert data to 1D array
   data_mask = new(dimsizes(data_1d),typeof(data),data@_FillValue)
;-- create 1D masked data array

   data_mask_1d = 0 ;--
initialize array to 0

   do n=0,dimsizes(ind_latlon)-1
      nn = ind_latlon(n) ;-- get
index of point we're checking
if(lat1d(nn).lt.minlat.or.lat1d(nn).gt.maxlat.or.lon1d(nn).lt.minlon.or.lon1d(nn).gt.maxlon)

         continue
      else
         data_mask_1d(nn) = 1 ;-- its
inside the sub-region
      end if
   end do

   if(num(data_mask_1d.eq.1).gt.0) then
      data_mask = where(data_mask_1d.eq.1,data_1d,data_1d@_FillValue)
;-- get data or _FillValue
   else
      print("WARNING: no points inside sub-region (Lat:
"+minlat+","+maxlat+" Lon: "+minlon+","+maxlon+")")
      exit
   end if

;-- for plotting purposes convert the masked data back to a 2D data array
   data_mask_2d = onedtond(data_mask,(/nlat,nlon/)) ;--
convert back to a 2d data array
   data_mask_2d!0 = "lat" ;-- set
named dimension lat
   data_mask_2d!1 = "lon" ;-- set
named dimension lon
   data_mask_2d&lat = lat ;-- use lat
   data_mask_2d&lon = lon ;-- use lon
   data_mask_2d&lat@units = "degrees_north" ;-- set
lat units
   data_mask_2d&lon@units = "degrees_east" ;-- set
lon units

   return(data_mask_2d) ;-- return 2D data mask array
end

Hope this helps,
Karin

Am 19.05.14 10:08, schrieb Sebastian Milinski:
> Dear all,
>
> I would like to mask a region in a tripolar grid. Outside this region the values should be 0. In the transition zone (e.g. 3°) I want the values to decay (linearly or exponentially) to 0.
> My current approach is to create a weight file with values from 0 to 1 on a regular grid. I will then interpolate it to the tripolar grid and multiply it with the data.
>
> Ho can I create this transition zone in NCL? For straight boundaries in N/S or E/W direction I could introduce some latitude/longitude dependence.
> But what about more complex shapes? E.g. a masked region along a coastline?
>
> I thought this should be a common problem when applying a regionally confined forcing in a model. I would very much appreciate any hints.
>
> Regards
> Sebastian
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>

-- 
Dipl. Geophys. Karin Meier-Fleischer
Application Support, Visualization
Deutsches Klimarechenzentrum GmbH      E-Mail:   meier-fleischer@dkrz.de
Bundesstrasse 45a                      Internet: http://www.dkrz.de/
20146 Hamburg                          Phone:    +49 (0)40 460094 126
Germany                                Fax:      +49 (0)40 460094 270
Geschäftsführer: Prof. Dr. Thomas Ludwig
Sitz der Gesellschaft: Hamburg
Amtsgericht Hamburg HRB 39784

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon May 19 04:56:04 2014