Automatic Mapping of the Lava Flows at Piton de la Fournaise Volcano, by Combining Thermal Data in Near and Visible Infrared

Z. Servadio1,2, N. Villeneuve1 and P. Bachelery3,4 1Laboratoire Geosciences Reunion, Universite de la Reunion, Institut de Physique du Globe de Paris, CNRS, UMR 7154, Geologie des Systemes Volcaniques, Saint Denis 2Institut de Recherche pour le Developpement, US 140, BP172, 97492 Sainte-Clotilde Cedex 3Clermont Universite, Universite Blaise Pascal, Laboratoire Magmas et Volcans, CNRS, UMR 6524, Observatoire de Physique du Globe de Clermont-Ferrand, BP 10448, F-63000 Clermont-Ferrand 4IRD, R 163, LMV, F-63038 Clermont-Ferrand France


Introduction
Knowing the eruptive history of a volcano is an essential key to the understanding of its functioning, and therefore of the evolution of the character of dangerousness of its eruptions. For an essentially effusive basaltic volcano such as the Piton de la Fournaise, the spatial and temporal distribution of lava flows allows to deduct numerous parameters of its activity, on a magmatic and a structural point of view. Satellite imaging brings more advantages than the methods used in aerial pictures studies, especially by supplying bigger temporal and spectral series. The revisiting of satellites over a region can allow the generation of dynamic mappings of the implementation of the lava flow, and also bring information on the phenomenology of the eruptions: Surface, volume, flow, spatial distribution… Furthermore, satellite images have the advantage of supplying data that grant a global visualization of the study area, and information on not easily accessible areas. The interpretation of these satellite data enables obtaining information on the surfaces and volumes of the lava field flows, but also on its nature and behavior. In a tropical environment such as La Reunion, where the climatological context presents a strong cloudiness, a satellite revisit is statistically necessary.
The optical satellite images have already been successfully used to realize mappings of lava flows. For example, in Nevado Subancaya in Peru (Legelay-Padovanie et al., 1997) or in Etna in Italy (Honda et al.,2002), the combination of spectral and morphological properties helped to elaborate surface lava flows mappings and also allowed to individualize the main The use of thermal satellite imaging to characterize the relative chronology of the implementation of lava flows was the aim of Kahle et al. (1998)'s work. Abrams et al. (1991, used it in association with optical satellite imaging to realize a chronological mapping of the lava flows in Hawaii. Other authors have realized mappings of lava flow's temperature (Hirn et al., 2005) and implementation maps by using a thermal camera (Harris et al., 2007;James et al., 2007;Lombardo et al., 2009).
At the Piton de la Fournaise, the first research that used photogrammetry to realize mappings of the outlines and thickness of lava flows took place on lava flows dated back to 1972 and 1976 (Lénat, 1987). Bonneville et al. (1989) mapped the main geological units by using SPOT1 images. After which, Despinoy (2000) realized a mapping of the lava flows above Les Grandes Pentes using a CASI hyperspectral sensor. Villeneuve (2000) realized outlines, volume calculations and a chronological follow up of the implementation thanks to stereophotogrammetry and the use of DGPS (Digital Global Positioning System). Lénat et al. (2001) associated RADAR and SPOT images to map the field of lava of the Enclos Fouqué's caldera. Recently, De Michele et Briole (2007) used a technic of correlation of images to extract lava flows which implemented between two series of aerial pictures. The study of the incoherencies in the interferograms helped realize a dynamic follow up and a mapping of the lava flows (Tinard, 2007 ;Froger et al., 2007 ;Froger et al., in press).
The aim of this article is to propose an original approach by combining thermal images with ones acquired in the visible and near infrared in order to extract independent outlines for each lava flow. The extraction is therefore independent of the operator's subjectivity. Only one automatic extraction is possible when associating thermal and optical images for the implementation of a lava flow. An automatic extraction of the outlines of lava flow is realized and then compared to a mapping of references realized by photo-interpretation. We can therefore estimate the precision of the automatic extraction.

