## Documentation Center |

On this page… |
---|

Raster geodata consists of georeferenced data grids and images
that in the MATLAB^{®} workspace are stored as matrices or objects.
While raster geodata looks like any other matrix of real numbers,
what sets it apart is that it is georeferenced, either to the globe
or to a specified map projection, so that each pixel of data occupies
a known patch of territory on the planet.

Whether a raster geodata set covers the entire planet or not, its placement and resolution must be specified. This additional information can be supplied in the form of a referencing object, a referencing matrix, or a referencing vector.

A referencing object is an instance of the `map.rasterref.GeographicRasterReference` class,
for raster data referenced to a geographic latitude-longitude system,
or the `map.rasterref.MapRasterReference` class, for raster data referenced to a planar
(projected) map coordinate system. A spatial referencing object encapsulates
the relationship between a geographic or planar coordinate system
and a system of intrinsic coordinates anchored to the columns and
rows of a 2-D spatially referenced raster grid or image. Unlike the
older referencing matrix and vector representations (described below),
a referencing object is self-documenting, providing a rich set of
properties to describe both the intrinsic and extrinsic geometry.
The use of referencing objects is preferred, but referencing matrices
and vectors continue to be supported for the purpose of compatibility.

A referencing matrix is a 3-by-2 matrix of doubles that describes
the scaling, orientation, and placement of the data grid on the globe.
For a given referencing matrix, `R`, one of the following
relations holds between rows and columns and coordinates (depending
on whether the grid is based on map coordinates or geographic coordinates,
respectively):

[x y] = [row col 1] * R, or [long lat] = [row col 1] * R

For additional details about and examples of using referencing
matrices, see the reference page for `makerefmat`.

In many instances (when the data grid or image is based on latitude
and longitude and is aligned with the geographic graticule), a referencing
matrix has more degrees of freedom than the data requires. In such
cases, you may encounter a more compact representation, a three-element *referencing
vector*. A referencing vector defines the pixel size and
northwest origin for a regular, rectangular data grid:

refvec = [cells-per-degree north-lat west-lon]

In MAT-files, this variable is often called `refvec` (or `maplegend`).
The first element, cells-per-degree, describes the angular extent
of each grid cell (e.g., if each cell covers five degrees of latitude
and longitude, cells-per-degree would be specified as `0.2`).
Note that if the latitude extent of cells differs from their longitude
extent you cannot use a referencing vector, and instead must specify
a referencing object or matrix. The second element, north-lat, specifies
the northern limit of the data grid (as a latitude), and the third
element, west-lon, specifies the western extent of the data grid (as
a longitude). In other words, north-lat, west-lon is the northwest
corner of the data grid. Note, however, that cell (1,1) is always
in the southwest corner of the grid. This need not be the case for
grids or images described by referencing objects or matrices.

All regular data grids require a referencing object, matrix, or vector, even if they cover the entire planet. Geolocated data grids do not, as they explicitly identify the geographic coordinates of all rows and columns. For details on geolocated grids, see Geolocated Data Grids.

Regular data grids are rectangular, non-sparse, matrices of
class `double`.

Imagine an extremely coarse map of the world in which each cell represents 60º. Such a map matrix would be 3-by-6.

miniZ = [1 2 3 4 5 6; 7 8 9 10 11 12; 13 14 15 16 17 18];

Now make a referencing object:

miniR = georasterref('RasterSize', size(miniZ), ... 'Latlim', [-90 90], 'Lonlim', [-180 180])

Your output appears like this:

miniR = GeographicCellsReference with properties: LatitudeLimits: [-90 90] LongitudeLimits: [-180 180] RasterSize: [3 6] RasterInterpretation: 'cells' ColumnsStartFrom: 'south' RowsStartFrom: 'west' CellExtentInLatitude: 60 CellExtentInLongitude: 60 RasterExtentInLatitude: 180 RasterExtentInLongitude: 360 XIntrinsicLimits: [0.5 6.5] YIntrinsicLimits: [0.5 3.5] CoordinateSystemType: 'geographic' AngleUnit: 'degree'

Set up an equidistant cylindrical map projection:

figure('Color','white') ax = axesm('MapProjection', 'eqdcylin'); axis off setm(ax,'GLineStyle','-', 'Grid','on','Frame','on')

Draw a graticule with parallel and meridian labels at 60º intervals:

setm(ax, 'MlabelLocation', 60, 'PlabelLocation',[-30 30],... 'MLabelParallel','north', 'MeridianLabel','on',... 'ParallelLabel','on','MlineLocation',60,... 'PlineLocation',[-30 30])

Map the data using

`geoshow`and display with a color ramp and legend:geoshow(miniZ, miniR, 'DisplayType', 'texturemap'); colormap('autumn') colorbar

Note that the first row of the matrix is displayed at the bottom of the map, while the last row is displayed at the top.

The latitude and longitude limits of a regular grid are among the most important properties of its referencing object. Given an older data grid that includes a referencing vector but not a referencing object, a good way to check the limits is to convert the referencing vector to a geographic raster reference object, as in the following simple exercise:

