I need to compare the values of two raster layers at a set of given points. The raster layers are at a national scale, split into tiles and are at 1m resolution (2+TB). The points number in the 100s of millions.
I have an implementation which works but is slow and consumes a lot of resources. Two virtual rasters of the layers are created and the points are split into ~1 million point chunks (arbitrarily selected). Rasterio is then used to sample the rasters at each point and the results (in vector format) are compared. The process takes about 20 mins / million points and tops out at about 10GB of RAM. I've parallelised the process but the limiting factor seems to be the RAM requirement. The CPU load is relatively low and the disk I/O does not appear to be a limiting factor.
The raster layers are already loaded in a local GeoServer. Using this to interface with the data should mitigate the RAM requirement and enable more parallel processes to run concurrently. However, I have been unable to find a way to sample a raster at a given point rather than downloading an area GeoTIFF.
I can go down the GeoTIFF route if that is the only option but it seem counterproductive to speeding up the process due to the fact that it must be written to disk first (even if its just a tempfile). So, I have two questions:
- What is the fastest way to sample a raster (of this size) at a given set of points
- Can it be done using a GeoServer WCS
NB This process is taking place on a scalable GCP VM. If required it can be scaled up to 96 vCPUs and/or 624GB of RAM but I'd prefer not to brute force it and find a more elegant way. The scaling could apply to the GeoServers specs too.