RECONSTRUCTION OF SKY ILLUMINATION DOMES FROM GROUND-BASED PANORAMAS

The knowledge of the sky illumination is important for radiometric corrections and for computer graphics applications such as relighting or augmented reality. We propose an approach to compute environment maps, representing the sky radiance, from a set of ground-based images acquired by a panoramic acquisition system, for instance a mobile-mapping system. These images can be affected by important radiometric artifacts, such as bloom or overexposure. A Perez radiance model is estimated with the blue sky pixels of the images, and used to compute additive corrections in order to reduce these radiometric artifacts. The sky pixels are then aggregated in an environment map, which still suffers from discontinuities on stitching edges. The influence of the quality of estimated sky radiance on the simulated light signal is measured quantitatively on a simple synthetic urban scene; in our case, the maximal error for the total sensor radiance is about 10%.


INTRODUCTION
1.1 Why does sky illumination matter?
The recent development of on-line viewers of street level images leads to the acquisition of a huge number of high-resolution, georeferenced, urban, terrestrial images.Panoramic images are also used for architectural or archaeological surveys.Understanding the physical process of the formation of digital images is a key to the development of many applications, for professionals as well as for the general public.For instance, relighting or augmented reality can be performed from an estimation of a reflectance map of the scene (Yu et al., 1999).While processing the radiometric information carried by the pixels values, it is of primary interest to estimate precisely the light sources of the scene.For the outdoor scenes, the single light source is the sun, and the medium through which its light reaches the scene is the atmosphere.The sky is a complex participating medium, that may be highly heterogeneous in the presence of clouds, and that combines multiple physical phenomena (molecular and aerosol scattering and attenuation): its modelization for radiative transfer is thus a difficult task.In this work, we propose an image-based approach to retrieve the downward sky radiance from the whole dome above the scene.The extraction and aggregation of sky pixels from the images highlight radiometric issues in the acquired images; these issues are tackled with a non-physical method.

Related work
A simplified model of sky radiance can be computed by using a radiative transfer code, such as 6S (Vermote et al., 1997).This method implies the knowledge of the atmospheric composition at the time of acquisition.This composition can be characterized by the optical depth of the different gases and aerosols, that can be obtained by meteorological measures.But as these measures are costly, they usually are not available .It can also be characterized by one simple observable parameter: the optical visibility, which is easier to measure.However, the radiative transfer computation is usually limited to atmosphere made of homogeneous parallel planes (Schanda, 1986).Though it is sufficient for low resolution applications such as meteorology or climatology, it can be a strong limitation for processing high resolution terrestrial images, in which local variations of the sky downward radiance can have a significant effect, as in the presence of clouds.In a non physically-based way, there exist empirical parametric models to represent the sky radiance; they are presented in (?).The model of (Perez et al., 1993) is the most popular, and its small number of parameters is valuable for inversion purpose.These models are valid for clear or continuously veiled skies, but cannot represent clouds.In computer graphics, (Debevec, 1998) proposes to retrieve an environment map by imaging a mirrored ball.It can be used to photo-realistically introduce objects in a scene.However, the estimated radiance with this method is valid only for the location of the mirrored ball, and cannot be used for other points of the scene in the generic case of non-lambertian materials and participating medium.That is why we prefer an estimation of sky radiance for the whole dome.
After a formulation of the problem (section 2), two approaches are presented for estimating a sky radiance map; they are then used jointly for attenuating the effect of radiometric artifacts (section 3).The environment maps resulting of these approaches are then presented for a urban dataset acquired by a mobile-mapping vehicle, and the influence of their quality on a simulation problem is then quantitatively measured (section 4).

Light sources model
Performing a complete simulation of the physical interaction between the sky and the solar light (scattering and attenuation) is a very costly computation using volume rendering (Pharr and Humphreys, 2004).Furthermore, it requires the knowledge of the composition of the atmosphere, that is usually not homogeneous, nor made of homogeneous parallel planes (mostly because of the presence of clouds).The main idea in the proposed paper International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XXXIX-B3, 2012 XXII ISPRS Congress, 25 August -01 September 2012, Melbourne, Australia is that if the sky is now considered as a light source rather than a part of the scene, its complete modelization and simulation is no longer required and it suffices to sample the downward radiance using sky-facing pixel values.The light source from the terrestrial point of view can thus be considered as the sum of a solar component and a sky component.
The sun component can then be modeled independently from the sky, using a point light source located at an infinite distance moving along the well-known solar trajectory.The sky illumination is however more complex : its is a time-varying 4D + T lightfield.Our main assumption is that every point of the scene receives the same radiance from the sky seen in a given direction This neglects time variations between acquired images and their parallax.These limitations are reasonable when the weather is good (e.g.there is no low-altitude fog or cloud), since the scene (typically a street) is relatively small and we can neglect the scattering in the low layers of the atmosphere.The sky component can thus be modeled as an environment map, i.e. a 2D map giving for each direction of the upper hemisphere the value of the downward radiance coming from the sky.

