New save capture fuction

[mdanilevicz]

Proposed new function to save processed/aligned captures in separate tifs (per spectral band) instead of stacked tif. This was useful for me because OpenDroneMap did not read the stacked tif as input to assemble the orthomosaic. This function is based on micasense save_capture_as_stack.

[poynting]

Hi @moDanilevicz thanks for the new functionality. Mostly it looks fine; I have a couple of comments.

1. It appears to support 5-band and 6-band data, but not 10-band data. I think a similar approach to that taken in Capture.save_capture_as_stack() is appropriate using the eo_list and lwir_list.

2. Saving thermal data in centi-Celcius to UINT16 will cause rollover for data below 0C. I've experienced this case flying myself, especially in the springtime when shadows are still very cold. I think we'll need to either save to INT16_t or FLOAT32 to prevent this. Alternatively we could convert back to centi-Kelvin to avoid rollovers.

[rrdpereira]

Thanks, i started to think in that solution!

[and-viceversa]

@mdanilevicz @poynting

To meet the criteria above, something like this?

    def save_capture_as_bands(self, out_file_name, photometric='MINISBLACK'):
        from osgeo.gdal import GetDriverByName, GDT_UInt16, GDT_Float32
        if self.__aligned_capture is None:
            raise RuntimeError("Call Capture.create_aligned_capture() prior to saving as stack.")

        # handle accidental .tif.tif values predictably
        if out_file_name.endswith('.tif'):
            out_file_path = out_file_name[:-4]
        else:
            out_file_path = out_file_name

        rows, cols, bands = self.__aligned_capture.shape
        driver = GetDriverByName('GTiff')

        for i in self.eo_indices():
            out_raster = driver.Create(out_file_path + f'_{i + 1}.tif', cols, rows, 1, GDT_UInt16,
                                       options=['INTERLEAVE=BAND', 'COMPRESS=DEFLATE', f'PHOTOMETRIC={photometric}'])
            out_band = out_raster.GetRasterBand(1)
            out_data = self.__aligned_capture[:, :, i]
            out_data[out_data < 0] = 0
            out_data[out_data > 2] = 2  # limit reflectance data to 200% to allow some specular reflections
            out_band.WriteArray(out_data * 32768)  # scale reflectance images so 100% = 32768
            out_band.FlushCache()

        for i in self.lw_indices():
            out_raster = driver.Create(out_file_path + f'_{i + 1}.tif', cols, rows, 1, GDT_Float32,
                                       options=['INTERLEAVE=BAND', 'COMPRESS=DEFLATE', f'PHOTOMETRIC={photometric}'])
            out_band = out_raster.GetRasterBand(1)
            out_data = (self.__aligned_capture[:, :, i])
            out_band.WriteArray(out_data)
            out_band.FlushCache()

I'm not sure whether C or cK is better... maybe an additional parameter to choose?

[poynting]

Yes that looks pretty good; I would like if there were a type option - UINT16 or FLOAT32, and based on that option it either wrote the files out as scaled UINT16 reflectance (32768=100%) and cK, or floating point reflectance (1.0=100%) and floating-point C. This way the file type will always be ensured to be consistent, which I think is important for downstream usability.

[and-viceversa]

Looking better with an output data type option. Defaults to GDT_UInt16 to match Capture.save_capture_as_stack() behavior.

I understand how this is working for reflectance images... but how do the units match up in case of radiance images?

For example in #77:

