Recovering High Dynamic Range Radiance Maps from Photographs

“Recovering High Dynamic Range Radiance Maps from Photographs” by Paul E. Debevec and Jitendra Malik of the University of California at Berkeley appeared at SIGGRAPH 1997, in the proceedings of the 24th Annual Conference on Computer Graphics and Interactive Techniques (pages 369–378). From a set of photographs of a static scene taken at known exposure times, it recovers the response curve of the whole imaging process, from film or sensor to the stored pixel values, up to a scale factor, and then merges the photographs into a radiance map whose values are proportional to scene radiance. No calibration chart is needed, and the curve is not given a parametric form. The method, the way it is used today, and a worked NumPy example of response recovery and merging are in high dynamic range imaging. This article covers the publication: the problem as it stood in 1997, the method and the demonstrations as printed, the places where the paper’s text, equations and code disagree, and what came after. Section, figure and footnote numbers below are from the authors’ copy, which matches a camera-ready SIGGRAPH paper but lacks the proceedings page numbers; the ACM Digital Library copy could not be read.

Problem

The paper starts from computer graphics, not vision. Scanned photographs were used as texture maps, image-based modelling and rendering used images as the basic modelling primitive, and rendered elements had to blend with photographs. All of this, the authors argue, needs to know how scene radiance became pixel values, and that mapping is rarely linear: if one pixel has twice the value of another, it probably did not see twice the radiance.

Their Fig. 1 traces the mapping through two pipelines. In film, exposure forms a latent image, development turns it into optical density, and a scanner converts density to voltages, digitizes them and remaps the result; printing adds its own curve. Film’s response to exposure is its characteristic curve, with a small response at zero exposure and saturation at high exposure. A digital camera’s CCD charge is proportional to irradiance, but most cameras remap it nonlinearly before storage, to mimic film, to anticipate the display’s nonlinearity, and often to reduce the 12-bit converter output to 8 bits. In both, the strongest nonlinearity is saturation, where every radiance above a level maps to the same maximum value. The composition of these steps is unknown in advance.

The paper’s motivation is limited dynamic range. Sunlit scenes and scenes with shiny surfaces or light sources cannot be captured without under-exposing some parts or saturating others; a series of exposures covers the range, but combining them requires the unknown curve. Section 1.1 lists what a radiance map would make possible:

• Image-based modelling and rendering. Existing systems assumed the same exposure settings and response for all images, which fails in a building with windows; recovering reflectance (up to full BRDFs) needs absolute radiance, including the values of specular highlights and light sources that a normal photograph clips.
• Image processing. Blurring, edge detection, colour correction and correspondence expect values proportional to radiance and misbehave near saturation; synthetic motion blur is the paper’s example.
• Compositing. Plates from a still camera, live action on another film stock and rendered elements could be combined as radiance and re-rendered through any one of the processes.
• A research tool. Measured radiance maps would allow quantitative comparison of rendering and global illumination with real scenes, and would let tone reproduction methods (the paper cites Tumblin and Rushmeier 1993, Schlick 1994, Ferwerda et al. 1996 and Ward, Rushmeier and Piatko 1997) be tested on real radiance maps as well as on synthetic ones.

Section 1.2 places the work. Photographers had long fitted scenes into film’s range by choosing shutter, aperture and development, with Ansel Adams’s system as one of the best-known examples, and scientific photography, in astronomy for instance, calibrated film’s nonlinearity with dedicated procedures. The authors wanted a self-calibrating procedure without charts or photometric instruments. Of the two prior works they cite, Madden (1993, a GRASP Laboratory technical report) combined several CCD integration times but had the sensor’s linear output and so faced no unknown response.

The other, Mann and Picard (1995), did address an unknown response. In their technical report they take two pictures whose exposures differ by a known factor $k$, so that a pixel recorded as $f(q)$ in one image is recorded as $f(kq)$ in the other; plotting one image’s values against the other’s (their “range-range” plot, in practice a joint histogram) and following a chain of pixels gives samples of $f$ at $q_0, kq_0, k^2 q_0, \ldots$ (Section 4). They then suggest restricting $f$ to the empirical film law $f(q) = \alpha + \beta q^{\gamma}$, whose exponent regression on the plot recovers, and combine the images with weights from “certainty functions”, the derivatives of the shifted response curves (Section 5). They describe the method as completely automatic. Debevec and Malik judged it “a rather limited method”: the parametric form is assumed arbitrarily, image noise is not treated satisfactorily, and the recovery uses only part of the data (Section 1.2).