Input data
The available data is a set of georeferenced terrestrial images acquired by a mobile-mapping vehicle in a urban area, with calibrated cameras, distributed so as to cover virtually the whole hemisphere.The images are georeferenced by processing the GPS and INS data.They are acquired with 12 bits per pixel, then are corrected of geometrical distortion and of vignetting (using a flat-field).The flat-field correction has to be performed carefully: the bit depth of the resulting image is set to 16 bits per pixel in order to avoid artificial overexposure in the corners (because of the division by a value < 1).Furthermore, rosy artifacts can appear in the areas where the blue channel is overexposed while the red one is not (an example can be seen on figure 1); that is due to a higher coefficient for the red channel for the white balance.For these reasons, the overexposed pixels are detected before the flat-field correction.
These images are affected by different undesirable effects due to inherent limitations of CCD sensors: overexposure in the solar  halo, flare effect, bloom effect, inaccuracy of the flat-field correction, and vertical smear (figure 2).The vertical smear, usually negligible in the absence of very high objects, can become very important when the sensor is worn-down (figure 3).These effects, hard to model and correct, affect the accuracy of the physical measure of the light's energy, and imply discontinuities while stitching several images into a panoramic.We use images taken in a very short period of time (about 1 minute) to generate one environment map, so the sky can be considered as static, even if all images are not taken from the same location.
Examples of images used for the results section (4) are shown on figure 4.

The problem: radiance estimation from images
Due to overexposure, the solar irradiance is not readily accessible in general purpose images.Its irradiance may be computed by knowing the optical thickness of the whole atmosphere, or by indirect estimation using pairs of points located on both sides of a shadow edge (Thomas et al., 2008).In this article, we only focus on sky radiance estimation.The problem is the computation of the sky radiance coming from all points of the upper hemisphere, using images acquired by a ground-based mobile-mapping vehicle, with an accurate geolocation, but affected with radiometric artifacts.The method aims at working for all weather conditions.This problem is ill-conditioned since the true signal due to sky radiance and the noise signal due to radiometric artifacts cannot be demixed.

SKY RADIANCE ESTIMATION
In this section, we propose two different approaches to estimate the radiance map from the images: aggregation of pixels detected as sky (including the clouds), and estimation parametric model with blue sky pixels (avoiding clouds).These two methods can be mixed in order to reduce the influence of artifacts and try to enhance the accuracy of the sky radiance estimation.The corresponding results are presented in section 4.1.

Extraction and aggregation of sky pixels
The first task is to extract the pixels that see the sky in the images, in order to set them apart from pixels of building or vegetation.We use simple thresholds on color values, set by experiments on the data set.These thresholds must lead to a good balance between accepting clouds and rejecting highly illuminated walls.The reflexions on cars are not a problem because pixels below the horizon are not used, but reflexions on windows are usually above the vehicle; a morphological filtering is used to remove most of the small areas detected among the buildings.The thresholds used are as follows, where the cmin threshold depends on the camera response and its bit depth : The main errors remaining are overdetection of large windows and of bright walls above the horizon, as well as a lack in detection of dark clouds.The sky areas affected by strong bloom are also badly detected, because of their having a high red component (figure 5).
The aggregation of the detected sky pixels is made by a simple bilinear interpolation on a regular grid, averaging the values from different images if they overlap.It leads to artifacts in the environment map: discontinuities are visible between the areas extracted from different images, for the reasons mentioned in section 2.2.Furthermore, there are holes in the environment map, for the directions where no sky is visible or detected in any image; this occurs principally for the areas close to the horizon.We limit the extent of these holes by using images taken while the vehicle goes through a crossroads, so we have a maximum solid angle of visible sky.Then, an interpolation is performed to fill the holes.This interpolation leads to more artifacts in the areas with no sky pixels extracted (figure 8).The discontinuities can be reduced by image processing techniques of blending, for instance by solving a Poisson equation, using only reliable pixels in the data term (Bhat et al., 2008).However, the discontinuities mainly affect the visual quality, and are not a problem in themselves for the light simulation.The problem is to retrieve an environment map with values closest as possible to the sky radiance at the time of acquisition.

Estimation of a parametric model
Another possible approach is to use a parametric model to reconstruct the environment map.The most widely used model is the Perez model (Perez et al., 1993), that uses 5 parameters a, b, c, d, e to describe the fall-off of the light around the sun position, plus 1 parameter Lz to set the zenithal radiance.However, this model is limited to the description of low spatial frequency phenomena, and therefore cannot model the clouds.The detection now avoids the clouds and extract only the background blue sky pixels, by using hue and saturation thresholds proposed by (Schmitt and Priese, 2009): The overexposed pixels, detected before the flat-field correction, are rejected.The parameters estimation is then performed by non-linear least square fitting, and gives coherent results as long as the sky pixels extracted are scattered homogeneously in the hemisphere.However, this estimation is affected by the radiometric effects mentioned in 2.2.Predominantly, the bloom effect in the images for which the sun is close to their field of view, leads to an overestimation of the values of the estimated Perez model (figure 8).

