2.3.8.4. vacumm.misc.grid.masking – Masking

Functions:

check_poly_islands() check_poly_straits() clip_shape() clip_shapes() convex_hull() create_polygon() d2m() envelop() erode_coast() get_avail_rate()

get_coast() get_coastal_indices() get_dist_to_coast() grid_envelop() grid_envelop_mask() mask2d() masked_polygon() merge_masks() plot_polygon()

polygon_mask() polygon_select() polygons() proj_shape() resol_mask() rsamp() t2uvmasks() uniq() zcompress()

Classes:

GetLakes indices() lakes() mask() ncells() ocean() plot() Lakes indices() lakes() mask() ncells() ocean() plot()

Inheritance diagram:  

Utilities to deal with masks and selections

See also

Tutorials: Les masques, Les polygones

GetLakes

alias of vacumm.misc.grid.masking.Lakes

class Lakes(mask, nmaxcells=None)

Bases: object

Find lakes in a 2D mask where 0 is water and 1 is land

Example:

>>> from vacumm.misc.grid import Lakes
>>> import numpy as N
>>> from pylab import pcolor,show,title
>>> # Build the mask
>>> mask=N.ones((500,500))
>>> mask[3:10,100:102]=0
>>> mask[103:110,200:210]=0
>>> mask[203:210,200:220]=0
>>> # Find the lakes
>>> lakes = Lakes(mask,nmaxcells=80)
>>> print 'Number of lakes:', len(lakes.indices())
Number of lakes: 2
>>> pcolor(lakes.mask(),shading=False)
>>> title('Lakes')
>>> show()
>>> first_two_lakes = lakes[:2]

indices(id=None, nbig=None)

Return a list of lake coordinates

• id: Select lake #<id>
• nbig: Consider the first nbig greatest lakes

lakes(id=None, nbig=None)

Return an array of same size as masks, but with points in lakes set to its lake id, and others to set to zero

• id: Select lake #<id>
• nbig: Consider the first nbig greatest lakes

mask(id=None, nbig=None, **kwargs)

Returns a boolean land/sea mask from lakes() (land is True)

• id: Select lake #<id>
• nmaxcells: Do not consider lakes with number of points greater than nmaxcells

Example:

>>> mask_lake2 = GetLakes(mask).mask(2)

ncells()

Number of cell point for each lake

ocean(mask=False)

Get the biggest lake or its mask (integer type)

• mask: If True, return the mask, not the lake [default: True]

plot(**kwargs)

Display the lakes

Keywords are passed to :func`~matplotlib.pyplot.pcolor`

check_poly_islands(mask, polys, offsetmin=0.85, offsetmax=1.5, dcell=2)

Check that islands as greater as a cell are taken into account

• mask: The initial mask. A cdms variable with X (longitude) and Y (latitude), or a tuple (lon, lat, mask).
• polys: Coastal polygons.
• offset: Islands whose area is > 1-offset and < 1+offset % of the mean cell area are checked [default: .95]
• dcell: dx and dy relative extension from the center of the cell in resolution units.

check_poly_straits(mask, polys, dcell=2, threshold=0.75)

Check that straits are opened by counting the number of polygon intersection in each cell and the area of water

• mask: The initial mask. A cdms variable with X (longitude) and Y (latitude).
• polys: Coastal polygons.
• dcell: dx and dy relative extension from the center of the cell in resolution units.
• threshold: Maximal fraction of cell area that must be covered of land (> .5) [default: .25]

clip_shape(shape, clip=None)

Clip a Point, a Polygon or a LineString shape using clipping polygon

Params:

• shape: A valid Point, a Polygon or a LineString instance.
• clip, optional: A clipping polygon.

Return:

A possible empty list of intersection shapes.

clip_shapes(shapes, clip=None)

Same as clip_shape() but applied to a list of shapes

convex_hull(xy, poly=False, method='delaunay')

Get the envelop of cloud of points

Params:

• xy: (x,y) or array of size (2,nxy) (see get_xy()).

• poly:

  • True: Return as Polygon instance.

  • False: Return two 1D arrays xpts,ypts

  • method:

    • "angles": Recursive scan of angles between points.
    • "delaunay": Use Delaunlay triangulation.

create_polygon(data, proj=False, mode='poly')

Create a simple Polygon instance using data

Param:

• data: Can be either a (N,2) or (2,N) array, a (xmin,ymin,xmax,ymax) list, tuple or array, or a Polygon.
• proj, optional: Geographical projection function.
• mode, optional: Output mode. "line" = LineString, "verts" = numpy vertices, else Polygon.

d2m(x, y)

envelop(*args, **kwargs)

Shortcut to convex_hull()

erode_coast(var, mask2d=None, copy=True, maxiter=10, corners=0.0, hardmask=None, **kwargs)

Erode coastal mask of var to fit to mask2d using 2D smoothing

Soothing is performed using shapiro2d(), in a convergence loop. Loop is ended if :

• the mask no longer changes,
• the number of iteration exceeds maxiter.

Warning

Must be used only for erodeing a thin border of the coast.

Params:

