2023-03-31
We consider the 1-dimensional case only. Let x be the location coordinate (where right is the positive direction) and t be the time variable. Then a wave is a function of the form \begin{aligned} u(x,t) = F(x-ct) + G(x+ct) \end{aligned} where F is a right-travelling function and G is a left-travelling function.
Next we consider reflection and transmission at the boundary of two media. Consider a right-travelling wave. Assume that u = A e^{i(kx - \omega t)} where A\in\mathbb{C}, k, \omega \in \mathbb{R}. | A | is the amplitude, the argument of A is the phase, k is the angular wavenumber, \omega is the angular frequency. We have the relation k=\omega / v where v is the speed of the wave in medium 1. We will generally focus on the imaginary part.
Assume the boundary of the two media is at x_{0}. Let v' be the speed of the wave and k' the angular wavenumber in medium 2. Of course, k' = \omega/ v'. Write u(x,t) = A e^{ikx_{0}} e^{i(k(x-x_{0}) - \omega t)}. Let u^{\mathrm{refl}} be the reflected wave and u^{\mathrm{trans}} the transmitted wave. It is worked out that \begin{aligned} u^{\mathrm{refl}} (x,t) &= \frac{v' - v}{v' + v} A e^{ikx_{0}} e^{i ( - k(x - x_{0}) - \omega t)} \\ u^{\mathrm{trans}} (x,t) &= \frac{2v'}{ v' + v} A e^{ikx_{0}} e^{i (k' (x-x_{0}) -\omega t)}. \end{aligned}
With this setup, we can apply this to the study of interference of membrane (consisting of medium 2). The transmitted wave is reflected by the back surface of the membrane, which is then transmitted back to medium 1. This wave would interfere with the reflected wave by the front surface of the membrane. Of course, the wave keeps reflecting between the front and back surfaces. We will ignore their contribution.
The main reference is from wikipedia 1.
Next we need to compute the irradiance of the total reflected wave. We now specialise to electromagnetic waves which have an electric part E and a magnetic part B. The electric wave and the magnetic wave are orthogonal to each other. We have the relation B=E/c, where c is the speed of light in the medium. Let \varepsilon be the permittivity of the medium. Let \mu be the permeativity of the medium. We have the relation c = \frac{1}{\sqrt{\varepsilon\mu}}.
It is worked out that the energy density is \begin{aligned} \varepsilon E(x,t)^{2}. \end{aligned} The power per unit area is \begin{aligned} c\varepsilon E(x,t)^{2}. \end{aligned} One basically converts the length to time, knowing the speed of light. What is important for perception is the irradiance, also called the intensity. This is the average power per unit area (with respect to t). We just need to average over one period of the wave. Thus the irradiance is \begin{aligned} \frac{c\varepsilon}{T} \int_{t - \frac{T}{2}} ^{t + \frac{T}{2}} E^{2}(x, s) ds. \end{aligned} Here T=2\pi/\omega is the period.
We note E will be the imaginary part of sums of Ae^{i(kx - \omega t)}’s, which is a sinusoidal wave. Let E(x,t) = A e^{i(kx - \omega t + \phi)} with A\in\mathbb{R}_{>0} and \phi the phase. The integral can be made explicit. We get the irradiance: \begin{aligned} \frac {1} {2}c \varepsilon A^{2}. \end{aligned} We can bring in the spectral irradiance distribution with respect to vaccuum wavelengths. Due to orthogonality, contributions to irradiance from different vaccuum wavelengths do not cancel and we find the total irradiance by integrating the contributions.
The main reference is this2.
By Section 1, we can compute the spectral irradiance distribution of the reflected light. For light of different wavelengths (in vaccuum), the IOR of the medium may be different. We will ignore this.
We can convert spectral power distribution (the viewer receives) to a ‘colour’ with respect to the 2^{\circ} CIE standard viewer. To display the colour on a device, we need to further convert it to a representation in P3/Display or sRGB etc.
First we have the CIE XYZ space which we denote by E. It is a 3-dimensional real vector space. Every point corresponds to a sensation of colour perceived by the standard (human) viewer, though due to human limitations, only a subdomain D of the XYZ space is used. The dimension corrresponds to the fact that humans are generally tristimulus. Mixing colours in this space is compatible with the additive structure of the vector space. The mapping of the space of spectral radiance distributions to the XYZ space is given by integration against the CIE’s color matching functions \bar{x}(\lambda), \bar{y}(\lambda), \bar{z}(\lambda). Note that radiance and irradiance are not the same. The emitted light needs to travel to the viewer. We need to compute how much power the viewer receives. We will assume a conversion rate of 1 from irradiance and the power the viewer receives. Let p(\lambda) denote the spectral power distribution. Let V denote the segment of visible vaccuum wavelengths. We have \begin{aligned} X &= \int_{V} p(\lambda) \bar{x}(\lambda) d\lambda, \\ Y &= \int_{V} p(\lambda) \bar{y}(\lambda) d\lambda, \\ Z &= \int_{V} p(\lambda) \bar{z}(\lambda) d\lambda . \end{aligned} See reference3. When we get a batch of XYZ coordinates, we normalise them so that the largest luminance is 1.
The quantity Y corresponds to luminance.
An CIE RGB colour space can be specified by fixing three points R,G,B in the subdomain D of the XYZ space E. We require that OR, OG and OB form a basis for E. Then every point in E can be written as a linear combination of OR, OG and OB. Only the points with in the linear combination with coefficients in [0,1] can be represented. The coefficients are then adjusted using a gamma curve and converted to an integer value between 0 and 255. We get the hex representation of the colour.
Set \begin{aligned} x &= \frac{X}{X+Y+Z}\\ y &= \frac{Y}{X+Y+Z}\\ z &= \frac{Z}{X+Y+Z}. \end{aligned}
This is the projection to the plane X+Y+Z=1 with respect to the origin. We focus on this affine plane and call it U. Every point can be described by using (x,y). Instead of specifying R,G,B, we fix three points r, g, b in U which are non-colinear and a white point w in the interior of the triangle RGB. They are assumed to lie in D. Let W be the point in the XYZ space, such that OW is parallel to Ow and Y_{W}=1. We can determine R,G,B that lifts r,g,b such that OW=OR+OG+OB. \label{eq:WRGB}\tag{eq:WRGB}
Let \begin{aligned} c = \begin{pmatrix} x_{c}\\ y_{c}\\z_{c} \end{pmatrix}, \quad \text{for $c=r,g,b$.} \end{aligned}
Let \begin{aligned} R = \alpha_{r} r, G=\alpha_{g} g, B=\alpha_{b} b \end{aligned} for \alpha_{c}\in\mathbb{R}_{>0}. Let \begin{aligned} M = \begin{pmatrix} r &g&b \end{pmatrix}. \end{aligned} Then [eq:WRGB] becomes \begin{aligned} W = M \begin{pmatrix} \alpha_{r} & & \\ & \alpha_{g} & \\ & & \alpha_{b} \end{pmatrix} \begin{pmatrix} 1\\1\\1 \end{pmatrix}. \end{aligned} Then clearly \begin{pmatrix} \alpha_{r} \\ \alpha_{g} \\\alpha_{b} \end{pmatrix} = M^{-1}W. This is the called wscale in the Python code.
Let \widetilde{M} = \begin{pmatrix} R & G & B \end{pmatrix}. Then the RGB component of a point v\in E is given by \begin{aligned} \widetilde{M}^{-1} v. \end{aligned} Denote the coefficients of v with respect to the basis R,G,B by v^{0}_{R}, v^{0}_{G},v^{0}_{B}. Ref4.
We can look up the r,g,b,w values from here5 and here6 for example.
Next we need to apply the gamma correction. For P3-Display and sRGB this is given by C =\begin{cases} 12.92C_{\text{linear}},&C_{\text{linear}}\leq 0.0031308\\ 1.055C_{\text{linear}}^{1/2.4}-0.055,&C_{\text{linear}}>0.0031308. \end{cases} We write C_{\text{linear}} for v^{0}_{R}, v^{0}_{G},v^{0}_{B} above. Now we have the gamma corrected values and call them v_{R}, v_{G},v_{B}. When we get negative values we may apply desaturation to land in the non-negative range. Then we may get values exceeding 1, I do not know the correct treatment for overexposure. My current treatment is to scale down by dividing out the largest value. Finally we convert them linearly to integer values between 0 and 255. We get the hex representation of the colour.