Points, Pixels, Grids, and Cells

Terms of Reference

Introduction

Satellite remote sensing has made possible the collection of data over large areas of the Earth. These data are often stored in "grids." Grids are an efficient means of storing data because the location of a value within the grid is implicit--it is not explicitly stored in the grid. The location is also constant, which makes it easy to compare data from different sensors or different time periods.

The purpose of this document is to give a brief introduction to mapping, gridding, and the associated terminology used in the mapx library at NSIDC. Along with comments in the source code, it also describes the format of the map and grid parameter files.

Maps

A map projection is a procedure or mathematical formula to transform a curved surface onto a plane. The curved surface is usually the surface of the Earth and the plane is what we call a "map." The derivation of a map transformation requires a mathematical representation of the surface of the Earth. The surface of constant geopotential referred to as mean sea level or the geoid is quite complicated and its mathematical form is correspondingly complex. When defining map projections then, it is common to approximate the surface of the earth with a simpler surface such as an ellipsoid or sphere.

The ellipsoid is flattened at the poles, relative to the sphere, by about 1 part in 300. For maps of large areas (continents or bigger) using the sphere introduces no significant error. When mapping smaller areas it is not uncommon to choose an ellipsoid specifically for each area; each area is mapped on an ellipsoid with a different radius and eccentricity.

Coordinates for points on the surface of the Earth are given in latitude and longitude. Geodetic latitude is the angle between a vertical ray at the point and the plane of the equator. Geocentric latitude is the angle between a ray from the center of the Earth to the point and the plane of the equator. On a sphere, the geodetic latitude is equivalent to the geocentric latitude. It's more common to use geodetic latitude. There is only one definition of longitude--the angle between a ray from the center of the Earth to the point and the plane of the Greenwich Meridian (an arbitrary reference point). Coordinates in the map plane (x,y) are defined in the usual manner with x horizontal and increasing to the right and y vertical and increasing upward.

Projections

Map transformations are referred to as projections because they are commonly visualized as geometric constructs. That is, we think of the map as a piece of paper wrapped around a globe with lines emanating from points on the globe and "projected" onto the paper. It's important to note that not all map projections can be constructed in this manner.

Map projections are specified as a set of equations giving x and y in terms of latitude and longitude. The manner in which we wrap the paper around the globe determines the form of the equations. Cylindrical projections wrap the paper around in a tube, conical projections use a cone, and azimuthal projections use a flat surface. Several parameters are used in these equations to specify exactly where the surface intersects the globe, the scale or size of the globe, and the translation and rotation of the map on the plane surface.

There is no one best map projection. Each projection has different properties and different "best" uses. Two of the most important characteristics of maps are whether they are conformal or equal-area. No map projection is both, and some are neither. If a map is conformal then angles within a small area are reproduced accurately. This means that shapes are preserved. A small circle on the globe will look like a small circle on the map. In formal terms, at any point on the map the scale h along a meridian of longitude is equal to the scale k along a parallel of latitude. If a map is equal-area then a small circle placed anywhere on the map will always cover the same amount of area on the globe and the product of h and k at any point is one.

For maps that are not equal-area, hk - 1 gives a measure of the areal distortion. For maps that are not conformal, k/h gives the aspect ratio, which is a measure of shape distortion. For example, on the Polar Stereographic map true at 70N (a conformal map) the areal distortion varies from -6% at the pole to +29% at 45N and goes up to +276% at the equator while, by definition, the aspect ratio remains 1:1 everywhere. The Azimuthal Equal-Area map over the pole has, of course, no areal distortion, while the aspect ratio varies from 1:1 at the pole to 1.17:1 at 45N and goes up to 2:1 at the equator.

Tissot's indicatrix combines the two previous measures. Tissot's indicatrix is a small ellipse which shows both the scale distortion and the maximum angular deformation at a particular point on the map. See Maling or Snyder for a complete discussion.

Grids

A grid is a rectangular array of points. Grids record regularly spaced samples over an area. When sampling over the surface of the earth a grid is determined by a map projection, a sampling interval, an origin, and the number of rows and columns. So called "lat/lon" or "equal angle" grids are equivalent to grids sampled on a Cylindrical Equidistant map projection.

