The 2014 – 2015 Lava Flow Field at Holuhraun , Iceland : Using Airborne Hyperspectral Remote Sensing for Discriminating the Lava Surface

The Holuhraun lava flow was the largest effusive eruption in Iceland for 230 years, with an estimated lava bulk volume of ~1.44 km3 and covering an area of ~84 km2. The six month long eruption at Holuhraun 2014–2015 generated a diverse surface environment. Therefore, the abundant data of airborne hyperspectral imagery above the lava field, calls for the use of time-efficient and accurate methods to unravel them. The hyperspectral data acquisition was acquired five months after the eruption finished, using an airborne FENIX-Hyperspectral sensor that was operated by the Natural Environment Research Council Airborne Research Facility (NERC-ARF). The data were atmospherically corrected using the Quick Atmospheric Correction (QUAC) algorithm. Here we used the Sequential Maximum Angle Convex Cone (SMACC) method to find spectral endmembers and their abundances throughout the airborne hyperspectral image. In total we estimated 15 endmembers, and we grouped these endmembers into six groups; (1) basalt; (2) hot material; (3) oxidized surface; (4) sulfate mineral; (5) water; and (6) noise. These groups were based on the similar shape of the endmembers; however, the amplitude varies due to illumination conditions, spectral variability, and topography. We, thus, obtained the respective abundances from each endmember group using fully constrained linear spectral mixture analysis (LSMA). The methods offer an optimum and a fast selection for volcanic products segregation. However, ground truth spectra are needed for further analysis.


