![[scrshot_images.jpg]](scrshot_images.jpg)
Imagine a photographic motif with exceptional contrast in brightness, say a whitewashed farm building in heavily backlight, with an open barn door in it, giving a view into the dark inside of the building. Naturally, the photographer has to decide here, by appropriate choice of exposure time and aperture, which part of the motif (whether sun or house wall or interior) he will get good contrast; the other parts will become heavily overexposed (white) or underexposed (black). Nevertheless, every part of the motif by itself could, by applying the appropriate exposure time, be moved into the working range of the camera and could be mapped rich in contrast.
An exposure bracketing would match the optimum exposure for every single part of the motif — unfortunately, with each optimum exposure a separate single image in a series of pictures. Now, the idea is to generate from an exposure bracketing one single image, where the "real" circumstances of brightness, the so to say physical radiance data are recorded. To be practical, we represent such data not by integers as in the case of ordinary digital pictures, but by floating point numbers, whose exponent makes it easy to express differences in magnitudes of many decimal powers. Such images are called HDR images, the initials "HDR" stand for "High Dynamic Range". A first quantitative orientation concerning the dynamic in brightness of an image provides the ratio of the brightest value to the (non-zero) darkest value. For an 8-bit integer image this ratio can be in a single color channel maximally 255 : 1, for 16-bit integers 65535 : 1, for a 32-bit floating point image about 10^32 : 10^-32 = 10^64.
Concerning contrast and detailedness, the brightness information coded in HDR images normally exceeds the dynamic range of usual reproduction media as monitor or paper. For instance, if the whitewashed wall is represented by the brightest, i.e. ungrayed white, then the brightest sun will appear just as white ("bright") — the white of the paper or screen, whereas, in terms of measurable radiance intensity the wall may have (in any arbitrary unit) a value of 1, the sun corresponds to a value of 1000! One of the advantages of HDR images is that the decision on which information will be of interest does not have to be made until the moment of looking at it (instead of during the shot). Another advantage is that HDR images provide a basis for improving the contrast at selected spots, or to generate more realistic images of high-contrast scenes by using algorithms modeling the automatic contrast compensation of our sense of sight.
The rest of this text consists of two parts, a theoretical one and a practical one. In order to work with the program it isn't necessary that you understand the theoretical chapter in detail, nevertheless a rough idea of the underlying theoretical concepts facilitates a reasonable evaluation of the results. If wished, you can jump straight to the practical part.
As far as I know, the method we use was first described by P. E. Debevec and J. Malik [1].
The radiation emanated or/and reflected by the objects of a scene goes through the aperture and lens system and causes a certain light intensity E at a pixel in the plane of the photo sensors. During the exposure the amount of E is supposed to be constant. Integrated over the exposure time Δt we get from E the light quantum X = E⋅Δt, which defines the exposure result. In a following digitization process X is transfered by an unknown non-linear function z=f(X) into an integer value z (e.g. z=0...255 for 8-bit converters) — values which we read out later as the digital RGB values. The "response function" z=f(X) is a device property, it characterizes the entire digitization process X --> z of a color channel including photo sensor characteristic and possible camera-internal correction algorithms. Even the intermediate step of a first chemical developed photo material, which has been digitized later, can be incorporated here. The task is, for given z-values of an exposure bracketing and known exposure times to reconstruct the function z=f(X)=f(E⋅Δt), and its inverse, and finally to calculate the values of E. The intensity values E constitute the content of a HDR image.
Let our bracketing series contain P exposures. Let us pick out N pixels with the corresponding (unknown) light intensities E_i (i=1,...,N). In the j-th exposure (j=1,...,P) the i-th pixel may give the output z_ij. Accordingly to z=f(X)=f(E⋅Δt) then has to be
z_ij = f(E_i⋅Δt_j) for all i=1...N, j=1...P
or inverted and logarithmised
g(z_ij) = ln(E_i) + ln(Δt_j) for all i=1...N, j=1...P,
where the abbreviation
g(z) := ln[f^{-1}(z)]
was introduced. Known are the times Δt_j, unknown are the ln(E_i) and the function g(z). The domain of the integer variable z is finite, let it be the M numbers z_1,..., z_M (e.g. 0,...,255 for 8-bit converters). Recovering the function g(z) means therefore finding out for the M possible arguments z_k the related values g(z_k). — The whole problem we now formulate as the task, to determine the M unknowns g(z_k) and the N unknowns ln(E_i) so, that the quadratic objective function
N P
Ω = ∑ ∑ [g(z_ij) - ln(E_i) - ln(Δt_j)]^2 + λ ∑ g''(z)^2
i j z
becomes a minimum. The first term secures the compliance of the above equations in a least square sense, the second one is a smoothing term, it contains the sum of squared second derivations and punishes curvatures of g(z). With the parameter λ>0 the balance between both terms can be adjusted, high λ takes curvatures more into account and promotes smoother curves at the cost of the accuracy of the first term.
Because the objective function Ω is quadratic in the unknowns g(z_k) and ln(E_i), minimizing of Ω leads to a linear least squares problem, which is solvable with straightforward standard methods like QR-decomposition or singular value decomposition (SVD).
After the function g(z) has been found in this way and so exp(g(z)) = f^{-1}(z), the z-data can be converted via f^{-1}(z) = E⋅Δt into their related intensities E: E_i = f^{-1}(z_ij) / Δt_j.
Three Remarks:
1.) The solutions for E and f^{-1}(z) are unique only so far, up to an arbitrary constant factor α>0, because the equations remain complied, if ln(E) is replaced by ln(E)+ln(α) and g(z) by g(z)+ln(α). To make it unique, the additional condition g(z_mid)=0, i.e. f^{-1}(z_mid)=1 is added, where z_mid means the middle z-value. I.e. it is stipulated that the middle z-value corresponds to the light quantum X=1. This is a priori a mere arbitrary stipulation about the absolute amount of the X- respectively of the E-values, which therefore can have only relative meaning. For absolute radiance data a rescaling would be necessary. Further this normalization happens for each color channel independently, so that a gray tone is postulated at this point. However, especially for camera raw data, a sensor output of (z_mid, z_mid, z_mid) doesn't have to represent a gray tone necessarily, and a subsequent color correction would be needed.
2.) Digital cameras often work with motif-dependent correction algorithms. If the same algorithm is used for all pictures of an exposure bracketing, this is not a problem, one will merely get (perhaps) different response functions depending on the motif. However, it would be problematic if different correction algorithms are used within an exposure bracketing (assumption 2 violated) or on the other hand, if zones within an image would be post-processed differently (assumption 3 violated). Though a result will be achieved then — in the end it is a minimization task — but perhaps not a good one. One should be on the safe side with camera raw data.
3.) The assumed monotonicity of the response function is not enforced in the numerical procedure, it can result also non-monotonic curves. This is often a hint at problematic ("against rules") input data.
The CinePaint tool for creation of HDR images, which we will introduce now, has become part of the official CinePaint distribution since version 0.20-0. The image editing program CinePaint was suggested as platform for an open source solution by the fact, that it can handle — unlike GIMP — floating point data (32-bit) inherently. In the following I describe some essentials of the new plug-in.
The tool is started in CinePaint via the <Toolbox>-menu "File" -> "New From" -> "Bracketing to HDR". It opens a separate main window, wherein under "File" -> "Open" the images of an exposure bracketing are loaded; all image formats supported by CinePaint are allowed. The images are sorted automatically by average brightness per pixel, the darkest one on top, and thus assumed to be the shortest exposed.
At the moment reading of EXIF-data (aperture, exposure time, sensitivity, ...) is not supported. However, at the moment anyhow only relative radiance values in the sense of remark 1 are created, no absolute ones. For this purpose it is sufficient to know the ratio of the light quanta with respect to the exposure times of the several images. For instance we can set the exposure time of the shortest exposed image arbitrarily =1 and refer the remaining times to it. With the "stops" choice menu one can generate exposure times increasing by a constant factor. Following usual photographic practice the factor F itself is not specified, but (in keeping with F = 2^stop) the number of stops of exposure difference between the exposures is entered (each stop representing a doubling/halving of the amount of light). This stop value corresponds just to the ratio of aperture values, an alteration by an aperture stop halves or doubles the aperture and so the light quantum X; in other words, doubling of the aperture and doubling of the exposure time are considered to be equivalent in this context. Thus, if the exposure bracketing was not accomplished by extending of the exposure time, but by successive increasing of the aperture (= decrementing of the aperture stop), it can specified here as well: steps by one aperture stop: stop=1, steps of a third: stops=1/3 etc. — If required, the exposure times can be edited also separately.
In the input field "grid points" the number N of pixels which are to be singled out from an image can be specified; values between 50 and 100 are a good choice. The larger N the closer meshed a response function is computed, but also the more time takes it. (To be more precise, N denotes not simply a number of pixels but, more in the abstract, of tuples, of which there is more in the next section.)
In the input field "smoothing" one can adjust the parameter λ mentioned above. The goal should be to get a monotonic, smooth response function with as small a λ as possible. Here one can experiment a little bit. I got fine results with λ=10...20 for input data I had confidence in. Understood that with a sufficient high λ>200...500 one can enforce smooth, monotonic response functions even for most heavily irregular input data (see section "Follow-up Curves", whose response curves oscillated wildly for λ=20), but of course the substance of such curves tends to zero.
On to the explanation of the diagrams I called "follow-up curves". Assume our exposure series contain three images. The numerical procedure expects then N triple [2] (N = grid points) of z-values as input information, say
(9,43,122), (15,55,135), ... (N exemplars)
with the following content: The z-value 9 in exposure number one corresponds to the z-value 43 in exposure number two and to the z-value 122 in exposure number three. The relationships of the increases with the exposure time is defined by the characteristic z=f(E⋅Δt), which we try to extract. Assuming that a certain pixel has in the three exposures the z-values (9, 43, 122). But in the first image there are altogether 1000 pixel with z=9. Theoretically, if all assumptions were fulfilled, all these 1000 points should have then in the second image the value 43 and in the third one the value 127. But de facto these "follow-up values" will scatter, perhaps because, in spite of assumption 4, some objects in the scene moved (leaves in the wind, clouds floating across the sky), perhaps because the camera positions were not exactly identical, not to mention fabrication tolerances of the photo sensors. Furthermore, the neighborhood of a pixel plays a role. A single z=9 value surrounded by z=0 and z=200 values might be the victim of blooming and make its follow-up values less certain than a z=9 value in similar surroundings. Obviously there is some creative leeway in the intelligent selection of the N tuples.
Currently it is done this way: At first the "follow-up curves" are acquired which means that: For all (z=0)-points in image 1 their corresponding values in image 2 are looked up and the average <z> of all those as well as their deviation Δz is determined. This gives for z=0 the triple (0, <z>, Δz). And so for all z-values 0, 1, ..., z_max in image 1. These averages <z> plotted over z gives the follow-up curve of image 2 to the reference image 1. Of course every other image could be used as a reference image as well. Such curves are shown in the tab "FollowUpCurves", the deviation is represented by vertical error bars. Follow-up curves are a good indicator of the quality of a data base. Ideally their deviation should be minimal and their slope uniform. Below on the left is an example, where the middle image (green curve) was the reference image. (The "pseudo follow-up curve" of the reference image itself is always a straight line and without deviation and is plotted only for orientation purposes.) On the right we see the same exposures, but the blue image was displaced sideways here by an arbitrary 20 pixels against both the others. The deviation of the blue curve is now considerable and the behavior non-monotonic. Clearly, a response function determined from such a data base will look strange and the related HDR image could have some odd color gradients. We also note this: in that situation, had the blue image be chosen as reference image, then both follow-up curves (green and red) would be abysmal instead of just the blue one. By the way, the black line shows the histogram: the number of points in the reference image with the corresponding z-value; what's missing is a scale at the right side.
![[scrshot_follpanA.jpg]](scrshot_follpanA.jpg)
|
![[scrshot_follpanA.jpg]](scrshot_follpanB.jpg)
|
From the follow-up curves, namely those for the currently chosen reference image, the N tuples mentioned above are selected internally (N vertical slices through a follow-up curve diagram). A valuation of the neighborhood of the pixel thus doesn't take place at present. Here remains room for alternate algorithms.
The module "Shifts" is planed to automatically correct small displacements of the input images, as happens in amateurish attempts without a tripod, camera put on a shovel handle and such, which the program author has had high hopes, tried, tried again, and then tried to correct and undo the damage. So far this module exists only in a very raw version, and it is mentioned here only in passing.
Let us start with a typical exposure series. By the way, the brightest parts of the filament lie at the saturation limit in even the darkest photo (i.e. out of the working range) and would certainly tolerate further darkening:
![[scrshot_reihe_klein.jpg]](scrshot_reihe_klein.jpg)
As soon as at least two images have been loaded, a HDR image can be created at any time by clicking the "HDR" button, using just the user supplied numerical parameters. It is better if one proceeds a bit more systematically, like this:
After loading the images and adjustment of the exposure data, take a first look at the quality of the input data in the "FollowUpCurves" diagram; play with different reference images! To avoid frequently switching between the two tabs cards "Images" and "FollowUpCurves" in the next steps, you can get the second panel in a separate window via the menu item "Analyze" -> "FollowUpCurves".
Now, to decide on the method of generating the N tuples. At present this is restricted to selecting a reference image, choosing one for which the follow-up curves take on a favorable shape. Later on, alternate algorithms could be written, but more probably this whole task will be automated, at least for the default choice.
![[scrshot_ccd.jpg not found]](scrshot_ccd.jpg)
Clicking the "Compute Response" button
starts the calculation of three response functions, one for each
color channel, which can be evaluated in a window which will open
automatically. Any time a non-monotonic characteristic is cause
for concern, making one check the choice of the reference image
(in the "FollowUpCurves" panel) or try a higher smoothing
coefficient λ and recompute the response functions. If
no reasonable characteristics can be found with moderate λ ≤ 50,
a problem with the input data probably exists. Currently there
are four parameters which affect the response functions for a
given image data set: 1) the exposure times, 2) the number N of
grid points, 3) the smoothing coefficient, 4) the reference image.
Because the calculation of response functions is the most time-consuming
step, re-calculations after changes in any one of the four are
not done automatically, but require clicking the "Compute
Response" button.
![[scrshot_hdr.jpg]](scrshot_hdr.jpg)
A click on the "HDR"
button computes the resultant HDR image using the last calculated
response functions; it appears in a new window as a CinePaint
image. The "32-bit IEEE Float" in the window title announces,
that indeed a floating point image has been created. With the
CinePaint color picker <image menu> -> "Tools"
-> "Color Picker" RGB values can be read out, and
you will notice that with floating point numbers where. By convention,
1.0 is the brightest representable value in each color channel,
in a HDR image, defined values can now occur which are greater
than 1.0.
![[scrshot_hdr_gamma.jpg]](scrshot_hdr_gamma.jpg)
![[scrshot_gammaexpose.jpg]](scrshot_gammaexpose.jpg)
The CinePaint tool called in the <image menu> under "Image"
-> "Colors" -> "Gamma-Expose" [3] allows simulation of various exposures with
its slide control "expose" and allows you a quasi review
of the light situations of the input images, but now continuously.
With the slide control "gamma" we can deviate from the
linearity (γ=1) of the brightness transmission profile.
Values over 1.0 emphasize the dark and the bright domains and
move differences in brightness, which had disappeared into the
black or white, into the visible, of course at the cost of the
contrast of the previous middles. This allows us to make visible
e.g. subtleties of the filament together with nuances from dark
regions in a single image, as demonstrated by the above example.
Of course, more sophisticated manipulations will avoid the graying
of the mid-tones, but this is beyond our subject here, where we
are generating the data base for such things.
Certain user actions make the buttons "Init" and "Compute Response" become yellow or blue colored and sometimes inactive (grayed out). This gives useful information about inner program states and whether the presented data is up-to-date.
In the program kernel there are two central instances, lets call them the "image container" and the "HDR calculator". The image container contains all images currently loaded and listed in the table. But often, even in experimental phases, one wishes to proceed with the HDR calculation (response function and image) using not all of these images, but only with some of them. To simplify this the green activation buttons were introduced, allowing single images to be deactivated. A deactivated image appears grayed out.
By default, newly loaded images are always activated.
The "HDR calculator" performs the real computations. It has to be initialized, and it is initialized in fact only with the activated images. Hitting the "Init" button forces such an initialization, with synchronization to the present state of the image container. Because every HDR computation needs at least two input images, the HDR calculator can only be initialized when at least two activated images are present in the image container. If an initialization is impossible in the current state of the program, the "Init" button appears inactive (grayed out). — After loading of new images the program always tries to (re)initialize the calculator automatically. If at first only one image was loaded into an empty image container, automatic initialization fails; the "Init" button appears after that grayed out — and blue. Blue means that no HDR calculator could be created at all. — If an HDR calculator exists already and an image gets de/activated, the "Init" button becomes yellow. Yellow indicates that the images used in the HDR calculator no longer coincide with the currently activated ones: if you want results for the present activated images a reinitialization would be necessary.
In summary: The "Init" button is blue if no HDR calculator exists, it is normal colored if an initialized HDR calculator exists and the used images coincide with the currently activated ones, it is yellow if a difference has occurred in activated images; and it is — in which color ever — inactive (grayed out), if a (re)initialization to the present set of activated images is impossible.
The display of colors of the "Compute Response" button acts quite similarly. It is blue if no characteristics have been calculated at all — or have been destroyed already again, as is the case after each initialization; it is normal colored, if calculated characteristics exist (always visible in the "Response functions" window) and belong to the current set of parameters; it is yellow if changes in parameters has been done which would lead to other characteristics than the present ones. And the button becomes inactive (grayed out) if the current program state doesn't allow a (re)calculation.
Hitting the "HDR" button causes an HDR image to be created, using the last determined response function, if such exist, in a straight computation, or the complete calculation, if no curves are available. The color of the adjacent "Compute Response" button indicates what will happen at any time: If it is gray, only the HDR image will be created; if it is yellow, ditto, though using out of date characteristics, which no longer represent the current parameter set; if it is blue, the complete computation will be done and accordingly it will take longer.
This is the first release of the program, so bugs and disasters are to be expected. Comments, suggestions, corrections concerning the program itself as well as this text are very welcome.
Thanks to Ted Miller, Elkhart, IN, USA, for revisions to the English translation.
Hartmut Sbosny <hartmut.sbosny@gmx.de>
Chemnitz, 2005
| [1] | P. E. Debevec und J. Malik: "Recovering High Dynamic Range Radiance Maps from Photographs", in SIGGRAPH 97; [pdf document, 1.4 MB]. Link to Debevec's HDR web site. (return to text) |
| [2] | Triple for 3 exposures, P-tuple for P exposures (return to text) |
| [3] | In version 0.20-0 an additional opportunity has become available under "View" -> "Expose". While with the previous "Gamma-Expose" the image data itself could be changed, the new "Expose" affects exclusively on the view of them. (return to text) |