Currently I'm using skimage for all but the LCH colorspace (which is implemented directly in rio-color)

Related to https://github.com/mapbox/rio-color/issues/28, we might want to make rio-color a more full-fledged colorspace module.

Let's identify which color spaces are useful for the purposes of histogram matching to guide that decision.

/cc @celoyd @virginiayung

Ran successfully on your test data 1 but not test data 2 (Andros Island, The Bahamas). I tried to run on my own satellite imagery and got the same error as test data 2: AttributeError: "DatasetReader' object has no attribute 'mask flags'.

Beautiful approach to histogram matching by the way! I am trying to integrate histogram matching into a personal project and am running into an error. from rio_hist.match import hist_match_worker hist_match_worker(mypath, myref, mypath_out, 1,creation_options, '1,2', 'LCH', None) What are the creation_options in this case? Are they used for writing and reading the profiles in rasterio? Are there any options to send to this command to make this work? I tried a similar approach to how you tested the cli but with no success. Thank you!

Here is a basic implementation of matching 16-bit RGB image to an 8-bit target image in LCH space. The image is then saved back in 16-bit RGB with the bit-depth corresponding to that of the sensor and a black point offset to allow further adjustmests without too much clipping of shadows or highlights. The RGB->LCH->RGB transform is lossy in colour for 8-bit outputs, 16-bit preserves the LCH space better and has less quantization. The colour space transfroms are all derived from sckit-image.

   from skimage.color import rgb2lab, lab2lch, lch2lab, lab2rgb
    with rasterio.open(test_image) as src:
        rgb_src = src.read()
        profile = src.profile
    
    with rasterio.open(sample_tif) as target:
        rgb_sample = target.read()[:3]
    
    float_rgb_src = rgb_src/(2**16)
    float_rgb_sample = rgb_sample/(2**8)
    rgb_src_image = reshape_as_image(float_rgb_src)
    rgb_sample_image = reshape_as_image(float_rgb_sample)
    
    lch_src = lab2lch(rgb2lab(rgb_src_image))
    lch_sample = lab2lch(rgb2lab(rgb_sample_image))
    
    lch_target = lch_src.copy()
    """Adjust Luminance"""
    lch_target[:,:,0] = hist_match(lch_src[:,:,0],lch_sample[:,:,0])
    """Adjust Chroma"""
    lch_target[:,:,1] = hist_match(lch_src[:,:,1],lch_sample[:,:,1])
    rgb_target_image = lab2rgb(lch2lab(lch_target))
    
    float_rgb_target = reshape_as_raster(rgb_target_image)
    rgb_target = float_rgb_target*(2**bit_depth)+ black_pt

Specifying --match-proportion 1.0 will do a full histogram match (the default)

Specifying --match-proportion 0.0 will effectively not change anything.

For every value in between, it moves the source closer to the reference by linear interpolation on the cumulative distribution.

Resolves #8

todo

• [x] try with some other scenes; i'm sure there are edge cases that don't work well.
• [x] try with other colorspaces. only RGB so far.
• [ ] unit and integration tests

cc @celoyd @virginiayung

The colorspace is called Lab, sic, and spelling it LAB is likely to trip up people who are familiar with it.

Because the colorspace string is downcased internally, this is entirely cosmetic.

Well, the results and the analyses looks so good. But I wonder how to matching the histogram or if there are any papers described this algorithm. Thanks!

Hue-Saturation-Value colorspace has numerous problems for analytical purposes which are well articulated here: https://en.wikipedia.org/wiki/HSL_and_HSV#Disadvantages

TL;DR; there are quantifiably better alternatives (Lab, Lch and Luv) so why bother with HSV

Dropping HSV means we can drop skimage too!

Currently rio-hist can handle nodata in RGB rasters.

We also need to support

• RGBA Alpha bands
• internal GDAL masks
• external GDAL masks (.msk sidecar)

Suppose a user inputs an array where masked values are not all unique (e.g. user is matching histograms for .tif red raster band with nodata value of 0. When rastertio.open() reads the .tif, machine epsilon may code some zeros a little bit off, like say 0.00000000001. Then the mask will cover 0 and 0.00000000001). Then, before the insertion of this edit, mask_index[0] would contain more than one value, and everything breaks. My edit replaces all masked values with np.nan up front to ensure np.unique (line 53) records only one masked value, thus preserving the mask_index[0] system below.

Thank you very much for this nice repo. I am also working on normalization. My problem is the reference image. If the content of the reference image is different with the source image, how can I normalize brightness/illumination using histogram matching?

usage of np.iinfo forces the user to use only arrays (images) with integer data type in cs_forward function in rio_hist/utils.py. So, a small change was added to the same function which will use np.finfo if the input arr is float and use np.iinfo otherwise