Introduction
Lava flow emplacement is an important constructive geological process that contributes to reshaping natural landscapes [1][2][3].To assess the hazards and long-term impacts posed by lava flows, it is vital to understand aspects such as the return period of effusive eruptions, to map the areas covered by eruptions in the past and to characterize the evolution of lava flow surfaces after emplacement [4,5].In high eruption frequency areas, lava flows often overlap each other.If the overlapping lava flows erupt within a short time span and have similar chemical and surface characteristics, discrimination will be further complicated by their similar spectral signatures.Spectral reflectance plays an important role in visible and shortwave infrared (VIS-SWIR) remote sensing.Each material absorbs and reflects the incoming radiation in a characteristic way.In the 400-2500 nm range, minerals display absorption features due to the interaction of light with cations (Fe, Mg, Al) and anions (OH, CO 3 ) [6].Reflectance spectra provide information about the specific material and their composition.They are used for different applications such as classification of remotely sensed data, identification of mineral features of rock, and environmental assessment [7,8].The interest in reflectance spectra of volcanic rocks has increased recently as they can play an important role as planetary analogues.In fact, these spectra can be used to identify compounds by data acquired by ongoing solar system exploration missions [9,10].
Characterization of surface spectral reflectance by satellite remote sensing is constrained by the spectral range and resolution (i.e., number of spectral bands) as well as by the spatial resolution of the imagery.Whereas multispectral imagery can be acquired at very high spatial resolution (e.g., WorldView [11,12]); the spatial resolution of hyperspectral satellite data remains low (e.g., EO-1 Hyperion with a ground resolution of 30 m x 30 m); and spectral mixing is thus a major issue [13].The spectral reflectance of lava of different compositions has also been documented using laboratory spectrometry with decimeter-size samples [14].For accessible volcanic terrains, field spectrometry offers a useful alternative approach for characterizing the spectral reflectance of contrasted lava surfaces and for documenting its spatial variation at different spatial scales [5,14].The great variety of morphologies observed in the 2014-2015 Holuhraun lava flows [1,15] encouraged a detailed study of their spectral characteristics, to obtain information about lava composition and detect possible differences in the spectra of the flow.In spectroscopy, the identification of the mineral constituents of major rock types is typically approached using spectral unmixing methods [5,16].Usually, in the visible and near-infrared spectral range, mafic rocks are characterized by very low reflectance due to the presence of large amounts of dark mafic minerals [14].The 2014-2015 lava flow at Holuhraun in NE Iceland offers an excellent diverse surface environment for investigating and characterizing lava deposits.Its intense volcanic activity [1,[17][18][19], geomorphological complexity [20], and well-documented flank eruptions [1] perplex the remote sensing monitoring of the bulk volcanic edifice.However, the detailed field mapping of lithologies is frequently obstructed by difficulties in accessibility, the scale of lava flow fields, topography, while remote sensing has become increasingly important in mapping volcanic terrains and specifically in mapping lava flows.Mapping individual lava flows using satellite remote sensing is challenging for at least three reasons: vegetation cover, spatial overlapping, and spectral similarity [3,4].Moreover, a high eruption frequency often leads to lava flows overlapping each other.If the overlapping lava flows are erupted within a short period and have similar chemical and surface characteristics, discrimination will be further complicated by their similar spectral signatures.
Hyperspectral remote sensing provides information on hundreds of distinct and contiguous channels of the electromagnetic spectrum, thus enabling the identification of multiple ground objects through their detailed spectral profiles.However, restrictions on the spatial resolution of hyperspectral data, the multiple scattering of the incident light between objects, and microscopic material mixing form the mixed pixel problem.Pixels are identified as mixed when they are composed of the spectral signatures of more than one ground object.Therefore, we adopted linear spectral mixture analysis (LSMA) techniques [8,21], which model the pixel spectra as a combination of pure components (endmembers) weighted by the fractions (abundances) that contribute to the total reflectance of the mixed pixel [22].Ideally, each selected endmember from the hyperspectral image under study has the maximum possible abundance of a single physical material present and minimum abundance of the rest of the physical materials.Spectral unmixing typically consists of two main substages: (a) endmember extraction; and (b) abundance estimation [22].In this paper, we focus on both endmember extraction and estimation of fractional abundances of the lava field products on 2014-2015 Holuhraun lava fields.For this purpose, an airborne hyperspectral image with an AisaFENIX sensor on board a NERC Airborne Research Facility (Natural Environment Research Council Airborne Research Facility) campaign was acquired at Holuhraun after the eruption and for the sub-pixel analysis we used the sequential maximum angle convex cone (SMACC) algorithm to identify the spectral image endmembers while the LSMA method was employed to retrieve the abundances.Our approach was narrowed to the eruptive fissure vent part since it is considered to have a more diverse surface.The resulting abundances from the LSMA method were both quantitatively and qualitatively compared with the spectral indices technique, aerial and field photographs, respectively.The objective was to retrieve the main lava surface type contributing to the signal recorded by airborne hyperspectral at the very top surface of Holuhraun.

The 2014-2015 Eruption at Holuhraun
The eruption took place in the tectonic fissure swarm between the Bárðarbunga-Veiðivötn and the Askja volcanic systems (Figure 1a).It lasted about six months (31 August 2014 to 27 February 2015) and produced a bulk volume ~1.44 km 3 of basaltic lava [1].Lava effusion rates during the eruption period range from 320 to 10 m 3 /s.Averaged values are ∼250, 100, and 50 m 3 /s during the initial (August-September 2014), intermediate (October-December 2014) and final phase (December 2014 to February 2015), respectively [1,17] (Figure 1b).The lava was emplaced on the sandur plains (glacial outwash sediment plains) north of the Vatnajökull/Dyngjujökull glacier, partially covering the previous two Holuhraun lava flow fields south of the Askja caldera [1].The area is gently sloping (average inclination <0.5%; i.e., ∼0.3 • ) to the east-northeast.The shallow gradient resulted in low topographic forcing of the flow and, therefore, rather slow lava flow advance.During its emplacement history, the lava field was initially dominated by channels and horizontal expansion.Then it transitioned to grow in volume primarily by inflation, tube-fed flow (i.e., transport of lava through roofed over partially or filled channels) and vertical stacking of lava-lobes.The 2014-2015 effusive eruption products originate from intense activity in the vent, in which high oxidation occurs in this area.The main lava channel shows significant inflation (5-10 m).Lava advancement rates were generally low ∼0.0167 m/s during the initial eruption phase [1] and dropped to ∼0.0017 m/s during the middle of November 2014 [23].The six-month-long effusive eruption features diverse surface structures and morphologies.The 2014-2015 lava flow at Holuhraun in NE Iceland offers an excellent diverse surface environment to investigate and characterize lava deposits.
Remote Sens. 2018, 10, x FOR PEER REVIEW 3 of 20 photographs, respectively.The objective was to retrieve the main lava surface type contributing to the signal recorded by airborne hyperspectral at the very top surface of Holuhraun.

The 2014-2015 Eruption at Holuhraun
The eruption took place in the tectonic fissure swarm between the Bárðarbunga-Veiðivötn and the Askja volcanic systems (Figure 1a).It lasted about six months (31 August 2014 to 27 February 2015) and produced a bulk volume ~1.44 km 3 of basaltic lava [1].Lava effusion rates during the eruption period range from 320 to 10 m 3 /s.Averaged values are ∼250, 100, and 50 m 3 /s during the initial (August-September 2014), intermediate (October-December 2014) and final phase (December 2014 to February 2015), respectively [1,17] (Figure 1b).The lava was emplaced on the sandur plains (glacial outwash sediment plains) north of the Vatnajökull/Dyngjujökull glacier, partially covering the previous two Holuhraun lava flow fields south of the Askja caldera [1].The area is gently sloping (average inclination <0.5%; i.e., ∼0.3°) to the east-northeast.The shallow gradient resulted in low topographic forcing of the flow and, therefore, rather slow lava flow advance.During its emplacement history, the lava field was initially dominated by channels and horizontal expansion.Then it transitioned to grow in volume primarily by inflation, tube-fed flow (i.e., transport of lava through roofed over partially or filled channels) and vertical stacking of lava-lobes.The 2014-2015 effusive eruption products originate from intense activity in the vent, in which high oxidation occurs in this area.The main lava channel shows significant inflation (5-10 m).Lava advancement rates were generally low ∼0.0167 m/s during the initial eruption phase [1] and dropped to ∼0.0017 m/s during the middle of November 2014 [23].The six-month-long effusive eruption features diverse surface structures and morphologies.The 2014-2015 lava flow at Holuhraun in NE Iceland offers an excellent diverse surface environment to investigate and characterize lava deposits.

Spectral Unmixing on Lava
Various spectroscopy studies [2,5,7,14,25] over the volcanic area have examined the mineralogical composition of the extensive lava fields.Usually, in the visible (VIS) and near-infrared (NIR) spectral range, mafic rocks are characterized by very low reflectance due to the presence of large amounts of dark mafic minerals [14].Spectral indices provide the first efficient way to emphasize subtle spectral variations at the surface [26].More elaborate methods have been developed to discriminate and quantify mixtures of mafic minerals.They have been used to derive composition maps of mafic minerals [27][28][29].However, some lava flows can have a similar chemical/mineralogical composition but dissimilar spectral behaviour due to the different grain size, surface texture, and presence of weathering [13,14].The main components of igneous rocks do not display any peculiar spectral features in the visible and near infrared spectral range.In the case of basalts, the only spectral feature commonly found is an absorption peak, due to iron, located around 1000 nm [26].However, in the case of hydrothermal alteration, hydroxyl bearing minerals show distinctive absorption features in the 2000-2500 nm spectral region [30].Because of the heterogeneity of the lava surface, mixed pixels are very common which is illustrated in Figure 2a,b [24], (b) coverage of the three main phases after Pedersen et al [1].

Spectral Unmixing on Lava
Various spectroscopy studies [2,5,7,14,25] over the volcanic area have examined the mineralogical composition of the extensive lava fields.Usually, in the visible (VIS) and near-infrared (NIR) spectral range, mafic rocks are characterized by very low reflectance due to the presence of large amounts of dark mafic minerals [14].Spectral indices provide the first efficient way to emphasize subtle spectral variations at the surface [26].More elaborate methods have been developed to discriminate and quantify mixtures of mafic minerals.They have been used to derive composition maps of mafic minerals [27][28][29].However, some lava flows can have a similar chemical/mineralogical composition but dissimilar spectral behaviour due to the different grain size, surface texture, and presence of weathering [13,14].The main components of igneous rocks do not display any peculiar spectral features in the visible and near infrared spectral range.In the case of basalts, the only spectral feature commonly found is an absorption peak, due to iron, located around 1000 nm [26].However, in the case of hydrothermal alteration, hydroxyl bearing minerals show distinctive absorption features in the 2000-2500 nm spectral region [30].Because of the heterogeneity of the lava surface, mixed pixels are very common which is illustrated in Figure 2a and Figure 2b.Spectral Mixing Analysis (SMA) has been specifically developed to account for mixtures [10].Analysis of the data sample can simply be performed on these abundance fractions rather than the sample itself.This method is well-suited for spectroscopic analysis because most of the spectral shapes are due to different materials.The signal detected by a sensor at a single pixel is frequently a combination of numerous disparate signals.Unmixing techniques were applied to the volcano of Nyamuragira for discriminating lava flows of different ages by Li et al. [5].The most recent study by Daskalopoulou et al. [16], used unmixing techniques to segregate lava flows and related products from the historical Mt.Etna.Nonetheless, there are no findings concerning lava flow delineation through unmixing in Iceland.

Airborne Hyperspectral Data Acquisitions
Airborne hyperspectral data were acquired on 4 September 2015 between 16.56 and 17.58 (local time) with an AisaFENIX sensor (Specim, Spectral Imaging Ltd, http://www.specim.fi)[31]   Spectral Mixing Analysis (SMA) has been specifically developed to account for mixtures [10].Analysis of the data sample can simply be performed on these abundance fractions rather than the sample itself.This method is well-suited for spectroscopic analysis because most of the spectral shapes are due to different materials.The signal detected by a sensor at a single pixel is frequently a combination of numerous disparate signals.Unmixing techniques were applied to the volcano of Nyamuragira for discriminating lava flows of different ages by Li et al. [5].The most recent study by Daskalopoulou et al. [16], used unmixing techniques to segregate lava flows and related products from the historical Mt.Etna.Nonetheless, there are no findings concerning lava flow delineation through unmixing in Iceland.

Airborne Hyperspectral Data Acquisitions
Airborne hyperspectral data were acquired on 4 September 2015 between 16.56 and 17.58 (local time) with an AisaFENIX sensor (Specim, Spectral Imaging Ltd, http://www.specim.fi)[31] on board a NERC Airborne Research Facility (Natural Environment Research Council Airborne Research Facility http://www.bas.ac.uk/nerc-arf) aircraft [32].Pushbroom VNIR and SWIR sensor, are two separate detectors with common fore-optics.The hyperspectral data contain 622 channels with spectral range from ~400 nm to 2500 nm (break at ~970 nm).The pixel size of this data is explained in Section 4.2.2.
In total, eight flights were acquired at the Holuhraun lava flow during this period with an average altitude of 2.4 km (Figure 3a).The data are delivered as level 1b ENVI BIL format files which means that radiometric calibration algorithms have been applied and navigation information has been synced to the image data (Figure 3b).In this study, we subset the data to focus on the area around the eruptive fissures vent (Figure 3c) which is thought to have a diverse surface and has field photographs.Very high-resolution aerial photographs of the lava field (0.5 m spatial resolution) from Loftmyndir ehf (http://www.loftmyndir.is/)[33] were used for comparison and validation of the unmixing results.
Remote Sens. 2018, 10, x FOR PEER REVIEW 5 of 20 separate detectors with common fore-optics.The hyperspectral data contain 622 channels with spectral range from ~400 nm to 2500 nm (break at ~970 nm).The pixel size of this data is explained in section 4.2.2.In total, eight flights were acquired at the Holuhraun lava flow during this period with an average altitude of 2.4 km (Figure 3a).The data are delivered as level 1b ENVI BIL format files which means that radiometric calibration algorithms have been applied and navigation information has been synced to the image data (Figure 3b).In this study, we subset the data to focus on the area around the eruptive fissures vent (Figure 3c) which is thought to have a diverse surface and has field photographs.Very high-resolution aerial photographs of the lava field (0.5 m spatial resolution) from Loftmyndir ehf (http://www.loftmyndir.is/)[33] were used for comparison and validation of the unmixing results.

Spectral Unmixing and Abundance Retrieval
The processing workflow towards unmixing and generating abundance consists of four steps: (1) Atmospheric correction to retrieve surface reflectance; (2) Data masking, geocorrection, reprojection, and resampling; (3) An endmember selection algorithm was adopted to select the endmembers; then a linear spectral mixing analysis method was employed to retrieve the abundance (Figure 4).

Atmospheric Correction
Remote-sensing applications require removing the atmospheric effect from the imagery, to retrieve the spectral reflectance of the surface materials.In this study, the data were atmospherically corrected using the quick atmospheric correction (QUAC) algorithm [34,35], since we had no prior

Spectral Unmixing and Abundance Retrieval
The processing workflow towards unmixing and generating abundance consists of four steps: (1) Atmospheric correction to retrieve surface reflectance; (2) Data masking, geocorrection, reprojection, and resampling; (3) An endmember selection algorithm was adopted to select the endmembers; then a linear spectral mixing analysis method was employed to retrieve the abundance (Figure 4).

Atmospheric Correction
Remote-sensing applications require removing the atmospheric effect from the imagery, to retrieve the spectral reflectance of the surface materials.In this study, the data were atmospherically corrected using the quick atmospheric correction (QUAC) algorithm [34,35], since we had no prior knowledge to perform empirical calibration [36,37].QUAC is an in-scene approach, requiring only an approximate specification of sensor band locations (i.e., central wavelengths) and their radiometric calibration; no additional metadata is required [35].QUAC does not involve first principles radiative transfer calculations, and therefore it is significantly faster than physics-based methods; however, it is also more approximate [35].

Data Masking, Geocorrection, Reprojection, and Resampling
In this study, we use the Airborne Processing Library (APL) software for processing the data [38].The first step of the APL processing is to apply the mask of bad channels to atmospherically corrected data, creating a new file with bad channels set to zero (Appendix A on Figure A1).The next step uses the navigation file, the view vector file, and the digital elevation file (DEM) to calculate the ground position for each pixel then change the projection to UTM (Universal Transverse Mercator) Zone 28N [38].We used satellite-based ASTER sensor for the DEM.In the final step we resampled output pixel size to ~3.5 m according to the height above ground level (AGL) that is given by the theoretical pixel size chart that can be found in Appendix A on Figure A2 (https://nerc-arf-dan.pml.ac.uk/trac/ wiki/Processing/PixelSize) [39].

Endmembers Selection
The conventional image-based endmember selection approach based on scatterplots of the image bands may not be effective in identifying a sufficient number of endmembers.In this paper, we employed the sequential maximum angle convex cone (SMACC) algorithm [34] to identify spectral image endmembers.Endmembers are spectra that represent pure surface materials in a spectral image.The extreme points were used to determine a convex cone, which defined the first endmember.A constrained oblique projection was applied to the existing cone to derive the next endmember.The cone was then increased to include a new endmember [8,40].This process was repeated until a projection derived an endmember that already existed within the convex cone, or until a specified number of endmembers was satisfied [21].When implemented with SMACC, the output endmember number was set as 5, 10, 15, 20, and 30 respectively.Better endmembers could be identified easily from the 15 endmembers output (more detail in Section 6.2).Then, we used the selected 15 endmembers for deriving the abundance.

Linear Spectral Mixture Analysis
The linear spectral mixture analysis (LSMA) approach was adopted to calculate the abundance of endmembers for each pixel.LSMA assumes that the spectrum measured by a sensor is a linear combination of the spectra of all components (endmembers) within the pixel, and the spectral proportions of the endmembers (i.e., their abundance) reflect the proportion of area covered by distinct features on the ground [8,21].The general equation for linear spectral mixing can be expressed as: where R ij,λ is the measured reflectance at wavelength λ for pixel ij, where i is the column pixel number, and j is the line pixel number; p ij,n is the fraction of endmembers n contributing to the image spectrum of pixel ij; N is the total number of endmembers; R n,λ is the reflectance of endmember n at wavelength λ; and E λ is the error at wavelength λ of the fit of N spectral endmembers.The fraction p ij,n can be solved using a least-square method with fully constrained unmixing.Fully constrained unmixing means that the sum of the endmember fractional (abundance) values for each pixel must equal unity, which requires a complete set of endmembers.Therefore, it should meet the following two conditions: In the majority of cases, the unmixing is only partially constrained because the extracted endmember set is incomplete for the image and only term (2) (i.e., Equation ( 2)) is satisfied.In this study, fully constrained LSMA were applied to the FENIX image to obtain the abundance result and both SMACC and LSMA were executed by ENVI 5.3 and IDL 8.5 language programming. , = 1 In the majority of cases, the unmixing is only partially constrained because the extracted endmember set is incomplete for the image and only term (2) (i.e., Equation ( 2)) is satisfied.In this study, fully constrained LSMA were applied to the FENIX image to obtain the abundance result and both SMACC and LSMA were executed by ENVI 5.3 and IDL 8.5 language programming.

Endmember Groups
The approximate locations of the 15 endmembers selected are shown in Figure 5a.SMACC first finds the brightest spectral in the image and defines it as the first endmember.In this study, the first endmember (endmember 1) represented saturated hot material.We grouped these 15 endmembers into six groups; (1) basalt; (2) hot material; (3) oxidized surface; (4) sulfate mineral; (5) water; and (6) noise (Figure 5b-g).These groups were based on the similar shape of the endmembers with the USGS spectral library; however, the amplitude of the endmembers within a group vary due to illumination conditions, spectral variability, and topography.We added up the abundances within the group to derive the abundance according to this endmembers group.

Endmember Groups
The approximate locations of the 15 endmembers selected are shown in Figure 5a.SMACC first finds the brightest spectral in the image and defines it as the first endmember.In this study, the first endmember (endmember 1) represented saturated hot material.We grouped these 15 endmembers into six groups; (1) basalt; (2) hot material; (3) oxidized surface; (4) sulfate mineral; (5) water; and (6) noise (Figure 5b-g).These groups were based on the similar shape of the endmembers with the USGS spectral library; however, the amplitude of the endmembers within a group vary due to illumination conditions, spectral variability, and topography.We added up the abundances within the group to derive the abundance according to this endmembers group.

Basalt Abundance
Figure 6a indicates the presence of the dominant basalt abundance pixel throughout the image.This abundance is associated with endmember 8 which is characterized by very low reflectance (Figure 5b) due to the presence of large amounts of dark mafic rock since the study area is dominated by basaltic lava (Figure 6b).This abundance is associated with endmember 8 which is characterized by very low reflectance (Figure 5b) due to the presence of large amounts of dark mafic rock since the study area is dominated by basaltic lava (Figure 6b).

Hot Material Abundance
As shown in Figure 7a, the hot material abundance map is very sparse.This abundance is described as blends of the endmember 1, 4, 5, and 13 which are characterized by very high reflectance in the SWIR due to the presence of hot material (Figure 5c). Figure 5a shows that endmembers 1, 4, 5, and 13 are located in the lower right corner and the upper part of the image, Figure 7b shows a false color (NIR-SWIR) image which agrees with the abundance map, i.e., some patches of hot material (red-yellow color) exist in the area.The false color image is created by stacking R: 2200 nm; G: 1600 nm, and B: 896 nm.This indicates that the lava field is still emitting hot material during the data acquisition.

Hot Material Abundance
As shown in Figure 7a, the hot material abundance map is very sparse.This abundance is described as blends of the endmember 1, 4, 5, and 13 which are characterized by very high reflectance in the SWIR due to the presence of hot material (Figure 5c). Figure 5a shows that endmembers 1, 4, 5, and 13 are located in the lower right corner and the upper part of the image, Figure 7b shows a false color (NIR-SWIR) image which agrees with the abundance map, i.e., some patches of hot material (red-yellow color) exist in the area.The false color image is created by stacking R: 2200 nm; G: 1600 nm, and B: 896 nm.This indicates that the lava field is still emitting hot material during the data acquisition.

Oxidized Surface Abundance
The oxidized surface endmembers (3, 6, and 12) have the highest abundance fraction at the vent as shown in Figure 8a.This agrees with a field observation shown in Figure 8b which highlights the matching dominant oxidized surface at the vent wall.

Oxidized Surface Abundance
The oxidized surface endmembers (3, 6, and 12) have the highest abundance fraction at the vent as shown in Figure 8a.This agrees with a field observation shown in Figure 8b which highlights the matching dominant oxidized surface at the vent wall.

Oxidized Surface Abundance
The oxidized surface endmembers (3, 6, and 12) have the highest abundance fraction at the vent as shown in Figure 8a.This agrees with a field observation shown in Figure 8b which highlights the matching dominant oxidized surface at the vent wall.

Sulfate Mineral Abundance
The sulfate mineral endmembers (2, 7, 10, 11, and 15) have the highest abundance fraction around the lava pond and there are four most prominent areas for the sulfate (Figure 9a).This surface mineral looked as if it had been dusted by snow (white color) commonly identified as thernadite (  ) [41].This can be directly seen from a true color image.This mineral formed as the flow cooled, a thin sublimate coating formed on the surface of the lava [41].Figure 9b, c shows the thernadite formed in surface lava at Holuhraun.

Sulfate Mineral Abundance
The sulfate mineral endmembers (2, 7, 10, 11, and 15) have the highest abundance fraction around the lava pond and there are four most prominent areas for the sulfate (Figure 9a).This surface mineral looked as if it had been dusted by snow (white color) commonly identified as thernadite (Na 2 SO 4 ) [41].This can be directly seen from a true color image.This mineral formed as the flow cooled, a thin sublimate coating formed on the surface of the lava [41].Figure 9b,c shows the thernadite formed in surface lava at Holuhraun.

Water Abundance
The water abundance (Figure 10a) has the highest abundance fraction at the location mainly recognized as a glacial river (Figure 10b).Endmember 14 represents water which is characterized by a relatively low reflectance and has the highest reflectance in the blue wavelength.Water has high absorption and virtually no reflectance in the NIR-SWIR wavelengths range (Figure 5f).

Water Abundance
The water abundance (Figure 10a) has the highest abundance fraction at the location mainly recognized as a glacial river (Figure 10b).Endmember 14 represents water which is characterized by a relatively low reflectance and has the highest reflectance in the blue wavelength.Water has high absorption and virtually no reflectance in the NIR-SWIR wavelengths range (Figure 5f).

Noise Abundance
Figure 11 shows the abundance map corresponding to endmember 9. We consider this endmember as representing noise due to an unrecognized spectral signature since this spectrum is characterized by saturated reflectance in channels ~2000 nm and ~2400 nm (Figure 5g).The saturated reflectance could be due to corrupted bands in some pixels.

Noise Abundance
Figure 11 shows the abundance map corresponding to endmember 9. We consider this endmember as representing noise due to an unrecognized spectral signature since this spectrum is characterized by saturated reflectance in channels ~2000 nm and ~2400 nm (Figure 5g).The saturated reflectance could be due to corrupted bands in some pixels.
Figure 11 shows the abundance map corresponding to endmember 9. We consider this endmember as representing noise due to an unrecognized spectral signature since this spectrum is characterized by saturated reflectance in channels ~2000 nm and ~2400 nm (Figure 5g).The saturated reflectance could be due to corrupted bands in some pixels.

False Color Abundance
The abundance results depicted as false color (R: Oxidized surface; G: Sulfate mineral; B: Basalt) images show that the majority of rocks or minerals in the study area are dominated by basalt as shown in the blue color in Figure 12a.The other colors such as magenta and yellow indicate a mixture.The mixture phenomenon is illustrated in Figure 12b, as the surface has 0.25 oxidized surface mix with 0.75 basalt resulting in the magenta color; and 0.25 oxidized surface mix with 0.75 sulfate mineral resulting in the yellow color pixel.

Validation
The very high-resolution aerial photograph was used for ground truth.The aerial photograph was classified into oxidized surface, sulfate, basalt, and water using visual image interpretation and used for validation of the unmixing results.We only validate three endmembers for basalt-oxidized, sulfate, and water-since the noise and hot material cannot be detected based on visual interpretation.We classified the endmembers that have fractional abundance > 0.5.Validation was based on 150 randomly generated point samples within each class.Table 1 show the validation results, with a resulting mean overall accuracy 79% and mean Kappa index of 0.73.This result shows that the abundances have moderate agreement with the sample points.

Validation
The very high-resolution aerial photograph was used for ground truth.The aerial photograph was classified into oxidized surface, sulfate, basalt, and water using visual image interpretation and used for validation of the unmixing results.We only validate three endmembers for basalt-oxidized, sulfate, and water-since the noise and hot material cannot be detected based on visual interpretation.We classified the endmembers that have fractional abundance > 0.5.Validation was based on 150 randomly generated point samples within each class.Table 1 show the validation results, with a resulting mean overall accuracy 79% and mean Kappa index of 0.73.This result shows that the abundances have moderate agreement with the sample points.The correlation between the spectral index images and the abundance image was analyzed.We only correlated the three endmembers since there are no reference spectral indices for sulfate mineral, hot material, and noise.Here we compared the basalt, oxidized, and water abundance images with the mafic, oxidized, and water index images proposed by Inzana et al.,  (Appendix B).We applied these indices to the hyperspectral image and compared them with the result from each abundance.Figure 13a-c shows the scatter plots results.The R 2 values were 0.46, 0.91, and 0.77 for the basalt, oxidized surface, and water, respectively.The oxidized surface and water indicate a good correlation with the indices (Figure 13b,c).This suggests that both oxidation and water generated from a spectral index are properly validated [2,44].Meanwhile, basalt shows a low correlation with the mafic index (Figure 13a) suggesting that the estimates of the basalt surface from the unmixing technique is an overestimation, since the basalt abundance shows the older lava flows as mafic with a relatively high fraction compared to the mafic index that only showed for fresh lava flow.This being due to a full spectrum of hyperspectral can easily differentiate between basalt surface and non-basalt.

Number of Endmembers
The determination of the number of endmembers is critical, since underestimation may result in a poor representation of the mixed pixels, whereas overestimation may result in an overly segregated area [16].Table 2 shows the relationship between the number of endmember and the number of pixels that have fractional abundance > 0.5 and the mean correlation with mafic, oxidized, and water index.We considered abundance >0.5 as high abundance.As the number of endmembers increase, the number of pixels also increases for an oxidized surface, sulfate mineral, water, and noise abundances, respectively.This is due to an increase of endmembers that is detected for each group.Meanwhile, the basalt abundance shows the opposite, as the endmembers increase the number of pixels with abundance >0.5 decreases.These results show that as more endmembers are considered the mixing of basalt with other endmembers increases resulting in a decrease of the fractional abundance of basalt.According to the results, we considered the 15 endmembers as an optimum number for this

Number of Endmembers
The determination of the number of endmembers is critical, since underestimation may result in a poor representation of the mixed pixels, whereas overestimation may result in an overly segregated area [16].Table 2 shows the relationship between the number of endmember and the number of pixels that have fractional abundance > 0.5 and the mean correlation with mafic, oxidized, and water index.We considered abundance >0.5 as high abundance.As the number of endmembers increase, the number of pixels also increases for an oxidized surface, sulfate mineral, water, and noise abundances, respectively.This is due to an increase of endmembers that is detected for each group.Meanwhile, the basalt abundance shows the opposite, as the endmembers increase the number of pixels with abundance >0.5 decreases.These results show that as more endmembers are considered the mixing of basalt with other endmembers increases resulting in a decrease of the fractional abundance of basalt.According to the results, we considered the 15 endmembers as an optimum number for this study since they have the highest mean correlation with mafic, oxidized, and water index.Clearly, the selection of appropriate endmembers in such a diverse volcanic environment, considering the particularities of the FENIX dataset, is of great importance in order to obtain accurate unmixing results.In addition, since only a small number of the available materials spectra are expected to be present in a single pixel, the abundance vectors are often sparse [45].

Size of Lava Field Area
As the methods were only tested on a subset area of the lava field vent, to apply the methods for the entire lava flow is challenging for several reasons.(1) The high spatial heterogeneity typically gives rise to mixed pixels containing multiple materials and it will increase the number of endmembers detected by SMACC [40].(2) Different illumination occurs within the different flight lines for the entire lava flow (Figure 3b) since the data acquisition time is acquired between 16.56 and 17.58 local times which results from the very low sun angle during the acquisition.This problem can be approached by collecting ground truth spectra, extensive calibration, and atmospheric correction using simultaneous and constrained calibration of multiple hyperspectral images through a new generalized empirical line model purposed by Kizel et al. [37].(3) The computation time to perform unmixing also must be considered for the entire lava field since the area is relatively large (84 km 2 ) and the hyperspectral data contains 622 channels with a 3.5-meter spatial resolution.In order to process the full set of data we need to consider using high performance computing (HPC) [46].

Using Full Optical Region for Mapping Recent Lava Flow (VIS-SWIR-TIR)
Hyperspectral VIS-SWIR image data is effective discrimination mafic, oxidation, sulfate etc.However, not all the minerals and surface type are always mapped uniquely with VIS-SWIR hyperspectral data.A typical surface such as rock forming minerals associated with unaltered rocks and alteration minerals associated with altered rocks can be identified with TIR (Thermal Infrared) data [47][48][49].Image processing methods that have become standard for hyperspectral VNIR/SWIR data analysis also work for hyperspectral TIR data [47].Vaughan et al [47] showed that pixel classification techniques based on spectral variability within the scene and mineral libraries for matching spectral emissivity features can be used for TIR-derived mineral maps using SEBASS hyperspectral TIR image data.Hyperspectral TIR instruments operational for airborne surveys are also available in the NERC Airborne Research Facility with a Specim AisaOWL sensor [48].A synergistic use of airborne data from both FENIX (VIS-SWIR) and OWL (TIR) allows great potential for lava discrimination in future study due to the complementary nature of the reflective (VIS-SWIR) and emissive (TIR) spectral regions.This might significantly improve our understanding of physical lava surface properties.Specifically, VIS-SWIR imaging spectrometers can discriminate surface materials and TIR data acquisitions can help to identify the thermal characteristics of different materials [47][48][49].For instance, combining emissivity spectra with reflectance spectra in a mixing model would improve discriminating lava from surfaces [50][51][52].

Conclusions
In this study, an application of potential spectral unmixing methods on 2014-2015 Holuhraun lava flow field was presented.In total, we acquired fifteen spectral endmembers and their abundances.The first endmember was chosen as the brightest pixel which represented saturated incandescent lava.We grouped these 15 endmembers into six groups (basalt, oxidized surface, sulfate mineral, hot material, water, and noise) based on the shape of the endmembers since the amplitude varies due to illumination conditions, spectral variability, and topography.The endmembers represent pure surface materials in a hyperspectral image.We concluded that the selection of appropriate endmembers in such a diverse volcanic environment, considering the particularities of the FENIX dataset, is of great importance in order to obtain accurate unmixing results.Combination of SMACC and LSMA methods offers an optimum and a fast selection for volcanic products segregation However, ground-truthing spectra are recommended for further analysis.A synergistic use of airborne data from both FENIX (VIS-SWIR) and OWL (TIR) gives a great potential for lava discrimination in future study due to the complementary nature of the reflective (VIS-SWIR) and emissive (TIR) spectral regions.This might significantly improve our understanding of physical lava surface properties.

Figure 1 .
Figure 1.Bárðarbunga volcano and the Holuhraun lava flow field.(a) geological setting by the Icelandic Meteorological Office (after modification) [24], (b) coverage of the three main phases after Pedersen et al [1].

Figure 2 .
Figure 2. Illustration of (a) the mixed pixel in the lava surface caused by the presence of small, sub pixel targets within the area; (b) variability of lava surfaces in Holuhraun lava field which include the oxidizing surface, sulfate mineral, and lava.

Figure 2 .
Figure 2. Illustration of (a) the mixed pixel in the lava surface caused by the presence of small, sub pixel targets within the area; (b) variability of lava surfaces in Holuhraun lava field which include the oxidizing surface, sulfate mineral, and lava.

Figure 3 .
Figure 3. (a) Map showing line acquisition of FENIX hyperspectral image in the Holuhraun lava field; (b) Image mosaic from eight FENIX lines collected during the campaign (red box shows the image subset location); (c) Image subset of the focusing study area in the eruptive fissure vent of Holuhraun.

Figure 3 .
Figure 3. (a) Map showing line acquisition of FENIX hyperspectral image in the Holuhraun lava field; (b) Image mosaic from eight FENIX lines collected during the campaign (red box shows the image subset location); (c) Image subset of the focusing study area in the eruptive fissure vent of Holuhraun.

Figure 4 .
Figure 4.The workflow processing to derive an abundance map from the FENIX hyperspectral data

Figure 4 .
Figure 4.The workflow processing to derive an abundance map from the FENIX hyperspectral data

Figure 5 .
Figure 5. (a) The spatial distribution of 15 endmembers extracted by SMACC; The numbers on the image indicate the approximate location of the pixels selected as the represented endmembers of (b) basalt; (c) hot material; (d) oxidized surface; (e) sulfate mineral; (f) water; and (g) noise, extracted by SMACC.

Figure
Figure6aindicates the presence of the dominant basalt abundance pixel throughout the image.This abundance is associated with endmember 8 which is characterized by very low reflectance (Figure5b) due to the presence of large amounts of dark mafic rock since the study area is dominated by basaltic lava (Figure6b).

Figure 6 .
Figure 6.(a) The abundance map for basalt endmember, yellow areas indicate the highest fraction of basalt meanwhile the black areas indicate the lowest fraction of basalt (the red box shows the approximate location of the field photo); (b) field photograph of basaltic lava field of the Holuhraun.

Figure 6 .
Figure 6.(a) The abundance map for basalt endmember, yellow areas indicate the highest fraction of basalt meanwhile the black areas indicate the lowest fraction of basalt (the red box shows the approximate location of the field photo); (b) field photograph of basaltic lava field of the Holuhraun.

Figure 7 .
Figure 7. (a) Abundance map for the hot material endmember, yellow areas indicate the highest fraction of hot material meanwhile the black areas indicate the lowest fraction of incandescent lava; (b) The false color (NIR-SWIR) image show that hot material (red-yellow color) exists in the area.

Figure 7 .
Figure 7. (a) Abundance map for the hot material endmember, yellow areas indicate the highest fraction of hot material meanwhile the black areas indicate the lowest fraction of incandescent lava; (b) The false color (NIR-SWIR) image show that hot material (red-yellow color) exists in the area.

Figure 7 .
Figure 7. (a) Abundance map for the hot material endmember, yellow areas indicate the highest fraction of hot material meanwhile the black areas indicate the lowest fraction of incandescent lava; (b) The false color (NIR-SWIR) image show that hot material (red-yellow color) exists in the area.

Figure 8 .
Figure 8.(a) Abundance map oxidized surface endmember; yellow areas indicate the highest fraction of oxidized surface meanwhile the black areas indicate the lowest fraction of oxidized surface; (b) Field photograph of an oxidized surface of the vent wall (red box and the line shows the approximate location of the field photograph)

Figure 8 .
Figure 8.(a) Abundance map oxidized surface endmember; yellow areas indicate the highest fraction of oxidized surface meanwhile the black areas indicate the lowest fraction of oxidized surface; (b) Field photograph of an oxidized surface of the vent wall (red box and the line shows the approximate location of the field photograph)

Figure 9 .
Figure 9. (a) The Abundance map for the sulfate mineral endmember, the yellow areas indicate the highest fraction of sulfate mineral meanwhile the black areas indicate the lowest fraction of the sulfate mineral; (b) Field photograph of sulfate mineral (white surface) formed on the surface of lava (the red boxes and lines show the approximate location of the field and aerial photo respectively); (c) aerial photograph of sulfate mineral (white surface) formed on the surface of lava (The numbers on the image indicate the approximate location of the sulfate for both the abundance and photograph).

Figure 10 .
Figure 10.(a) The abundance map for water endmember, the highest abundance fraction indicated by a yellow color, and the lowest abundance fraction indicated by a black red box shows the approximate location of the aerial photograph); (b) aerial photograph of the glacial river.

Figure 10 .
Figure 10.(a) The abundance map for water endmember, the highest abundance fraction indicated by a yellow color, and the lowest abundance fraction indicated by a black red box shows the approximate location of the aerial photograph); (b) aerial photograph of the glacial river.

Figure 11 .
Figure 11.The abundance map for the noise endmember, the highest abundance fraction is indicated by the yellow color, and the lowest abundance fraction is indicated by the black color.

Figure 11 .
Figure 11.The abundance map for the noise endmember, the highest abundance fraction is indicated by the yellow color, and the lowest abundance fraction is indicated by the black color.

5. 8 . 20 Figure 12 .
Figure 12.(a) False color of abundance highlighting for R: oxidized surface; G: sulfate mineral; and B: Basalt; (b) Illustration of the mixed pixels in the area, 0.25 oxidized surface mix with 0.75 basalt resulting in the magenta color; and 0.25 oxidized surface mix with 0.75 sulfate mineral resulting in the yellow color pixel.

Figure 12 .
Figure 12.(a) False color of abundance highlighting for R: oxidized surface; G: sulfate mineral; and B: Basalt; (b) Illustration of the mixed pixels in the area, 0.25 oxidized surface mix with 0.75 basalt resulting in the magenta color; and 0.25 oxidized surface mix with 0.75 sulfate mineral resulting in the yellow color pixel.
Remote Sens. 2018, 10, x FOR PEER REVIEW 14 of 20 abundance shows the older lava flows as mafic with a relatively high fraction compared to the mafic index that only showed for fresh lava flow.This being due to a full spectrum of hyperspectral can easily differentiate between basalt surface and non-basalt.

Figure 13 .
Figure 13.Linear regression analysis between the spectral index images and the (a) basalt abundance; (b) oxidized surface abundance; and (c) water index. .

Table 1 .
Validation of the endmembers that have abundance > 0.5.

Table 1 .
Validation of the endmembers that have abundance > 0.5.

Table 2 .
Comparison number of pixels that have abundance >0.5, R 2 and number of endmembers.