Tone curves or HDR tone mapping curves handled with simple 1D transfer functions applied separately on the 3 RGB channels are inherentily bad because they can affect the color saturation differently in highlights or shadows, and shift hues. This video shows it all :
The solution pointed in the conference is a ratio-preserving curve, using a well-choosen norm as the ratio denominator. In darktable’s filmic, I use the max(RGB), but they show how it can darken colors. The problem is their yellow-weighted norm is kind of black magic, and different norms yield different results, each of them having different drawbacks and bottlenecks.
My claim is a good norm should be as close as possible as the physical energy of the light spectrum. (To be checked).
The spectrum is a 2D graph linking wavelengths and their intensities, like this daylight spectrum :
(Source)
This energy could be expressed like:
E = \int_\Omega I(\lambda)^2 d\lambda
with I the distribution of intensities among wavelengths, where \lambda is the wavelength and \Omega is the domain of visible wavelengths where \lambda \in [370;750] nm.
Unfortunately, in a photo software editor, what we have is not a spectrum but a tristimulus vector (RGB). So I(\lambda) and E have no trivial expressions.
But a camera sensor is 3 colorimeters stacked at once doing essentially:
\begin{cases}
R &= \int_\Omega I(\lambda) \bar R(\lambda) d\lambda \\
G &= \int_\Omega I(\lambda) \bar G(\lambda) d\lambda \\
B &= \int_\Omega I(\lambda) \bar B(\lambda) d\lambda \\
\end{cases}
where R, G, B are 3 scalars representing the pixel components, and \bar R, \bar G, \bar B are the transmittances of the color filters. These are specific to each sensor and are usually not known, even if they can be profiled. However, RGB spaces have base vectors defined in XYZ space (known), and \bar X, \bar Y, \bar Z transmittance functions are known too :
So we can translate RGB chromaticities to XYZ space using the primaries of the camera sensor RGB space, then :
\begin{cases}
X &= \int_\Omega I(\lambda) \bar X(\lambda) d\lambda \\
Y &= \int_\Omega I(\lambda) \bar Y(\lambda) d\lambda \\
Z &= \int_\Omega I(\lambda) \bar Z(\lambda) d\lambda \\
\end{cases}
Now, we know the transmittances functions, but we still don’t know the I(\lambda) functions, so we have no way to isolate \int_\Omega I(\lambda)^2 d\lambda = E in the above equations (with an integration by parts, for example). The most we can do is using the Cauchy-Schwarz inequality and express an upper bound:
X \leq \left(\int_\Omega I(\lambda)^2 d\lambda \right)^{1/2} \left(\int_\Omega\bar X(\lambda)^2 d\lambda\right)^{1/2} \Rightarrow \dfrac{X^2}{\int_\Omega\bar X(\lambda)^2 d\lambda} \leq E.
Propagating this result, we can use each tristimulus value to get the lower bound of E :
E \geq max\left(\dfrac{X^2}{\int_\Omega\bar X(\lambda)^2 d\lambda} ; \dfrac{Y^2}{\int_\Omega\bar Y(\lambda)^2 d\lambda} ; \dfrac{Z^2}{\int_\Omega\bar Z(\lambda)^2 d\lambda}\right) .
Now, we need to find an upper bound of E ? Someone has an idea ? What do you think so far ?
[Full disclosure : I wrote that trapped at home with a bad cold and fever]
All the best. Sometimes rest is all you need to make new discoveries. 
