If we specify extents contained entirely within a DEM tile and require no transformations (e.g. do not remove geoid and keep pixel center or area center coordinate system of the original tiles), we want dem-stitcher to return a subset of that tile that contains the extent. This is not the case. We will post a notebook to demonstrate this issue is resolved as well as an integration test.
In original version, we (well, the error was mine) used rasterio.merge (link) with the bounds keyword argument. I believed this function to be obtaining a window around the bound extents; however, as the rasterio source code reveals, there is resampling that is done to ensure that the extents fit nicely into the CRS, which may be beneficial in many, if not most, GIS applications.
More precisely, to go from Pixel center to Area center coordinates, we need to translate the original rasters and then perform any necessary resampling. If resampling is done before a translation of the geo-transform, then things go wrong (majorly).
Further, the window operations were also resampling under the hood. All these issues need to be fixed.
In a general sense, if we have two affine transformations T_r and T_t representing "resampling" and "translation" resepectively, then T_r * T_t is not the same as T_t * T_r (it's worth checking this on simple examples i.e. seeing that both the geo-transform/affine transformaion and the resulting rasters are different).
This was further complicated when we wanted to change the resolution at higher latitudes because glo-30 shrinks latitude spacing as documented here. Thus, more resampling was done.
Which order is correct? We will take the approach T_r * T_t i.e.:
First translating profiles (if required) on raw data
Resampling (if required) at the end of all transformations.
To Reproduce
To be filled in later.
Additional context
The issue #31 points out is due the fact this resampling is done so a shift was selected based on investigating sites. However, this is not a totally correct approach.
Work to be done
We need tests to make sure that the above bug is removed.
The bug
Simple demonstration of bug:
If we specify extents contained entirely within a DEM tile and require no transformations (e.g. do not remove geoid and keep pixel center or area center coordinate system of the original tiles), we want
dem-stitcherto return a subset of that tile that contains the extent. This is not the case. We will post a notebook to demonstrate this issue is resolved as well as an integration test.In original version, we (well, the error was mine) used
rasterio.merge(link) with theboundskeyword argument. I believed this function to be obtaining a window around the bound extents; however, as the rasterio source code reveals, there is resampling that is done to ensure that the extents fit nicely into the CRS, which may be beneficial in many, if not most, GIS applications.More precisely, to go from Pixel center to Area center coordinates, we need to translate the original rasters and then perform any necessary resampling. If resampling is done before a translation of the geo-transform, then things go wrong (majorly).
Further, the window operations were also resampling under the hood. All these issues need to be fixed.
In a general sense, if we have two affine transformations
T_randT_trepresenting "resampling" and "translation" resepectively, thenT_r * T_tis not the same asT_t * T_r(it's worth checking this on simple examples i.e. seeing that both the geo-transform/affine transformaion and the resulting rasters are different).This was further complicated when we wanted to change the resolution at higher latitudes because
glo-30shrinks latitude spacing as documented here. Thus, more resampling was done.Which order is correct? We will take the approach
T_r * T_ti.e.:To Reproduce
To be filled in later.
Additional context
The issue #31 points out is due the fact this resampling is done so a shift was selected based on investigating sites. However, this is not a totally correct approach.
Work to be done
We need tests to make sure that the above bug is removed.