AUTOMATIC CLOUD DETECTION FROM MULTI-TEMPORAL SATELLITE IMAGES : TOWARDS THE USE OF PLÉIADES TIME SERIES

Contrary to aerial images, satellite images are often affected by the presence of clouds. Identifying and removing these clouds is one of the primary steps to perform when processing satellite images, as they may alter subsequent procedures such as atmospheric corrections, DSM production or land cover classification. The main goal of this paper is to present the cloud detection approach, developed at the French Mapping agency. Our approach is based on the availability of multi-temporal satellite images (i.e. time series that generally contain between 5 and 10 images) and is based on a region-growing procedure. Seeds (corresponding to clouds) are firstly extracted through a pixel-to-pixel comparison between the images contained in time series (the presence of a cloud is here assumed to be related to a high variation of reflectance between two images). Clouds are then delineated finely using a dedicated region-growing algorithm. The method, originally designed for panchromatic SPOT5-HRS images, is tested in this paper using time series with 9 multi-temporal satellite images. Our preliminary experiments show the good performances of our method. In a near future, the method will be applied to Pléiades images, acquired during the in-flight commissioning phase of the satellite (launched at the end of 2011). In that context, this is a particular goal of this paper to show to which extent and in which way our method can be adapted to this kind of imagery.


INTRODUCTION
Detecting clouds is one of the primary steps to perform when processing satellite images.In an industrial context, more particularly in the context of the French mapping agency (IGN-F), clouds may affect production lines in three different manners:

•
During the acquisition of satellite data.At the wavelengths of the visible light, clouds are opaque and also hide the ground surface from Earth observation satellites.As a work-around, it is possible to acquire a scene in several passes: that results in time series containing a high number of satellite images (in general between 5 and 10) acquired at different epochs.

•
During the computation of DSM.The presence of clouds in satellite stereoviews is a classical cause of failure of stereo-matching procedures (Grün, 2000) (Eckert et al., 2005) • During the computation of large cloud-free orthomosaics from input satellite images.The presence of clouds in input satellite images has here several consequences.As mentioned in the first item, that implies the acquisition of time series i.e. the acquisition of multitemporal satellite images.Therefore, computing large seamless orthomosaics from these multi-temporal images involves identifying cloud pixels in input data so that only cloud-free pixels are used for generating the final product (Din-Chang et al., 2008).In addition, due to the fact that images are not acquired simultaneously, there are radiometric heterogeneities between them.These differences may be related to variations in the Bidirectional Reflectance Distribution Function (BRDF) of ground surfaces or to variations in the solar illumination.These differences must also be considered when generating orthomosaics through specific radiometric equalization procedures (Chandelier and Martinoty, 2009) (Falala et al., 2008).Eventually, the cloud shadows may cause some problems in the procedure, more particularly during the radiometric equalization.In particular, they may cause a "leopard skin" effect in the final orthomosaics (Soille, 2008).As a consequence, these shadows must also be detected for a subsequent de-shadowing procedure that aims at enhancing the display of the corresponding areas in the final orthomosaics (Simpson and Stitt, 1998) (Richter and Müller, 2005).Guyane, Martinique, Guadeloupe) where the cloud cover is known to be important all year round.To limit this problem, providers of satellite data generally wait for a favourable weather window so that the cloud cover is minimal.However, this solution is not operable in areas (e.g.Guyane) that are almost always nebulous.Here, the solution, already mentioned above, consists in acquiring images regardless the weather conditions.That results in obtaining time series i.e. a pile of satellite images with a high number of images, as depicted in Figure 2 for our Hispaniola test area.If this solution involves difficult pos-procedures for generating a virtual cloud-free satellite image from those contained in time series, it also presents the advantage of maximizing the acquisition capability of the satellite.This is the main reason why this acquisition configuration is considered by an increasing number of satellite data providers.In our project, we assumed to have such multitemporal satellite images and we built our cloud detection system upon this hypothesis.Until now, the analysis of satellite images for cloud identification and removal has been carried out manually.That entails a long, fastidious and costly work.As a consequence, it becomes necessary to introduce a certain degree of automation in the cloud detection procedure i.e. to develop expert systems that are able to focus the operator's attention on potential cloud areas.