• var: Atleast-2D A MV2 variable with mask.
• mask2d, optional: Reference 2D mask.
• maxiter, optional: Maximal number of iteration for convergence loop. If equal to None, no check is performed.
• corners, optional: Weights of corners when calling shapiro2d().
• copy, optional: Modify the variable in place or copy it.
• hardmask, optional: Mask always applied after an erosion step.
• Other keywords are passed to shapiro2d().

Return:

• A MV2 variable

Tutorial:

Rognage de la côte

get_avail_rate(data, num=True, mode='rate')

Get the availability rate of a masked array along its forst dimension

Params:

• data: Mask or masked array.

• num, optional: Get result as pure numpy array or formatted as input array when possible (MV2 array).

• mode, optional:

  • “rate”: Output is between 0 and 1.
  • “pct”: Same but in percents from 0 to 100.
  • “count”: Count valid data.

get_coast(mask, land=True, b=True, borders=True, corners=True)

Get a mask integer of ocean coastal points from a 2D mask

• mask: Input mask with 0 on water and 1 on land.
• land: If True, coast is on land [default: True]
• corners: If True, consider borders as coastal points [default: True]
• borders: If True, consider corners as coastal points [default: True]

At each point, return 0 if not coast, else an interger ranging from 1 to (2**8-1) to describe the coast point:

128 4  64
8   +   2
16  1  32

get_coastal_indices(mask, coast=None, asiijj=False, **kwargs)

Get indices of ocean coastal points from a 2D mask

Params:

• mask: Input mask with 0 on water and 1 on land.
• boundary: If True, consider outside boundary as land [default: True]
• asiijj: Get results as II,JJ instead of [(j0,i0),(j1,i1)…]

get_dist_to_coast(grid, mask=None, proj=True)

Get the distance to coast on grid

Params:

• grid: (x,y), grid or variable with a grid.
• mask: Land/sea mask or variable with a mask. If not provided, gets mask from grid if it is a masked variable, or estimates it using polygon_mask() and a GHSSC shoreline resolution given by gshhs_autores().
• proj: Distance are computed by converting coordinates from degrees to meters.

Example:

>>> dist = get_dist_to_coast(sst.getGrid(), sst.mask)
>>> dist = get_dist_to_coast(sst)

See also:

get_coastal_indices()

grid_envelop(gg, centers=False, poly=True)

Return a polygon that encloses a grid

• gg: A cdms grid or a tuple of (lat,lon)

• centers:

  • True: The polygon is at the cell center.
  • False: The polygon is at the cell corners.

• poly:

  • True: Return as Polygon instance.
  • False: Return two 1D arrays xpts,ypts

grid_envelop_mask(ggi, ggo, poly=True, **kwargs)

Create a mask on output grid from the bounds of an input grid: points of output grid that are not within bounds of input grid are masked.

Params:

• ggi: Input grid or (lon,lat) or cdms variable.
• ggo: Output grid or (lon,lat) or cdms variable.
• Other keywords are passed to grid_envelop()

Returns:

A mask on output grid, where True means “outside bounds of input grid”.

mask2d(mask, method='and')

Convert a 3(or more)-D mask to 2D by performing compression on first axes

Note

Mask is False on ocean

• mask: At least 2D mask

• method: Compression method

  • "and": Only one point must be unmask to get the compressed dimension unmasked.
  • "or": All points must be unmask to get the compressed dimension unmasked.

masked_polygon(vv, polys, copy=0, reverse=False, **kwargs)

Mask a variable according to a polygon list

Params:

• vv: The MV2 array.
• polys, optional: Polygons (or something like that, see polygons()).
• copy, optional: Copy data?.
• reverse, optional: Reverse the mask?
• Other keywords are passed to polygon_mask().

merge_masks(varlist, copy=False, mode='max')

Merge the mask of a list of variable

plot_polygon(poly, ax=None, **kwargs)

Simply plot a polygon on the current plot using matplotlib.pyplot.plot()

Params:

• poly: A valid Polygon instance.
• Extra keyword are passed to plot function.

polygon_mask(gg, polys, mode='intersect', thresholds=[0.5, 0.75], ocean=False, fractions=0, yclean=True, premask=None, proj=False)

Create a mask on a regular or curvilinear grid according to a polygon list

Params:

• gg: A cdms grid or variable or a tuple of (xx,yy).

• polys: A list of polygons or GMT resolutions or Shapes instance like shorelines.

• mode, optional: Way to decide if a grid point is masked. Possible values are:

  • intersect, 1, area (default): Masked if land area fraction is > thresholds[0]. If more than one intersections, land area fraction must be > thresholds[1] to prevent masking straits.
  • else: Masked if grid point inside polygon.

• thresholds, optional: See intersect mode [default: [.5, .75]]

• ocean, optional: If True, keep only the ocean (= biggest lake).

• fractions: If True or 1, return the total fraction of cells covered by the polygons; if 2, returns the total area for each cell.

• proj, optional: Geographical projection function. See get_proj(). Use False to not convert coordinates to meters, and speed up the routine.

Result:

A numpy array of boolean (mask) with False on land.

polygon_select(xx, yy, polys, zz=None, mask=False)

Select unstructered points that are inside polygons

Params:

• xx: Positions along X.
• yy: Positions along Y.
• polys: Polygons.
• zz, optional: Values at theses positions.
• mask, optional: Return masked positions and values of True or 1,

Examples:

>>> xs, ys = polygon_select(x, y, [poly1, poly2])
>>> xs, ys, zs = polygon_select(x, y, [poly1, poly2], z)

>>> xmasked, ymasked, zmasked = polygon_select(x, y, polys, z, mask=True)
>>> print N.allclose(xs, xmasked.compressed())

>>> valid = polygon_select(x, y, polys, mask=2)
>>> print N.allclose(xs, x[valid])

Outputs:

Depends on mask

• False or 0: selected x,y or x,y,z
• True/1: the same but as masked arrays
• 2: boolean array of valid points (the inverse of the mask)

polygons(polys, proj=None, clip=None, shapetype=2, **kwargs)

Return a list of Polygon instances

Params:

• polys: A tuple or list of polygon proxies (see examples).
• shapetype: 1 = Polygons, 2=Polylines(=LineStrings) [default: 2]
• proj: Geographical projection to convert positions. No projection by default.
• clip, optional: Clip all polygons with this one.

Example:

>>> import numpy as N
>>> X = [0,1,1,0]
>>> Y = N.array([0,0,1,1])
>>> polygons( [(X,Y)]] )                                # from like grid bounds
>>> polygons( [zip(X,Y), [X,Y], N.array([X,Y])] )       # cloud
>>> polygons( [polygons(([X,Y],), polygons(([X,Y],)])   # from polygins
>>> polygons( [[min(X),min(Y),max(X),max(Y)],] )        # from bounds

proj_shape(shape, proj)

Project a Point, a Polygon or a LineString shape using a geographical projection

Params:

• shape: A Point, Polygon, or a LineString instance with coordinates in degrees.
• proj: A projection function of instance. See get_proj().

Return:

A similar instance with its coordinates converted to meters

resol_mask(grid, res=None, xres=None, yres=None, xrelres=None, yrelres=None, relres=None, scaler=None, compact=False, relmargin=0.05, nauto=20)

Create a mask based on resolution criteria for undersampling 2D data

Params:

• grid: A cdms array with a grid, a cdms grid or a tuple of axes.

• res, optional: Horizontal resolution of arrows

  (in both directions) for undersampling [default: None]. If 'auto', resolution is computed so as to have at max nauto arrow in along an axis. If it is a complex type, its imaginary part set the nauto parameter and res is set to 'auto'.

  Warning

  Use of res makes the supposition that X and y units are consistent.

• xres, optional: Same along X [default: res]

• yres, optional: Same along Y [default: res]

• relres, optional: Relative resolution

  (in both directions).

  • If > 0, = median(res)*relres.
  • If < -1, =``min(res)*abs(relres)``.
  • If < 0 and > -1, =max(res)*abs(relres)

• xrelres, optional: Same along X [default: relres]

• yrelres, optional: Same along Y [default: relres]

• scaler, optional: A callable object that transform X and Y coordinates. Transformed coordinates are used instead of original coordinates if resolutions are negatives. A typical scaler is for example a geographic projection that convert degrees to meters (like a Basemap instance).

• compact, optional: If no unsersampling is efective, returns False.

Note

A resolution value set to False implies no undersampling.

Example:

>>> mask = resol_mask((x2d, y2d), relres=2.5) # sampling=2.5
>>> mask = resol_mask(var.getGrid(), xres=2., scaler=mymap, yres=-100.) # 100m
>>> mask = resol_mask(var.getGrid(), res=20j) # at max 20 arrows per axis
>>> mask = resol_mask(var.getGrid(), res='auto', nauto=20) # equivalent

Note

It is usually better to regrid with a convenient method instead of unsersample because of aliasing.

rsamp(x, y, r, z=None, rmean=0.7, proj=False, rblock=3, getmask=False)

Undersample data points using a radius of proximity

• x: 1D x array.
• y: 1D Y array.
• r: Radius of proximity.
• z: 1D Z array.
• proj: Geographic projection instance to convert coordinates.
• rmean: Radius of averaging relative to r (0<rmean<1)
• rblock: Size of blocks relative to r for computation by blocks.

t2uvmasks(tmask, getu=True, getv=True)

Compoute masks at U and V points from a mask at T points on a C grid (True is land)

• tmask: A 2D mask.
• getu: Get the mask at U-points
• getv: Get the mask at V-points

Return umask,vmask OR umask OR vmask depending on getu and getv.

uniq(data)

Remove duplicates in data

Example:

>>> import numpy as N
>>> a = N.arange(20.).reshape((10,2))
>>> a[5] = a[1]
>>> b = uniq(a)

zcompress(z, *xy, **kwargs)

Compress 1D arrays according to the mask of the first one

• z: Reference (masked) array
• xy: Additional arrays to be compressed
• numpify: Force convertion to numpy array

Example:

>>> z, x, y = zcompress(z, x, y, numpify=True)