Contribution

As the paper presents it:

• Response recovery from ordinary photographs. The response of any imaging process, film or digital, is recovered up to scale from photographs taken from one viewpoint at known exposure durations, assuming only that the curve is smooth and monotonic. The fit is a linear least-squares problem, small enough to solve in seconds, and the complete MATLAB routine, gsolve.m, is printed in the appendix.
• A high dynamic range radiance map, a weighted combination of all the exposures at every pixel, stored in floating point or in the RADIANCE format.
• Practical guidance: how many photographs are needed, how to extend the range of single film scans, how to obtain absolute radiance, how to handle colour, and how to render “virtual photographs” from a radiance map through any response curve.
• Demonstrations on a digital camera and on scanned colour print film, with renderings by a tone reproduction operator and a motion blur experiment that compares synthetic and real blur.

Method

Response recovery

The method rests on reciprocity: film density, or a sensor’s charge, depends only on the exposure $X = E\,\Delta t$, the product of irradiance and exposure time. The paper states that in typical print films reciprocity holds to within 1/3 stop for exposure times from 10 seconds to 1/10,000 of a second, without citing a source for the figure; for CCDs it holds when each site counts all photons absorbed during the integration time (Section 2.1). The input is a set of photographs taken from the same vantage point at known exposure durations $\Delta t_j$, of a static scene under unchanging light. With $Z_{ij}$ the value of pixel $i$ in photograph $j$, $Z_{ij} = f(E_i\,\Delta t_j)$ (Eq. 1), and with $g = \ln f^{-1}$, Eq. 2 reads

$$g(Z_{ij}) = \ln E_i + \ln \Delta t_j.$$

Since pixel values are integers from $Z_{\min}$ to $Z_{\max}$, $g$ is a finite table, and the paper solves for its values and the $N$ log irradiances by minimizing the objective of Eq. 3,

$$\mathcal{O} = \sum_{i=1}^{N} \sum_{j=1}^{P} \left[ g(Z_{ij}) - \ln E_i - \ln \Delta t_j \right]^2 + \lambda \sum_{z=Z_{\min}+1}^{Z_{\max}-1} g''(z)^2,$$

with $g''(z) = g(z-1) - 2g(z) + g(z+1)$. The second term couples neighbouring values of $g$; $\lambda$ is to be chosen for the noise expected in the pixel values. Fig. 2 gives the intuition: the values of one pixel across five exposures trace a segment of the curve whose height depends on the unknown $\ln E_i$, and the fit slides the segments of different pixels up and down until they line up into one smooth curve.

Three refinements complete the algorithm (Section 2.1). The scale is fixed by adding the equation $g(Z_{\text{mid}}) = 0$, so a pixel midway between $Z_{\min}$ and $Z_{\max}$ is assigned unit exposure. A hat-shaped weight $w(z)$ (Eq. 4), rising linearly from $Z_{\min}$ and falling to $Z_{\max}$, multiplies both terms, because the curve is steep and fits poorly near its ends. And only a sample of pixels is used: with $N$ pixels and $P$ photographs, the authors want $N(P-1) > Z_{\max} - Z_{\min}$ for a sufficiently overdetermined system, and for a pixel-value range of 255 and 11 photographs they judge about 50 pixels more than enough. Pixels should spread evenly over the range of values and over the image, and come from regions of low intensity variance, so that radiance is constant across the pixel and optical blur matters little; the authors picked them by hand, noting that this could easily be automated. Monotonicity is not enforced. The authors note it could be, by turning the problem into a non-negative least-squares problem, but they found the smoothness term enough. Running times were a few seconds. The weighted objective, the hat function and the merge are written out in high dynamic range imaging.

Radiance map, storage and number of images

With $g$ known, each exposure gives $\ln E_i = g(Z_{ij}) - \ln \Delta t_j$ (Eq. 5), and the radiance map averages these over all exposures with the same hat weights (Eq. 6). The curve also converts any later photograph taken through the same process, not only those used to recover it. The paper credits the averaging with reducing noise and film grain, and notes that, since the weights ignore saturated values, blooming (charge or light spilling from saturated sites into neighbours, footnote 5) has little effect (Section 2.2).