Related works
In the last few years, many researches have been carried out to automate the cloud detection in satellite images.Some works can be found in the meteorology community, for example in the paper by Derrien and Le Gléau (2010) who propose to combine a temporal analysis to a region-growing technique in order to improve the detection of low clouds at twilight.If this work appears to be interesting, the meteorology context appears to be very different from the context generally considered in the geoscience community.Thus, the images used in the above-mentioned project are delivered by the SEVIRI system that offers 11 spectral bands (including thermal bands), a high temporal resolution (one image every 15 minutes) but a low spatial resolution (3km at nadir).By contrast, the images used in the geoscience community generally offer fewer spectral channels but a better spatial resolution, with a Ground Sample Distance (GSD) of e.g.5m for SPOT5-HRS, 70cm for Pléiades and 50cm for Worldview.Most of the cloud detection approaches found in the geoscience literature are dedicated to single-date satellite images (Irish, 2000) (Irish et al., 2006) (Le Hégarat andAndré, 2009).By contrast, only a few of them make a full use of multi-temporal satellite images.In this category, the method by Sedano et al. (2011), originally designed for single-date High-Resolution images, uses information from multi-temporal Low-Resolution MODIS images in order to perform the cloud detection.Thus, the difference between MODIS and input images is analysed to extract the seed points that correspond to clouds, which are then used in a subsequent region-growing procedure in order to delineate clouds finely.The main drawback of the method is related to the low resolution of MODIS images (featured by a GSD of 250m) that makes the detection of small clouds impossible.Multi-temporal LANDSAT and FORMOSAT-2 images are used in (Hagolle et al., 2010).The extraction of cloud pixels is here based on the assumption that a sudden increase of reflectance in the blue band between two images is due to the presence of a cloud.This leads to the computation of a preliminary cloud mask that is then refined using additional radiometric tests, performed on the remaining bands.If this method is shown to give good results, it is limited to the availability of several bands and can also not be used in every context, especially when using panchromatic SPOT5-HRS images.

Presentation
The main goal of this paper is to present the method developed at the French Mapping Agency for identifying and delineating clouds in High Resolution multi-temporal satellite images.We here propose to couple the ideas of the papers by Sedano et al. (2011) and by Hagolle et al. (2010) to build our cloud detection system.Following the design of the method by Sedano et al. (2011), our method is composed of two steps: seeds are firstly extracted and secondly extended during a region-growing procedure in order to delineate clouds finely.However, our method differs from this latter method in the way to extract seeds.We here follow the assumption made by Hagolle et al.
(2010) and we consider that a high variation of reflectance between two images (in time series) is related to the presence of a cloud in the scene.
Our method was originally designed for SPOT5-HRS images, mainly because these images are still used in our production lines, in particular for the production of the Reference3D® database (Bouillon et al., 2006)  The production of large area seamless orthomosaics (Falala et al., 2008) •

Cloud detection
In that context, this is also a particular goal of this paper to show to which extent and in which way the method presented here can be adapted to Pléiades images.
The remaining of the paper is organized as follows.Section 2 presents input data.Section 3 presents our cloud detection method.Section 4 shows the preliminary results that we obtained, starting from SPOT5-HRS images.Eventually, Section 5 highlights the perspectives related to the use of Pléiades time series for detecting clouds.

INPUT DATA
Two different kinds of input data are also considered in our project: panchromatic SPOT5-HRS images and (in a near future) Pléiades images.
SPOT5-HRS data.As shown in Figure 3 and as detailed in (Bouillon et al., 2006), the SPOT5-HRS instrument is composed of two telescopes with a viewing angle (along the track) of 20° forward and 20° aft.This configuration allows the acquisition of single pass stereopairs, with a time delay between 2 images of 90 seconds and a corresponding Base-to-Height (B/H) ratio of about 0.8.The system is featured by a swath of 120 km, a pixel GSD of 5 m along the track and 10 m across the track.For our 1 http://smsc.cnes.fr/PLEIADES/index.htm.Last visited: 2012/4/16 project, multi-temporal panchromatic images were acquired.That resulted in the availability of time series i.e. a pile of images, acquired at different epochs, with an average overlapping of 6 images, as depicted in Figure 2.
Figure 3. Acquisition configuration of the SPOT5-HRS instrument (Bouillon et al., 2006) Pléiades data.Contrary to SPOT5-HRS images, Pléiades images are not only panchromatic but have 4 channels: Red, Green, Blue and Near Infrared (NIR).In addition, they are featured by a ground pixel of 70cm.In our project, we also assume to have time series, in a similar way as for SPOT5-HRS data.

METHOD
As introduced in Section 1.2 and as justified and detailed in this section, our method is composed of 3 steps: A region growing step

