Description:
MATLAB's Mapping Toolbox only provides a
"geotiffread" function, but it does not have a
"geotiffwrite" function. This is the MATLAB program to write a 2D or
3D array to a single or multi-band GeoTIFF file, where data can be either 1-bit
monochrome data (i.e. binary or bilevel), 8-bit unsigned integer, 16-bit signed
/ unsigned integer, 32-bit signed integer, or 32-bit floating point.
This program is based on GeoTIFF Format Specification under:
http://www.remotesensing.org/geotiff/spec/geotiffhome.html
or http://download.osgeo.org/geotiff/spec
It does not need MATLAB's Mapping Toolbox, or any other
library.
The motivation to create such a program is because I want to
convert the free satelite data from NASA & USGS to the GARMIN 3D topo map
for my Garmin GPS receiver. DEM2TOPO (http://people.uleth.ca/~brad.gom/dem2topo)
is a wonderful program that can nicely extract contours from the elevation
data. However, only 1201x1201 pixels data are supported. For USGS24K DEM data,
1-arc-second SRTM1 data, etc, DEM2TOPO cannot process them directly, which must
be converted into GeoTIFF format first (see "Known Issues" in the
above site). Although there is a free Windows program "3DEM" that can
convert USGS24K DEM & SRTM1 data into GeoTIFF format, sometimes you will
still need to interpolate those DEM data for void or missing elevation points.
If you just want to write DEM data to a GeoTIFF file, all
you have to provide is a bounding box (2 latitude & 2 longitude) and a 2D
or 3D array. Otherwise, you can provide proper GeoTIFF fields to the option
argument. To make your life easier, a GUI program “make_option” is also
provided. So you can select the proper GeoTIFF fields to generate the option
argument.
I recommend that you use "make_option" program to
generate option argument instead of constructing it manually.
Usage:
[option bbox] =
make_option([bbox]);
geotiffwrite(filename,
bbox, image, [bit_depth, option]);
Argument |
Value |
|
filename |
GeoTIFF file name to be written. |
|
bbox |
Bounding box with West/East of longitude and South/North
of latitude. The latitude & longitude must be in Decimal Degree (DD)
format, and they must be in a 2x2 array with the order of: [ longitude_West,
latitude_South; longitude_East,
latitude_North ] Note: If GTModelTypeGeoKey is not ModelTypeGeographic, bbox will be reset to empty. |
|
image |
2D or 3D array. Although any orientation can be specified
through option.Orientation, it is simpler to just use the following default
orientation: First row at North edge; last
row at South edge; First column at West edge; last
column at East edge; Note: “geotiffwrite” does not support full-color RGB image, and the 3rd dimension only indicates different bands. Therefore, when the 3rd dimension of a 3D array equals to 3, it only represents a 3-band GeoTIFF image, not a full-color RGB image. |
|
bit_depth |
(optional) Number of bits to represent a data point (i.e.
bits per sample): 1 for 1-bit monochrome data
(i.e. binary or bilevel) 8 for 8 bit unsigned integer 16 for 16 bit signed integer 32 for 32 bit floating point
(default) -16 for 16-bit unsigned integer -32 for 32-bit signed integer Note: If you set values to option.ColorMap, 16 will be 16-bit unsigned integer, and you cannot write in 32-bit depth or 1-bit depth. |
|
option |
(optional) A structure of GeoTIFF fields: First, please run: [option bbox] = make_option( [bbox] );
This will include most option fields, except the following: option.NaN By default, NaN value (void data point) in the image variable will be kept untouched. However, if you decide to replace it with other value (e.g. -32768), you can define it here. option.ColorMap
[Nx3 array] Where N = 2^bit_depth, and bit_depth can only be 8 or 16. Values in option.ColorMap range from 0 to 65535, and value 0 in image points to the first row in the ColorMap. Columns in option.ColorMap are [R G B]. Black is represented by [0 0 0], and white is represented by [65535 65535 65535]. option.Orientation = [1] | 2 | 3 | 4 | 5 | 6 | 7 | 8 1: Row from Top, Col from Left;
(default) 2: Row from Top, Col from Right; 3: Row from Bottom, Col from
Right; 4: Row from Bottom, Col from
Left; 5: Row from Left, Col from Top; 6: Row from Right, Col from Top; 7: Row from Right, Col from
Bottom; 8: Row from Left, Col from
Bottom; option.ModelTransformationTag [4x4 array] In most case, ModelPixelScaleTag
and ModelTiepointTag from "make_option.m" program is sufficient to
transform from raster to model space. If raster image requires rotation or
shearing, you can use this ModelTransformationTag to define a 4x4 affine
transformation matrix. |
Examples:
- Although not illustrated, it would be easier to use
"make_option" for option argument
Data from UTM projected coordinate systems:
ftp://ftp.remotesensing.org/pub/geotiff/samples/spot/chicago/UTM2GTIF.TIF
or http://download.osgeo.org/geotiff/samples/spot/chicago/UTM2GTIF.TIF
img =
imread('UTM2GTIF.TIF');
clear option
option.GTModelTypeGeoKey
= 1;
option.ModelPixelScaleTag
= [10;10;0];
option.ModelTiepointTag
= [0;0;0;444650;4640510;0];
option.GeogEllipsoidGeoKey
= 7008;
option.ProjectedCSTypeGeoKey
= 26716;
option.GTCitationGeoKey
= 'UTM Zone 16N NAD27"';
option.GeogCitationGeoKey
= 'Clarke, 1866 by Default';
geotiffwrite('UTM2GTIFX.TIF',[],img,8,option)
Data from colormap indexed scan image:
ftp://ftp.remotesensing.org/pub/geotiff/samples/usgs/o41078a5.tif
or http://download.osgeo.org/geotiff/samples/usgs/o41078a5.tif
[img cmap] =
imread('o41078a5.tif');
clear option
option.GTModelTypeGeoKey
= 1;
option.ModelPixelScaleTag
= [2.4384;2.4384;0];
option.ModelTiepointTag
= [0;0;0;698753.304798;4556059.506392;0];
option.ProjectedCSTypeGeoKey
= 32617;
option.PCSCitationGeoKey
= 'UTM Zone 17 N with WGS84';
option.ColorMap =
65535*cmap;
geotiffwrite('o41078a5X.tif',[],img,8,option)
DEM data:
ftp://ftp.remotesensing.org/pub/geotiff/samples/usgs/mlatlon.tif
or http://download.osgeo.org/geotiff/samples/usgs/mlatlon.tif
img =
imread('mlatlon.tif');
clear option
option.GeographicTypeGeoKey
= 4267;
option.GTCitationGeoKey
= 'Geographic Model/GeoTIFF 1.0';
option.GeogCitationGeoKey
= 'NAD 27 Datum';
geotiffwrite('mlatlonX.tif',[-122.6261,37.4531;-122.0777,38],img,16,option)
Even simpler:
http://dds.cr.usgs.gov/srtm/version2_1/SRTM3/North_America/N46W082.hgt.zip
srtm =
load_srtm('N46W082.hgt');
geotiffwrite('N46W082.tif',[-82,46;-81,47],flipud(srtm.dat))
6-band GeoTIFF data is available under:
http://burnseverity.cr.usgs.gov/fires/2009/yuch/padd09a_postm.zip
Versions:
1. 3D array
is supported to write Multi-band GeoTIFF data.
2. “make_option”
program is created to provide a GUI window to generate option argument for
“geotiffwrite” program.
Download:
Please download “geotiffwrite” and “make_option” proggrams
under:
Feedback:
Please kindly provide your feedbacks, ratings, or comments
under:
http://www.mathworks.com/matlabcentral/fileexchange/27959#feedback