READHGT: Import/download NASA SRTM data files (.HGT)

Imports or download SRTM worldwide digital elevation model data files.

Description

This function imports .HGT "height" binary data files from NASA SRTM global digital elevation model of Earth land, corresponding to 1x1 degree tiles of 3-arc seconds resolution (SRTM3, around 90 m) and 1-arc second (SRTM1, around 30 m) for USA territory, and returns coordinates vectors latitude and longitude, and a matrix of elevation values.

The function includes also an automatic download of data from the USGS SRTM webserver, so indicating latitude and longitude is sufficient to get the data and instant map plot anywhere in the World.

Examples:

   readhgt(46:47,-123:-122,'merge')

plots a map of Seattle and surrounding volcanoes,

   X = readhgt(48,2,'plot');

plots a map of Paris (France) and returns downloaded SRTM data in structure X.

Type "doc readhgt" for syntax, examples and use.

• Required Products: MATLAB
• MATLAB release: MATLAB 6.5 (R13)

[Jeff]

Regarding my previous comment, I just now noticed you recently added the 'crop' argument update which is exactly what I needed. You are way ahead of me... Amazing what a little bit of reading can do :) Thanks.

[François Beauducel]

A user asked me how to extract an elevation profile from SRTM data... here is a little example (Mount Etna, Sicily):

X = readhgt(37,14:15,'merge','plot');
[x,y] = ginput(2)

then left-click on 2 different points: you will get 2 coordinates in x and y vectors. Then interpolate to create a profile vector:

xi = linspace(x(1),x(2));
yi = linspace(y(1),y(2));
zi = interp2(X.lon,X.lat,double(X.z),xi,yi);

note you need to use the DOUBLE command since SRTM elevation data are int16 class. Plot the profile on the DEM:

hold on,plot(xi,yi,'r'),hold off

Since x and y are longitudes and latitudes, the linear profile is an equirectangular approximation (rhumb line or loxodrome) of the great circle... To compute distance along this path:

di = sqrt(((xi-x(1)).*cosd(yi)).^2 + (yi-y(1)).^2)*6370*pi/180;

where cosd(yi) reduces longitude distance versus latitude, and last factor is latitude degree length in km. Then plot the profile:

figure
plot(di,zi)

that's all. If you know the profile coordinates, you don't need to use GINPUT and you directly define xi and yi with LINSPACE.

[Martin Richard]

 

[François Beauducel]

Thanks Warwick for your kind feedback. I'll add a 2x2 tile option, in few days.
François.

[Warwick]

thank you. this is terrific, and surprisingly fast. A suggestion for a small increase in quick and easy usability would be an example of a 2by2 tile.