Context
The Piton de la Fournaise (Reunion Island, Indian Ocean) is a basaltic volcano whose functioning is connected to the activity of a hot spot (Courtillot et al. 1986). The eruptions take place inside the Enclos Fouqué's caldera and have contributed to the creation of a volcanic cone whose summit is occupied by two craters: the Bory on the east side and the Dolomieu on the west side (figure 1).
During the past three decades and besides the 1986 and 1998 exceptions, the lava flows have been taking place inside de Enclos Fouqué's caldera ( figure 1 and 4) Bachèlery (1981).
intrusions situated along the rift zones (Bachèlery, 1981). The intrusions are interconnected at the same level as the central cone. Eruptions are of three kinds: summit zone eruptions, proximal eruptions, and distal eruptions (Peltier et al. 2009). It is clear that in the recent history of the Piton de La Fournaise, summit zone eruptions are the most frequent (Villeneuve et Bachèlery, 2006;Peltier et al., 2009). Most of these eruptions take place according to leveled cracks that progressively migrate by furthering from the central cone (Bachèlery, 1981).
Degassing at the event during eruptions can generate lava fountains which cause pyroclastic deposits and cones of several meters high. Two morphologies of lava flow are observed at Piton de la Fournaise: 'a'a type lava flows and pahoehoe type lava flows.
The occidental part of the Enclos Fouqué is largely recovered by a vast field of lava that Lénat et al. 2001 name the CLEF (Champ de Lave de l'Enclos Fouqué). This field of lava, essentially formed by a pahoehoe lava type flow, may have been constituted from slow emissions from the volcano's summit zone, between 1750 and 1794.
At the Piton de la Fournaise, the flow of lava emitted by eruptive cracks represent long (several meters to several kilometers), thin, (about ten meters), shallow (one meter in the slopes to 5 meters in flat zones (Letourneur et al., 2008)) lava flows, which shows the poor viscosity of the emitted magma. The juxtaposition of several individual flows during a same phase will contribute to the constitution of fields of lava, and particularly for long time eruptions (more than one month). In this case, in regards to the initial lava flow, the new income contributes essentially to the thickening and widening of the lava field flow. In this case, the lava field will be considerated here as a same unit.
The spectral properties of the lava flows differ according to the type of the surface (mainly 'a'a type at Piton de la Fournaise), but also according to the age, either because of a chemical transformation of the rock on the surface, or because of the implantation (always very fast in Reunion Island) of a vegetation cover (lichen, moss, shrub…). Lava flows which implementations were separated by several years can therefore be distinguished by their spectral properties. For the summit zone of the volcano, all the more in the Dolomieu crater, where the rocks are superimposed with only a few years or a few months of interval, the spectral properties of the diverse lava flows can then be very similar.
The oldest known eruption at the Piton de la Fournaise dates back to 1644. About 200 events have been counted since that date thanks to archives, 95% of them took place in the caldera (Lacroix, 1936;Stieltjes et al., 1989;Peltier et al. 2009). Never the less, this database is incomplete, particularly in the case of short time and low scale eruptions, before 1980. The mean average magma emission at Piton de la Fournaise estimated over a century, is 0,01 km 3 .an -1 (Lénat et Bachèlery, 1987), or 0.3 m 3 .s -1 . The debit estimations show a temporal estimation. For example, Stieltjes et al. (1989), calculate a mean debit of 0.3 m 3 .s -1 over 54 years , but obtain 0.78 m 3 .s -1 for a period of 25 years . These variations are partly due to the existence of long periods of inactivity. For example, no eruption took place during 1992 and 1998; witch is to say 6 years of inactivity. Also, another inactivity as long was observed between 1966 and 1972 (Villeneuve, 2000). Peltier et al. (2009) illustrate these debits variation and show a more important activity since 1998. Between January 1990 and January 2010, 61 eruptions have been registered with a total volume of emitted lava estimated at 473 Mm 3 (figure 2), and 33 eruptions between 1998 and 2010, with a total volume of emitted lava of 313 Mm 3 (Peltier et al. 2009, OVPF 20092010). From these observations, we have calculated a mean debit estimated between 0.45 m³.s -1 and de 0.82 m 3 .s -1 , from 1980 to 2010. These estimations are superior to those obtained by Stieltjes et al. (1989), on former periods.

Method
The originality of this research states in the use of thermal data as an analyze mask. In spite of its low spatial resolution (90m), thermal data brings essential information in our automatic mapping method. It allows determining with certitude the zone where the newly implemented lava flow is localized. The automatic extraction of the outline can be realized in this analyze mask. Also, its utilization enables treating the lava flows separately from one another because for one thermal image, only o n e l a v a f l o w i s a s s o c i a t e d i n t h i s methodology. This is particularly adapted in the case of the constitution of a lava field flow.

Optical and thermal data
The automatic extraction method has been realized by the combination of thermal and optical data. SPOT and ASTER data have been used. SPOT data have a wavelength from the visible to mid infrared, and a spatial resolution from 2,5 to 20 meters. ASTER data have a spectrum from visible to thermal infrared and a spatial resolution that varies from 15m (VNIR) to 90m (TIR).
The ASTER TIR thermal data have to be acquired at the end of the eruption or very little time afterwards, in order for the thermal anomaly to be clearly visible on the entire zone. The maximum post eruption delay of acquisition is variable and depends on the thickness of the lava flow and therefore on its speed of cooling. In most of the cases, it is less than a month. ASTER VNIR and SPOT data can be acquired long after the lava flow's implementation. The principal is not having a new lava flow implementing on the same zone. The recent lava flows present low reflectances between the visible and mid infrared wavelength. The basalt spectrum, in the visible and short wavelength of the infrared (0.4-2.4 µm), is dominated by the presence of iron, which, at different levels of crescent oxidation, increases the reflectance (Despinoy, 2000). In the same way, the presence of lichen that grows on the lava flow increases the reflectance. In near and medium infrared, the presence of chlorophyll in the vegetation induces a strong signal (Kahle et al. 1995), permitting to discriminate precise outlines in zones with vegetation cover, especially near the Grandes Pentes area as for the Piton de la Fournaise.
The KALIDEOS project from the CNES (Centre National d'Etudes Spatiales, and GEO Grid (AIST/METI), and the NASA, grant free satellite data in the case of a research program. The data used in this article are from SPOT data from the KALIDEOS program, and from the GEO Grid program for the ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) data. The climatological conditions must be optimal because the presence of clouds masks, deforms or reduces the thermal anomaly, which in a tropical zone and especially in the volcano's zone is frequent. Only seven ASTER TIR data (the ASTER satellite was released in orbit in December 1999), acquired at the end of the eruption don't present these types of issues, and have therefore been used for the treatments (Table 1 and figure 3). Optical data acquired after the eruption have been associated to the former images.

Error calculation and outlines precision of photo-interpretation
The precision of the outlines extracted by automatic method leans entirely on their comparison with a base of outlines realized from very high-resolution satellite data' photo interpretation, and aerial pictures from IGN. This base is considerate as a reference.  The outlines realized after photo interpretation, also known as reference, are not without errors. They contain errors linked to the operator's subjectivity and from the resolution of the data used to extract the outlines. One has to take into consideration these imprecisions before interpreting the results. This will be approached in the photo interpretation chapter, where outlines extracted by satellite data' photo interpretation are tested and compared on the base of those extracted from IGN's aerial pictures.
The method chosen lies on the constitution of an error matrix, because it permits to know the precision of all or part of the classification. It's expressed in percent or area.

Error matrix
The error matrix enables comparing two thermal maps. These matrixes are constructed according to a methodology developed by the Remote sensing Center of Canada (http://ccrs.nrcan.gc.ca/glossary/index_e.php?id=3124).
In our study case, an outline is compared to another that will be considered as reference (table 2).
In the matrix, numbers contained by cells correspond to areas. The total of the areas by column represent the total area of the class obtained by automatic classification. The total of the areas of each line corresponds to what has been correctly classified. The area values present in the diagonal correspond to the ones correctly classified. The other cells represent the areas which were classified wrong, either by omission, or by commission (table 2). Table 2. An error matrix for two classes (http://ccrs.nrcan.gc.ca/glossary/index_e.php?id=3124 (modifiée)).
It is therefore possible to calculate a global precision (in percent) of the classification, which is to say the area of what has been classified correctly in regards to the total area of the zone that has been classified Pg ∑ Well class area total area of the study * It is also interesting to look at the mean precision of each class of classification:

Pm
Pmi i The precision of each class is different. Some pixels can be attributed to classes that don't suit them, if the spectral properties are similar. The exactitude producer (table 2) bases itself on the area of the pixels correctly classified in regards to the area of the lava flow considered as reference. As for the exactitude user (table 2), it bases on the area of the pixel correctly classified in regard to the area of the lava flow obtained by extraction.
The exactitude producer of the "lava flow" class allows knowing the accuracy of the outline obtained by automatic extraction:

Pm lava
Well class area Real lava flow area * The mean distance separating outlines vectors is calculated from the area which was not correctly classified, divided by the smallest of the two compared lava flow outline's perimeter: Error β ∑ Not correctly class area Perimeter *

Photo interpretation
The « referenced » outlines have been realized on a period from 1980 to 2010, by vectorization on the base of IGN aerial data from 1997 and 2003. As for the eruptions post August 2003, the research was done on a base of SPOT THR, or SPOT 4 AND 5 satellite data ( Figure 4). The vectorization of the lava flows by using aerial photos is one of the most precise (Paparoditis et al. 2006), but it is also one of the most binding to realize. It could ask for several worked days. For the entire French territories, the temporal recurrence in the IGN's acquisitions is 5 years.
It is a problem when dealing with superimposed lava flows. Satellite imaging is then a complementary tool because, even if the spatial resolution of SPOT data is lower (2,5m to 20m with optical SPOT data, or inferior or equal to 1 meter with aerial photos delivered by IGN for the BDORTHOs), the temporal resolution grants data between eruptions and defines outlines for the superimposed zones. Unfortunately, the resolution of the outlines is www.intechopen.com Automatic Mapping of the Lava Flows at Piton de la Fournaise Volcano, by Combining Thermal Data in Near and Visible Infrared 209 then less precise. The SPOT satellite was put into function in 1986, the outline data base was completed by the mappings of the OVPF and the BRGM (Stieltjes et al., 1985(Stieltjes et al., et 1989Billard, 1974 ;Bussière, 1967 ;Mc Dougall, 1971).
Our methodology leans on different types of satellite data, it is important to know the influence of the spatial resolution on the extracted outlines. A comparison between the outlines obtained by photo-interpretation of different data and the referenced outlines is then done by using an error matrix and a mean distance between the outlines (table 3). According to the used satellite data, the awaited error on the outline is about the same size of the pixel (table 3). The classification's precision can tell that a mapping by satellite data' photo-interpretation is 85% more reliable for satellite data with a 10m to 20 m resolution, and 95% reliable for THR SPOT data (table 3). The other error due to the referenced outlines extraction can be the consequence of the operator's subjectivity. It presents an error from 2% to 5%. The same test was run from aerial photos, and the errors didn't exceed 2%.

Automatic process: TIR-VNIR method
The lava flows which implement in zones where spectral properties are different, such as a vegetation zone or a soil that presents high spectral reflectances, are the easiest to identify (figure 5). To automatically extract an outline, different methods can be used: classifications, threshold, or automatic detection of change (Inglada et al., 2003;Habib et al., 2007). The distance between the referenced outlines and the outlines automatically extracted is then proportional to the pixel's size, and the lava flow's classification precision is between 95% and 99%. On the other hand, lava flows which implement in low spectral reflectance zones, such as the central cone of the Piton de la Fournaise or the upperstream part of the Grandes Pentes, ask for more complex treatments. A data treatment methodology is then put together by using ENVI software ( Figure 6). The visible data Principal Composant Analysis (PCA) is applied in order to maximize the data's anti-correlation. The thermal bands, near www.intechopen.com infrared and the PCA are grouped in a multi-band data. A multi-level binary classification named decisional tree classification is then applied by using the following steps:  By using thermal data, distinguish « hot » from « cold » zones.  Extract and keep only « hot » zones.  Discretization of the bare soil and covered vegetation zones by using the Near Infra-red band.  Classify the pixels with low anti-correlation values. The seven lava flows for which ASTER data were available (table 1) were tested. For half of them, the tested eruptions are on a weak slope zone, with a substrate composed of lava flows with similar spectral properties; the other half shows various substrates and slopes (figure 7).   Table 3. Error matrices of outlines obtained by photo-interpretation, comparing those from the photo-interpretation of aerial photographs with those of the different types of satellite data.

www.intechopen.com
Automatic Mapping of the Lava Flows at Piton de la Fournaise Volcano, by Combining Thermal Data in Near and Visible Infrared

213
The classification achieved, it is then possible to export the « lava flow » class as a vector that represents the outline of the lava flow (figure 8). The extraction is realized in less than an hour, once the data collected.

Results and interpretations
In order to validate the outlines automatic extraction method, each extraction result has been compared to the referenced outline. The error matrixes (table 4) represent the tests run on the seven tested objects localized on figure 7.
The error matrixes enable calculating global errors of classifications between 77% and 96% (table 4). We can observe a disparity between the summit zone eruptions and those implemented on sloppy substrates with various spectral properties. The first ones show a global precision between 77% and 88%, whereas the others variate from 91% to 96%. This is partly due to low reflectances observed for the substrate at the summit of the volcano. For example, the outlines of a lava flow newly implemented are hardly distinguished from the intra Dolomieu lava effusion zone that presents similar ages. The mean precisions show the same disparity, with values included between 77% and 84% at the summit zone, and included between 86% and 90% for the lava flows situated on the cone's flanks and on the slopes. The area of the lava flows also play a role in the classification precision's difference,  Table 4. Error matrices for the seven outlines obtained by automatic extraction, by comparing those from the photo-interpretation of aerial photographs with those from automatic extraction.
the summit zone lava flows have less superficy than the ones observed along the slopes. An error on the outline will echo as importantly as the aria of the mapped lava flow is low. As to say, the more the mapped lava flow's area is important, the less the error will echo on the precision of the classification. The exactitude of the producer of the "lava flow" class is the most important in the validation of an automatic lava flow outline extraction method. The later varies from 74% and 87% on the entire tested lava flow. There is no observed difference between the eruptions that occur in the Dolomieu crater and the ones that occur on the flanks and the Grandes Pentes of the volcano. On another hand, the morphology of the lava flows plays a role on this precision (figure 8). The more the lava flow's superficy is important and the more compact it is, like a cercle shape form or very large, the more precise the "lava flow" class will be.
The distances between the automatic extracted lava flow outlines, in the case of superimposed lava flows, and those considered as reference is about 30m (table 4 and figure  9). Now, the satellite data used in the automatic extraction have a visible, near and medium infrared spatial resolution from 10m to 20m. If the extraction was realized without any possible confusion, we should have mean distance of about the size of these pixels. Nevertheless, in some low reflectances zones, only the thermal imaging is able to bound the outline's localization zone (figure 8 and 9). It therefore considerably reduces the precision because this data has a 90m pixel. We notice that in zones where the lava regularly effuses, the reflectance's difference cannot allow a free from error extraction of the outline, and it is therefore deducted from the thermal mask. This mask is essential because it allows the definition of a zone of interest inside witch the outline can be determinated. There is no possible confusion: to one thermal data, only one lava flow can be associated.
In vegetalized zones, the extraction can be considered as a simple threshold of the infrared, the outline is then about the size of the pixel. In some zones, especially the Dolomieu crater or certain parts of the cone, the differences of reflectances of the lava flows don't allow to differentiate them. The outline is then obtained thanks to the thermal data. It can represent up to 40% of the summit zone lava flows outline, but only 5% to 0% of the proximal and distal eruption's outlines. The mean error observed on the outline's extraction can then reach twice the size of the optical data's pixel. Some light diffusion and wedging effects can interact on the precision.  (1), a low reflectance with few or no thermal diffusion zone (2), and a various reflectance with few or no thermal diffusion zone (3). The weakest precision is for the first zone, for the second zone, the low reflectance is due to the substrate's nature and to luminosity issues in the zone, because the shadow projected by the rampart can interfere. A luminosity parameter is to be taken into consideration when choosing and time acquiring SPOT data.
As for the extraction of the third zone, the thermal infrared band essentially obtains it, because the lava flow implemented in a vegetalized zone with a high spectral signature difference compared to the lava flow. It is the zone where the extraction is the most precise because only based on the SPOT data.

Discussion and conclusion
The vectors of the outlines obtained by automatic extraction properly match with the referenced outlines since we obtain a mean exactitude producer of 80% for the "lava flow" class. If we compare our automatic extraction results to those obtained by DGPS, they are less precise. The DGPS has a precision of about one centimeter at its antenna, now it's the operator who transports it. The error is due to the positioning of the antenna regards to the outlines. The error then is about one meter, for pluri-meter in our extractions. Nevertheless, the effusion zones are not all accessible, and cannot let realizing the outlines in their integrality, especially in the Grandes Pentes. Computer assisted drawing methods are on the other hand less reliable than our methodology for the data' distortion are not taken into account and can provoke hundred of meters error. The results of our classification are close to those obtained in other contexts' literature. (Azerzaq et al., 1997 ;Messar et Messar, 1997 ;Yüksel et al., 2008).
The vectorization time by photo-interpretation can take several days, whereas the automatic method enables obtaining this vector in less than one hour. There is therefore a considerable gain of time. The gaps of observations due to clouds, or zone with high thermal diffusion and low reflectances, could be filled by using RADAR imaging in treatment sequences (Weisseil et al.2004). However, by adding Bi data will increase the treatment time, especially by using coherence data.
The association of thermal and optical data has already been realized in other automatic classifications with outlines extractions contexts: glaciology (Raciviteanu et al. 2008), the canopy (Joshi et al., 2006), agriculture (Kasdan, 1979;Saito et al., 2001). For similar spectral resolution data, the error matrix results are comparable to those obtained by our methodology.
This methodology was developed in order to automate lava flow outline's extraction and therefore ensure a fast update of the Piton de la Fournaise's database. The lava flow map was updated thanks to photo-interpretation and automatic extraction (figure 4). It allowed us to test the reliability of the outlines extracted according to each methodology used in this article, and to know their precision.
The errors measured by the matrixes give us the extracted surface's error, by comparing the automatic lava flow's area and the referenced one. We saw that: 1/ the difference of interpretation between two operators for the same data can be of 2% for the aerial photos and it varies between 2% and 5% with satellite data. 2/if the outlines extracted from satellite data's photo interpretation are compared to those extracted from aerial data, the exactitude producer of the "lava flow" class can reduce to 80%, but is at an average of 90%. 3/ Those obtained with the same type of data varies from 74% to 87%, which is to say an average of less than 10%. It influences the obtained area in less than 10%, exception made for low scale summit zone eruptions which represent an error of 20%, which remains modest and without major consequences on the volume and production estimation rates for the volcano.

Acknowledgements
These reasearch were financed by the "Region Reunion" and the ministry of "l'Enseignement Supérieur et de la Recherche" program. The authors would like to thank the CNES (Centre National d'Etudes Spatial) for the free access to the SPOT data via the KALIDEOS program (http://kalideos.cnes.fr).
This piece of work also used free ASTER data obtained thanks to AIST GEO Grid (Ministry of Economy, Trade and Industry) programm.
Thanks to Rebecca Roger for the translation.

Note
All the lava flow outlines' vectors will be available online on the Laboratoire de Geosciences Reunion website.