Using the parametric model to correct the environment map
In order to reduce the influence of radiometric artifacts, we propose to use a correction based on an estimated Perez model.We select images that do not have the sun in the hemisphere in front of the camera, so they are only little affected by radiometric artifacts, thus their radiometry is reliable.The parameters Lz, a, b, c, d, e are estimated with these selected images.We then compute the difference between the blue sky pixels of the original images and the corresponding reconstructed values with the Perez model.This difference, computed independently for each channel, does not fit a simple parametric model (see figure 6).Indeed, this difference is due to radiometric artifacts (noisy and hard to model in a physical way), but also to the misfit of the Perez model in the presence of heterogeneous veil.We assume that images classified a priori as having artifacts only suffer from a veil (i.e. an additive constant, estimated as their mean error).They are thus corrected by subtracting their mean error.The complete processing algorithm is summarized in figure 7.

Environment maps
Using 13 images, the environment maps computed with the different techniques presented in 3 are shown on figure 8.The Perez model is estimated with only 7 reliable images.The per-channel corrections computed from the estimated Perez model slightly attenuates the discontinuities in the stitching process, mostly by correcting the red values, that tends to be increased by artifacts.Though the Perez model is visually nice, it is probably more accurate for the later simulated light to take into accounts the effect of clouds.We can notice that this dataset is particularly tricky, because the sun itself is masked by clouds; the overexposed areas inside these clouds cannot be filled by the values of the Perez model.
The influence of the quality of this sky radiance estimation on the light signal is discussed in the next section (4.2).
Figure 7: The steps of the proposed algorithm.

Influence on the simulated light
The purpose of the sky radiance reconstruction is to use physicallybased methods for extracting information form the pixels of the images.The estimation of the reflectance of the materials of the scene is an important step for these methods, and is usually performed by minimizing the difference between acquired images and simulated images with unknown reflectances (Coubard et al., 2011).So the influence of the environment map on the simulated light is very important to know.That is why we compare the simulated signal in a simple urban scene (a street lined with buildings), with several estimated environment maps, to quantify the impact of their inaccuracy.The light is simulated by raytracing, with a code based on LuxRender (Pharr and Humphreys, 2004).The maximum number of interreflections between objects is set to 2, and the simulation is made for the red channel.The radiance computed with the three environment maps presented in section 4.1 is given for two points of the scene (figure 9).
The Perez model clearly underestimates the total light scattered by the sun, because it doe not taken into account the clouds, that are usually brighter than blue sky.This underestimation is about 20% with respect to the environment map computed by aggregation of all sky pixels.In a scene of urban canyon, the scattered radiance is the main contribution to total radiance in the shadows (e.g. point B on figure 9), with the same order of magnitude than reflected radiance (Martinoty, 2005).This would lead in our case to a maximal error with an order of magnitude of 10% on the total sensor signal.

CONCLUSION AND FUTURE WORK
In this paper, we propose an approach for retrieving the sky radiance from every direction of the upper hemisphere of an outdoor scene, using panoramic images.The radiometric issues inherent to CCD sensors are hard to modelize and correct.A method using an estimated parametric model of blue sky for computing per-image corrections is proposed and reduces slightly the visual artifacts.Anyway, the influence of the errors of the environment map are small on the total simulated signal.
The detection of sky pixels with colorimetric thresholds is not perfect; it could be enhanced by using LIDAR data when they are available with the acquired images, or simply by using a 3D model and a tree detection.Another interesting improvement could come from exposure bracketing: the overexposed areas would be reduced, and the artifacts could me more easily avoided as they depend on the exposure.

Figure 1 :
Figure 1: Example of rosy artifacts in the case of a flat-field correction without prior rejection of overexposed pixels.

Figure 2 :
Figure 2: Extracts from images from the urban dataset with visible radiometric artifacts: bloom and overexposure.

Figure 3 :
Figure 3: Image affected by a strong vertical smear, because of a worn-out sensor.

Figure 4 :
Figure 4: Images used for the radiance retrieval, from left and right cameras, for 6 vehicle locations (all images displayed with the same settings).

Figure 5 :
Figure 5: Result of sky+clouds (a) and blue sky (b) detection.

Figure 6 :
Figure 6: Original image (top) and difference between the estimated Perez model and the blue sky pixels (homogeneous white areas are those not detected as blue sky; the display settings are chosen for enhancing contrast).

Figure 8 :Figure 9 :
Figure 8: Estimated environment maps: a) Aggregation of sky pixels, b) Perez model and c) Aggregation of corrected sky pixels.The aggregation results are cut for angles close to the horizon, where there is no extracted data.d) Difference between aggregated environment maps with or without correction (different display parameters).The red deviation of artifacts is highlighted.