Save geotif

docs/changelog.md

@@ -2,6 +2,10 @@
 
 All notable changes to this project will be documented below.
 
+#### geospatial.save_geotif
+
+* v0.1dev: **geospatial.save_geotif**(*out_img, in_img, out_path, nodata=None*)
+
 #### geospatial.read_geotif
 
 * v0.1dev: spectral = **geospatial.read_geotif**(*filename, bands="B,G,R", cropto=None*)

docs/save_geotif.md

@@ -0,0 +1,35 @@
+## Save Geo-tif Data
+
+Save data from segmentation or modification of a geotif. The resulting geotif file will have the same georeferenced coordinates as the original image. 
+
+**plantcv.geospatial.save_geotif**(*out_img, in_img, out_path, nodata=None*)
+
+- **Parameters:**
+    - out_img - Image to be saved, could be a segmented mask generated by PlantCV.
+    - in_img - Spectral_data object created by read_geotif. Used to access metadata.
+    - out_path - Path to save the output geotif.
+    - nodata - Value to assign as "no data" during saving out. Default is "None". 
+
+- **Example use:**
+    - below
+
+
+```python
+import plantcv.geospatial as geo
+import plantcv.plantcv as pcv
+
+# Read geotif in
+img = geo.read_geotif(filename="./rgb.tif", bands="R,G,B")
+# Threshold using PlantCV
+gray_rgb = pcv.rgb2gray_hsv(rgb_img=img.pseudo_rgb, channel="s")
+thresh_rgb = pcv.threshold.binary(gray_img=gray_rgb, threshold=60, object_type="light")
+# Save mask as a geotif
+geo.save_geotif(thresh_rgb, img, "./masked_rgb.tif")
+
+```
+
+![Screenshot](documentation_images/rgb_tif.png)
+![Screenshot](documentation_images/gray_tif.png)
+![Screenshot](documentation_images/masked_tif.png)
+
+**Source Code:** [Here](https://github.com/danforthcenter/plantcv-geospatial/blob/main/plantcv/geospatial/save_geotif.py)

mkdocs.yml

@@ -18,6 +18,7 @@ nav:
   - 'PlantCV Namespace':
     - 'Geopspatial Tools':
         - Read Geo-tif Data: read_geotif.md
+        - Save Geo-tif Data: save_geotif.md
         - Transform coordinate points: transform_points.md
         - Transform coordinate polygons: transform_polygons.md
 markdown_extensions:

plantcv/geospatial/__init__.py

@@ -2,12 +2,14 @@
 from plantcv.geospatial.read_geotif import read_geotif
 from plantcv.geospatial.transform_points import transform_points
 from plantcv.geospatial.transform_polygons import transform_polygons
+from plantcv.geospatial.save_geotif import save_geotif
 
 # Auto versioning
 __version__ = version("plantcv-geospatial")
 
 __all__ = [
     "read_geotif",
     "transform_points",
-    "transform_polygons"
+    "transform_polygons",
+    "save_geotif"
     ]

plantcv/geospatial/save_geotif.py

@@ -0,0 +1,55 @@
+# Save georeferenced TIF files from plantcv segmented masks
+
+import rasterio
+
+
+def save_geotif(out_img, in_img, out_path, nodata=None):
+    """Save Georeferenced TIF file from image.
+
+    Parameters
+    ----------
+    out_img : ndarray
+        Image to be saved, could be a segmented mask generated by PlantCV.
+    in_img : str
+        Spectral_data object created by read_geotif.
+    out_path : str
+        Path to save the output geotif.
+    nodata : int or float, optional
+        Value to assign as "no data" during saving out.
+    """
+    metadata = in_img.metadata.copy()
+    # Update metadata - change count in case input image has mask layer
+    metadata.update({
+        "driver": "GTiff",
+        "height": out_img.shape[0],
+        "width": out_img.shape[1],
+        "dtype": "uint8",
+        "count": 3,
+        "nodata": nodata
+    })
+    # Read in channels to know if reordering is necessary
+    waves = in_img.wavelength_dict
+
+    # For color outputs
+    if len(out_img.shape) == 3:
+        # Reorder dimensions
+        out_img = out_img.transpose(2, 0, 1)
+        # Reorder channels if needed
+        chord = [waves[480], waves[560], waves[650]]
+        # For multispec
+        if out_img.shape[0] > 3:
+            metadata.update({
+                "count": out_img.shape[0]
+            })
+            for i in range(out_img.shape[0]-3):
+                chord.append(i+3)
+        out_img = out_img[chord, :, :]
+        with rasterio.open(out_path, 'w', **metadata) as dest:
+            dest.write(out_img)
+    # For binary outputs
+    elif len(out_img.shape) == 2:
+        metadata.update({
+            "count": 1
+        })
+        with rasterio.open(out_path, 'w', **metadata) as dest:
+            dest.write(out_img, 1)

tests/test_geospatial_read_geotif.py

@@ -1,4 +1,4 @@
-"""Tests for geospatial.readd_geotif."""
+"""Tests for geospatial.read_geotif."""
 
 import pytest
 from plantcv.geospatial import read_geotif

tests/test_geospatial_save_geotif.py

@@ -0,0 +1,31 @@
+"""Tests for geospatial.save_geotif."""
+
+import os
+import plantcv.plantcv as pcv
+from plantcv.geospatial import save_geotif, read_geotif
+
+
+def test_save_geotif_color(test_data, tmpdir):
+    """Test for PlantCV."""
+    # Create a test tmp directory
+    cache_dir = tmpdir.mkdir("cache")
+    # Create test image
+    img = read_geotif(filename=test_data.cropped_tif, bands="B,G,R,RE,N")
+    filename = os.path.join(cache_dir, 'test.tif')
+    save_geotif(out_img=img.array_data, in_img=img, out_path=filename)
+    # Assert that the file was created
+    assert os.path.exists(filename)
+
+
+def test_save_geotif_binary(test_data, tmpdir):
+    """Test for PlantCV."""
+    # Create a test tmp directory
+    cache_dir = tmpdir.mkdir("cache")
+    # Create test image
+    img = read_geotif(filename=test_data.rgb_tif, bands="R,G,B")
+    gray = pcv.rgb2gray_hsv(rgb_img=img.pseudo_rgb, channel="s")
+    thresh = pcv.threshold.binary(gray_img=gray, threshold=60, object_type="light")
+    filename = os.path.join(cache_dir, 'test_binary.tif')
+    save_geotif(out_img=thresh, in_img=img, out_path=filename)
+    # Assert that the file was created
+    assert os.path.exists(filename)