How are the Mallett2019 basis functions calculated?

[NogginBops]

In the paper the authors suggest using the CIE 2006 LMS 10° Cone fundamentals (but I guess their XYZ color matching function equivalents).

The documentation however seems to imply that CIE 1931 2 Degree Standard Observer is used?
https://colour.readthedocs.io/en/develop/generated/colour.recovery.spectral_primary_decomposition_Mallett2019.html

I'm trying to verify that "lossless" RGB -> sd -> RGB is possible with the method, but I'm struggling to get those results

illuminant = 'D65'
observer = 'CIE 1931 2 Degree Standard Observer'
colorspace = 'ITU-R BT.709'
shape = color.SpectralShape(380, 780, 5)

illuminant_ccs = color.CCS_ILLUMINANTS[observer][illuminant]

# / 100 is to normalize D65, there might be a better way to do this?
illuminant = color.SDS_ILLUMINANTS[illuminant].copy().align(shape) / 100
observer = color.MSDS_CMFS[observer].copy().align(shape)
colorspace = color.RGB_COLOURSPACES[colorspace]

red   = [1, 0, 0]

red_sd = color.recovery.RGB_to_sd_Mallett2019(red).align(shape)
red_xyz = color.colorimetry.sd_to_XYZ_integration(red_sd, cmfs=observer, illuminant=illuminant, k=1)
red_rgb = color.XYZ_to_RGB(red_xyz, colorspace, illuminant=illuminant_ccs, apply_cctf_encoding=False)
red_srgb = color.RGB_to_RGB(red_rgb, colorspace, 'sRGB', chromatic_adaptation_transform=None, apply_cctf_decoding=False, apply_cctf_encoding=True)

np.set_printoptions(suppress=True)
print(red, red_rgb, red_srgb, red_rgb[0], red_rgb[1], red_rgb[2])

[1, 0, 0] [ 1.00645053 0. -0. ] [ 1.00278025 0.00136576 0.00031839] 1.00645053412 5.18912857128e-10 -3.36285367608e-10

This result seems pretty close as both green and blue are on the order of e-10 but the red value has an error on the order of e-4 which is much higher than what is expected.

I'm not sure what is going wrong, but my guess is that it either has something to do with my normalization of D65 or my observer doesn't match the observer used to generate the Mallett2019 basis functions.

Maybe it would be a good idea to document more clearly how the values where generated if they are not generated as recommended in the paper?

Any help is appreciated! 🙂

[KelSolaar]

Hello @NogginBops,

This is a really good question and we should have written the details. I'm fairly confident that those were computed using the Matlab code in the supplemental.zip file associated with the paper: https://diglib.eg.org/handle/10.2312/sr20191216

@enneract : Would it be please possible to confirm?

Cheers,

Thomas

[NogginBops]

Does my code for checking this look reasonable?
I'm not super that with color science as a field so I'm not sure if some of the operations I'm doing are dumb or not.

Like, can I just do D65 / 100? It seems the spectrum has values >100 so maybe I should divide by D65.values.max()? (maybe the important part is matching what the matlab code is doing)
Is it ok to call .align() on "standard" spectra like D65 or should calculations be done in their defined shape?

[KelSolaar]

Given that we are recovering reflectance (as an illuminant is passed to the integration part and we are assuming the values are under D65), you should let k=None and divide the integration result by 100 to convert it from percentage to float scalar:

import colour as color
import numpy as np

illuminant = "D65"
observer = "CIE 1931 2 Degree Standard Observer"
colorspace = "ITU-R BT.709"
shape = color.SpectralShape(380, 780, 5)

illuminant_ccs = color.CCS_ILLUMINANTS[observer][illuminant]

illuminant = color.SDS_ILLUMINANTS[illuminant].copy().align(shape)
observer = color.MSDS_CMFS[observer].copy().align(shape)
colorspace = color.RGB_COLOURSPACES[colorspace]

sample = [1, 1, 1]

sample_sd = color.recovery.RGB_to_sd_Mallett2019(sample).align(shape)
sample_xyz = color.colorimetry.sd_to_XYZ_integration(
    sample_sd, cmfs=observer, illuminant=illuminant
) / 100
sample_rgb = color.XYZ_to_RGB(sample_xyz, colorspace)

np.set_printoptions(suppress=True)
print(sample, sample_rgb)

At this point, I get this:

[1, 1, 1] [ 1.00004314  1.00001476  0.99972668]

Not perfect but relatively good!