Load the Korea 5-arc-minute elevation grid and inspect the referencing vector,

`refvec`:korea = load('korea','map','refvec')

korea = map: [180x240 double] refvec: [12 45 115]

The referencing vector,

`korea.refvec`, indicates that there are 12 cells per angular degree. This horizontal resolution is 5 times finer than that of the topo grid, which is one cell per degree. It also indicates that the northwest corner is at 45 degrees North, 115 degrees East, but it needs to be combined with the size of the data grid before the southern and eastern limits can be determined.Convert the referencing vector to a geographic raster reference object (providing a size vector):

`korea.R = refvecToGeoRasterReference(korea.refvec, ... size(korea.map))`

korea = map: [180x240 double] refvec: [12 45 115] R: [1x1 map.rasterref.GeographicCellsReference]

Examine the

`Latlim`and`Lonlim`properties:korea.R.LatitudeLimits korea.R.LongitudeLimits

ans = 30 45 ans = 115 135

You can access and manipulate gridded geodata using either geographic
or intrinsic raster coordinates. Use the `russia.mat` file
to explore this. As was demonstrated above, the north, south, east,
and west limits of the mapped area can be determined as follows:

russia = load('russia','map','refvec'); R = refvecToGeoRasterReference(russia.refvec, size(russia.map)); R.LatitudeLimits R.LongitudeLimits

ans = 35 80 ans = 15 190

Display a map of Russia:

figure('Color','white') worldmap(R.LatitudeLimits,R.LongitudeLimits) cmap = jet(4); geoshow(russia.map,cmap,R)

The `map.rasterref.GeographicCellsReference.intrinsicToGeographic` method
can be used to retrieve the geographic coordinates at the center of
a given grid cell. For example, consider the cell in row 23, column
79. In intrinsic raster coordinates, the center of this cell is located
at:

xIntrinsic = 79; yIntrinsic = 23;

This corresponds to the following location in latitude-longitude,
obtained via the `intrinsicToGeographic` method:

[lat, lon] = intrinsicToGeographic(R, xIntrinsic, yIntrinsic)

Your output appears like this:

lat = 39.5000 lon = 30.7000

The `geographicToIntrinsic` method does the
reverse, converting from latitude-longitude to intrinsic *x* and *y*:

[xIntrinsic, yIntrinsic] = geographicToIntrinsic(R, lat, lon)

Your output appears as follows:

xIntrinsic = 79 yIntrinsic = 23

Finally, if you know the latitude and longitude limits of a region, you can calculate the required matrix size and an appropriate referencing object for any desired map resolution and scale. However, before making a large, memory-taxing data grid, you should first determine what its size will be. For a map of the continental U.S. at a scale of 10 cells per degree, do the following:

Specify latitude limits of 25ºN to 50ºN and longitudes from 60ºW to 130ºW:

latlim = [ 25 50]; lonlim = [-130 -60];

Specify 10 cells per degree, and compute the number of rows and columns that this cell density implies:

cellsPerDegree = 10; numRows = cellsPerDegree * diff(latlim); numCols = cellsPerDegree * diff(lonlim);

Construct a referencing object and verify that its size is a reasonable 250 by 700 cells:

R = georasterref('RasterSize',[numRows numCols], ... 'Latlim', latlim, 'Lonlim', lonlim); R.RasterSize

Your output appears like this:

ans = 250 700

In addition to regular data grids, the toolbox provides another
format for geodata: *geolocated data grids*. These
multivariate data sets can be displayed, and their values and coordinates
can be queried, but unfortunately much of the functionality supporting
regular data grids is not available for geolocated data grids.

The examples thus far have shown maps that covered simple, regular quadrangles, that is, geographically rectangular and aligned with parallels and meridians. Geolocated data grids, in addition to these rectangular orientations, can have other shapes as well.

To define a geolocated data grid, you must define three variables:

A matrix of indices or values associated with the mapped region

A matrix giving cell-by-cell latitude coordinates

A matrix giving cell-by-cell longitude coordinates

The following exercise demonstrates this data representation:

Load the MAT-file example of an irregularly shaped geolocated data grid called

`mapmtx`:load mapmtx whos Name Size Bytes Class Attributes description 1x54 108 char lg1 50x50 20000 double lg2 50x50 20000 double lt1 50x50 20000 double lt2 50x50 20000 double map1 50x50 20000 double map2 50x50 20000 double source 1x43 86 char

Two geolocated data grids are in this data set, each requiring three variables. The values contained in

`map1`correspond to the latitude and longitude coordinates, respectively, in`lt1`and`lg1`. Notice that all three matrices are the same size. Similarly,`map2`,`lt2`, and`lg2`together form a second geolocated data grid. These data sets were extracted from the`topo`data grid shown in previous examples. Neither of these maps is regular, because their columns do not run north to south.To see their geography, display the grids one after another:

close all axesm mercator gridm on framem on h1=surfm(lt1,lg1,map1); h2=surfm(lt2,lg2,map2);