A grid coordinate system is defined in the map plane with axes parallel to the rows and columns of the grid and units equal to the sampling interval. The grid sample locations are at the whole integer grid coordinate points. To conform to mathematical array conventions, the grid coordinates (r,s) start at (0,0) in the upper left corner with r increasing to the right, and s increasing downward (this also conforms to digital image processing conventions). A grid cell (i,j) is defined as the area between grid coordinates i-.5 and i+.5, and j-.5 and j+.5. To conform to rounding conventions, the lower bound is included in the cell while the upper bound is not (round up at .5). Note that the r coordinate corresponds to the grid column number j, and the s coordinate corresponds to the grid row number i.

Think of the grid as a map drawn on graph paper. To find a grid cell location, count graph lines over from the left and down from the top, then read the latitude and longitude off the map.

Binning

The process of sampling data on a grid is sometimes referred to as binning, especially when the sampling method is to average all data that falls into the grid cell. Other methods used to sample data are: take the closest value to the grid point (nearest neighbor) or interpolate between surrounding grid points. Also, binning doesn't need to be the average of multiple data points; it could be the minimum, the maximum, the median, or the latest value in the grid cell. The characteristics of the data should determine the sampling method.

Pictures

Pixels are analogous to grid cells. To display gridded data as a picture, define a transformation from grid coordinates to pixel coordinates and then sample the gridded data at the whole integer pixel coordinate points. The most common technique is to map grid cells onto pixels one for one. More sophisticated techniques, often referred to as "resampling", allow for scrolling, zooming, and rotating.

Parameters

Many programs we use at NSIDC depend on the "mapx" module to handle map projections. These include: GISMO, PSQ, regrid, gridomatic, and the EASE-Grid processing software. Most of these tools provide access to predefined maps, but, it is relatively easy to define your own map with a map projection parameters (.mpp) file.

Map Projection Parameters (.mpp) File Format

Each parameter is described by a "keyword: value" pair on a single line. Parameters can be followed by comments, which begin with a semi-colon. Parameters can appear in any order. The default value is used when the corresponding keyword is not found. Unrecognized keywords are silently ignored. This means that if you misspell a keyword you'll get the default value without warning. The following parameters define the map projection.

Keyword Format Default Description
Map Projection string required field projection name (see list below)
Map Reference Latitude decimal degrees required field reference latitude for map projection
Map Reference Longitude decimal degrees required field reference longitude for map projection
Map Second Reference Latitude decimal degrees none used by some projections
Map Rotation decimal degrees 0.0 rotation counter-clockwise
Map Scale float 1.0 map scale factor (radius units per map unit)
Map ISin Justify int 2=even justify flag for ISin map
Map Origin Latitude decimal degrees refrence latitude translated map origin latitude
Map Origin Longitude decimal degrees reference longitude translated map origin longitude
Map False Easting float 0.0 offset map x origin in map units
Map False Northing float 0.0 offset map y origin in map units
Map Eccentricity float Clark 1866 eccentricity of ellipsoid
Map Equatorial Radius float Clark 1866 ellipsoid or equivalent authallic sphere radius of ellipsoid or sphere

The following parameters are used only by programs that actually draw maps.

Keyword Format Default Description
Map Southern Bound decimal degrees 90.00S bottom of map, starting point for latitude graticule
Map Northern Bound decimal degrees 90.00N top of map
Map Western Bound decimal degrees 180.00W left side of map, starting point for longitude graticule
Map Eastern Bound decimal degrees 180.00E right side of map
Map Graticule Latitude Interval decimal degrees 30. graticule spacing
Map Graticule Longitude Interval decimal degrees 30. graticule spacing
Map Graticule Label Latitude decimal degrees 0.00N where to label meridians
Map Graticule Label Longitude decimal degrees 0.00E where to label parallels
Map CIL Detail Level integer 1 level of detail for coastlines, islands, and lakes
Map BDY Detail Level integer 0 level of detail for political boundaries
Map RIV Detail Level integer 0 level of detail for rivers

note: For Integerized Sinusoidal, set scale to size of each row.

Map Projection Names

• Albers Conic Equal-Area
• Azimuthal Equal-Area
• Azimuthal Equal-Area (ellipsoid)
• Cylindrical Equal-Area
• Cylindrical Equal-Area (ellipsoid)
• Cylindrical Equidistant
• Integerized Sinusoidal
• Interupted Homolosine Equal-Area
• Lambert Conic Conformal (ellipsoid)
• Mercator
• Mollweide
• Orthographic
• Polar Stereographic
• Sinusoidal

