5.2 The SampledSpectrum Class

SampledSpectrum uses the CoefficientSpectrum infrastructure to represent an SPD with uniformly spaced samples between a starting and an ending wavelength. The wavelength range covers from 400 nm to 700 nm—the range of the visual spectrum where the human visual system is most sensitive. The number of samples, 60, is generally more than enough to accurately represent complex SPDs for rendering. (See the “Further Reading” section for background on sampling rate requirements for SPDs.) Thus, the first sample represents the wavelength range , the second represents , and so forth. These values can easily be changed here as needed.

By inheriting from the CoefficientSpectrum class, SampledSpectrum automatically has all of the basic spectrum arithmetic operations defined earlier. The main methods left to define for it are those that initialize it from spectral data and convert the SPD it represents to other spectral representations (such as RGB). The constructor for initializing it with a constant SPD is straightforward.

We will often be provided spectral data as a set of samples, where the th sample has some value at wavelength . In general, the samples may have an irregular spacing and there may be more or fewer of them than the SampledSpectrum stores. (See the directory scenes/spds in the pbrt distribution for a variety of SPDs for use with pbrt, many of them with irregular sample spacing.)

The FromSampled() method takes arrays of SPD sample values v at given wavelengths lambda and uses them to define a piecewise linear function to represent the SPD. For each SPD sample in the SampledSpectrum, it uses the AverageSpectrumSamples() utility function, defined below, to compute the average of the piecewise linear function over the range of wavelengths that each SPD sample is responsible for.

The AverageSpectrumSamples() function requires that the values be sorted by wavelength. The SpectrumSamplesSorted() function checks that they are; if not, SortSpectrumSamples() sorts them. Note that we allocate new storage for the sorted samples and do not modify the values passed in by the caller; in general, doing so would likely be unexpected behavior for a user of this function (who shouldn’t need to worry about these requirements of its specific implementation). We won’t include the implementations of either of these two functions here, as they are straightforward.

In order to compute the value for the th spectral sample, we compute the range of wavelengths that it’s responsible for—lambda0 to lambda1—and use the AverageSpectrumSamples() function to compute the average value of the given piecewise linear SPD over that range. This is a 1D instance of sampling and reconstruction, a topic that will be discussed in more detail in Chapter 7.

Figure 5.2 shows the basic approach taken by AverageSpectrumSamples(): it iterates over each of the linear segments between samples that are partially or fully within the range of wavelengths, lambdaStart to lambdaEnd. For each such segment, it computes the average value over its range, scales the average by the wavelength range the segment covers, and accumulates a sum of these values. The final average value is this sum divided by the total wavelength range.

The function starts by checking for and handling the edge cases where the range of wavelengths to average over is outside the range of provided wavelengths or the case where there is only a single sample, in which case the average value is trivially computed. We assume that the SPD has a constant value (the values at the two respective endpoints) outside of the provided sample range; if this isn’t a reasonable assumption for a particular set of data, the data provided should have explicit values of (for example) 0 at the endpoints.

Having handled these cases, the next step is to check to see if part of the range to average over goes beyond the first and/or last sample value. If so, we accumulate the contribution of the constant segment(s), scaled by the out-of-bounds wavelength range.

And now we advance to the first index i where the starting wavelength of the interpolation range overlaps the segment from to . A more efficient implementation would use a binary search rather than a linear search here. ^† However, this code is currently only called at scene initialization time, so the lack of these optimizations doesn’t currently impact rendering performance.

The loop below iterates over each of the linear segments that the averaging range overlaps. For each one, it computes the average value over the wavelength range segLambdaStart to segLambdaEnd by averaging the values of the function at those two points. The values in turn are computed by interp(), a lambda function that linearly interpolates between the two endpoints at the given wavelength.

The std::min() and std::max() calls below compute the wavelength range to average over within the segment; note that they naturally handle the cases where lambdaStart, lambdaEnd, or both of them are within the current segment.

5.2.1 XYZ Color

