colour.recovery.spectral_primary_decomposition_Mallett2019#

colour.recovery.spectral_primary_decomposition_Mallett2019(colourspace: RGB_Colourspace, cmfs: MultiSpectralDistributions | None = None, illuminant: SpectralDistribution | None = None, metric: Callable = np.linalg.norm, metric_args: tuple = (), optimisation_kwargs: dict | None = None) MultiSpectralDistributions[source]#

Perform the spectral primary decomposition as described in Mallett and Yuksel (2019) for given RGB colourspace.

Parameters:
  • colourspace (RGB_Colourspace) – RGB colourspace.

  • cmfs (MultiSpectralDistributions | None) – Standard observer colour matching functions, default to the CIE 1931 2 Degree Standard Observer.

  • illuminant (SpectralDistribution | None) – Illuminant spectral distribution, default to CIE Standard Illuminant D65.

  • metric (Callable) –

    Function to be minimised, i.e., the objective function.

    metric(basis, *metric_args) -> float

    where basis is three reflectances concatenated together, each with a shape matching shape.

  • metric_args (tuple) – Additional arguments passed to metric.

  • optimisation_kwargs (dict | None) – Parameters for scipy.optimize.minimize() definition.

Returns:

Basis functions for given RGB colourspace.

Return type:

colour.MultiSpectralDistributions

References

[MY19]

Notes

  • In-addition to the BT.709 primaries used by the sRGB colourspace, [MY19] tried BT.2020, P3 D65, Adobe RGB 1998, NTSC (1987), Pal/Secam, ProPhoto RGB, and Adobe Wide Gamut RGB primaries, every one of which encompasses a larger (albeit not-always-enveloping) set of CIE L*a*b* colours than BT.709. Of these, only Pal/Secam produces a feasible basis, which is relatively unsurprising since it is very similar to BT.709, whereas the others are significantly larger.

Examples

>>> from colour import MSDS_CMFS, SDS_ILLUMINANTS, SpectralShape
>>> from colour.models import RGB_COLOURSPACE_PAL_SECAM
>>> from colour.utilities import numpy_print_options
>>> cmfs = (
...     MSDS_CMFS["CIE 1931 2 Degree Standard Observer"]
...     .copy()
...     .align(SpectralShape(360, 780, 10))
... )
>>> illuminant = SDS_ILLUMINANTS["D65"].copy().align(cmfs.shape)
>>> msds = spectral_primary_decomposition_Mallett2019(
...     RGB_COLOURSPACE_PAL_SECAM,
...     cmfs,
...     illuminant,
...     optimisation_kwargs={"options": {"ftol": 1e-5}},
... )
>>> with numpy_print_options(suppress=True):
...     print(msds)  
[[ 360.            0.3395134...    0.3400214...    0.3204650...]
 [ 370.            0.3355246...    0.3338028...    0.3306724...]
 [ 380.            0.3376707...    0.3185578...    0.3437715...]
 [ 390.            0.3178866...    0.3351754...    0.3469378...]
 [ 400.            0.3045154...    0.3248376...    0.3706469...]
 [ 410.            0.2935652...    0.2919463...    0.4144884...]
 [ 420.            0.1875740...    0.1853729...    0.6270530...]
 [ 430.            0.0167983...    0.054483 ...    0.9287186...]
 [ 440.            0.       ...    0.       ...    1.       ...]
 [ 450.            0.       ...    0.       ...    1.       ...]
 [ 460.            0.       ...    0.       ...    1.       ...]
 [ 470.            0.       ...    0.0458044...    0.9541955...]
 [ 480.            0.       ...    0.2960917...    0.7039082...]
 [ 490.            0.       ...    0.5042592...    0.4957407...]
 [ 500.            0.       ...    0.6655795...    0.3344204...]
 [ 510.            0.       ...    0.8607541...    0.1392458...]
 [ 520.            0.       ...    0.9999998...    0.0000001...]
 [ 530.            0.       ...    1.       ...    0.       ...]
 [ 540.            0.       ...    1.       ...    0.       ...]
 [ 550.            0.       ...    1.       ...    0.       ...]
 [ 560.            0.       ...    0.9924229...    0.       ...]
 [ 570.            0.       ...    0.9970703...    0.0025673...]
 [ 580.            0.0396002...    0.9028231...    0.0575766...]
 [ 590.            0.7058973...    0.2941026...    0.       ...]
 [ 600.            1.       ...    0.       ...    0.       ...]
 [ 610.            1.       ...    0.       ...    0.       ...]
 [ 620.            1.       ...    0.       ...    0.       ...]
 [ 630.            1.       ...    0.       ...    0.       ...]
 [ 640.            0.9835925...    0.0100166...    0.0063908...]
 [ 650.            0.7878949...    0.1265097...    0.0855953...]
 [ 660.            0.5987994...    0.2051062...    0.1960942...]
 [ 670.            0.4724493...    0.2649623...    0.2625883...]
 [ 680.            0.3989806...    0.3007488...    0.3002704...]
 [ 690.            0.3666586...    0.3164003...    0.3169410...]
 [ 700.            0.3497806...    0.3242863...    0.3259329...]
 [ 710.            0.3563736...    0.3232441...    0.3203822...]
 [ 720.            0.3362624...    0.3326209...    0.3311165...]
 [ 730.            0.3245015...    0.3365982...    0.3389002...]
 [ 740.            0.3335520...    0.3320670...    0.3343808...]
 [ 750.            0.3441287...    0.3291168...    0.3267544...]
 [ 760.            0.3343705...    0.3330132...    0.3326162...]
 [ 770.            0.3274633...    0.3305704...    0.3419662...]
 [ 780.            0.3475263...    0.3262331...    0.3262404...]]