diff --git a/satpy/tests/writer_tests/test_mitiff.py b/satpy/tests/writer_tests/test_mitiff.py index 9c0ef86a39..a54782889e 100644 --- a/satpy/tests/writer_tests/test_mitiff.py +++ b/satpy/tests/writer_tests/test_mitiff.py @@ -691,6 +691,93 @@ def test_convert_proj4_string(self): proj4_string = w._add_proj4_string(ds1, ds1) self.assertEqual(proj4_string, check['proj4']) + def test_save_dataset_palette(self): + """Test writer operation as palette.""" + import os + import numpy as np + from libtiff import TIFF + from satpy.writers.mitiff import MITIFFWriter + + expected = np.full((100, 200), 0) + + exp_c = ([0, 3, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], + [1, 4, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], + [2, 5, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]) + + color_map = [[0, 3], [1, 4], [2, 5]] + pal_desc = ['test', 'test2'] + unit = "Test" + + dataset = self._get_test_one_dataset() + palette = {'palette': True, + 'palette_color_map': color_map, + 'palette_description': pal_desc, + 'palette_unit': unit, + 'palette_channel_name': dataset.attrs['name']} + + w = MITIFFWriter(base_dir=self.base_dir) + w.save_dataset(dataset, **palette) + filename = "{:s}_{:%Y%m%d_%H%M%S}.mitiff".format(dataset.attrs['name'], + dataset.attrs['start_time']) + tif = TIFF.open(os.path.join(self.base_dir, filename)) + # Need to check PHOTOMETRIC is 3, ie palette + self.assertEqual(tif.GetField('PHOTOMETRIC'), 3) + colormap = tif.GetField('COLORMAP') + # Check the colormap of the palette image + self.assertEqual(colormap, exp_c) + IMAGEDESCRIPTION = 270 + imgdesc = (tif.GetField(IMAGEDESCRIPTION)).decode('utf-8').split('\n') + found_color_info = False + unit_name_found = False + name_length_found = False + name_length = 0 + names = [] + unit_name = None + for key in imgdesc: + if name_length_found and name_length > len(names): + names.append(key) + continue + elif unit_name_found: + name_length = int(key) + name_length_found = True + unit_name_found = False + elif found_color_info: + unit_name = key + unit_name_found = True + found_color_info = False + elif 'COLOR INFO:' in key: + found_color_info = True + # Check the name of the palette description + self.assertEqual(name_length, 2) + # Check the name and unit name of the palette + self.assertEqual(unit_name, ' Test') + # Check the palette description of the palette + self.assertEqual(names, [' test', ' test2']) + for image in tif.iter_images(): + np.testing.assert_allclose(image, expected, atol=1.e-6, rtol=0) + def suite(): """The test suite for this writer's tests. diff --git a/satpy/writers/mitiff.py b/satpy/writers/mitiff.py index b720eb5edc..df0c91924b 100644 --- a/satpy/writers/mitiff.py +++ b/satpy/writers/mitiff.py @@ -1,6 +1,6 @@ #!/usr/bin/env python # -*- coding: utf-8 -*- -# Copyright (c) 2018 Satpy developers +# Copyright (c) 2018, 2019 Satpy developers # # This file is part of satpy. # @@ -20,7 +20,6 @@ """ import logging - import numpy as np from satpy.writers import ImageWriter @@ -51,6 +50,7 @@ def __init__(self, name=None, tags=None, **kwargs): self.mitiff_config = {} self.translate_channel_name = {} self.channel_order = {} + self.palette = False def save_image(self): raise NotImplementedError("save_image mitiff is not implemented.") @@ -61,6 +61,8 @@ def save_dataset(self, dataset, filename=None, fill_value=None, def _delayed_create(create_opts, dataset): try: + if 'palette' in kwargs: + self.palette = kwargs['palette'] if 'platform_name' not in kwargs: kwargs['platform_name'] = dataset.attrs['platform_name'] if 'name' not in kwargs: @@ -160,6 +162,11 @@ def _make_channel_list(self, datasets, **kwargs): channels.append( ds.attrs['prerequisites'][ch][0]) break + elif self.palette: + if 'palette_channel_name' in kwargs: + channels.append(kwargs['palette_channel_name'].upper()) + else: + LOG.error("Is palette but can not find palette_channel_name to name the dataset") else: for ch, ds in enumerate(datasets): channels.append(ch + 1) @@ -244,8 +251,8 @@ def _add_proj4_string(self, datasets, first_dataset): if '+a=6378137.0 +b=6356752.31414' in proj4_string: proj4_string = proj4_string.replace("+a=6378137.0 +b=6356752.31414", "+ellps=WGS84") - if '+units=m' in proj4_string: - proj4_string = proj4_string.replace("+units=m", "+units=km") + if '+units=m' in proj4_string: + proj4_string = proj4_string.replace("+units=m", "+units=km") if not any(datum in proj4_string for datum in ['datum', 'towgs84']): proj4_string += ' +towgs84=0,0,0' @@ -267,7 +274,20 @@ def _add_proj4_string(self, datasets, first_dataset): proj4_string += ' +y_0=%.6f' % ( (-datasets.attrs['area'].area_extent[1] + datasets.attrs['area'].pixel_size_y) + y_0) - + elif '+x_0=0' in proj4_string and '+y_0=0' in proj4_string and isinstance(datasets, list): + proj4_string = proj4_string.replace("+x_0=0", '+x_0=%.6f' % ( + (-first_dataset.attrs['area'].area_extent[0] + + first_dataset.attrs['area'].pixel_size_x) + x_0)) + proj4_string = proj4_string.replace("+y_0=0", '+y_0=%.6f' % ( + (-first_dataset.attrs['area'].area_extent[1] + + first_dataset.attrs['area'].pixel_size_y) + y_0)) + elif '+x_0=0' in proj4_string and '+y_0=0' in proj4_string: + proj4_string = proj4_string.replace("+x_0=0", '+x_0=%.6f' % ( + (-datasets.attrs['area'].area_extent[0] + + datasets.attrs['area'].pixel_size_x) + x_0)) + proj4_string = proj4_string.replace("+y_0=0", '+y_0=%.6f' % ( + (-datasets.attrs['area'].area_extent[1] + + datasets.attrs['area'].pixel_size_y) + y_0)) LOG.debug("proj4_string: %s", proj4_string) proj4_string += '\n' @@ -362,6 +382,17 @@ def _add_calibration_datasets(self, ch, datasets, reverse_offset, reverse_scale, # http://stackoverflow.com/questions/1598579/rounding-decimals-with-new-python-format-function return skip_calibration, _table_calibration, _reverse_offset, _reverse_scale, _decimals + def _add_palette_info(self, datasets, palette_unit, palette_description, **kwargs): + # mitiff key word for palette interpretion + _palette = '\n COLOR INFO:\n' + # mitiff info for the unit of the interpretion + _palette += ' {}\n'.format(palette_unit) + # The length of the palette description as needed by mitiff in DIANA + _palette += ' {}\n'.format(len(palette_description)) + for desc in palette_description: + _palette += ' {}\n'.format(desc) + return _palette + def _add_calibration(self, channels, cns, datasets, **kwargs): _table_calibration = "" @@ -547,7 +578,11 @@ def _make_image_description(self, datasets, **kwargs): else: LOG.debug("Area extent: %s", datasets.attrs['area'].area_extent) - _image_description += self._add_calibration(channels, cns, datasets, **kwargs) + if self.palette: + LOG.debug("Doing palette image") + _image_description += self._add_palette_info(datasets, **kwargs) + else: + _image_description += self._add_calibration(channels, cns, datasets, **kwargs) return _image_description @@ -565,6 +600,48 @@ def _calibrate_data(self, dataset, calibration, min_val, max_val): (float(max_val) - float(min_val))) * 255. return _data.clip(0, 255) + def _save_as_palette(self, tif, datasets, **kwargs): + # MITIFF palette has only one data channel + if len(datasets.dims) == 2: + LOG.debug("Palette ok with only 2 dimensions. ie only x and y") + # 3 = Palette color. In this model, a color is described with a single component. + # The value of the component is used as an index into the red, green and blue curves + # in the ColorMap field to retrieve an RGB triplet that defines the color. When + # PhotometricInterpretation=3 is used, ColorMap must be present and SamplesPerPixel must be 1. + tif.SetField('PHOTOMETRIC', 3) + + # As write_image can not save tiff image as palette, this has to be done basicly + # ie. all needed tags needs to be set. + tif.SetField('IMAGEWIDTH', datasets.sizes['x']) + tif.SetField('IMAGELENGTH', datasets.sizes['y']) + tif.SetField('BITSPERSAMPLE', 8) + tif.SetField('COMPRESSION', tif.get_tag_define('deflate')) + + if 'palette_color_map' in kwargs: + tif.SetField('COLORMAP', kwargs['palette_color_map']) + else: + LOG.ERROR("In a mitiff palette image a color map must be provided: palette_color_map is missing.") + + data_type = np.uint8 + # Looks like we need to pass the data to writeencodedstrip as ctypes + tif.WriteEncodedStrip(0, np.ascontiguousarray(datasets.data.astype(data_type), data_type).ctypes.data, + datasets.sizes['x'] * datasets.sizes['y']) + tif.WriteDirectory() + + def _save_as_enhanced(self, tif, datasets, **kwargs): + """Save datasets as an enhanced RGB image""" + img = get_enhanced_image(datasets.squeeze(), enhance=self.enhancer) + if 'bands' in img.data.sizes and 'bands' not in datasets.sizes: + LOG.debug("Datasets without 'bands' become image with 'bands' due to enhancement.") + LOG.debug("Needs to regenerate mitiff image description") + image_description = self._make_image_description(img.data, **kwargs) + tif.SetField(IMAGEDESCRIPTION, (image_description).encode('utf-8')) + for i, band in enumerate(img.data['bands']): + chn = img.data.sel(bands=band) + data = chn.values.clip(0, 1) * 254. + 1 + data = data.clip(0, 255) + tif.write_image(data.astype(np.uint8), compression='deflate') + def _save_datasets_as_mitiff(self, datasets, image_description, gen_filename, **kwargs): """Put all togehter and save as a tiff file with the special tag @@ -572,7 +649,7 @@ def _save_datasets_as_mitiff(self, datasets, image_description, """ from libtiff import TIFF - tif = TIFF.open(gen_filename, mode='w') + tif = TIFF.open(gen_filename, mode='wb') tif.SetField(IMAGEDESCRIPTION, (image_description).encode('utf-8')) @@ -618,18 +695,10 @@ def _save_datasets_as_mitiff(self, datasets, image_description, tif.write_image(data.astype(np.uint8), compression='deflate') break + elif self.palette: + LOG.debug("Saving dataset as palette.") + self._save_as_palette(tif, datasets, **kwargs) else: LOG.debug("Saving datasets as enhanced image") - img = get_enhanced_image(datasets.squeeze(), enhance=self.enhancer) - if 'bands' in img.data.sizes and 'bands' not in datasets.sizes: - LOG.debug("Datasets without 'bands' become image with 'bands' due to enhancement.") - LOG.debug("Needs to regenerate mitiff image description") - image_description = self._make_image_description(img.data, **kwargs) - tif.SetField(IMAGEDESCRIPTION, (image_description).encode('utf-8')) - for i, band in enumerate(img.data['bands']): - chn = img.data.sel(bands=band) - data = chn.values.clip(0, 1) * 254. + 1 - data = data.clip(0, 255) - tif.write_image(data.astype(np.uint8), compression='deflate') - + self._save_as_enhanced(tif, datasets, **kwargs) del tif