.travis.yml

@@ -8,16 +8,16 @@ env:
         - PYTHON_VERSION=$TRAVIS_PYTHON_VERSION
         - NUMPY_VERSION=stable
         - MAIN_CMD='python setup.py'
-        - CONDA_DEPENDENCIES='pillow gdal xarray dask toolz mock coverage coveralls codecov'
+        - CONDA_DEPENDENCIES='pillow gdal xarray dask toolz mock coverage coveralls codecov rasterio'
         - SETUP_XVFB=False
         - EVENT_TYPE='push pull_request'
         - SETUP_CMD='test'
         - CONDA_CHANNELS='conda-forge'
+        - CONDA_CHANNEL_PRIORITY=True
+
 install:
-# Get rasterio 1.0
     - git clone --depth 1 git://github.com/astropy/ci-helpers.git
     - source ci-helpers/travis/setup_conda.sh
-    - conda install -c conda-forge/label/dev rasterio
 script: coverage run --source=trollimage setup.py test
 after_success:
     - if [[ $PYTHON_VERSION == 3.6 ]]; then coveralls; fi

CHANGELOG.md

-## Version 1.6.3 (2018/12/20)
+## Version 1.7.0 (2019/02/28)
+
+### Issues Closed
 
+* [Issue 27](https://github.com/pytroll/trollimage/issues/27) - Add "overviews" to save options
+* [Issue 5](https://github.com/pytroll/trollimage/issues/5) - Add alpha channel to Colormaps
+
+In this release 2 issues were closed.
+
+### Pull Requests Merged
+
+#### Bugs fixed
+
+* [PR 42](https://github.com/pytroll/trollimage/pull/42) - Fix stretch_linear to be dask serializable
+* [PR 41](https://github.com/pytroll/trollimage/pull/41) - Refactor XRImage pil_save to be serializable
+
+#### Features added
+
+* [PR 44](https://github.com/pytroll/trollimage/pull/44) - Add support for adding overviews to rasterio-managed files
+* [PR 43](https://github.com/pytroll/trollimage/pull/43) - Add support for jpeg2000 writing
+* [PR 40](https://github.com/pytroll/trollimage/pull/40) - Modify colorize routine to allow colorizing using colormaps with alpha channel
+* [PR 39](https://github.com/pytroll/trollimage/pull/39) - Add 'keep_palette' keyword argument 'XRImage.save' to prevent P -> RGB conversion on save
+* [PR 36](https://github.com/pytroll/trollimage/pull/36) - Add support for saving gcps
+
+In this release 7 pull requests were closed.
+
+
+## Version 1.6.3 (2018/12/20)
 
 ### Pull Requests Merged
 
@@ -12,7 +38,6 @@ In this release 1 pull request was closed.
 
 ## Version 1.6.2 (2018/12/20)
 
-
 ### Pull Requests Merged
 
 #### Bugs fixed
@@ -25,7 +50,6 @@ In this release 2 pull requests were closed.
 
 ## Version 1.6.1 (2018/12/19)
 
-
 ### Pull Requests Merged
 
 #### Bugs fixed

appveyor.yml

@@ -3,8 +3,9 @@ environment:
     PYTHON: "C:\\conda"
     MINICONDA_VERSION: "latest"
     CMD_IN_ENV: "cmd /E:ON /V:ON /C .\\ci-helpers\\appveyor\\windows_sdk.cmd"
-    CONDA_DEPENDENCIES: "pillow gdal xarray dask toolz mock coverage coveralls codecov"
+    CONDA_DEPENDENCIES: "pillow gdal xarray dask toolz mock coverage coveralls codecov rasterio"
     CONDA_CHANNELS: "conda-forge"
+    CONDA_CHANNEL_PRIORITY: "True"
 
   matrix:
     - PYTHON: "C:\\Python27_64"

debian/changelog

+trollimage (1.7.0-1) unstable; urgency=medium
+
+  * New upstream release.
+  * debian/patches
+    - refresh all patches
+
+ -- Antonio Valentino <antonio.valentino@tiscali.it>  Fri, 01 Mar 2019 06:54:52 +0000
+
 trollimage (1.6.3-1) unstable; urgency=medium
 
   * Initial version (Closes: #916550)

debian/patches/0001-No-display.patch

@@ -8,10 +8,10 @@ Skip tests that require display.
  1 file changed, 1 insertion(+)
 
 diff --git a/trollimage/tests/test_image.py b/trollimage/tests/test_image.py
-index fb72c6c..60ce4a8 100644
+index 816d4f2..a949c8f 100644
 --- a/trollimage/tests/test_image.py
 +++ b/trollimage/tests/test_image.py
-@@ -1243,6 +1243,7 @@ class TestXRImage(unittest.TestCase):
+@@ -1790,6 +1790,7 @@ class TestXRImage(unittest.TestCase):
      def test_putalpha(self):
          pass
  

trollimage/colormap.py

@@ -26,6 +26,7 @@
 import numpy as np
 from trollimage.colorspaces import rgb2hcl, hcl2rgb
 
+
 def colorize(arr, colors, values):
     """Colorize a monochromatic array *arr*, based *colors* given for
     *values*. Interpolation is used. *values* must be in ascending order.
@@ -52,6 +53,7 @@ def colorize(arr, colors, values):
     except AttributeError:
         return channels
 
+
 def palettize(arr, colors, values):
     """From start *values* apply *colors* to *data*.
     """
@@ -124,6 +126,14 @@ class Colormap(object):
                         (self.values.max() - self.values.min()))
                        * (max_val - min_val) + min_val)
 
+    def to_rio(self):
+        """Converts the colormap to a rasterio colormap.
+        """
+        self.colors = (((self.colors * 1.0 - self.colors.min()) /
+                        (self.colors.max() - self.colors.min())) * 255)
+        return dict(zip(self.values, tuple(map(tuple, self.colors))))
+
+
 # matlab jet "#00007F", "blue", "#007FFF", "cyan", "#7FFF7F", "yellow",
 # "#FF7F00", "red", "#7F0000"
 
@@ -160,7 +170,6 @@ purples = Colormap((0.0, (252 / 255.0, 251 / 255.0, 253 / 255.0)),
 reds = Colormap((0.0, (1.0, 245 / 255.0, 240 / 255.0)),
                 (1.0, (103 / 255.0, 0.0, 13 / 255.0)))
 
-
 # * Multihue *
 
 # BuGn
@@ -231,7 +240,6 @@ sequential_colormaps = [blues, greens, greys, oranges, purples, reds,
                         bugn, bupu, gnbu, orrd, pubu, pubugn, purd, rdpu,
                         ylgn, ylgnbu, ylorbr, ylorrd]
 
-
 # * Diverging *
 
 brbg = Colormap((0.0, (84 / 255.0, 48 / 255.0, 5 / 255.0)),
@@ -246,7 +254,6 @@ brbg = Colormap((0.0, (84 / 255.0, 48 / 255.0, 5 / 255.0)),
                 (0.9, (1 / 255.0, 102 / 255.0, 94 / 255.0)),
                 (1.0, (0 / 255.0, 60 / 255.0, 48 / 255.0)))
 
-
 piyg = Colormap((0.0, (142 / 255.0, 1 / 255.0, 82 / 255.0)),
                 (0.1, (197 / 255.0, 27 / 255.0, 125 / 255.0)),
                 (0.2, (222 / 255.0, 119 / 255.0, 174 / 255.0)),
@@ -259,8 +266,6 @@ piyg = Colormap((0.0, (142 / 255.0, 1 / 255.0, 82 / 255.0)),
                 (0.9, (77 / 255.0, 146 / 255.0, 33 / 255.0)),
                 (1.0, (39 / 255.0, 100 / 255.0, 25 / 255.0)))
 
-
-
 prgn = Colormap((0.0, (64 / 255.0, 0 / 255.0, 75 / 255.0)),
                 (0.1, (118 / 255.0, 42 / 255.0, 131 / 255.0)),
                 (0.2, (153 / 255.0, 112 / 255.0, 171 / 255.0)),
@@ -273,7 +278,6 @@ prgn = Colormap((0.0, (64 / 255.0, 0 / 255.0, 75 / 255.0)),
                 (0.9, (27 / 255.0, 120 / 255.0, 55 / 255.0)),
                 (1.0, (0 / 255.0, 68 / 255.0, 27 / 255.0)))
 
-
 puor = Colormap((0.0, (127 / 255.0, 59 / 255.0, 8 / 255.0)),
                 (0.1, (179 / 255.0, 88 / 255.0, 6 / 255.0)),
                 (0.2, (224 / 255.0, 130 / 255.0, 20 / 255.0)),
@@ -286,7 +290,6 @@ puor = Colormap((0.0, (127 / 255.0, 59 / 255.0, 8 / 255.0)),
                 (0.9, (84 / 255.0, 39 / 255.0, 136 / 255.0)),
                 (1.0, (45 / 255.0, 0 / 255.0, 75 / 255.0)))
 
-
 rdbu = Colormap((0.0, (103 / 255.0, 0 / 255.0, 31 / 255.0)),
                 (0.1, (178 / 255.0, 24 / 255.0, 43 / 255.0)),
                 (0.2, (214 / 255.0, 96 / 255.0, 77 / 255.0)),
@@ -299,7 +302,6 @@ rdbu = Colormap((0.0, (103 / 255.0, 0 / 255.0, 31 / 255.0)),
                 (0.9, (33 / 255.0, 102 / 255.0, 172 / 255.0)),
                 (1.0, (5 / 255.0, 48 / 255.0, 97 / 255.0)))
 
-
 rdgy = Colormap((0.0, (103 / 255.0, 0 / 255.0, 31 / 255.0)),
                 (0.1, (178 / 255.0, 24 / 255.0, 43 / 255.0)),
                 (0.2, (214 / 255.0, 96 / 255.0, 77 / 255.0)),
@@ -312,7 +314,6 @@ rdgy = Colormap((0.0, (103 / 255.0, 0 / 255.0, 31 / 255.0)),
                 (0.9, (77 / 255.0, 77 / 255.0, 77 / 255.0)),
                 (1.0, (26 / 255.0, 26 / 255.0, 26 / 255.0)))
 
-
 rdylbu = Colormap((0.0, (165 / 255.0, 0 / 255.0, 38 / 255.0)),
                   (0.1, (215 / 255.0, 48 / 255.0, 39 / 255.0)),
                   (0.2, (244 / 255.0, 109 / 255.0, 67 / 255.0)),
@@ -325,7 +326,6 @@ rdylbu = Colormap((0.0, (165 / 255.0, 0 / 255.0, 38 / 255.0)),
                   (0.9, (69 / 255.0, 117 / 255.0, 180 / 255.0)),
                   (1.0, (49 / 255.0, 54 / 255.0, 149 / 255.0)))
 
-
 rdylgn = Colormap((0.0, (165 / 255.0, 0 / 255.0, 38 / 255.0)),
                   (0.1, (215 / 255.0, 48 / 255.0, 39 / 255.0)),
                   (0.2, (244 / 255.0, 109 / 255.0, 67 / 255.0)),
@@ -350,7 +350,6 @@ spectral = Colormap((0.0, (158 / 255.0, 1 / 255.0, 66 / 255.0)),
                     (0.9, (50 / 255.0, 136 / 255.0, 189 / 255.0)),
                     (1.0, (94 / 255.0, 79 / 255.0, 162 / 255.0)))
 
-
 diverging_colormaps = [brbg, piyg, prgn, puor, rdbu, rdgy, rdylbu, rdylgn,
                        spectral]
 
@@ -442,6 +441,7 @@ qualitative_colormaps = [set1, set2, set3,
                          paired, accent, dark2,
                          pastel1, pastel2]
 
+
 def colorbar(height, length, colormap):
     """Return the channels of a colorbar.
     """
@@ -451,6 +451,7 @@ def colorbar(height, length, colormap):
 
     return colormap.colorize(cbar)
 
+
 def palettebar(height, length, colormap):
     """Return the channels of a palettebar.
     """

trollimage/image.py

@@ -76,6 +76,7 @@ def check_image_format(fformat):
              "png": "png",
              "xbm": "xbm",
              "xpm": "xpm",
+             'jp2': 'jp2',
              }
     fformat = fformat.lower()
     try:

trollimage/tests/test_colormap.py

@@ -23,14 +23,15 @@
 """Test colormap.py
 """
 
-
 import unittest
 from trollimage import colormap
 import numpy as np
 
+
 class TestColormapClass(unittest.TestCase):
     """Test case for the colormap object.
     """
+
     def test_colorize(self):
         """Test colorize
         """
@@ -79,8 +80,6 @@ class TestColormapClass(unittest.TestCase):
         self.assertTrue(np.allclose(colors, cm_.colors))
         self.assertTrue(all(channels == [0, 0, 1, 2, 3, 3]))
 
-
-
     def test_set_range(self):
         """Test set_range
         """
@@ -161,6 +160,20 @@ class TestColormapClass(unittest.TestCase):
         self.assertTrue(np.allclose(channel, np.arange(4)))
         self.assertTrue(np.allclose(palette, cm_.colors))
 
+    def test_to_rio(self):
+        """Test conversion to rasterio colormap
+        """
+        cm_ = colormap.Colormap((1, (1, 1, 0)),
+                                (2, (0, 1, 1)),
+                                (3, (1, 1, 1)),
+                                (4, (0, 0, 0)))
+
+        d = cm_.to_rio()
+        exp = {1: (255, 255, 0), 2: (0, 255, 255),
+               3: (255, 255, 255), 4: (0, 0, 0)}
+
+        self.assertEqual(d, exp)
+
 
 def suite():
     """The suite for test_colormap.

trollimage/tests/test_image.py

@@ -58,7 +58,6 @@ class CustomScheduler(object):
 
 
 class TestEmptyImage(unittest.TestCase):
-
     """Class for testing the mpop.imageo.image module
     """
 
@@ -251,7 +250,6 @@ class TestEmptyImage(unittest.TestCase):
 
 
 class TestImageCreation(unittest.TestCase):
-
     """Class for testing the mpop.imageo.image module
     """
 
@@ -340,7 +338,6 @@ class TestImageCreation(unittest.TestCase):
 
 
 class TestRegularImage(unittest.TestCase):
-
     """Class for testing the mpop.imageo.image module
     """
 
@@ -696,7 +693,6 @@ class TestRegularImage(unittest.TestCase):
 
 
 class TestFlatImage(unittest.TestCase):
-
     """Test a flat image, ie an image where min == max.
     """
 
@@ -728,7 +724,6 @@ class TestFlatImage(unittest.TestCase):
 
 
 class TestNoDataImage(unittest.TestCase):
-
     """Test an image filled with no data.
     """
 
@@ -805,15 +800,22 @@ class TestXRImage(unittest.TestCase):
         import dask.array as da
         from dask.delayed import Delayed
         from trollimage import xrimage
+        from trollimage.colormap import brbg, Colormap
 
-        data = xr.DataArray(np.arange(75).reshape(5, 5, 3) / 75., dims=[
+        # RGBA colormap
+        bw = Colormap(
+            (0.0, (1.0, 1.0, 1.0, 1.0)),
+            (1.0, (0.0, 0.0, 0.0, 0.5)),
+        )
+
+        data = xr.DataArray(np.arange(75).reshape(5, 5, 3) / 74., dims=[
             'y', 'x', 'bands'], coords={'bands': ['R', 'G', 'B']})
         img = xrimage.XRImage(data)
         with NamedTemporaryFile(suffix='.png') as tmp:
             img.save(tmp.name)
 
         # Single band image
-        data = xr.DataArray(np.arange(75).reshape(15, 5, 1) / 75., dims=[
+        data = xr.DataArray(np.arange(75).reshape(15, 5, 1) / 74., dims=[
             'y', 'x', 'bands'], coords={'bands': ['L']})
         # Single band image to JPEG
         img = xrimage.XRImage(data)
@@ -824,7 +826,21 @@ class TestXRImage(unittest.TestCase):
         with NamedTemporaryFile(suffix='.png') as tmp:
             img.save(tmp.name)
 
-        data = xr.DataArray(da.from_array(np.arange(75).reshape(5, 5, 3) / 75.,
+        # Single band image palettized
+        data = xr.DataArray(np.arange(75).reshape(15, 5, 1) / 74., dims=[
+            'y', 'x', 'bands'], coords={'bands': ['L']})
+        # Single band image to JPEG
+        img = xrimage.XRImage(data)
+        img.palettize(brbg)
+        with NamedTemporaryFile(suffix='.png') as tmp:
+            img.save(tmp.name)
+        # RGBA colormap
+        img = xrimage.XRImage(data)
+        img.palettize(bw)
+        with NamedTemporaryFile(suffix='.png') as tmp:
+            img.save(tmp.name)
+
+        data = xr.DataArray(da.from_array(np.arange(75).reshape(5, 5, 3) / 74.,
                                           chunks=5),
                             dims=['y', 'x', 'bands'],
                             coords={'bands': ['R', 'G', 'B']})
@@ -832,7 +848,7 @@ class TestXRImage(unittest.TestCase):
         with NamedTemporaryFile(suffix='.png') as tmp:
             img.save(tmp.name)
 
-        data = data.where(data > (10 / 75.0))
+        data = data.where(data > (10 / 74.0))
         img = xrimage.XRImage(data)
         with NamedTemporaryFile(suffix='.png') as tmp:
             img.save(tmp.name)
@@ -852,8 +868,9 @@ class TestXRImage(unittest.TestCase):
         import rasterio as rio
 
         # numpy array image - scale to 0 to 1 first
-        data = xr.DataArray(np.arange(75.).reshape(5, 5, 3) / 75., dims=[
-            'y', 'x', 'bands'], coords={'bands': ['R', 'G', 'B']})
+        data = xr.DataArray(np.arange(75).reshape(5, 5, 3) / 75.,
+                            dims=['y', 'x', 'bands'],
+                            coords={'bands': ['R', 'G', 'B']})
         img = xrimage.XRImage(data)
         with NamedTemporaryFile(suffix='.tif') as tmp:
             img.save(tmp.name)
@@ -1001,6 +1018,7 @@ class TestXRImage(unittest.TestCase):
         import dask.array as da
         from trollimage import xrimage
         import rasterio as rio
+        from rasterio.control import GroundControlPoint
 
         # numpy array image
         data = xr.DataArray(np.arange(75).reshape(5, 5, 3), dims=[
@@ -1044,6 +1062,90 @@ class TestXRImage(unittest.TestCase):
             da.store(*delay)
             delay[1].close()
 
+        # GCPs
+        class FakeArea():
+            def __init__(self, lons, lats):
+                self.lons = lons
+                self.lats = lats
+
+        gcps = [GroundControlPoint(1, 1, 100.0, 1000.0, z=0.0),
+                GroundControlPoint(2, 3, 400.0, 2000.0, z=0.0)]
+        crs = 'epsg:4326'
+
+        lons = xr.DataArray(da.from_array(np.arange(25).reshape(5, 5), chunks=5),
+                            dims=['y', 'x'],
+                            attrs={'gcps': gcps,
+                                   'crs': crs})
+
+        lats = xr.DataArray(da.from_array(np.arange(25).reshape(5, 5), chunks=5),
+                            dims=['y', 'x'],
+                            attrs={'gcps': gcps,
+                                   'crs': crs})
+
+        data = xr.DataArray(da.from_array(np.arange(75).reshape(5, 5, 3), chunks=5),
+                            dims=['y', 'x', 'bands'],
+                            coords={'bands': ['R', 'G', 'B']},
+                            attrs={'area': FakeArea(lons, lats)})
+        img = xrimage.XRImage(data)
+        with NamedTemporaryFile(suffix='.tif') as tmp:
+            img.save(tmp.name)
+            with rio.open(tmp.name) as f:
+                fgcps, fcrs = f.gcps
+            for ref, val in zip(gcps, fgcps):
+                self.assertEqual(ref.col, val.col)
+                self.assertEqual(ref.row, val.row)
+                self.assertEqual(ref.x, val.x)
+                self.assertEqual(ref.y, val.y)
+                self.assertEqual(ref.z, val.z)
+            self.assertEqual(crs, fcrs)
+
+        # with rasterio colormap provided
+        exp_cmap = {i: (i, 255 - i, i, 255) for i in range(256)}
+        data = xr.DataArray(da.from_array(np.arange(81).reshape(9, 9, 1), chunks=9),
+                            dims=['y', 'x', 'bands'],
+                            coords={'bands': ['P']})
+        img = xrimage.XRImage(data)
+        with NamedTemporaryFile(suffix='.tif') as tmp:
+            img.save(tmp.name, keep_palette=True, cmap=exp_cmap)
+            with rio.open(tmp.name) as f:
+                file_data = f.read()
+                cmap = f.colormap(1)
+            self.assertEqual(file_data.shape, (1, 9, 9))  # no alpha band
+            exp = np.arange(81).reshape(9, 9, 1)
+            np.testing.assert_allclose(file_data[0], exp[:, :, 0])
+            self.assertEqual(cmap, exp_cmap)
+
+        # with trollimage colormap provided
+        from trollimage.colormap import Colormap
+        t_cmap = Colormap(*tuple((i, (i, i, i)) for i in range(20)))
+        exp_cmap = {i: (int(i * 255 / 19), int(i * 255 / 19), int(i * 255 / 19), 255) for i in range(20)}
+        exp_cmap.update({i: (0, 0, 0, 255) for i in range(20, 256)})
+        data = xr.DataArray(da.from_array(np.arange(81).reshape(9, 9, 1), chunks=9),
+                            dims=['y', 'x', 'bands'],
+                            coords={'bands': ['P']})
+        img = xrimage.XRImage(data)
+        with NamedTemporaryFile(suffix='.tif') as tmp:
+            img.save(tmp.name, keep_palette=True, cmap=t_cmap)
+            with rio.open(tmp.name) as f:
+                file_data = f.read()
+                cmap = f.colormap(1)
+            self.assertEqual(file_data.shape, (1, 9, 9))  # no alpha band
+            exp = np.arange(81).reshape(9, 9, 1)
+            np.testing.assert_allclose(file_data[0], exp[:, :, 0])
+            self.assertEqual(cmap, exp_cmap)
+
+        # with bad colormap provided
+        bad_cmap = [[i, [i, i, i]] for i in range(256)]
+        data = xr.DataArray(da.from_array(np.arange(81).reshape(9, 9, 1), chunks=9),
+                            dims=['y', 'x', 'bands'],
+                            coords={'bands': ['P']})
+        img = xrimage.XRImage(data)
+        with NamedTemporaryFile(suffix='.tif') as tmp:
+            self.assertRaises(ValueError, img.save, tmp.name,
+                              keep_palette=True, cmap=bad_cmap)
+            self.assertRaises(ValueError, img.save, tmp.name,
+                              keep_palette=True, cmap=t_cmap, dtype='uint16')
+
         # with input fill value
         data = np.arange(75).reshape(5, 5, 3)
         # second pixel is all bad
@@ -1083,6 +1185,45 @@ class TestXRImage(unittest.TestCase):
             np.testing.assert_allclose(file_data[2], exp[:, :, 2])
             np.testing.assert_allclose(file_data[3], exp_alpha)
 
+    @unittest.skipIf(sys.platform.startswith('win'), "'NamedTemporaryFile' not supported on Windows")
+    def test_save_jp2_int(self):
+        """Test saving jp2000 when input data is int."""
+        import xarray as xr
+        from trollimage import xrimage
+        import rasterio as rio
+
+        # numpy array image
+        data = xr.DataArray(np.arange(75).reshape(5, 5, 3), dims=[
+            'y', 'x', 'bands'], coords={'bands': ['R', 'G', 'B']})
+        img = xrimage.XRImage(data)
+        self.assertTrue(np.issubdtype(img.data.dtype, np.integer))
+        with NamedTemporaryFile(suffix='.jp2') as tmp:
+            img.save(tmp.name, quality=100, reversible=True)
+            with rio.open(tmp.name) as f:
+                file_data = f.read()
+            self.assertEqual(file_data.shape, (4, 5, 5))  # alpha band added
+            exp = np.arange(75).reshape(5, 5, 3)
+            np.testing.assert_allclose(file_data[0], exp[:, :, 0])
+            np.testing.assert_allclose(file_data[1], exp[:, :, 1])
+            np.testing.assert_allclose(file_data[2], exp[:, :, 2])
+            np.testing.assert_allclose(file_data[3], 255)
+
+    @unittest.skipIf(sys.platform.startswith('win'), "'NamedTemporaryFile' not supported on Windows")
+    def test_save_overviews(self):
+        """Test saving geotiffs with overviews."""
+        import xarray as xr
+        from trollimage import xrimage
+        import rasterio as rio
+
+        # numpy array image
+        data = xr.DataArray(np.arange(75).reshape(5, 5, 3), dims=[
+            'y', 'x', 'bands'], coords={'bands': ['R', 'G', 'B']})
+        img = xrimage.XRImage(data)
+        self.assertTrue(np.issubdtype(img.data.dtype, np.integer))
+        with NamedTemporaryFile(suffix='.tif') as tmp:
+            img.save(tmp.name, overviews=[2, 4])
+            with rio.open(tmp.name) as f:
+                self.assertEqual(len(f.overviews(1)), 2)
 
     def test_gamma(self):
         """Test gamma correction."""
@@ -1335,7 +1476,13 @@ class TestXRImage(unittest.TestCase):
         import dask
         import xarray as xr
         from trollimage import xrimage
-        from trollimage.colormap import brbg
+        from trollimage.colormap import brbg, Colormap
+
+        # RGBA colormap
+        bw = Colormap(
+            (0.0, (1.0, 1.0, 1.0, 1.0)),
+            (1.0, (0.0, 0.0, 0.0, 0.5)),
+        )
 
         arr1 = np.arange(150).reshape(1, 15, 10) / 150.
         arr2 = np.append(arr1, np.ones(150).reshape(arr1.shape)).reshape(2, 15, 10)
@@ -1431,10 +1578,10 @@ class TestXRImage(unittest.TestCase):
             img.palettize(brbg)
             pal = img.palette
 
-            img = img.convert('RGBA')
-            self.assertTrue(np.issubdtype(img.data.dtype, np.floating))
-            self.assertTrue(img.mode == 'RGBA')
-            self.assertTrue(len(img.data.coords['bands']) == 4)
+            img2 = img.convert('RGBA')
+            self.assertTrue(np.issubdtype(img2.data.dtype, np.floating))
+            self.assertTrue(img2.mode == 'RGBA')
+            self.assertTrue(len(img2.data.coords['bands']) == 4)
 
         # PA -> RGB (float)
         img = xrimage.XRImage(dataset3)
@@ -1447,7 +1594,23 @@ class TestXRImage(unittest.TestCase):
 
         self.assertRaises(ValueError, img.convert, 'A')
 
+        # L -> palettize -> RGBA (float) with RGBA colormap
+        with dask.config.set(scheduler=CustomScheduler(max_computes=0)):
+            img = xrimage.XRImage(dataset1)
+            img.palettize(bw)
+
+            img2 = img.convert('RGBA')
+            self.assertTrue(np.issubdtype(img2.data.dtype, np.floating))
+            self.assertTrue(img2.mode == 'RGBA')
+            self.assertTrue(len(img2.data.coords['bands']) == 4)
+            # convert to RGB, use RGBA from colormap regardless
+            img2 = img.convert('RGB')
+            self.assertTrue(np.issubdtype(img2.data.dtype, np.floating))
+            self.assertTrue(img2.mode == 'RGBA')
+            self.assertTrue(len(img2.data.coords['bands']) == 4)
+
     def test_colorize(self):
+        """Test colorize with an RGB colormap."""
         import xarray as xr
         from trollimage import xrimage
         from trollimage.colormap import brbg
@@ -1553,7 +1716,29 @@ class TestXRImage(unittest.TestCase):
                                    alpha.reshape((1,) + alpha.shape)))
         np.testing.assert_allclose(values, expected)
 
+    def test_colorize_rgba(self):
+        """Test colorize with an RGBA colormap."""
+        import xarray as xr
+        from trollimage import xrimage
+        from trollimage.colormap import Colormap
+
+        # RGBA colormap
+        bw = Colormap(
+            (0.0, (1.0, 1.0, 1.0, 1.0)),
+            (1.0, (0.0, 0.0, 0.0, 0.5)),
+        )
+
+        arr = np.arange(75).reshape(5, 15) / 74.
+        data = xr.DataArray(arr.copy(), dims=['y', 'x'])
+        img = xrimage.XRImage(data)
+        img.colorize(bw)
+        values = img.data.compute()
+        self.assertTupleEqual((4, 5, 15), values.shape)
+        np.testing.assert_allclose(values[:, 0, 0], [1.0, 1.0, 1.0, 1.0], rtol=1e-03)
+        np.testing.assert_allclose(values[:, -1, -1], [0.0, 0.0, 0.0, 0.5])
+
     def test_palettize(self):
+        """Test palettize with an RGB colormap."""
         import xarray as xr
         from trollimage import xrimage
         from trollimage.colormap import brbg
@@ -1572,6 +1757,27 @@ class TestXRImage(unittest.TestCase):
             [8, 8, 8, 8, 8, 8, 8, 9, 9, 9, 9, 9, 9, 9, 10]]])
         np.testing.assert_allclose(values, expected)
 
+    def test_palettize_rgba(self):
+        """Test palettize with an RGBA colormap."""
+        import xarray as xr
+        from trollimage import xrimage
+        from trollimage.colormap import Colormap
+
+        # RGBA colormap
+        bw = Colormap(
+            (0.0, (1.0, 1.0, 1.0, 1.0)),
+            (1.0, (0.0, 0.0, 0.0, 0.5)),
+        )
+
+        arr = np.arange(75).reshape(5, 15) / 74.
+        data = xr.DataArray(arr.copy(), dims=['y', 'x'])
+        img = xrimage.XRImage(data)
+        img.palettize(bw)
+
+        values = img.data.values
+        self.assertTupleEqual((1, 5, 15), values.shape)
+        self.assertTupleEqual((2, 4), bw.colors.shape)
+
     def test_merge(self):
         pass
 

trollimage/version.py

@@ -23,9 +23,9 @@ def get_keywords():
     # setup.py/versioneer.py will grep for the variable names, so they must
     # each be defined on a line of their own. _version.py will just call
     # get_keywords().
-    git_refnames = " (HEAD -> master, tag: v1.6.3)"
-    git_full = "69942415f7bfbd9160189534ca60782eaa2e4b59"
-    git_date = "2018-12-20 15:10:36 -0600"
+    git_refnames = " (HEAD -> master, tag: v1.7.0)"
+    git_full = "d35a7665ad475ff230e457085523e21f2cd3f454"
+    git_date = "2019-02-28 12:49:51 -0600"
     keywords = {"refnames": git_refnames, "full": git_full, "date": git_date}
     return keywords
 

trollimage/xrimage.py

@@ -47,7 +47,6 @@ try:
 except ImportError:
     rasterio = None
 
-
 try:
     # rasterio 1.0+
     from rasterio.windows import Window
@@ -58,7 +57,6 @@ except ImportError:
         """Replace the missing Window object in rasterio < 1.0."""
         return (y_off, y_off + y_size), (x_off, x_off + x_size)
 
-
 logger = logging.getLogger(__name__)
 
 
@@ -72,6 +70,7 @@ class RIOFile(object):
         self.kwargs = kwargs
         self.rfile = None
         self._closed = True
+        self.overviews = kwargs.pop('overviews', None)
 
     def __setitem__(self, key, item):
         """Put the data chunk in the image."""
@@ -104,6 +103,9 @@ class RIOFile(object):
 
     def close(self):
         if not self._closed:
+            if self.overviews:
+                logger.debug('Building overviews %s', str(self.overviews))
+                self.rfile.build_overviews(self.overviews)
             self.rfile.close()
             self._closed = True
 
@@ -164,7 +166,11 @@ def color_interp(data):
 
 
 class XRImage(object):
-    """Image class using an :class:`xarray.DataArray` as internal storage."""
+    """Image class using an :class:`xarray.DataArray` as internal storage.
+
+    It can be saved to a variety of image formats, but if Rasterio is installed,
+    it can save to geotiff and jpeg2000 with geographical information.
+    """
 
     def __init__(self, data):
         """Initialize the image with a :class:`~xarray.DataArray`."""
@@ -219,7 +225,7 @@ class XRImage(object):
         return ''.join(self.data['bands'].values)
 
     def save(self, filename, fformat=None, fill_value=None, compute=True,
-             **format_kwargs):
+             keep_palette=False, cmap=None, **format_kwargs):
         """Save the image to the given *filename*.
 
         Args:
@@ -229,6 +235,9 @@ class XRImage(object):
                            `rasterio` or `PIL` libraries ('jpg', 'png',
                            'tif'). By default this is determined by the
                            extension of the provided filename.
+                           If the format allows, geographical information will
+                           be saved to the ouput file, in the form of grid
+                           mapping or ground control points.
             fill_value (float): Replace invalid data values with this value
                                 and do not produce an Alpha band. Default
                                 behavior is to create an alpha band.
@@ -237,6 +246,11 @@ class XRImage(object):
                             a `dask.Delayed` object or a tuple of
                             ``(source, target)`` to be passed to
                             `dask.array.store`.
+            keep_palette (bool): Saves the palettized version of the image if
+                                 set to True. False by default.
+            cmap (Colormap or dict): Colormap to be applied to the image when
+                                     saving with rasterio, used with
+                                     keep_palette=True. Should be uint8.
             format_kwargs: Additional format options to pass to `rasterio`
                            or `PIL` saving methods.
 
@@ -250,33 +264,45 @@ class XRImage(object):
 
         """
         fformat = fformat or os.path.splitext(filename)[1][1:4]
-        if fformat == 'tif' and rasterio:
+        if fformat in ('tif', 'jp2') and rasterio:
             return self.rio_save(filename, fformat=fformat,
                                  fill_value=fill_value, compute=compute,
+                                 keep_palette=keep_palette, cmap=cmap,
                                  **format_kwargs)
         else:
             return self.pil_save(filename, fformat, fill_value,
                                  compute=compute, **format_kwargs)
 
     def rio_save(self, filename, fformat=None, fill_value=None,
-                 dtype=np.uint8, compute=True, tags=None, **format_kwargs):
-        """Save the image using rasterio."""
+                 dtype=np.uint8, compute=True, tags=None,
+                 keep_palette=False, cmap=None,
+                 **format_kwargs):
+        """Save the image using rasterio.
+
+        Overviews can be added to the file using the `overviews` kwarg, eg::
+
+          img.rio_save('myfile.tif', overviews=[2, 4, 8, 16])
+
+        """
         fformat = fformat or os.path.splitext(filename)[1][1:4]
         drivers = {'jpg': 'JPEG',
                    'png': 'PNG',
-                   'tif': 'GTiff'}
+                   'tif': 'GTiff',
+                   'jp2': 'JP2OpenJPEG'}
         driver = drivers.get(fformat, fformat)
 
         if tags is None:
             tags = {}
 
-        data, mode = self.finalize(fill_value, dtype=dtype)
+        data, mode = self.finalize(fill_value, dtype=dtype,
+                                   keep_palette=keep_palette, cmap=cmap)
         data = data.transpose('bands', 'y', 'x')
         data.attrs = self.data.attrs
 
         crs = None
+        gcps = None
         transform = None
-        if driver == 'GTiff':
+        if driver in ['GTiff', 'JP2OpenJPEG']:
             if not np.issubdtype(data.dtype, np.floating):
                 format_kwargs.setdefault('compress', 'DEFLATE')
             photometric_map = {
@@ -298,13 +324,20 @@ class XRImage(object):
                 transform = rasterio.transform.from_bounds(west, south,
                                                            east, north,
                                                            width, height)
+
+            except KeyError:  # No area
+                logger.info("Couldn't create geotransform")
+            except AttributeError:
+                try:
+                    gcps = data.attrs['area'].lons.attrs['gcps']
+                    crs = data.attrs['area'].lons.attrs['crs']
+                except KeyError:
+                    logger.info("Couldn't create geotransform")
+
             if "start_time" in data.attrs:
                 stime = data.attrs['start_time']
                 stime_str = stime.strftime("%Y:%m:%d %H:%M:%S")
                 tags.setdefault('TIFFTAG_DATETIME', stime_str)
-
-            except (KeyError, AttributeError):
-                logger.info("Couldn't create geotransform")
         elif driver == 'JPEG' and 'A' in mode:
             raise ValueError('JPEG does not support alpha')
 
@@ -314,11 +347,25 @@ class XRImage(object):
                          count=data.sizes['bands'],
                          dtype=dtype,
                          nodata=fill_value,
-                         crs=crs, transform=transform, **format_kwargs)
+                         crs=crs,
+                         transform=transform,
+                         gcps=gcps,
+                         **format_kwargs)
         r_file.open()
+        if not keep_palette:
             r_file.colorinterp = color_interp(data)
         r_file.rfile.update_tags(**tags)
 
+        if keep_palette and cmap is not None:
+            if data.dtype != 'uint8':
+                raise ValueError('Rasterio only supports 8-bit colormaps')
+            try:
+                from trollimage.colormap import Colormap
+                cmap = cmap.to_rio() if isinstance(cmap, Colormap) else cmap
+                r_file.rfile.write_colormap(1, cmap)
+            except AttributeError:
+                raise ValueError("Colormap is not formatted correctly")
+
         if compute:
             # write data to the file now
             res = da.store(data.data, r_file)
@@ -343,11 +390,8 @@ class XRImage(object):
             # Take care of GeoImage.tags (if any).
             format_kwargs['pnginfo'] = self._pngmeta()
 
-        def _create_save_image(fill_value, filename, fformat, format_kwargs):
-            img = self.pil_image(fill_value)
-            img.save(filename, fformat, **format_kwargs)
-        delay = dask.delayed(_create_save_image)(
-            fill_value, filename, fformat, format_kwargs)
+        img = self.pil_image(fill_value, compute=False)
+        delay = img.save(filename, fformat, **format_kwargs)
         if compute:
             return delay.compute()
         return delay
@@ -449,28 +493,36 @@ class XRImage(object):
         """Convert the image from P or PA to RGB or RGBA."""
         self._check_modes(("P", "PA"))
 
-        if self.mode.endswith("A"):
-            alpha = self.data.sel(bands=["A"]).data
-            mode = mode + "A" if not mode.endswith("A") else mode
-        else:
-            alpha = None
-
         if not self.palette:
             raise RuntimeError("Can't convert palettized image, missing palette.")
-
         pal = np.array(self.palette)
         pal = da.from_array(pal, chunks=pal.shape)
-        flat_indexes = self.data.data[0].ravel().astype('int64')
-        new_shape = (3,) + self.data.shape[1:3]
+
+        if pal.shape[1] == 4:
+            # colormap's alpha overrides data alpha
+            mode = "RGBA"
+            alpha = None
+        elif self.mode.endswith("A"):
+            # add a new/fake 'bands' dimension to the end
+            alpha = self.data.sel(bands="A").data[..., None]
+            mode = mode + "A" if not mode.endswith("A") else mode
+        else:
+            alpha = None
+
+        flat_indexes = self.data.sel(bands='P').data.ravel().astype('int64')
+        dim_sizes = ((key, val) for key, val in self.data.sizes.items() if key != 'bands')
+        dims, new_shape = zip(*dim_sizes)
+        dims = dims + ('bands',)
+        new_shape = new_shape + (pal.shape[1],)
         new_data = pal[flat_indexes].reshape(new_shape)
         coords = dict(self.data.coords)
         coords["bands"] = list(mode)
 
         if alpha is not None:
-            new_arr = da.concatenate((new_data, alpha), axis=0)
-            data = xr.DataArray(new_arr, coords=coords, attrs=self.data.attrs, dims=self.data.dims)
+            new_arr = da.concatenate((new_data, alpha), axis=-1)
+            data = xr.DataArray(new_arr, coords=coords, attrs=self.data.attrs, dims=dims)
         else:
-            data = xr.DataArray(new_data, coords=coords, attrs=self.data.attrs, dims=self.data.dims)
+            data = xr.DataArray(new_data, coords=coords, attrs=self.data.attrs, dims=dims)
 
         return data
 
@@ -531,14 +583,14 @@ class XRImage(object):
             new_img.palette = self.palette
         return new_img
 
-    def _finalize(self, fill_value=None, dtype=np.uint8):
+    def _finalize(self, fill_value=None, dtype=np.uint8, keep_palette=False, cmap=None):
         """Wrapper around 'finalize' method for backwards compatibility."""
         import warnings
         warnings.warn("'_finalize' is deprecated, use 'finalize' instead.",
                       DeprecationWarning)
-        return self.finalize(fill_value, dtype)
+        return self.finalize(fill_value, dtype, keep_palette, cmap)
 
-    def finalize(self, fill_value=None, dtype=np.uint8):
+    def finalize(self, fill_value=None, dtype=np.uint8, keep_palette=False, cmap=None):
         """Finalize the image to be written to an output file.
 
         This adds an alpha band or fills data with a fill_value (if specified).
@@ -550,10 +602,16 @@ class XRImage(object):
         in the ``DataArray`` ``.attrs`` dictionary.
 
         """
+        if keep_palette and not self.mode.startswith('P'):
+            keep_palette = False
+
+        if not keep_palette:
             if self.mode == "P":
-            return self.convert("RGB").finalize(fill_value=fill_value, dtype=dtype)
+                return self.convert("RGB").finalize(fill_value=fill_value, dtype=dtype,
+                                                    keep_palette=keep_palette, cmap=cmap)
             if self.mode == "PA":
-            return self.convert("RGBA").finalize(fill_value=fill_value, dtype=dtype)
+                return self.convert("RGBA").finalize(fill_value=fill_value, dtype=dtype,
+                                                     keep_palette=keep_palette, cmap=cmap)
 
         if np.issubdtype(dtype, np.floating) and fill_value is None:
             logger.warning("Image with floats cannot be transparent, so "
@@ -563,6 +621,7 @@ class XRImage(object):
         final_data = self.data
         # if the data are integers then this fill value will be used to check for invalid values
         ifill = final_data.attrs.get('_FillValue') if np.issubdtype(final_data, np.integer) else None
+        if not keep_palette:
             if fill_value is None and not self.mode.endswith('A'):
                 # We don't have a fill value or an alpha, let's add an alpha
                 alpha = self._create_alpha(final_data, fill_value=ifill)
@@ -586,11 +645,24 @@ class XRImage(object):
         final_data.attrs = self.data.attrs
         return final_data, ''.join(final_data['bands'].values)
 
-    def pil_image(self, fill_value=None):
-        """Return a PIL image from the current image."""
+    def pil_image(self, fill_value=None, compute=True):
+        """Return a PIL image from the current image.
+
+        Args:
+            fill_value (int or float): Value to use for NaN null values.
+                See :meth:`~trollimage.xrimage.XRImage.finalize` for more
+                info.
+            compute (bool): Whether to return a fully computed PIL.Image
+                object (True) or return a dask Delayed object representing
+                the Image (False). This is True by default.
+
+        """
         channels, mode = self.finalize(fill_value)
-        res = np.asanyarray(channels.transpose('y', 'x', 'bands').values)
-        return PILImage.fromarray(np.squeeze(res), mode)
+        res = channels.transpose('y', 'x', 'bands')
+        img = dask.delayed(PILImage.fromarray)(np.squeeze(res.data), mode)
+        if compute:
+            img = img.compute()
+        return img
 
     def xrify_tuples(self, tup):
         """Make xarray.DataArray from tuple."""
@@ -657,6 +729,25 @@ class XRImage(object):
         else:
             raise TypeError("Stretch parameter must be a string or a tuple.")
 
+    @staticmethod
+    def _compute_quantile(data, dims, cutoffs):
+        """Helper method for stretch_linear.
+
+        Dask delayed functions need to be non-internal functions (created
+        inside a function) to be serializable on a multi-process scheduler.
+
+        Quantile requires the data to be loaded since it not supported on
+        dask arrays yet.
+
+        """
+        # numpy doesn't get a 'quantile' function until 1.15
+        # for better backwards compatibility we use xarray's version
+        data_arr = xr.DataArray(data, dims=dims)
+        # delayed will provide us the fully computed xarray with ndarray
+        left, right = data_arr.quantile([cutoffs[0], 1. - cutoffs[1]], dim=['x', 'y'])
+        logger.debug("Interval: left=%s, right=%s", str(left), str(right))
+        return left.data, right.data
+
     def stretch_linear(self, cutoffs=(0.005, 0.005)):
         """Stretch linearly the contrast of the current image.
 
@@ -668,21 +759,13 @@ class XRImage(object):
         logger.debug("Left and right quantiles: " +
                      str(cutoffs[0]) + " " + str(cutoffs[1]))
 
-        # Quantile requires the data to be loaded, not supported on dask arrays
-        def _compute_quantile(data, cutoffs):
-            # delayed will provide us the fully computed xarray with ndarray
-            left, right = data.quantile([cutoffs[0], 1. - cutoffs[1]],
-                                        dim=['x', 'y'])
-            logger.debug("Interval: left=%s, right=%s", str(left), str(right))
-            return left.data, right.data
-
         cutoff_type = np.float64
         # numpy percentile (which quantile calls) returns 64-bit floats
         # unless the value is a higher order float
         if np.issubdtype(self.data.dtype, np.floating) and \
                 np.dtype(self.data.dtype).itemsize > 8:
             cutoff_type = self.data.dtype
-        left, right = dask.delayed(_compute_quantile, nout=2)(self.data, cutoffs)
+        left, right = dask.delayed(self._compute_quantile, nout=2)(self.data.data, self.data.dims, cutoffs)
         left_data = da.from_delayed(left,
                                     shape=(self.data.sizes['bands'],),
                                     dtype=cutoff_type)
@@ -880,9 +963,12 @@ class XRImage(object):
             return np.concatenate(channels, axis=0)
 
         new_data = l_data.data.map_blocks(_colorize, colormap,
-                                          chunks=(3,) + l_data.data.chunks[1:], dtype=np.float64)
+                                          chunks=(colormap.colors.shape[1],) + l_data.data.chunks[1:],
+                                          dtype=np.float64)
 
-        if alpha is not None:
+        if colormap.colors.shape[1] == 4:
+            mode = "RGBA"
+        elif alpha is not None:
             new_data = da.concatenate([new_data, alpha.data], axis=0)
             mode = "RGBA"
         else: