Can't determine if CRS is correct
sagimann opened this issue · 3 comments
I have cropped a tif file using latest elevation lib, e.g.:
eio clip -o dem.tif --bounds 35.079732 30.901217 35.845468 31.870474
Then inspected the CRS of the file using rasterio
, which returned "EPSG:4326", plus some point inside the area (lat,lng), which returned 69:
import rasterio
with rasterio.open('my.tif') as dataset:
print(f'CRS: {dataset.crs}') # output EPSG:4326
row, col = dataset.index(lng, lat)
elevation = dataset.read(1, window=((row, row+1), (col, col+1)))
print(elevation[0, 0]) # output 69
I'm new to this field, but looking up EPSG:4326 says it is "WGS 84, latitude/longitude coordinate system based on the Earth's center of mass, used by the Global Positioning System among others."
Does this mean the generated tif file contains heights above the WGS84 ellipsoid? I'm confused because we have an independent and highly-accurate measuring system that explicitly reports 69m above MSL and NOT above the ellipsoid. And that system also tells us they are ~20m apart at that point, so they can't both be right..
What am I missing here?
4326 is a 2-dimensional coordinate system so it does not make sense for elevation of any kind of be expressed in it.
The real question is where dem.tif
came from, what the CRS is of that, and what the documentation/metadata says. Clipping should not really change that.
In WGS84 there are two concepts:
- height above ellipsoid (HAE) which is purely geometric and unrelated to gravity. Unless you are doing surveying or geodesy you probably don't want this
- WGS84 orthometric height: The distance above the geoid model that is part of WGS84. Today, that means EGM2008, but people often use EMG96 anyway, even though it is formally obsolete. People also use truncated representations of the model, leading to errors in the output (compared to the full model which is by definition correct for WGS84 orthometric height).
Also, "above MSL" is a fuzzy term and I am surprised if anything "highly accurate" is using it. Around the earth the average sea level varies in height relative to the geoid. Ask them what they really mean. Perhaps it is WGS84 orthometric height -- but it is difficult to measure that precisely. Perhaps they mean ITRF HAE converted via EGM2008, which is almost the same thing. Perhaps they mean some national height datum. Most of those are very close to WGS84 orthometric height, at the meter or two level, and to each other. There is one in Europe (Holland?) which is ~2m different.
I should have said earlier that a DEM that was in HAE, while conceptually sensible, seems very odd to me.
For a much longer treatise on height, see https://opencommons.uconn.edu/cgi/viewcontent.cgi?article=1000&context=nrme_monos&httpsredir=1&referer=
The real question is where
dem.tif
came from, what the CRS is of that, and what the documentation/metadata says. Clipping should not really change that.
As mentioned, dem.tif is the output of the eio
command line provided by this project. If my understanding is correct, it comes from SRTM which provides data in both ellipsoid and EGM96 frames of references. The problem is gdalinfo
says the file's CRS is WGS84 but when querying the elevation values, I get EGM96 values instead. So either I don't understand the header value or eio
gets the EGM96 variant from NASA but puts a WGS84 header on it...