The last thing will be to update capture.save_capture_as_stack to allow for radiance images and in this case change the stack output type to floating point, since radiance images won't have a fixed output range relative to reflectance images and won't intuitively fit into 16b (a single raw image can fit in 16b, but a group of raw images converted to radiance, can have a much larger dynamic range. For example one image may have gain 1 exposure 1ms, and another could have gain 8 exposure 10ms, and then changing 1 bit will represent a substantially different change in radiance.

Therefore if img_type == radiance do we need to force GDT_Float32 to write the correct units?

    def save_capture_as_bands(self, out_file_name, out_data_type='GDT_UInt16', photometric='MINISBLACK'):
        """
        Output the Images in the Capture object as separate GTiffs.
        :param out_file_name: str system file path without file extension
        :param out_data_type: str GDT_Float32 or GDT_UInt16
            Default: GDT_UInt16 will write images as scaled reflectance values. EO (32769=100%)
            and LWIR in centi-Kelvin (0-65535).
            GDT_Float32 will write images as floating point reflectance values. EO (1.0=100%)
            and LWIR in floating point Celsius.
            https://gdal.org/api/raster_c_api.html#_CPPv412GDALDataType
        :param photometric: str GDAL argument for GTiff color matching
        :return: None
        """
        from osgeo.gdal import GetDriverByName, GDT_UInt16, GDT_Float32
        if self.__aligned_capture is None:
            raise RuntimeError("Call Capture.create_aligned_capture() prior to saving as stack.")

        # handle incorrect datatype. default to GDT_UInt16.
        if out_data_type.strip() == 'GDT_UInt16':
            gdal_type = GDT_UInt16
        elif out_data_type.strip() == 'GDT_Float32':
            gdal_type = GDT_Float32
        else:
            warnings.warn(message='Output data type in Capture.save_capture_as_bands() was called as {}. '
                                  'Must use "GDT_UInt16" or "GDT_Float32". Defaulting to GDT_UInt16...'
                          .format(out_data_type),
                          category=UserWarning)
            gdal_type = GDT_UInt16

        # predictably handle accidental .tif.tif values
        if out_file_name.endswith('.tif'):
            out_file_path = out_file_name[:-4]
        else:
            out_file_path = out_file_name

        rows, cols, bands = self.__aligned_capture.shape
        driver = GetDriverByName('GTiff')

        for i in self.eo_indices():
            out_raster = driver.Create(out_file_path + f'_{i + 1}.tif', cols, rows, 1, gdal_type,
                                       options=['INTERLEAVE=BAND', 'COMPRESS=DEFLATE', f'PHOTOMETRIC={photometric}'])
            out_band = out_raster.GetRasterBand(1)
            out_data = self.__aligned_capture[:, :, i]
            out_data[out_data < 0] = 0
            out_data[out_data > 2] = 2  # limit reflectance data to 200% to allow some specular reflections
            # if GDT_UInt16, scale reflectance images so 100% = 32768. GDT_UInt16 resolves to 2.
            out_band.WriteArray(out_data * 32768 if gdal_type == 2 else out_data)
            out_band.FlushCache()

        for i in self.lw_indices():
            out_raster = driver.Create(out_file_path + f'_{i + 1}.tif', cols, rows, 1, gdal_type,
                                       options=['INTERLEAVE=BAND', 'COMPRESS=DEFLATE', f'PHOTOMETRIC={photometric}'])
            out_band = out_raster.GetRasterBand(1)
            # if GDT_UInt16, scale data from float degC to back to centi-Kelvin to fit into UInt16.
            out_data = (self.__aligned_capture[:, :, i] + 273.15) * 100 if gdal_type == 2 \
                else self.__aligned_capture[:, :, i]
            out_band.WriteArray(out_data)
            out_band.FlushCache()

[and-viceversa]

The more I look at this, it seems like a Capture should have an img_type object variable that gets assigned during Capture.create_aligned_capture(). That way, it is simple to enforce the correct dtype in the various save methods.

[poynting]

The more I look at this, it seems like a Capture should have an img_type object variable that gets assigned during Capture.create_aligned_capture(). That way, it is simple to enforce the correct dtype in the various save methods.

Yes, I think this is a good option, and yes, I think if the type is radiance, it's necessary to force to floating point. It creates an opportunity to confuse these files with floating point reflectance, but I think that's ok.