The map is computed in single-precision floating point and can be stored in the image format of Ward’s RADIANCE rendering system, which gives each pixel an 8-bit mantissa per colour channel and one shared 8-bit exponent, so that a high dynamic range colour image needs only a third more storage than a conventional RGB image (Section 2.2.1).

Section 2.3 counts photographs. Recovering the curve needs at least two, whose exposures are close enough that some pixels fall in the film’s working range (between the toe and the shoulder) in both; more reduce sensitivity to noise. Covering a scene needs at least $\lceil R/F \rceil$ photographs, with $R$ the scene’s range of radiances and $F$ the film’s working range. To use few photographs, the authors suggest recovering the curve once from three or four exposures one or two stops apart, reading off the working range (as many as five or six stops for the processes they had seen) and then shooting each scene with as many exposures as its range needs, sometimes a single one for a diffuse interior.

Film scans, absolute radiance, colour and virtual photographs

Sections 2.4 to 2.7 extend the method, mostly in a paragraph each:

• Single film scans (Section 2.4, done with Gregory Ward Larson). Many scanners, and Kodak’s PhotoCD process, delivered 8 bits per channel and thus only part of the range of print film. Scanning the negative with slide-film settings makes most scanners record its whole detectable range, at the cost of coarser quantization; scanning it twice at two density settings and combining the scans, each through its own recovered curve, avoids that cost.
• Absolute radiance (Section 2.5). The scale factor can be approximated from the film’s ASA rating, the shutter speeds and formulas in the IES Lighting Handbook, which the authors consider good enough for simulating glare and predicting scotopic vision, or obtained exactly by photographing a calibration luminaire of known radiance.
• Colour (Section 2.6). Each channel gets its own curve and therefore its own unknown scale. By default, $g(Z_{\text{mid}}) = 0$ in every channel makes the pixel $(Z_{\text{mid}}, Z_{\text{mid}}, Z_{\text{mid}})$ achromatic, which is right only if the channels respond equally to achromatic light near the middle of the range. Films are balanced for a particular illuminant such as sunlight, so the paper recommends scaling the channels to that illuminant’s colour, or calibrating against a luminaire of known colour.
• Virtual photographs (Section 2.7). Passing radiance back through a response curve (Eq. 1) renders the map as if photographed by the recovered process or by another one. Footnote 9 notes the assumptions: equal spectral responses, and no model of grain, chromatic aberration, blooming or the modulation transfer function.

Footnotes record practical details. In footnote 4 the authors checked the shutter of their Canon EOS Elan by recording its sound with a Macintosh: the times were accurate to within a thousandth of a second and were exact powers of two (1/64, 1/32, …, 16, 32 s) rather than the rounded values shown on the camera (1/60, 1/30, …, 15, 30). The same footnote discourages varying the aperture instead of the time, because of vignetting. Irradiance is proportional to radiance at each pixel, but with a factor that varies across the sensor, the $\cos^4$ falloff given by Horn, which, the authors write, most modern lenses compensate for at f/8 and smaller apertures (footnote 3).

Worked Example: How Many Pixels, and Which

The high dynamic range imaging example already varies $\lambda$ and compares merging weights. This one tests the paper’s sampling advice: the bound $N(P-1) > Z_{\max} - Z_{\min}$, and the rule to take pixels from regions of low intensity variance. A rendered scene spans more than five decades of radiance: a ramp of 4.5 decades from left to right, a gentler gradient from top to bottom, and horizontal bars that change radiance by a factor of four. A simulated camera with a film-like S-shaped response records it at the paper’s DCS460 schedule of 11 exposures, one stop apart, with 1% noise and 8-bit output. The response is then recovered from $N$ random pixels, with $\lambda = 100$. The paper’s MATLAB code multiplies the smoothness rows by its parameter directly, which therefore acts as $\sqrt{\lambda}$ (see Inconsistencies below); this code multiplies them by $\sqrt{\lambda}$, so its $\lambda$ is the $\lambda$ of Eq. 3. A second experiment shifts each frame by a random offset of about half a pixel, a stand-in for imperfect registration, which the paper’s PhotoCD scans needed, and compares random pixels with pixels whose 5 × 5 neighbourhood is smooth.

