Adding folium_helper.py

.gitignore

@@ -2,6 +2,9 @@
 __pycache__/
 *.py[cod]
 *$py.class
+.DS_Store
+*.hdf
+*.gzip
 
 # C extensions
 *.so

src/folium_helper.py

@@ -0,0 +1,129 @@
+import folium
+import glob
+import numpy as np
+import os
+import rasterio as rio
+import tempfile
+
+from rasterio.warp import calculate_default_transform, reproject, Resampling
+from pyproj import Transformer 
+
+# -----------------------------------------------------------------------------
+# Uses rasterio to open a raster, get the metadata and crs
+# associated with it and get all the subdatasets in the file.
+# This is very useful for hdf files such as MODIS hdfs.
+# -----------------------------------------------------------------------------
+def print_subdatasets(filename):
+    bands_to_return = []
+    with rio.open(filename) as dataset:
+        meta_data = dataset.meta
+        crs = dataset.read_crs()
+
+        print([name for name in dataset.subdatasets if search_term in name])
+
+# -----------------------------------------------------------------------------
+# Gets a tiff that has the correct metadata for that tile, gets the metadata
+# from the source tif and copies to a destination tiff. 
+# -----------------------------------------------------------------------------     
+def add_metadata_to_annual_product(filepath, model_type, year, tile):
+    metadata_pull_src = [fv for fv in glob.glob(os.path.join(filepath, "{}-1*-{}-MOD-*.tif".format(year, tile)))][0]
+    with rio.open(metadata_pull_src) as src:
+        src_meta = src.meta
+    dst_tiffs = [os.path.join(filepath, fn) for fn in os.listdir(filepath) if "{0}-{1}".format(year, tile) in os.path.basename(fn)]
+    [copy_meta(dst_tiff, src_meta, metadata_pull_src) for dst_tiff in dst_tiffs]
+
+# -----------------------------------------------------------------------------
+# Given a path to a tiff with no metadata, assign the metadata given to that
+# tiff.
+# -----------------------------------------------------------------------------     
+def copy_meta(dst_path, src_meta, src_name):
+    print('Copying metadata from {} to {}'.format(src_name, dst_path))
+    with rio.open(dst_path, 'r+') as dst:
+        dst.crs = src_meta['crs']
+        dst.transform = src_meta['transform']        
+
+# -----------------------------------------------------------------------------
+# Given a tiff file as input, open the tiff and get the transform needed to
+# reproject from the tiff's source crs to the one we want (EPSG:3857).
+# For each band in the tiff, open then reproject it into the desired crs
+# then write to a temporary file. Return the path to the temp file.
+# -----------------------------------------------------------------------------
+def reproject_to_3857(input_tiff):
+    # Set desitnation CRS
+    dst_crs = f"EPSG:3857"
+
+    # set out path
+    out_path_rproj = os.path.join(tempfile.gettempdir(), input_tiff.split('/')[-1].replace('.tif','-3857.tif'))
+
+    with rio.open(input_tiff) as src:
+        # get src bounds and transform
+        transform, width, height = calculate_default_transform(src.crs, dst_crs, src.width, src.height, *src.bounds)
+        print('Transform: {}'.format(transform))
+        print('Width: {} Height: {}'.format(width, height))
+        kwargs = src.meta.copy()
+        kwargs.update({'crs': dst_crs,
+                   'transform': transform,
+                   'width': width,
+                   'height': height})
+
+        # reproject and write to file
+        with rio.open(out_path_rproj, 'w', **kwargs) as dst:
+            for i in range(1, src.count + 1):
+                reproject(source=rio.band(src, i),
+                      destination=rio.band(dst, i),
+                      src_transform=src.transform,
+                      src_crs=src.crs,
+                      dst_transform=transform,
+                      dst_crs=dst_crs,
+                      resampling=Resampling.nearest)
+    return out_path_rproj
+
+# -----------------------------------------------------------------------------
+# In order for folium to work properly we need to pass it the bounding box
+# of the tiff in the form of lat and lon. This is done by using rasterio.
+# -----------------------------------------------------------------------------
+def get_bounds(tiff_3857):
+    with rio.open(tiff_3857) as src:
+        src_crs = src.crs['init'].upper()
+        min_lon, min_lat, max_lon, max_lat = src.bounds
+    bounds_orig = [[min_lat, min_lon], [max_lat, max_lon]]
+    bounds = []
+    dst_crs = 'EPSG:4326'
+    for item in bounds_orig:   
+        #converting to lat/lon
+        lat = item[0]
+        lon = item[1]
+        proj = Transformer.from_crs(int(src_crs.split(":")[1]), int(dst_crs.split(":")[1]), always_xy=True)
+        lon_n, lat_n = proj.transform(lon, lat)
+        bounds.append([lat_n, lon_n])
+    center_lon = bounds[0][1] + (bounds[1][1] - bounds[0][1])/2
+    center_lat = bounds[0][0] + (bounds[1][0] - bounds[0][0])/2
+    return {'bounds': bounds, 'center': (center_lon, center_lat)}
+
+# -----------------------------------------------------------------------------
+# Use rasterio to open and read in the desired band name as a nd-array.
+# -----------------------------------------------------------------------------
+def open_and_get_band(file_name, band_num=1):
+    with rio.open(file_name) as data:
+        b = data.read(band_num)
+    return b
+
+# -----------------------------------------------------------------------------
+# Given an nd-array (band) and the bounds in lat lon of the nd-array, return
+# a folium layer. To add on the map.
+# -----------------------------------------------------------------------------
+def get_overlay(band, meta_dict, name, opacity=1.0, show=True):
+    return folium.raster_layers.ImageOverlay(band, 
+                                             bounds=meta_dict['bounds'], 
+                                             name=name, 
+                                             opacity=opacity, 
+                                             show=show)
+
+# -----------------------------------------------------------------------------
+# We don't need to keep those temp files we made for the reprojections around.
+# -----------------------------------------------------------------------------
+def cleanup(filename):
+    if os.path.exists(filename):
+        os.remove(filename)
+    else:
+        print('No file: {} exists.'.format(filename))