Name
Geo::Elevation::HGT - Elevation service with terrain data provided by Mapzen and Amazon AWS S3
Version
Version 0.08
Synopsis
use Geo::Elevation::HGT;
my ($lat, $lon) = (45.8325, 6.86444444444444); # MontBlanc
my $geh = Geo::Elevation::HGT->new();
print $geh->get_elevation_hgt ($lat, $lon)."\n";
# 4790.99999999998
Description
This module implements an elevation service with terrain data provided by Mapzen and Amazon AWS S3.
You provide the latitude and longitude in decimal degrees, with south latitude and west longitude being negative.
The return is the elevation for this latitude and longitude in meters. Bilinear interpolation is applied to the elevations at the four grid points adjacent to the latitude plus longitude pair. To use bicubic interpolation using the sixteen grid points adjacent to the latitude plus longitude pair set the bicubic
parameter to 1.
You can also use your own terrain tiles by providing the corresponding path, see below. A good source for Europe that I am using was compiled by Sonny -- many thanks to him; found at https://data.opendataportal.at/dataset/dtm-europe
In addition you can specify a cache folder for subsequent use of downloaded tiles, see below.
There are only core dependencies
Carp
IO::Uncompress::AnyUncompress
HTTP::Tiny
POSIX
File::Find
List::Util
Notice
Geo::Elevation::HGT loads the required terrain tiles (from .HGT format files, see below) into the returned object, i.e. into memory.
Any following query requiring the same tile will be much faster since it only involves memory access, instead of a download from the internet.
In a typical application of getting elevations for a gpx track of an outdoor activity, all track points are normally on one tile, maybe two.
To get the elevations of a few thousand gpx track points is therefore normally quite fast.
Here is a benchmark I did with bilinear interpolation on my 2015 NUC5i3RYK with Intel 5010U dual-core processor
- 4.5 s for the first elevation with download of the terrain tile from Amazon AWS S3
- 0.5 s for the first elevation with the terrain tile stored on my NAS
- 50,000 elevations per second with the terrain tile in memory
It is the user's responsibility to respect the license and terms of use for the data provided by Mapzen and Amazon AWS S3.
Example
Get elevation in meters of any latitude plus longitude pair by the 'get_elevation_hgt' method.
Pass latitude and longitude in decimal notation, i.e. 45.8325, 6.86444444444444 for MontBlanc.
use Geo::Elevation::HGT;
# Latitude, longitude, elevation of some famous mountains (Wikipedia)
my @mountains = ( {Denali => {lat => 63+( 4*60+ 7)/3600, lon => -(151+( 0*60+28)/3600), ele_wiki => 6190 } },
{Everest => {lat => 27+(59*60+16)/3600, lon => 86+(55*60+29)/3600 , ele_wiki => 8848 } },
{MontBlanc => {lat => 45+(49*60+57)/3600, lon => 6+(51*60+52)/3600 , ele_wiki => 4810 } },
{Aconcagua => {lat => -(32+(39*60+13)/3600), lon => -(70+(00*60+40)/3600), ele_wiki => 6960.8} } );
# Get elevations of these famous mountains
my $geh = Geo::Elevation::HGT->new();
for my $mountain ( @mountains ) {
my ($name) = %$mountain;
my ($lat, $lon, $ele_wiki) = ($mountain->{$name}{'lat'}, $mountain->{$name}{'lon'}, $mountain->{$name}{'ele_wiki'});
my $ele_geh = $geh->get_elevation_hgt ($lat, $lon);
print join( ' ', $name, $lat, $lon, $ele_wiki, $ele_geh, $ele_geh-$ele_wiki, "\n");
}
# returns
# Denali 63.0686111111111 -151.007777777778 6190 6126.9999999999 -63.0000000000955
# Everest 27.9877777777778 86.9247222222222 8848 8713.0000000001 -134.999999999898
# MontBlanc 45.8325 6.86444444444444 4810 4791.00000000002 -18.9999999999764
# Aconcagua -32.6536111111111 -70.0111111111111 6960.8 6923.00000000011 -37.7999999998856
The elevations returned for these famous mountains are lower than their effective ones since the underlying data, which originate from SRTM (Shuttle Radar Topography Mission) and other sources, smooth peaks and thereby show them lower than commonly known.
Nevertheless the elevations returned should be useful to obtain a good representation of the elevation profile of a gpx track.
According to my experience these are much better than elevations recorded by a gps logger on a mobile phone. Such elevation data include a lot of noise with erratic jumps by 100 m.
Methods
new
$geh = Geo::Elevation::HGT->new( %parameters )
Constructor, returns a new Geo::Elevation::HGT object.
Valid parameters, all optional:
* folder
- the path to a folder where the terrain tiles to use are stored; no default
* url
- the url of the terrain tiles to use; default https://elevation-tiles-prod.s3.amazonaws.com/skadi
* cache_folder
path to an existing folder where the terrain tiles downloaded from $geh->{url}
will be stored for subsequent use; no default
Note that cache will not expire and will have to be cleared from outside of Geo::Elevation::HGT
. Thinking is that terrain data will not change very frequently.
* debug
- set to 1 to get some debug output to STDERR; default 0
* bicubic
- set to 1 to perform bicubic interpolation; default 0 = bilinear interpolation
Note that bicubic interpolation involves more multiplications (40 per elevation value instead of 3), reducing performance on my NUC to about 50%. Bicubic interpolation has an advantage when determining the elevation of a sharp peak in that it also takes into account the elevation gradient around the peak. This will allow the calculated elevation to be above that of the adjacent grid points. The effect is not very prominent however and will also depend on how much the gradient is reflected in the underlying data.
get_elevation_hgt
$ele_geh = $geh->get_elevation_hgt ($lat, $lon)
Arguments: latitude and longitude in decimal degrees, with south latitude and west longitude being negative
Returns the elevation for this latitude and longitude in meters
Flags
* $geh->{switch}
- is set to 1 in case a required tile is not found under the $geh->{folder}
path in which case the data source is switched to $geh->{url}
; otherwise 0
* $geh->{cache}
- is set to 1 in case a required tile is used from $geh->{cache_folder}
; otherwise 0
* $geh->{fail}
- is set to 1 in case access to a required tile under $geh->{url}
fails in which case all corresponding elevations will be set to 0; otherwise 0.
Status Description
* $geh->{status_descr}
- is set to a string describing how the terrain tile was found, using the following key words, or a combination thereof
- Memory
- the tile was in memory from a previous call
- Folder
- the user provided a valid folder path, if nothing follows the tile was found and used
- Switch
- the tile was not found under the $geh->{folder}
path
- Cached
- the tile was found under the $geh->{cache_folder}
path
- Url
- the tile was downloaded from the $geh->{url}
path, unless subsequent Failed
indicates failure
- Failed
- access to the tile under $geh->{url}
failed
- Saved
- the downloaded tile was saved to cache under the $geh->{cache_folder}
path
get_elevation_batch_hgt
$ele_geh = $geh->get_elevation_batch_hgt ($latlon)
Argument: an array reference with arrays of latitude-longitude pairs, i.e. [[lat1, lon1], [lat2, lon2], ...]
Returns an array reference with the associated elevations [ele1, ele2, ...]
Provided for user convenience. Calls method get_elevation_hgt
in turn for every latitude-longitude pair.
Flags and Status Description reflect the state after the call with the last latitude-longitude pair.
HGT DEM (digital elevation model) files
The names of individual DEM files refer to the latitude and longitude of the lower-left (south-west) corner of the tile.
e.g. N37W105 has its lower left corner at 37 degrees north latitude and 105 degrees west longitude.
The DEM is provided as 16-bit signed integer data in a simple binary raster.
There are no header or trailer bytes embedded in the file.
Each file is a series of signed 16-bit integers (two bytes) giving the height of each cell in meters arranged from west to east and then north to south.
Elevations are in meters referenced to the WGS84/EGM96 geoid.
Byte order is Motorola ("big-endian") standard with the most significant byte first.
Grid size is 3601x3601 for 1-minute DEMs or 7201x7201 for 0.5-minute DEMs or 1201x1201 for 3-minute DEMs.
The rows at the north and south edges as well as the columns at the east and west edges of each cell overlap and are identical to the edge rows and columns in the adjacent cell.
HGT directory tree
https://elevation-tiles-prod.s3.amazonaws.com/skadi
stores HGT files in subdirectories by their latitude
...
├── N45
│ ├── N45E000.hgt.gz
│ ├── N45E001.hgt.gz
│ ├── N45E002.hgt.gz
...
├── N46
│ ├── N46E000.hgt.gz
│ ├── N46E001.hgt.gz
│ ├── N46E002.hgt.gz
...
├── S46
│ ├── S46W000.hgt.gz
│ ├── S46W001.hgt.gz
│ ├── S46W002.hgt.gz
...
$geh->{folder}
will work without following the above storage pattern, as long as a file with the correct name is found somewhere under that path.
Similarly, $geh->{cache_folder}
will be built to the above storage pattern, but will also work if files are stored in any other way under that path.
HGT files need to be compressed to GNU zip (gzip) or ZIP (zip) compression format with .hgt.gz or .zip file extension, respectively.
Author
Ulrich Buck, <ulibuck at cpan.org>
Bugs
Please report any bugs or feature requests to bug-geo-elevation-hgt at rt.cpan.org
, or through the web interface at https://rt.cpan.org/NoAuth/ReportBug.html?Queue=Geo-Elevation-HGT. I will be notified, and then you'll automatically be notified of progress on your bug as I make changes.
Support
You can find documentation for this module with the perldoc command.
perldoc Geo::Elevation::HGT
You can also look for information at:
RT: CPAN's request tracker (report bugs here)
AnnoCPAN: Annotated CPAN documentation
CPAN Ratings
Search CPAN
Acknowledgements
Inspired by
Using DEMs to get GPX elevation profiles
plus others
License and Copyright
This software is Copyright (c) 2020 by Ulrich Buck.
This is free software, licensed under:
The Artistic License 2.0 (GPL Compatible)