A remarkable property of the human visual system makes it possible to represent colors for human perception with just three floating-point numbers. The tristimulus theory of color perception says that all visible SPDs can be accurately represented for human observers with three values, , , and . Given an emissive SPD , these values are computed by integrating its product with the spectral matching curves , , and :

These curves were determined by the Commission Internationale de l’Éclairage (CIE) standards body after a series of experiments with human test subjects and are graphed in Figure 5.3. It is believed that these matching curves are generally similar to the responses of the three types of color-sensitive cones in the human retina. Remarkably, SPDs with substantially different distributions may have very similar , , and values. To the human observer, such SPDs actually appear the same visually. Pairs of such spectra are called metamers.

This brings us to a subtle point about representations of spectral power distributions. Most color spaces attempt to model colors that are visible to humans and therefore use only three coefficients, exploiting the tristimulus theory of color perception. Although XYZ works well to represent a given SPD to be displayed for a human observer, it is not a particularly good set of basis functions for spectral computation. For example, although XYZ values would work well to describe the perceived color of lemon skin or a fluorescent light individually (recall Figure 5.1), the product of their respective XYZ values is likely to give a noticeably different XYZ color than the XYZ value computed by multiplying more accurate representations of their SPDs and then computing the XYZ value.

pbrt provides the values of the standard , , and response curves sampled at 1-nm increments from 360 nm to 830 nm. The wavelengths of the th sample in the arrays below are given by the th element of CIE_lambda; having the wavelengths of the samples explicitly represented in this way makes it easy to pass the XYZ samples into functions like AverageSpectrumSamples() that take such an array of wavelengths as a parameter.

SampledSpectrum uses these samples to compute the XYZ matching curves in its spectral representation (i.e., themselves as SampledSpectrums).

The SampledSpectrum XYZ matching curves are computed in the SampledSpectrum::Init() method, which is called at system startup time by the pbrtInit() function defined in Section A.2.

Given the wavelength range and number of samples for SampledSpectrum, computing the values of the matching functions for each sample is just a matter of computing the sample’s wavelength range and using the AverageSpectrumSamples() routine.

All Spectrum implementations in pbrt must provide a method that converts their SPD to coefficients. This method is called, for example, in the process of updating pixels in the image. When a Spectrum representing the light carried along a ray from the camera is provided to the Film, the Film converts the SPD into XYZ coefficients as a first step in the process of finally turning them into RGB values used for storage and/or display.

To compute XYZ coefficients, SampledSpectrum computes the integrals from Equation (5.1) with a Riemann sum:

and so forth.

The coefficient of XYZ color is closely related to luminance, which measures the perceived brightness of a color. Luminance is described in more detail in Section 5.4.3. We provide a method to compute alone in a separate method as often only the luminance of a spectrum is desired. (For example, some of the light transport algorithms in Chapters 14, 15, and 16 use luminance as a measure of the relative importance of light-carrying paths through the scene.)

5.2.2 RGB Color

When we display an RGB color on a display, the spectrum that is actually displayed is basically determined by the weighted sum of three spectral response curves, one for each of red, green, and blue, as emitted by the display’s phosphors, LED or LCD elements, or plasma cells. ^† Figure 5.4 plots the red, green, and blue distributions emitted by a LED display and a LCD display; note that they are remarkably different.

Figure 5.5 in turn shows the SPDs that result from displaying the RGB color on those displays. Not surprisingly, the resulting SPDs are quite different as well. This example illustrates that using RGB values provided by the user to describe a particular color is actually only meaningful given knowledge of the characteristics of the display they were using when they selected the RGB values.

Given an representation of an SPD, we can convert it to corresponding RGB coefficients, given the choice of a particular set of SPDs that define red, green, and blue for a display of interest. Given the spectral response curves , , and , for a particular display, RGB coefficients can be computed by integrating the response curves with the SPD and using the tristimulus theory of color perception:

The integrals of the products of and so forth can be precomputed for given response curves, making it possible to express the full conversion as a matrix:

The conversion routines implemented in pbrt are based on a standard set of these RGB spectra that has been defined for high-definition television.

The inverse of this matrix gives coefficients to convert given RGB values expressed with respect to a particular set of RGB response curves to coefficients.