Pre-processing
The pre-processing step aims at making the subsequent analysis/comparison easier.It consists on the one hand in a geometric correction and on the other hand in a radiometric correction.
Geometric pre-processing.Because our method is based on a pixel-to-pixel comparison between the several satellite images contained in the time series (See Section 3.2 for more details), input satellite images must be co-registered.In our project, input satellite images are orthorectified using an in-house DTM (Reference3D®).Note that we only use (in the rest of the algorithm) the orthophotos generated at this step.It should be noted here that the method did not appear to be sensitive to the accuracy of the DTM used for orthorectification.Thus, the outcomes produced with orthophotos computed with GDEM (Global Digital Elevation Model, derived from ASTER images through stereo-matching algorithms) are not significantly different from the outcomes based on the Reference3D® or SRTM (Shuttle Radar Topography Mission) DTMs, even though these two latter DTMs are known to have a better accuracy (in altimetry) than GDEM.
Radiometric pre-processing.In addition to this geometric correction and because the comparison (detailed in Section 3.2) is made on a radiometrical basis, input data must be radiometrically corrected.For that purpose, we followed the recommendations found in (Lillesand et al., 2008) In the two following sections, we propose to describe the two steps involved in our cloud detection approach: the seed extraction step (Section 3.2) and the region-growing step (Section 3.3).

Seed extraction
This first step also aims at achieving the following goal: To extract cloud seeds in one orthophoto, we proceed in the following way.Firstly, we take one orthophoto as a reference (this Reference Orthophoto RO corresponds to the orthophoto we want to detect clouds in).Secondly, we select the other orthophotos of the pile that intersect this Reference Orthophoto RO and we crop them according to their overlapping on the Reference Orthophoto (Note that these cropped orthophotos are called Pile Orthophotos PO in the rest of this paper).From there, a pixel-to-pixel analysis is carried out between the Reference Orthophoto RO and each PO (Figure 4).If the radiometric difference is better than a given threshold T R (in our project, 8 bits images are used and T R =170), the pixel of the Reference Orthophoto receives a positive vote.Eventually (after this process is carried out for each PO), a pixel of RO is labelled "cloud" if a 2/3 majority votes for it.
At the end of this step, pixels corresponding to clouds are also identified.However, some cloud pixels are missing.This is the case for example for pixels covered with a cloud in all the orthophotos of the pile (in that case, the pixel-to-pixel difference is never higher than the T R threshold).This issue is solved by operating a region-growing procedure, starting from these seed points, as detailed in the next section.

Region growing
This step also aims at delineating clouds finely.It is based on the two following hypotheses, found e.g. in (Le Hégarat and André, 2009): • Clouds are brighter than the underlying landscape • Clouds are connex objects The region-growing algorithm used in our system is based on two criteria: a radiometric criterion and a homogeneity criterion.
To be aggregated to a cloud region, the radiometry of the pixel must also lie within a range of given values (in our project, 8 bits images are used and range = [T RG ; 255] with T RG =200).In addition, all the neighbours of this pixel must satisfy this first criterion.Note that the size of the neighbourhood structuring element S N can be specified to the algorithm and is set by default to S N =1.Some results, obtained with panchromatic satellite images, are presented in Figures 5 and 6.

RESULTS AND DISCUSSION
The method was tested in the Hispaniola (including the Dominican Republic and Haïti) test site that contains 52 SPOT5-HRS images, with an average overlapping of 6 images.
In this paper, we present the preliminary outcomes delivered by our method for a pile of 9 images.These results were evaluated visually and are discussed in this section.
As can be seen in Figures 5 and 6, we observe that most of clouds are correctly detected by our method.In particular, the method is effective in detecting cumulus clouds that appear to be fleecy in the scene.This is due to the fact that these clouds have sharp outlines and are brighter than the surrounding (ground) objects.In that sense, they perfectly fit the two assumptions used in the region-growing step: this explains why our procedure is so successful.By contrast, in a similar way as the other methods found in literature (Sedano et al., 2011), the detection of mists appears to be more difficult.This is due to the fact that these clouds are generally thinner and less bright than e.g.cumulus clouds.
Then, the spectral variation of the corresponding pixels is not significantly high (it lies in the normal range) and the first assumption used for extracting seeds is not verified.
In addition, we notice that confusions appear with some bright objects like white sand beaches.This is related to the fact that the corresponding regions in input SPOT5-HRS images may be saturated because of the poor radiometric resolution of input images.As a consequence, these regions are wrongly (and systematically) alerted as clouds at the end of the second step.
To remove them, a spatial analysis is necessary.It consists in comparing the number of pixels of the resulting region (computed at the end of the region-growing procedure) and the number of initial seeds.Our experiments show that, in the case of a saturated region, this ratio is high.By contrast, in the case of a cloud, this ratio is always low.This criterion is also used (with a threshold T S =6, applied to this ratio) to remove misclassified regions from the outcomes.