Keywords must be entered verbatim. Capitalization and white space matter. Unrecognized (or mispelled) keywords are silently ignored. Map names are not case sensitive, white space and dashes are ignored, and most word permutations are acceptable. For example, "Equal-Area Cylindrical" is equivalent to "CYLINDRICAL EQUALAREA".

The reference latitude and longitude specify the original location and orientation of the map projection. The origin of the rectangular map coordinate system can be rotated and translated with the rotation and origin keywords.

For the Cylindrical Equal-Area and Polar Stereographic projections, the second reference latitude specifies the latitude of "true" scale. Conceptually, this is where the projection plane intersects the surface. For the Albers Conic Equal-Area and Lambert Conic Conformal projections, the conic plane is secant to the surface between the reference latitude and the second reference latitude.

The equatorial radius and eccentricity specify the surface to be mapped (usually the Earth's surface). The scale and the equatorial radius determine the map units. For example, if the radius is in kilometers (as is the case with the default) and the scale is 1.0, then the map units will be kilometers. PostScript units are inches, so for PostScript output you would specify the scale in kilometers per inch. For a grid you can define the scale in terms of kilometers per grid cell. The gridding module "grids" also has parameters that specify grid cells per map unit. This allows the same .mpp file to be used for multiple nested grids.

The remaining parameters are used by programs that actually draw maps (e.g. psmap, mapenum). The southern, northern, western, and eastern bounds of the map are the starting point for the graticule and are used to speed up the search in the map outline data base. They cannot be counted on to clip the map accurately. All longitudes should be in the range -180 to +360. West to east should not span more than 360 degrees. West specifies the left side of the map and east the right, not necessarily the minimum and maximum longitudes.

The latitude and longitude graticule intervals specify the spacing between graticule lines. The graticule is the overlay of latitude and longitude lines on the map. By default, a parallel is drawn every 30 degrees starting at the bottom of the map (90.00S) and a meridian is drawn every 30 degrees starting at the left side of the map (180.00W). Label latitude and longitude specify the parallel and meridian along which to draw the graticule line labels.

The last three numbers are used by the database search routine to specify the level of detail for map lines. The level refers to the number and size of features that will be selected, as opposed to resolution of each feature. The higher the number, the more detail will be included. For example, level 1 will include only the largest features, such as whole continents, while level 2 would also include some smaller features like large lakes. In either case the continental outline will appear with the same amount of detail (the same jaggedness).

Grid Parameter Definition (.gpd) File Format

Grids are defined in grid parameter definition (.gpd) files.

Keyword Format Default Description
Grid Width integer required field number of columns in grid
Grid Height integer required field number of rows in grid
Grid Map Origin Column float 0.0 column coordinate of map origin
Grid Map Origin Row float 0.0 row coordinate of map origin
Grid Cells per Map Unit float 1.0 sets both column and row scale
Grid Columns per Map Unit float 1.0 sets column scale if not already set
Grid Rows per Map Unit float 1.0 sets row scale if not already set
Grid Map Units per Cell float 1.0 alternative method to set both column and row scale
Grid Map Units per Column float 1.0 sets column scale if not already set
Grid Map Units per Row float 1.0 sets row scale if not already set
Grid MPP File string none name of .mpp file, map projection parameters can also be specified in the same file with the grid definition parameters

Limitations

The mapx library is very flexible but it does impose certain restrictions on the kinds of grids that can be defined. All grids must be rectangular arrays. Rectangular arrays are easy to use with existing scientific visualization tools, like IDL, but they may be inefficient in terms of storage space for non-rectangular data sets. If some flag value (for example zero) is stored in the unused cells, then this wasted space can be compacted by any good data compression scheme. This same invalid data flag can be used within the grid as well, leading to further storage savings. Another restriction is that each grid must be based on a map projection defined in mapx. Many of the most useful projections have already been implemented and fortunately it is easy, given a mathematical formulation of the forward and inverse transformations, to add a new projection to mapx (see source code for instructions). The mathematical formulations of most common map projections can be found in the references below.

References

American Society of Photogrammetry. 1983. Manual of Remote Sensing. Robert N. Colwell, editor. Second edition. Falls Church, VA.

Gonzalez, Rafael C. and Paul Wintz. 1987.Digital Image Processing.Second edition. Reading, MA: Addison-Wesley.

Maling, D. H. 1992. Coordinate Systems and Map Projections. Second edition. Elmsford, N.Y.: Pergamon Press.

Snyder, John P. 1987.Map Projections, A Working Manual. U. S. Geological Survey Professional Paper 1395. Department of the Interior. Washington, D. C.