Given these functions, a SampledSpectrum can convert to RGB coefficients by first converting to XYZ and then using the XYZToRGB() utility function.

An RGBSpectrum can also be created easily, using the ToRGBSpectrum() method.

Going the other way and converting from RGB or XYZ values to a SPD isn’t as easy: the problem is highly under-constrained. Recall that an infinite number of different SPDs have the same (and thus RGB) coefficients. Thus, given an RGB or value, there are an infinite number of possible SPDs that could be chosen for it. There are a number of desirable criteria that we’d like a conversion function to have:

• If all of the RGB coefficients have the same value, the resulting SPD should be constant.
• In general, it’s desirable that the computed SPD be smooth. Most real-world objects have relatively smooth spectra. (The main source of spiky spectra is light sources, especially fluorescents. Fortunately, actual spectral data are more commonly available for illuminants than they are for reflectances.)

The smoothness goal is one of the reasons why constructing an SPD as a weighted sum of a display’s , , and SPDs is not a good solution: as shown in Figure 5.4, those functions are generally irregular and spiky, and a weighted sum of them will thus not be a very smooth SPD. Although the result will be a metamer of the given RGB values, it’s likely not an accurate representation of the SPD of the actual object.

Here we implement a method for converting RGBs to SPDs suggested by Smits (1999) that tries to achieve the goals above. This approach is based on the observation that a good start is to compute individual SPDs for red, green, and blue that are smooth and such that computing the weighted sum of them with the given RGB coefficients and then converting back to RGB give a result that is close to the original RGB coefficients. He found such spectra through a numerical optimization procedure.

Smits observed that two additional improvements could be made to this basic approach. First, rather than representing constant spectra by the sums of the computed red, green, and blue SPDs, the sum of which isn’t exactly constant, it’s better to represent constant spectra with constant SPDs. Second, mixtures of colors like yellow (a mixture of red and green) that are a mixture of two of the primaries are better represented by their own precomputed smooth SPDs rather than the sum of SPDs for the two corresponding primaries.

The following arrays store SPDs that meet these criteria, with their samples’ wavelengths in RGB2SpectLambda[] (these data were generated by Karl vom Berge).

If a given RGB color describes illumination from a light source, better results are achieved if the conversion tables are computed using the spectral power distribution of a representative illumination source to define “white” rather than using a constant spectrum as they are for the tables above that are used for reflectances. The RGBIllum2Spect arrays use the D65 spectral power distribution, which has been standardized by the CIE to represent midday sunlight. (The D65 illuminant will be discussed more in Section 12.1.2.)

The fragment <<Compute RGB to spectrum functions for SampledSpectrum>>, which is called from SampledSpectrum::Init(), isn’t included here; it initializes the following SampledSpectrum values by resampling the RGBRefl2Spect and RGBIllum2Spect distributions using the AverageSpectrumSamples() function.

The SampledSpectrum::FromRGB() method converts from the given RGB values to a full SPD. In addition to the RGB values, it takes an enumeration value that denotes whether the RGB value represents surface reflectance or an illuminant; the corresponding rgbIllum2Spect or rgbRefl2Spect values are used for the conversion.

Here we’ll show the conversion process for reflectances. The computation for illuminants is the same, just using the different conversion values. First, the implementation determines whether the red, green, or blue channel is the smallest.

Here is the code for the case of a red component being the smallest. (The cases for green and blue are analogous and not included in the book here.) If red is the smallest, we know that green and blue have greater values than red. As such, we can start to convert the final SPD to return by assigning to it the value of the red component times the white spectrum in rgbRefl2SpectWhite. Having done this, the remaining RGB value left to process is . The code in turn determines which of the remaining two components is the smallest. This value, times the cyan (green and blue) spectrum, is added to the result and we’re left with either or . Based on whether the green or blue channel is non-zero, the green or blue SPD is scaled by the remainder and the conversion is complete.

Given the method to convert from RGB, converting from XYZ color is easy. We first convert from XYZ to RGB and then use the FromRGB() method.

Finally, we provide a constructor that converts from an instance of the RGBSpectrum class, again using the infrastructure above.