import numpy as np
from scipy.ndimage import shift, uniform_filter

# --- A rendered scene: a 4.5-decade ramp with high-contrast bars ------------
H, W = 96, 128
yy, xx = np.mgrid[0:H, 0:W]
L = 10 ** (-1.5 + 4.5 * xx / W) * (0.6 + 0.4 * yy / H)  # radiance 0.02 to ~900
L *= np.where((yy // 12) % 2 == 0, 1.0, 0.25)            # horizontal bars, 4:1

# --- Camera with a known film-like response, 8-bit output ------------------
def f(X):                          # exposure X -> pixel value in [0, 255]
    return 255 / (1 + np.exp(-1.6 * (np.log10(np.maximum(X, 1e-12)) - 0.5)))

def capture(t, rng, jitter):       # jitter: misregistration in pixels (std)
    Ls = shift(L, rng.normal(0, jitter, 2), order=1, mode="nearest") if jitter else L
    X = Ls * t * (1 + 0.01 * rng.standard_normal(L.shape))  # 1% noise
    return np.clip(np.round(f(X)), 0, 255).astype(int)

times = 2.0 ** np.arange(-5, 6)    # 11 exposures, 1/32 s to 32 s, one stop apart
B, P = np.log(times), len(times)

Xs = np.logspace(-6, 6, 400001)    # true g = ln f^-1, scaled so g(128) = 0
g_true = np.log(np.interp(np.arange(256), f(Xs), Xs))
g_true -= g_true[128]
w = np.minimum(np.arange(256), 255 - np.arange(256)).astype(float)  # Eq. 4 hat

def gsolve(Z, lam):                # weighted Eq. 3, smoothness rows times sqrt(lam)
    N, n = Z.shape[0], 256
    A = np.zeros((N * P + n - 1, n + N)); b = np.zeros(A.shape[0])
    k = np.arange(N * P); i, j = k // P, k % P
    A[k, Z[i, j]], A[k, n + i], b[k] = w[Z[i, j]], -w[Z[i, j]], w[Z[i, j]] * B[j]
    A[N * P, 128] = 1                                       # g(128) = 0
    for z in range(1, n - 1):
        A[N * P + z, z - 1:z + 2] = np.sqrt(lam) * w[z] * np.array([1, -2, 1])
    return np.linalg.lstsq(A, b, rcond=None)[0][:n]

lv = uniform_filter(np.log(L) ** 2, 5) - uniform_filter(np.log(L), 5) ** 2
smooth = np.flatnonzero(lv.ravel() < 0.05)   # 5x5 windows clear of the bar edges

def trial(seed, N, jitter, pool):
    rng = np.random.default_rng(seed)
    Z = np.stack([capture(t, rng, jitter).ravel() for t in times], 1)
    idx = rng.choice(pool, N, replace=False)
    err = np.abs(gsolve(Z[idx], 100) - g_true)[5:251]
    return np.sqrt(np.mean(err ** 2)), err.max()

def row(label, runs):
    r = np.array(runs)
    print(f"  {label:<28} rms {np.median(r[:,0]):.3f} [{r[:,0].min():.3f}-"
          f"{r[:,0].max():.3f}]   max {np.median(r[:,1]):.3f} [{r[:,1].min():.3f}-"
          f"{r[:,1].max():.3f}]")

print(f"{len(smooth)} of {H*W} pixels lie in smooth 5x5 windows")
print("error of g in natural-log units, codes 5-250, lambda = 100;"
      " median [min-max] over seeds 0-9")
print("registered frames, random pixels (N(P-1) > 255 needs N >= 26):")
for N in [5, 15, 26, 50, 200]:
    row(f"N = {N:<3} N(P-1) = {N*(P-1)}", [trial(s, N, 0, H * W) for s in range(10)])
print("frames misregistered by 0.5 px (std), N = 50:")
row("random pixels", [trial(s, 50, 0.5, H * W) for s in range(10)])
row("pixels from smooth regions", [trial(s, 50, 0.5, smooth) for s in range(10)])

Output (NumPy 2.3.5, SciPy 1.16.3; identical on two runs; about 2 s):

8704 of 12288 pixels lie in smooth 5x5 windows
error of g in natural-log units, codes 5-250, lambda = 100; median [min-max] over seeds 0-9
registered frames, random pixels (N(P-1) > 255 needs N >= 26):
  N = 5   N(P-1) = 50          rms 0.031 [0.014-0.108]   max 0.317 [0.059-1.009]
  N = 15  N(P-1) = 150         rms 0.011 [0.009-0.029]   max 0.064 [0.045-0.327]
  N = 26  N(P-1) = 260         rms 0.010 [0.008-0.013]   max 0.049 [0.035-0.104]
  N = 50  N(P-1) = 500         rms 0.006 [0.005-0.009]   max 0.039 [0.024-0.052]
  N = 200 N(P-1) = 2000        rms 0.005 [0.003-0.006]   max 0.029 [0.021-0.040]
frames misregistered by 0.5 px (std), N = 50:
  random pixels                rms 0.050 [0.037-0.128]   max 0.141 [0.105-0.235]
  pixels from smooth regions   rms 0.023 [0.011-0.038]   max 0.074 [0.022-0.111]

The bound counts equations against unknowns: $NP$ data equations against $N$ log irradiances and 255 free values of $g$ once the scale is fixed, so it asks the data alone to overdetermine the system (our reading). It is not a solvability condition. The smoothness rows tie every code to its neighbours, so the system has a solution for any $N$ (our reasoning), and even five pixels give a median root-mean-square error of about 3% in exposure, because each pixel sweeps through many codes as the exposure doubles. Below the bound, though, some codes are never observed and the smoothness term has to interpolate across them, and draws that miss whole stretches of the curve go wrong there: with five pixels the largest error along the curve has a median near 0.3 in log units over the ten seeds and reaches 1.0, a factor of $e$. Once $N$ reaches the bound, the worst error over the ten seeds stays near 0.1 or below, and adding pixels beyond the paper’s 50 improves the root-mean-square error only slightly. This matches the paper’s choice of about 50 pixels for 11 photographs.

Misregistration matters more than the number of pixels. A pixel near a bar edge sees different radiances in different frames, which breaks Eq. 2’s assumption of one $E_i$ per pixel; at 50 random pixels it raises the median error about eightfold. Taking the 50 pixels from smooth neighbourhoods, as the paper advises, halves the median error but does not remove it, since the ramp still varies within each neighbourhood. This is our reading of a simulation; the paper itself gives no such measurement.

Results

The results are demonstrations, shown in figures; there is no table and no comparison with independently measured radiance (our reading of the whole paper).

Kodak DCS460 (Section 3, Figs. 3–5). Eleven grayscale photographs of an indoor scene, filtered down to 765 × 509 pixels, were taken at f/8 with exposure times from 1/30 s to 30 s, each twice the previous one. The curve recovered from 45 pixel locations (Fig. 4, natural logarithm) shows that the camera remaps the CCD’s linear output nonlinearly, presumably, the authors write, to mimic film; outliers in the underlying data come from sensor artifacts, light horizontal bands across some of the darker images. The radiance map is displayed through the logarithm of its values (Fig. 5), with nothing saturated and shadow detail visible. The relative radiances printed at five marked pixels run from 1.0 to 15116.0, which the authors describe as over four orders of magnitude of useful dynamic range.

A church on print film (Figs. 6–8). Sixteen photographs were taken with a Canon 35 mm SLR on Fuji 100 ASA colour print film through a 15 mm fisheye lens at f/8, from 30 s to 1/1000 s in one-stop steps. The film was developed professionally and scanned on a Kodak PhotoCD scanner with its automatic brightness and contrast adjustment (“Scene Balance Adjustment”, footnote 10) turned off, so that every frame went through the same response. PhotoCD does not scan exactly the same area of each negative, so the frames were registered to sub-pixel accuracy by normalized correlation; footnote 11 notes that cinema film, scanned by its sprocket holes, avoids the problem. Curves for red, green and blue were recovered from 28 pixel locations (Fig. 7). Red and green agree closely, while blue rises above them at low exposures, which reveals a slight blue cast in dark regions; because the curves model it, it does not affect the relative radiances.

Fig. 8 shows the church radiance map six ways: (a) one 2-second photograph, lacking detail in its dark regions and saturated in the two rightmost windows; (b, c) linear mappings to the display with scale factors a thousand times apart (c maps the lowest 0.1% of the range), both with useful detail; (d) a false-colour image that, per the caption, shows over five orders of magnitude, with the highest listed radiance nearly 250,000 times the lowest, according to the text; and (e, f) renderings by the tone reproduction operator of Ward, Rushmeier and Piatko, (e) with adaptive histogram compression and (f) also simulating glare, contrast sensitivity and scotopic vision. Images (d–f) are credited to Gregory Ward Larson. The renderings show detail both in the rafters and in the rightmost window, which is more than a thousand times brighter, with the sun directly behind it.

Motion blur and virtual photography (Fig. 9). The authors blurred (a) a 37 × 1 pixel box filter over an ordinary photograph and (b) the same filter over the radiance map, then re-photographed (b) through the recovered curves at the photograph’s exposure time; (c) is a real photograph blurred by rotating the camera on its tripod during an exposure of the same length. In (a) the windows’ clipped values smear into a muddy, unsaturated blur; in (b) the true radiances stay above the saturation point over most of the blurred region and leave crisply defined saturated regions, which match (c) qualitatively in the bright regions. The real blur path was curved, which the demonstration did not model.

Inconsistencies in the Paper

A close reading of the authors’ copy finds several places where the text, the equations and the appendix code disagree. None changes the method, but readers porting the code should know them.

1. Equation reference. Section 2.1 says the appendix routine minimizes “Equation 5”, which is the per-exposure radiance formula; the appendix correctly names the objective of Eq. 3, whose weighted form is unnumbered.
2. $\lambda$ against l. The weighted objective has $\lambda \sum [w(z)\,g''(z)]^2$, but gsolve.m builds the smoothness rows as l*w(i+1)*[1 -2 1], so the squared residual carries l squared. The code’s l is $\sqrt{\lambda}$; a value of $\lambda$ taken from the equation and passed to the code smooths as if it were $\lambda^2$.
3. “SVD”. The text says the system is solved with the singular value decomposition, and the code comment says the same, but the code solves it with MATLAB’s backslash operator, x = A\b. For a full-rank system both give the same least-squares solution (our observation).
4. Argument documentation. The code header documents B(j), the log exposure time of image j, but the code reads B(i,j), one entry per pixel and image.
5. Scale constraint. The text sets $g(Z_{\text{mid}}) = 0$ with $Z_{\text{mid}} = \tfrac{1}{2}(Z_{\min} + Z_{\max})$, which is 127.5 for 8-bit data; the code fixes A(k,129) = 1, that is, $g(128) = 0$. The text also counts $Z_{\max} - Z_{\min} + 1$ values of $g$ in one paragraph and $Z_{\max} - Z_{\min}$ in another. The code allocates $n + 1$ rows beyond the data rows but fills only one scale row and $n - 2$ smoothness rows, which leaves two zero rows that do not affect the solution.
6. $\lceil R/F \rceil$. If $R$ and $F$ are ratios of maximum to minimum radiance, the number of photographs needed to tile the scene’s range is about $\lceil \log R / \log F \rceil$, not $\lceil R/F \rceil$: a scene of $10^5$ and a film working range of $10^2$ would need about three photographs, not a thousand. The formula is right if both ranges are measured in stops, which the text does not say (our reading; an ambiguity rather than a confirmed error).
7. Mann and Picard’s venue. The reference list gives the “IS&T 46th annual conference”; their technical report’s header gives IS&T’s 48th Annual Conference, Washington, D.C., 7–11 May 1995, with the same pages, 422–428. Kirk and Andersen (2006) repeat “46th”, and Mitsunaga and Nayar (1999) print “48th”.

Impact

The paper is widely cited in both graphics and vision; the works below show how it was used and judged.

• Image-based lighting. A year later, Debevec (1998) rendered synthetic objects into real scenes lit by measured light. A radiance map of a mirrored ball, assembled with the 1997 technique from digital photographs in one-stop steps from 1/4 to 1/10,000 s, served as a light probe; renderings lit by a single ordinary photograph instead, whose brightest regions were truncated to about 2% of their true values, came out visibly wrong for rough and diffuse materials (Section 3). Lighting virtual objects with measured light in this way is one of the techniques described under augmented reality.
• Response recovery in vision. Mitsunaga and Nayar (1999, Section 2) describe the 1997 algorithm as using about 11 high-quality photographs with precisely known exposures and no restrictive model of the response, well suited when the images are not noisy and the exposure times are precise. Their own method fits a polynomial and estimates the exposure ratios from rough guesses, such as the f-number readings of an inexpensive lens. Grossberg and Nayar (2003, Section 1) note that every such method makes some prior assumption about the response, and place the 1997 method at one extreme, imposing smoothness without any particular form, with Mann and Picard’s gamma-curve form at the other and Mitsunaga and Nayar’s polynomials in between.
• Merging weights. The hat weights became a reference point. Mitsunaga and Nayar (1999, Section 6.3) call both the hat and Mann and Picard’s gradient weights “somewhat ad-hoc”, the latter less so, and propose signal-to-noise weights instead. Kirk and Andersen (2006, Section 2.2) classify Eq. 6 as log-domain weighted averaging and note that the hat was chosen with the expected shape of the response’s derivative in mind. Their review also describes Robertson, Borman and Stevenson (2003), who start from inverse variances and weight by certainty times squared exposure time. Kirk and Andersen derive inverse-variance weights from a camera noise model; in their predicted signal-to-noise curves for a linear camera, these and signal-to-noise weights come out best, while the hat gives a less jagged curve because it down-weights images gradually before they saturate (Section 4). They also verified their camera’s nominal shutter times because of the paper’s report of a camera whose times followed powers of two (Section 3.1, footnote 1).
• The baseline for later HDR work. Twenty years later, Kalantari and Ramamoorthi (2017, Section 1) open by calling the merging of differently exposed images, citing Debevec and Malik, the most common way to capture HDR images, and they built the ground truth of their training data with a triangle weighting they describe as similar to the paper’s (Section 4).

Limitations

Some limitations are stated in the paper:

• Static, registered, known exposures. The input must come from one viewpoint, with known exposure times, a static scene and unchanging light (Section 2.1). The church scans needed registration, and footnote 4 shows how far a camera’s displayed shutter times can be from the true ones: 1/60 against 1/64 s, or 15 against 16 s, about 6–7% (our arithmetic).
• Reciprocity fails at very long and very short exposures (Section 2.1).
• Manual sampling and unenforced monotonicity. Pixels were chosen by hand, and monotonicity was left to the smoothness term (Section 2.1).
• Colour balance depends on three unknown scale factors, which need an assumption or a calibration (Section 2.6).
• Spectral response. The quantity recovered is irradiance weighted by each channel’s spectral response, neither a radiometric nor a photometric quantity in the strict sense (Section 2.1), and virtual photographs assume equal spectral responses and ignore grain, chromatic aberration, blooming and the MTF (footnote 9).

Later work identified others:

• Dependence on precise exposures and clean images (Mitsunaga and Nayar, 1999, Section 2), which their method relaxes.
• Monotonicity in practice. On three exposures of a scene taken with a Nikon 990 Coolpix at exposures $e$, $2e$ and $4e$, Grossberg and Nayar (2003, Section 9, Fig. 8) recovered inverse responses with the 1997 method at $\lambda = 8$, 32 and 128 and with the polynomial method; all were reasonable fits, but only their own constrained model was monotonic, and it came closest to chart measurements. They do not say whether their $\lambda$ is that of Eq. 3 or the code’s parameter. On this three-exposure set, then, the smoothness term did not keep the curve monotonic, as it had in the authors’ experience.
• Noise. The hat weights a value by its code alone, not by its noise; Kirk and Andersen’s (2006) predictions favour variance-based weights, and the example in high dynamic range imaging shows the cost on a simulated camera with shot noise.
• Motion. Merging produces ghosting when the scene moves or the camera is hand-held (Kalantari and Ramamoorthi, 2017, Section 1), outside the paper’s assumptions.

What Came After

Two changes in cameras narrowed the problem the paper solved. Many consumer cameras can store linear RAW data, and most machine-vision cameras can run in a linear mode, which removes the uncertainty of the response curve (Kirk and Andersen, 2006, Section 1); and sensors that extend the range within one frame avoid bracketing altogether. The work moved to what the 1997 paper assumed away: alignment and ghost removal for moving scenes, learned merging, single-image reconstruction of clipped regions, and the display of radiance maps by tone mapping. High dynamic range imaging describes each of these, and computational photography places multi-exposure capture among burst and multi-frame methods. The response-recovery problem itself survives wherever only processed images are available, which is the case the paper was written for.