Showing coastlines will help to orient you to these skewed grids:

load coast plotm(lat,long,'r')

Notice that neither

`topo`matrix is a regular rectangle. One looks like a diamond geographically, the other like a trapezoid. The trapezoid is displayed in two pieces because it crosses the edge of the map. These shapes can be thought of as the geographic organization of the data, just as rectangles are for regular data grids. But, just as for regular data grids, this organizational logic does not mean that displays of these maps are necessarily a specific shape.Now change the view to a polyconic projection with an origin at 0ºN, 90ºE:

setm(gca,'MapProjection','polycon', 'Origin',[0 90 0])

As the polyconic projection is limited to a 150º range in longitude, those portions of the maps outside this region are automatically trimmed.

Mapping Toolbox™ software supports three different interpretations of geolocated data grids:

First, a map matrix having the same number of rows and columns as the latitude and longitude coordinate matrices represents the values of the map data at the corresponding geographic points (centers of data cells).

Next, a map matrix having one fewer row and one fewer column than the geographic coordinate matrices represents the values of the map data within the area formed by the four adjacent latitudes and longitudes.

Finally, if the latitude and longitude matrices have smaller dimensions than the map matrix, you can interpret them as describing a coarser

*graticule*, or mesh of latitude and longitude cells, into which the blocks of map data are warped.

This section discusses the first two interpretations of geolocated data grids. For more information on the use of graticules, see The Map Grid.

**Type 1: Values associated with the upper left grid coordinate. **As an example of the first interpretation, consider a 4-by-4
map matrix whose cell size is 30-by-30 degrees, along with its corresponding
4-by-4 latitude and longitude matrices:

Z = [ ... 1 2 3 4; ... 5 6 7 8; ... 9 10 11 12; ... 13 14 15 16]; lat = [ ... 30 30 30 30; ... 0 0 0 0; ... -30 -30 -30 -30; ... -60 -60 -60 -60]; lon = [ ... 0 30 60 90;... 0 30 60 90;... 0 30 60 90;... 0 30 60 90];

Display the geolocated data grid with the values of `map` shown
at the associated latitudes and longitudes:

figure('Color','white','Colormap',autumn(64)) set(gcf,'Renderer','zbuffer') axesm('pcarree','Grid','on','Frame','on',... 'PLineLocation',30,'PLabelLocation',30) mlabel; plabel; axis off; tightmap h = geoshow(lat,lon,Z,'DisplayType','surface'); set(h,'ZData',zeros(size(Z))) ht = textm(lat(:),lon(:),num2str(Z(:)), ... 'Color','blue','FontSize',14); colorbar('southoutside')

Notice that only 9 of the 16 total cells are displayed. The
value displayed for each cell is the value at the upper left corner
of that cell, whose coordinates are given by the corresponding `lat` and `lon` elements.
By convention, the last row and column of the map matrix are not displayed,
although they exist in the `CData` property of
the surface object.

**Type 2: Values centered within four adjacent coordinates. **For the second interpretation, consider a 3-by-3 map matrix
with the same `lat` and `lon` variables:

delete(h) delete(ht) Z3by3 = [ ... 1 2 3; ... 4 5 6; ... 7 8 9]; h = geoshow(lat,lon,Z3by3,'DisplayType','texturemap'); tlat = [ ... 15 15 15; ... -15 -15 -15; ... -45 -45 -45]; tlon = [ ... 15 45 75; ... 15 45 75; ... 15 45 74]; textm(tlat(:),tlon(:),num2str(Z3by3(:)), ... 'Color','blue','FontSize',14)

Display a surface plot of the map matrix, with the values of `map` shown
at the center of the associated cells:

All the map data is displayed for this geolocated data grid. The value of each cell is the value at the center of the cell, and the latitudes and longitudes in the coordinate matrices are the boundaries for the cells.

**Ordering of Cells. **You may have noticed that the first row of the matrix is displayed
as the top of the map, whereas for a regular data grid, the opposite
was true: the first row corresponded to the bottom of the map. This
difference is entirely due to how the `lat` and `lon` matrices
are ordered. In a geolocated data grid, the order of values in the
two coordinate matrices determines the arrangement of the displayed
values.

**Transforming Regular to Geolocated Grids. **When required, a regular data grid can be transformed into a
geolocated data grid. This simply requires that a pair of coordinates
matrices be computed at the desired spatial resolution from the regular
grid. Do this with the `meshgrat` function,
as follows:

load topo [lat,lon] = meshgrat(topo,topolegend); Name Size Bytes Class Attributes lat 180x360 518400 double lon 180x360 518400 double topo 180x360 518400 double topolegend 1x3 24 double topomap1 64x3 1536 double topomap2 128x3 3072 double

**Transforming Geolocated to Regular Grids. **Conversely, a regular data grid can also be constructed from
a geolocated data grid. The coordinates and values can be embedded
in a new regular data grid. The function that performs this conversion
is `geoloc2grid`; it takes a geolocated data grid
and a cell size as inputs.

Was this topic helpful?