Skip to content

Commit

Permalink
Merge pull request #777 from colour-science/feature/rosch-macadam
Browse files Browse the repository at this point in the history
PR: Implement support to generate the "Rösch-MacAdam" colour solid hue lines.
  • Loading branch information
KelSolaar committed Feb 7, 2021
2 parents 362ccee + d1cca65 commit eb837b4
Show file tree
Hide file tree
Showing 4 changed files with 206 additions and 47 deletions.
5 changes: 3 additions & 2 deletions colour/volume/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@
from .mesh import is_within_mesh_volume
from .pointer_gamut import is_within_pointer_gamut
from .spectrum import (generate_pulse_waves, XYZ_outer_surface,
is_within_visible_spectrum)
solid_RoschMacAdam, is_within_visible_spectrum)
from .rgb import (RGB_colourspace_limits, RGB_colourspace_volume_MonteCarlo,
RGB_colourspace_volume_coverage_MonteCarlo,
RGB_colourspace_pointer_gamut_coverage_MonteCarlo,
Expand All @@ -18,7 +18,8 @@
__all__ += ['is_within_mesh_volume']
__all__ += ['is_within_pointer_gamut']
__all__ += [
'generate_pulse_waves', 'XYZ_outer_surface', 'is_within_visible_spectrum'
'generate_pulse_waves', 'XYZ_outer_surface', 'solid_RoschMacAdam',
'is_within_visible_spectrum'
]
__all__ += [
'RGB_colourspace_limits', 'RGB_colourspace_volume_MonteCarlo',
Expand Down
151 changes: 130 additions & 21 deletions colour/volume/spectrum.py
Original file line number Diff line number Diff line change
@@ -1,9 +1,10 @@
# -*- coding: utf-8 -*-
"""
Visible Spectrum Volume Computations
====================================
Rösch-MacAdam colour solid - Visible Spectrum Volume Computations
=================================================================
Defines objects related to visible spectrum volume computations.
Defines objects related to *Rösch-MacAdam* colour solid, visible spectrum
volume computations.
References
----------
Expand All @@ -13,11 +14,16 @@
- :cite:`Mansencal2018` : Mansencal, T. (2018). How is the visible gamut
bounded? Retrieved August 19, 2018, from
https://stackoverflow.com/a/48396021/931625
- :cite:`Martinez-Verdu2007` : Martínez-Verdú, F., Perales, E., Chorro, E.,
de Fez, D., Viqueira, V., & Gilabert, E. (2007). Computation and
visualization of the MacAdam limits for any lightness, hue angle, and light
source. Journal of the Optical Society of America A, 24(6), 1501.
doi:10.1364/JOSAA.24.001501
"""

import numpy as np

from colour.colorimetry import (MSDS_CMFS, msds_to_XYZ, SpectralShape, sd_ones)
from colour.colorimetry import MSDS_CMFS, msds_to_XYZ, SpectralShape, sd_ones
from colour.constants import DEFAULT_FLOAT_DTYPE
from colour.volume import is_within_mesh_volume
from colour.utilities import zeros
Expand All @@ -31,7 +37,7 @@

__all__ = [
'SPECTRAL_SHAPE_OUTER_SURFACE_XYZ', 'generate_pulse_waves',
'XYZ_outer_surface', 'is_within_visible_spectrum'
'XYZ_outer_surface', 'solid_RoschMacAdam', 'is_within_visible_spectrum'
]

SPECTRAL_SHAPE_OUTER_SURFACE_XYZ = SpectralShape(360, 780, 5)
Expand All @@ -46,10 +52,11 @@
_CACHE_OUTER_SURFACE_XYZ_POINTS = {}


def generate_pulse_waves(bins):
def generate_pulse_waves(bins, pulse_order='Bins', filter_jagged_pulses=False):
"""
Generates the pulse waves of given number of bins necessary to totally
stimulate the colour matching functions.
stimulate the colour matching functions and produce the *Rösch-MacAdam*
colour solid.
Assuming 5 bins, a first set of SPDs would be as follows::
Expand Down Expand Up @@ -81,6 +88,29 @@ def generate_pulse_waves(bins):
----------
bins : int
Number of bins of the pulse waves.
pulse_order : unicode, optional
**{'Bins', 'Pulse Wave Width'}**,
Method for ordering the pulse waves. *Bins* is the default order, with
*Pulse Wave Width* ordering, instead of iterating over the pulse wave
widths first, iteration occurs over the bins, producing blocks of pulse
waves with increasing width.
filter_jagged_pulses : bool, optional
Whether to filter jagged pulses. When ``pulse_order`` is set to
*Pulse Wave Width*, the pulses are ordered by increasing width. Because
of the discrete nature of the underlying signal, the resulting pulses
will be jagged. For example assuming 5 bins, the center block with
the two extreme values added would be as follows::
0 0 0 0 0
0 0 1 0 0
0 0 1 1 0 <--
0 1 1 1 0
0 1 1 1 1 <--
1 1 1 1 1
Setting the ``filter_jagged_pulses`` parameter to `True` will result
in the removal of the two marked pulses above which avoid jagged lines
when plotting and having to resort to excessive ``bins`` values.
Returns
-------
Expand All @@ -89,7 +119,7 @@ def generate_pulse_waves(bins):
References
----------
:cite:`Lindbloom2015`, :cite:`Mansencal2018`
:cite:`Lindbloom2015`, :cite:`Mansencal2018`, :cite:`Martinez-Verdu2007`
Examples
--------
Expand All @@ -116,14 +146,59 @@ def generate_pulse_waves(bins):
[ 1., 1., 0., 1., 1.],
[ 1., 1., 1., 0., 1.],
[ 1., 1., 1., 1., 1.]])
>>> generate_pulse_waves(5, 'Pulse Wave Width')
array([[ 0., 0., 0., 0., 0.],
[ 1., 0., 0., 0., 0.],
[ 1., 1., 0., 0., 0.],
[ 1., 1., 0., 0., 1.],
[ 1., 1., 1., 0., 1.],
[ 0., 1., 0., 0., 0.],
[ 0., 1., 1., 0., 0.],
[ 1., 1., 1., 0., 0.],
[ 1., 1., 1., 1., 0.],
[ 0., 0., 1., 0., 0.],
[ 0., 0., 1., 1., 0.],
[ 0., 1., 1., 1., 0.],
[ 0., 1., 1., 1., 1.],
[ 0., 0., 0., 1., 0.],
[ 0., 0., 0., 1., 1.],
[ 0., 0., 1., 1., 1.],
[ 1., 0., 1., 1., 1.],
[ 0., 0., 0., 0., 1.],
[ 1., 0., 0., 0., 1.],
[ 1., 0., 0., 1., 1.],
[ 1., 1., 0., 1., 1.],
[ 1., 1., 1., 1., 1.]])
>>> generate_pulse_waves(5, 'Pulse Wave Width', True)
array([[ 0., 0., 0., 0., 0.],
[ 1., 0., 0., 0., 0.],
[ 1., 1., 0., 0., 1.],
[ 0., 1., 0., 0., 0.],
[ 1., 1., 1., 0., 0.],
[ 0., 0., 1., 0., 0.],
[ 0., 1., 1., 1., 0.],
[ 0., 0., 0., 1., 0.],
[ 0., 0., 1., 1., 1.],
[ 0., 0., 0., 0., 1.],
[ 1., 0., 0., 1., 1.],
[ 1., 1., 1., 1., 1.]])
"""

square_waves = []
square_waves_basis = np.tril(
np.ones((bins, bins), dtype=DEFAULT_FLOAT_DTYPE))[0:-1, :]
for square_wave_basis in square_waves_basis:

if pulse_order.lower() == 'bins':
for square_wave_basis in square_waves_basis:
for i in range(bins):
square_waves.append(np.roll(square_wave_basis, i))
else:
for i in range(bins):
square_waves.append(np.roll(square_wave_basis, i))
for j, square_wave_basis in enumerate(square_waves_basis):
square_waves.append(np.roll(square_wave_basis, i - j // 2))

if filter_jagged_pulses:
square_waves = square_waves[::2]

return np.vstack([
zeros(bins),
Expand All @@ -135,18 +210,44 @@ def generate_pulse_waves(bins):
def XYZ_outer_surface(cmfs=MSDS_CMFS['CIE 1931 2 Degree Standard Observer']
.copy().align(SPECTRAL_SHAPE_OUTER_SURFACE_XYZ),
illuminant=sd_ones(SPECTRAL_SHAPE_OUTER_SURFACE_XYZ),
point_order='Bins',
filter_jagged_points=False,
**kwargs):
"""
Generates the *CIE XYZ* colourspace outer surface for given colour matching
functions using multi-spectral conversion of pulse waves to *CIE XYZ*
tristimulus values.
Generates the *Rösch-MacAdam* colour solid, i.e. *CIE XYZ* colourspace
outer surface, for given colour matching functions using multi-spectral
conversion of pulse waves to *CIE XYZ* tristimulus values.
Parameters
----------
cmfs : XYZ_ColourMatchingFunctions, optional
Standard observer colour matching functions.
illuminant : SpectralDistribution, optional
Illuminant spectral distribution.
point_order : unicode, optional
**{'Bins', 'Pulse Wave Width'}**,
Method for ordering the underlying pulse waves used to generate the
*Rösch-MacAdam* colour solid. *Bins* is the default order, with
*Pulse Wave Width* ordering, instead of iterating over the pulse wave
widths first, iteration occurs over the bins, producing blocks of pulse
waves with increasing width.
filter_jagged_points : bool, optional
Whether to filter the underlying jagged pulses. When ``point_order`` is
set to *Pulse Wave Width*, the pulses are ordered by increasing width.
Because of the discrete nature of the underlying signal, the resulting
pulses will be jagged. For example assuming 5 bins, the center block
with the two extreme values added would be as follows::
0 0 0 0 0
0 0 1 0 0
0 0 1 1 0 <--
0 1 1 1 0
0 1 1 1 1 <--
1 1 1 1 1
Setting the ``filter_jagged_points`` parameter to `True` will result
in the removal of the two marked pulses above which avoid jagged lines
when plotting and having to resort to excessive ``bins`` values.
Other Parameters
----------------
Expand All @@ -157,11 +258,12 @@ def XYZ_outer_surface(cmfs=MSDS_CMFS['CIE 1931 2 Degree Standard Observer']
Returns
-------
ndarray
Outer surface *CIE XYZ* tristimulus values.
*Rösch-MacAdam* colour solid, *CIE XYZ* outer surface tristimulus
values.
References
----------
:cite:`Lindbloom2015`, :cite:`Mansencal2018`
:cite:`Lindbloom2015`, :cite:`Mansencal2018`, :cite:`Martinez-Verdu2007`
Examples
--------
Expand Down Expand Up @@ -207,18 +309,23 @@ def XYZ_outer_surface(cmfs=MSDS_CMFS['CIE 1931 2 Degree Standard Observer']
settings = {'method': 'Integration', 'shape': cmfs.shape}
settings.update(kwargs)

key = (hash(cmfs), hash(illuminant), str(settings))
key = (hash(cmfs), hash(illuminant), point_order, filter_jagged_points,
str(settings))
XYZ = _CACHE_OUTER_SURFACE_XYZ.get(key)

if XYZ is None:
pulse_waves = generate_pulse_waves(len(cmfs.wavelengths))
pulse_waves = generate_pulse_waves(
len(cmfs.wavelengths), point_order, filter_jagged_points)
XYZ = msds_to_XYZ(pulse_waves, cmfs, illuminant, **settings) / 100

_CACHE_OUTER_SURFACE_XYZ[key] = XYZ

return XYZ


solid_RoschMacAdam = XYZ_outer_surface


def is_within_visible_spectrum(
XYZ,
cmfs=MSDS_CMFS['CIE 1931 2 Degree Standard Observer']
Expand All @@ -227,8 +334,9 @@ def is_within_visible_spectrum(
tolerance=None,
**kwargs):
"""
Returns if given *CIE XYZ* tristimulus values are within visible spectrum
volume / given colour matching functions volume.
Returns if given *CIE XYZ* tristimulus values are within the visible
spectrum volume, i.e. *Rösch-MacAdam* colour solid, for given colour
matching functions and illuminant.
Parameters
----------
Expand All @@ -250,7 +358,8 @@ def is_within_visible_spectrum(
Returns
-------
bool
Is within visible spectrum.
Are *CIE XYZ* tristimulus values within the visible spectrum volume,
i.e. *Rösch-MacAdam* colour solid.
Notes
-----
Expand All @@ -276,7 +385,7 @@ def is_within_visible_spectrum(
vertices = _CACHE_OUTER_SURFACE_XYZ_POINTS.get(key)

if vertices is None:
_CACHE_OUTER_SURFACE_XYZ_POINTS[key] = vertices = (XYZ_outer_surface(
_CACHE_OUTER_SURFACE_XYZ_POINTS[key] = vertices = (solid_RoschMacAdam(
cmfs, illuminant, **kwargs))

return is_within_mesh_volume(XYZ, vertices, tolerance)
92 changes: 70 additions & 22 deletions colour/volume/tests/test_spectrum.py
Original file line number Diff line number Diff line change
Expand Up @@ -41,30 +41,78 @@ def test_generate_pulse_waves(self):
np.testing.assert_array_equal(
generate_pulse_waves(5),
np.array([
[0.00000000, 0.00000000, 0.00000000, 0.00000000, 0.00000000],
[1.00000000, 0.00000000, 0.00000000, 0.00000000, 0.00000000],
[0.00000000, 1.00000000, 0.00000000, 0.00000000, 0.00000000],
[0.00000000, 0.00000000, 1.00000000, 0.00000000, 0.00000000],
[0.00000000, 0.00000000, 0.00000000, 1.00000000, 0.00000000],
[0.00000000, 0.00000000, 0.00000000, 0.00000000, 1.00000000],
[1.00000000, 1.00000000, 0.00000000, 0.00000000, 0.00000000],
[0.00000000, 1.00000000, 1.00000000, 0.00000000, 0.00000000],
[0.00000000, 0.00000000, 1.00000000, 1.00000000, 0.00000000],
[0.00000000, 0.00000000, 0.00000000, 1.00000000, 1.00000000],
[1.00000000, 0.00000000, 0.00000000, 0.00000000, 1.00000000],
[1.00000000, 1.00000000, 1.00000000, 0.00000000, 0.00000000],
[0.00000000, 1.00000000, 1.00000000, 1.00000000, 0.00000000],
[0.00000000, 0.00000000, 1.00000000, 1.00000000, 1.00000000],
[1.00000000, 0.00000000, 0.00000000, 1.00000000, 1.00000000],
[1.00000000, 1.00000000, 0.00000000, 0.00000000, 1.00000000],
[1.00000000, 1.00000000, 1.00000000, 1.00000000, 0.00000000],
[0.00000000, 1.00000000, 1.00000000, 1.00000000, 1.00000000],
[1.00000000, 0.00000000, 1.00000000, 1.00000000, 1.00000000],
[1.00000000, 1.00000000, 0.00000000, 1.00000000, 1.00000000],
[1.00000000, 1.00000000, 1.00000000, 0.00000000, 1.00000000],
[1.00000000, 1.00000000, 1.00000000, 1.00000000, 1.00000000],
[0.0, 0.0, 0.0, 0.0, 0.0],
[1.0, 0.0, 0.0, 0.0, 0.0],
[0.0, 1.0, 0.0, 0.0, 0.0],
[0.0, 0.0, 1.0, 0.0, 0.0],
[0.0, 0.0, 0.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 0.0, 1.0],
[1.0, 1.0, 0.0, 0.0, 0.0],
[0.0, 1.0, 1.0, 0.0, 0.0],
[0.0, 0.0, 1.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 1.0, 1.0],
[1.0, 0.0, 0.0, 0.0, 1.0],
[1.0, 1.0, 1.0, 0.0, 0.0],
[0.0, 1.0, 1.0, 1.0, 0.0],
[0.0, 0.0, 1.0, 1.0, 1.0],
[1.0, 0.0, 0.0, 1.0, 1.0],
[1.0, 1.0, 0.0, 0.0, 1.0],
[1.0, 1.0, 1.0, 1.0, 0.0],
[0.0, 1.0, 1.0, 1.0, 1.0],
[1.0, 0.0, 1.0, 1.0, 1.0],
[1.0, 1.0, 0.0, 1.0, 1.0],
[1.0, 1.0, 1.0, 0.0, 1.0],
[1.0, 1.0, 1.0, 1.0, 1.0],
]))

np.testing.assert_array_equal(
generate_pulse_waves(5, 'Pulse Wave Width'),
np.array([
[0.0, 0.0, 0.0, 0.0, 0.0],
[1.0, 0.0, 0.0, 0.0, 0.0],
[1.0, 1.0, 0.0, 0.0, 0.0],
[1.0, 1.0, 0.0, 0.0, 1.0],
[1.0, 1.0, 1.0, 0.0, 1.0],
[0.0, 1.0, 0.0, 0.0, 0.0],
[0.0, 1.0, 1.0, 0.0, 0.0],
[1.0, 1.0, 1.0, 0.0, 0.0],
[1.0, 1.0, 1.0, 1.0, 0.0],
[0.0, 0.0, 1.0, 0.0, 0.0],
[0.0, 0.0, 1.0, 1.0, 0.0],
[0.0, 1.0, 1.0, 1.0, 0.0],
[0.0, 1.0, 1.0, 1.0, 1.0],
[0.0, 0.0, 0.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 1.0, 1.0],
[0.0, 0.0, 1.0, 1.0, 1.0],
[1.0, 0.0, 1.0, 1.0, 1.0],
[0.0, 0.0, 0.0, 0.0, 1.0],
[1.0, 0.0, 0.0, 0.0, 1.0],
[1.0, 0.0, 0.0, 1.0, 1.0],
[1.0, 1.0, 0.0, 1.0, 1.0],
[1.0, 1.0, 1.0, 1.0, 1.0],
]))

np.testing.assert_equal(
np.sort(generate_pulse_waves(5), axis=0),
np.sort(generate_pulse_waves(5, 'Pulse Wave Width'), axis=0))

np.testing.assert_array_equal(
generate_pulse_waves(5, 'Pulse Wave Width', True),
np.array([
[0.0, 0.0, 0.0, 0.0, 0.0],
[1.0, 0.0, 0.0, 0.0, 0.0],
[1.0, 1.0, 0.0, 0.0, 1.0],
[0.0, 1.0, 0.0, 0.0, 0.0],
[1.0, 1.0, 1.0, 0.0, 0.0],
[0.0, 0.0, 1.0, 0.0, 0.0],
[0.0, 1.0, 1.0, 1.0, 0.0],
[0.0, 0.0, 0.0, 1.0, 0.0],
[0.0, 0.0, 1.0, 1.0, 1.0],
[0.0, 0.0, 0.0, 0.0, 1.0],
[1.0, 0.0, 0.0, 1.0, 1.0],
[1.0, 1.0, 1.0, 1.0, 1.0],
]))


class TestXYZOuterSurface(unittest.TestCase):
"""
Expand Down
Loading

0 comments on commit eb837b4

Please sign in to comment.