PERSPECTIVES
In a near future, we plan to test our cloud detection method with multi-spectral Pléiades images.For that purpose, and because our method is designed for images with one band, we will use the panchromatic channel of these images.The major drawbacks -encountered and presented in the previous sectionshould not appear in these new experiments.Thus, regarding the errors related to saturated regions in input data, they should be discarded because 12-bits Pléiades images are shown to have a (much) higher dynamic range than the SPOT5-HRS images.In  addition, regarding the false negative detection related to mists, the information given by the NIR channel [750-950 nm] of Pléiades images might be very useful.Some tests are being made in that direction.
Eventually, the detection of shadows, already mentioned in the introduction and that appears to be difficult to carry out with panchromatic SPOT5-HRS images, becomes possible with multi-spectral Pléiades images.In that context, a cloud shadow method is being considered.In a similar way as the cloud detection method presented in this paper, it is based on a regiongrowing algorithm.It however differs in the way to extract seeds.This step is here carried out by selecting the pixels of the Reference Orthophoto RO that appear to be "darker" than the corresponding pixels of the other orthophotos PO contained in the pile.We moreover use the information given by the NIR channel that appears to better discriminate shadow pixels (Richter and Mueller, 2005).The detection of shadows is all the more important as it enables -in a similar way as (Le Hégarat-Mascle and André, 2009) -a joint detection of clouds and corresponding shadows.In our framework, clouds and shadows will be detected separately.We will then compare their respective size and shape (given the incidence angle of the sun, the local terrain model and an approximate altitude for clouds), which should lead to a mutual validation (or invalidation) of their detection.

CONCLUSION
In this paper, we have presented a method for the automatic identification and delineation of clouds in High-Resolution satellite images.Our method, originally designed for panchromatic SPOT5-HRS multi-temporal images, is based on a region-growing algorithm.It is easy to implement and requires only 4 parameters (T R , T RG , S N , T S ) that can be easily set to default.As an advantage to other methods found in literature, our method does not require thermal bands and works on 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 images with one panchromatic channel only.In that sense, the method appears to be sensor-independent and can be directly applied to e.g.multispectral Pléiades images, by using the panchromatic channel of these images.It should be noted that our method requires the availability of a pile of several images, acquired at different epochs.If this constraint prima facie appears to be strong, it actually corresponds to a normal configuration for satellite data acquisition, especially in areas where the cloud cover is important all year round.Our experiments were carried out in a difficult context i.e. with images containing a lot of radiometrically saturated pixels.That caused large false positive cloud detections and necessitated developing a specific post-procedure.These drawbacks will not appear with Pléiades images, as these images have a higher dynamic range.In addition, the availability of a NIR channel will allow us to build more robust procedures including a joint detection of clouds and shadows.

Figure 1 .
Figure 1.The first images acquired by the Pléiades system highlight one of the primary problems that satellite acquisition generally copes with: the presence of clouds in the scene.A view over the historical part of Paris (2011/12/22).© CNES

•
2D change detection(Champion et al., 2010) (Le- Bris andChehata, 2011) for updating building databases•The update of Land Cover / Land Use databases(Hermosilla et al., 2011) • 3D change detection by comparison of Digital Surface Models, computed from satellite images acquired at two different dates(Guérin et al., 2012) •The 3D reconstruction of buildings(Durupt and Taillandier, 2006) (Lafarge et al., 2008) • Figure 4. Seed Extraction.Illustration of the pixel-to-pixel comparison between the Reference Orthophoto RO and the other orthophotos PO contained in the pile (time series)

International
Figure 5. (top) Initial panchromatic SPOT5-HRS orthophoto.(bottom) The outcomes of our cloud detection algorithm.Pixels labelled "clouds" are highlighted in blue.

Figure 6 .
Figure 6.Zoom on the area highlighted in red in Figure 5. (top) Initial images.(bottom) Clouds detected by the method are highlighted in blue.
. However, with the advent of High-Resolution satellite images such as Geoeye, Worldview and Pléiades, it becomes necessary to adapt (at least, to test) the existing procedures to these new sensors.With this objective in mind, IGN-F is taking part in 2012 in the user's thematic commissioning that is conducted by the French Space Agency (CNES) and that aims at validating the future products and services based on Pléiades.More particularly, we plan to use Pléiades images to replay some studies, already carried out during the ORFEO 1 accompaniment program.These experiments will concern: