geomaskdocumentation
Thegeomaskreturns true for locations within a specified geographic region.
Back to Climate Data Tools Contents
Contents
Syntax
mask = geomask(lat,lon,latv,lonv) mask = geomask(lat,lon,latv,lonv,'inclusive') [mask,coords] = geomask(...)
Description
mask = geomask(lat,lon,latv,lonv)returns a mask the size oflat和lonthat istruefor all points within the bounds given bylatv,lonv.
- scalarlatv,lonv: Iflatv,lonvare scalar values, the outputmaskwill betruefor only the pixel closest tolatv,lonv.
- two-element arrays: Iflatv,lonvare two-element arrays (e.g.,[40 50],[110 120]) the outputmaskwill be true for alllat,lonvalues within the geographic quadrangle bounded bylatv,lonv.
- polygon defined bylatv,lonv: Iflatv,lonvcontain more than two elements, the outputmaskwill be true for all elements of lat,lon within the polygon defined bylatv,lonv.
- polygons in cell format: Iflatv,lonvare cell arrays (as is common for multiple areas in a shapefile), the outputmaskis true for all elements within any of the polygons inlatv,lonv.
mask = geomask(lat,lon,latv,lonv,'inclusive')includeslat,lonpoints that are on the boundary or boundaries defined bylatv,lonv.
[mask,coords] = geomask(...)Iflatv,lonv标量,是可选的坐标outpu吗t is a structure containing the coordinates of the pixel in the mask. The coords structure includescoords.row和coords.colwhich are the row and column of of thelat,lon电网,coords.lat和coords.lon, which are the geographic location of the output grid cell.
Example 1: Geographic quadrangle
Sometimes you're interested in values that are withing a geographic quadrangle. One such quadrangle is the Nino 3.4 box is defined as (5°S to 5°N, 170°W to 120°W). Here's some sample sea surface temperature data
% Load some sample data:loadpacific_sst.mat% Plot the sample data:imagescn(lon,lat,mean(sst,3)) axisxyimagecmoceanthermalxlabel'longitude'ylabel'latitude'% Define the Nino 3.4 box:latv = [-5 -5 5 5 -5]; lonv = [-170 -120 -120 -170 -170];% Plot the Nino 3.4 box:holdonplot(lonv,latv,'k-','linewidth',2)
For a simple geographic quadrangle the easiest way to get a mask is to define the limits of the bounding box (similar toingeoquadif you have Matlab's Mapping Toolbox).
Start by converting thelat,lonarrays into a grid, then get the mask by specifing the limits of the region interest:
% Convert lat,lon into a grid:[Lon,Lat] = meshgrid(lon,lat);% Get a mask of grid cells within the quadrangle:mask = geomask(Lat,Lon,[-5 5],[-170 -120]);% Plot the mask:figure imagescn(lon,lat,mask) xlabel'longitude'ylabel'latitude'title'yellow indicates *true* grid cells (the Nino 3.4 box)'
Example 2: Nearest grid cell
Sometimes you need a time series near a single point of interest. For example, you might want to plot the SST time series near Honolulu (21.3 N, 157.8 W). Let's find the grid cell closest to Honolulu:
[honolulu,coords] = geomask(Lat,Lon,21.3,-157.8); figure imagescn(lon,lat,honolulu) xlabel'longitude'ylabel'latitude'holdonborders% plots national boundaries for context
A single-pixel mask may be of limited use, but for such masks,geomaskalso offers the row, column, and geographic locations of the nearest pixel, see:
coords
coords = struct with fields: row: 20 col: 12 lat: 21.5000 lon: -157.5000
As proof, below we plot the mean SST map with a red star over our desired Honolulu location and a blue square over the grid cell closest to Honolulu:
% Plot the mean SST for context:figure imagescn(lon,lat,mean(sst,3)) axisxycmoceanthermalxlabel'longitude'ylabel'latitude'% Plot Honolulu as a red star and its nearest grid cell as a blue square:holdonplot(-157.8,21.3,'rp') plot(coords.lon,coords.lat,'bs') text(-157.8,21.3,' Honolulu','color','r','vert','top') text(coords.lon,coords.lat,' nearest grid cell','color','b','vert','bottom')
And with the row and column information in thecoordsoutput, it's easy to get a time series for the Honolulu grid cell. Note thatsstis a 3D matrix, so we have to use thesqueezecommand to remove singleton dimensions. For clarity, I
sst_honolulu = sst (coords.row、coords.col:);% remove singleton dimensions:sst_honolulu = squeeze(sst_honolulu); figure plot(t,sst_honolulu) axistightboxoffylabel('Honolulu SST') datetick('x','keeplimits')
Example 3: Arbitrarily-shaped polygon(s)
If you want to make a mask of an arbitrary polygon or multiple polygons, you can use the standard Matlab functioninpolygon, or you can letgeomaskdo it for you. Withgeomask,latv和lonvcan be numeric arrays just like you would use withinpolygon, orlatv和lonvcan be cell arrays with polygons in each cell.
In this example, we have a global grid and we want to know which grid cells are within the borders of Latin American countries. CDT comes with sampleborderdata, and by going through the country names, you can manually pick out which countries correspond to Latin America. Below I've picked 20 Latin American countries which we'll use in to make the mask.
This example may take a couple of seconds to compute because the national outlines inborderdataare somewhat high resolution and it takes time forgeomaskto compare the grid cells to the outline of 41,500 border data vertices.
% For this sample 1 degree resolution grid:[Lat,Lon] = geogrid;% Load some national border data:B = load('borderdata.mat');% Use only Latin American countries:latv = B.lat([8 17 21 33 38 39 41 48 49 55 75 78 79 120 159 161 162 165 211 214]); lonv = B.lon([8 17 21 33 38 39 41 48 49 55 75 78 79 120 159 161 162 165 211 214]);% Find grid cells corresponding to Latin American countries:mask = geomask(Lat,Lon,latv,lonv);% Plot the mask:figure imagescn(Lon,Lat,mask) xlabel'longitude'ylabel'latitude'
With themaskfor Latin American countries, it is now possible to uselocalto get a time series, say, of area-averaged land surface temperatures in Latin American countries.
Author Info
Thegeomaskfunction and supporting documentation were written byChad A. Greeneof the University of Texas at Austin, Institute for Geophysics (UTIG).