Support multiple GeoTIFF files (#109)

* Support GeoTIFF as well as NetCDF multifile raster

* Cache grid file for MultiFileRasters

* Use rioxarray merge_arrays rather than xarray.merge

Author: ianthomas23

example_multiband/create_multiband_netcdf.py

@@ -35,6 +35,7 @@
 
         x = np.linspace(lon*scale + offset + dx/2, (lon + 1)*scale + offset - dx/2, nx)
         y = np.linspace(lat*scale + offset + dy/2, (lat + 1)*scale + offset - dy/2, ny)
+        y = y[::-1]
         filename = os.path.join("example_multiband", f"dummy_{lon}_{lat}.nc")
 
         red = rng.random(size=(ny, nx), dtype=np.float32)

mapshader/core.py

@@ -19,7 +19,7 @@
 
 from mapshader.mercator import MercatorTileDefinition
 from mapshader.sources import MapSource
-from .multifile import MultiFileNetCDF
+from .multifile import MultiFileRaster
 
 import spatialpandas
 
@@ -76,8 +76,9 @@ def create_agg(source: MapSource,
     agg_func = source.agg_func
     geometry_type = source.geometry_type
 
-    if isinstance(source.data, MultiFileNetCDF):
+    if isinstance(source.data, MultiFileRaster):
         dataset = source.data.load_overview(z, source.band)
+        # Note this is really an xr.DataArray.
         if dataset is None:
             dataset = source.data.load_bounds(xmin, ymin, xmax, ymax, source.band,
                                               source.transforms)

mapshader/io.py

@@ -1,4 +1,4 @@
-from os.path import expanduser
+from os.path import expanduser, splitext
 
 import geopandas as gpd
 import numpy as np
@@ -34,31 +34,37 @@ def load_raster(file_path, transforms, force_recreate_overviews,
         The loaded data.
     """
 
-    if file_path.endswith('.tif'):
+    file_extension = splitext(file_path)[1]
+    if not file_extension:
+        raise RuntimeError(f"file_path does not have a file extension: {file_path}")
 
-        arr = xr.open_rasterio(expanduser(file_path),
-                               chunks={'y': 512, 'x': 512})
+    arr = None
 
-        if hasattr(arr, 'nodatavals'):
+    if '*' in file_path:
+        # Multiple files.
+        if file_extension in ['.nc', '.tif']:
+            arr = SharedMultiFile.get(file_path, transforms, force_recreate_overviews)
 
-            if np.issubdtype(arr.data.dtype, np.integer):
-                arr.data = arr.data.astype('f8')
+    else:
+        # Single file.
+        if file_extension == '.tif':
+            arr = xr.open_rasterio(expanduser(file_path), chunks={'y': 512, 'x': 512})
 
-            for val in arr.nodatavals:
-                arr.data[arr.data == val] = np.nan
+            if hasattr(arr, 'nodatavals'):
+                if np.issubdtype(arr.data.dtype, np.integer):
+                    arr.data = arr.data.astype('f8')
 
-        arr.name = file_path
+                for val in arr.nodatavals:
+                    arr.data[arr.data == val] = np.nan
 
-    elif file_path.endswith('.nc'):
+            arr.name = file_path
 
-        if '*' in file_path:
-            arr = SharedMultiFile.get(file_path, transforms, force_recreate_overviews)
-        else:
+        elif file_path.endswith('.nc'):
             # TODO: add chunk parameter to config
             arr = xr.open_dataset(file_path, chunks={'x': 512, 'y': 512})[layername]
             arr['name'] = file_path
 
-    else:
+    if arr is None:
         raise TypeError(f"Unable to load raster {file_path}")
 
     return arr