Calculating of Taftan volcano displacement using PSI technique and Sentinel 1 images
الموضوعات : فصلنامه علمی پژوهشی سنجش از دور راداری و نوری و سیستم اطلاعات جغرافیایی
Mahdieh Shirmohamadi
1
(Department of Remote sensing and Geographic Information Systems, Faculty of Geographical Sciences, Kharazmi university, TEHRAN, IRAN)
Seyed Mohammad Tavakkoli Sabour
2
(Department of Remote sensing and Geographic Information Systems, Faculty of Geographical Sciences, Kharazmi university, TEHRAN, IRAN)
Parviz Zeaiean Firooz Abadi
3
(Associated professor of remote sensing and GIS, Faculty of Geography, Kharazmi University)
Javad Sadidi
4
(Department of Remote sensing and Geographic Information Systems,
Faculty of Geographical Sciences, Kharazmi university ,Tehran,)
الکلمات المفتاحية: GPS, Taftan volcano, Sentinel 1, PSI,
ملخص المقالة :
Taftan is a semi-active volcano located in southeast Iran with many craters. The main objective of this study is to investigate whether subsidence or uplift in Taftan Peak. 58 Sentinel 1-A images acquired from January 2015 to December 2020 in the ascending orbit mode, and 102 Sentinel 1-A Sentinel 1-B images acquired from October 2014 to June 2020 in the descending orbit mode were preprocessed for this purpose. The interferograms with the permanent scatterer interferometry (PSI) method were created using SARPROZ and StaMPS softwares, in which atmospheric corrections were made automatically, and then the surface displacement of Taftan volcano estimated.The Line of Sight (LOS) displacement corresponding to the uplift was 0.5 mm to 1 mm yr-1 for the ascending orbit and 1 mm yr-1for the descending orbit. Because no GPS station was close to Taftan volcano, the GPS measurements of one station located in the study area (Saravan station) was used to check the accuracy of PSI method. The GPS station of SARAVAN has situated inside the town ,and it is appropriate to use and analyze PSI technique in this station. As a result, the researchers found that the PSI method is in accordance with the GPS data.
| Available online at https://sanad.iau.ir/Journal/jrors Journal Of Radar and Optical Remote Sensing and GIS JRORS 3 (2023) 26–37 |
|
Calculating of Taftan volcano displacement using PSI technique and Sentinel 1
Mahdieh Shirmohamadia 1, Seyed Mohammad Tavakkoli Sabourb, Parviz Zeaiean Firooz Abadic, Javad Sadidid
a Department of Remote sensing and Geographic Information Systems, Faculty of Geographical Sciences, Kharazmi university, Tehran, Iran
b Department of Remote sensing and Geographic Information Systems, Faculty of Geographical Sciences, Kharazmi university, Tehran, Iran
c Associated professor of remote sensing and GIS, Faculty of Geography, Kharazmi University
d Department of Remote sensing and Geographic Information Systems, Faculty of Geographical Sciences, Kharazmi university ,Tehran, Iran
Received 25 January 2023, Revised 07 June 2023, Accepted 20 September 2023
Abstract
Taftan is a semi-active volcano located in southeast of Iran with a number of craters. The main objective of this study is to investigate whether subsidence or uplift in Taftan peak. A total number of 58 images of Sentinel 1-A acquired between January 2015 to December 2020 in the ascending orbit mode and 102 images of both Sentinel 1-A, Sentinel 1-B acquired between October 2014 to June 2020 in the descending orbit mode were preprocessed for this purpose. The interferograms with the permanent scatterer interferometry (PSI) method were created using SARPROZ and StaMPS softwares in which atmospheric corrections were made automatically and following that get surface displacement of Taftan volcano. The results of the Line Of Sight (LOS) displacement corresponding to the uplift was observed to be 0.5 mm to 1 mm yr-1for the ascending orbit and 1 mm yr-1for the descending orbit. Because of GPS station lack close to Taftan volcano, the GPS measurements of one station located in the study area (Saravan station) was used to check the accuracy of PSI method, the GPS station of SARAVAN has been located inside of town and it is appropriate to use and analyze PSI technique in this station. As a result, it was found that the PSI method is in good agreement with the GPS data.
Keywords: Sentinel 1, PSI, Taftan volcano, GPS
1. Introduction and background
Nowadays, it is very common to investigate natural disasters such as floods, earthquakes, volcanoes, storms, etc. using satellite data. Determining the appropriate strategy to prevent and reduce the damage of natural disasters is one most important application of remote sensing. Remote sensing technology including active and passive remote sensing has already proved its potentials to be utilized in many different situations before, during and after such disasters. With the advances in active remote sensing technologies such as Sentinel 1 satellite from European Space Agency (ESA), surface dynamic monitoring became a reality, making it possible to estimate ground deformation and instability to be measured on time and with high accuracy. There are many developed algorithms to find surface displacements from satellite remote sensing data especially microwave images. Amongst them, Permanent Scatterer Interferometry method (PSI) is a widespread method for small ground surface displacements. This method is one of the most advanced differential interferometry techniques. Moreover, the number of processed images is more than the traditional method, and the reason is that it is possible to identify the points of permanent scatterers in several images. In this technique, unlike the traditional methods, the most important and required corrections are performed in the form of a formulated mathematical model.
Several investigators have reported surface displacements calculation using radar interferometric methods. Massonnet et al., identified the subsidence of Mount Etna located in the island of Sicily, Italy. They indicated that, in practice, it is possible to measure these changes from space without ground control points (Massonnet et al., 1995). In another study, Hooper used 22 ERS-1 and ERS-2 images acquired between 1992 to 2000 from Long Valley volcano, Eastern California and validated it with GPS data and reported an uplift about 3mm yr-1 (Hooper, 2006). Radar interferometric technique was used to estimate land subsidence/uplift for Campi Flegrei volcano near the populous city of Naples, southern Italy (Manconi et al., 2010). They found a maximum land uplift of 2 cm yr-1. Bathke et al., applied the same technique on two image datasets for two periods before and after eruption of Laima volcano Chile (Bathke et al., 2011). Shirzaei et al., monitored the surface deformation of the Damavand volcano located in the central part of the Alborz Mountain in the south of the Mazandaran Sea, Iran by Envisat satellite data during 2003 to 2008. They reached subsidence about 5 mm per year at the peak of the Damavand volcano in 2011. They depicted the velocity vector obtained from the GPS data of two stationary stations, Polour and Absard to the geometry of the satellite view and they achieved a good agreement between the amount of displacement obtained from radar interferometric and GPS observations (Shirzaei et al., 2011). Gorzi Zanjani et al., monitored the surface deformation of Sabalan volcano located in Ardabil in the northwest of Iran by Envisat satellite data during 2003 to 2010 and they reached subsidence about 3 mm per year on 2018 (Gorzi Zanjaniet al, 2019). Surface subsidence in Taftan and Bazman volcanoes were studied using Envisat satellite images acquired between 2003 to 2010 and 2003 to 2005 by Saadatmand. He indicated a subsidence about 5 mm per year for Taftan and 6 mm per year for Bazman volcanoes (Saadatmand, 2018).
1.1. DInSAR details
Briefly speaking, radar deferential interferometry (DInSAR), as an effective SAR data processing technique, makes use of a combination of two SAR images of the same area acquired from slightly different positions (Mao et al., 2021; 2023; Mao et al., 2023). As a result, a new image known as interferogram, whose phase component is created by the following terms:
𝛥𝛷= 𝛥𝛷𝑓𝑙𝑎𝑡+𝛥𝛷ℎ𝑒𝑖𝑔ℎ𝑡+𝛥𝛷𝑑𝑖𝑠𝑝lacement+𝛥𝛷𝑎𝑡𝑚𝑜sphere+ 𝛥𝛷noise (1)
Where 𝛥𝛷𝑓𝑙𝑎𝑡 is the flat earth phase that can be estimated from orbital data and easily removed, 𝛥𝛷ℎ𝑒𝑖𝑔ℎ𝑡 is the topographic phase error which is due to the inaccuracy of the external DEM, 𝛥𝛷𝑑𝑖𝑠placement is displacement phase, 𝛥𝛷𝑎𝑡𝑚𝑜sphere is atmospheric delay, 𝛥𝛷noise is residual noise. The goal of DInSAR is to remove all the interferometric phases except displacement. Conventional DInSAR is able to remove all of them except atmosphere noise, it can be resolved when employing advanced techniques, such as the permanent scatterer interferometry (PSI).
1.2. PSI Technique
The limitations of traditional radar interferometric method led to the creation of a technique called permanent scatterers. It utilizes ground features called permanent scatterers having relatively stable backscatter characteristics during several months and even several years. Image pixels containing permanent scatterers have a uniform and stable phase behavior over time. One of the images is selected as the master image based on its maximum cumulative coherence value, and the other images as the slave images. The selection of permanent scatter points is done in two steps, the first step is based on the backscatter amplitude and the second step is based on the corrected phase of the scatter points.
Equation 2 is the general equation for estimating phase differences in the permanent scatterer points having k+1 SAR images, k interferograms and H points of the permanent scatterer. If the pixels coordinate of the SAR image in the azimuth and range directions are x, η, then it is possible to calculate the following items:
- the amount of the phase difference in the correction of topographic effect (ΦΔ) in the set of permanent scatterer points (H) and
- The interferogram phase related to K interferograms produced in a matrix relation to the dimensions K*H equation 2.
ΦΔ=a1T + pxxT + pηηT + BΔqT + TvT + E (2)
Where a1T + pxxT + pηηT is atmospheric phase, BΔqT is topographic phase, TvT displacement phase, E is residual vector. In this non-linear equation, K*H independent equations and 3K+2H unknown variables are exist. In the first stage of solving equation 2, the displacement value is modeled and the topography is approximated. In this way, The residual phase value includes atmospheric delay, unmodeled displacement and other residual errors in the model. In order to remove the atmospheric effects, and at the beginning, the atmospheric contribution in the master image is approximated, then assuming that atmospheric conditions is not correlated with time but the modeled displacement is correlated, the contribution of the atmosphere is separated from the unmodeled displacement. After this separation, in the candidate permanent scatterer points, this value is spatially modeled and predicted for the entire image in each interferogram using the Kriging and variogram interpolation method. In the next step, the effects of noise are removed from the atmospheric signal. After estimating the atmospheric phase in each pixel, this signal is subtracted from the initial phase. Once, all the interferograms were corrected from atmospheric effects, the phase value and desired parameters are approximated at the points of the second order permanent scatterer. After selecting the final set of permanent scatterers, the displacement parameters, the atmospheric phase and the residual phase are determined (Ferretti et al., 2001).
2. Materials and Methods
2.1. Study area and data used
Taftan volcanic crater located in in 42 km northwest of Khash in Sistan and Baluchestan province between 61° 00‘E to 61° 15‘E longitudes and 28°30‘N to 28°40‘N latitudes (Figure 1). Taftan is a young and semi-active volcano with an area of 398.98 square kilometres. Its highest crater called Anjark has a height of 3943 meters. It is also 4050 meters above mean sea level.
In this study, 58 Sentinel 1-A images acquired between January 2015 to December 2020 in the ascending orbit mode and 102 Sentinel 1-A, Sentinel 1-B images acquired between October 2014 to June 2020 in the descending orbit mode have been used. For the validation of PSI method in Saravan station, 41 images of Sentinel 1-A acquired in 2014 to 2018 along with GPS data have been used. Satellite images have been analyzed with SARPROZ, Stamps, SNAP software packages.
Figure1. Study Area, Taftan Volcano, southeast Iran
2.2. Methodology
2.2.1 PSI implementation
58 Sentinel 1-A images in Taftan area were analyzed with SARPROZ. To speed up the process, the smaller area of interest (AOI) was defined. The master images were extracted based on maximum the cumulative coherence of all the images. The time of radar image acquisition and temporal and spatial baselines have been displayed as star graphs in figure 2. The central red point represents the master image in figure 2. The image of ascending orbit acquired on 01/09/2017 for the Taftan peak was selected as master image in figure a and the image acquired on 16/05/2016 for SARAVAN station was selected as master image in figure b. Then the images were co-registered on the basis of precision orbits and the digital elevation model (DEM).
In SARPROZ Amplitude Stability Index (ASI) is applied to select the permanent scatterers (PS):
ASI = 1- DA = 1 – (sA/mA) (3)
where DA represents the amplitude dispersion, mA is the mean deviation of amplitude in time, and sA is the standard deviation of amplitude in time. (Ferretti et al., 2001) recommended a suitable threshold to select the first PSS when ASI > 0.75 (DA <0.25). After choosing the first order PS, it is necessary to estimate a reference network by connecting the PSS by Delaunay triangulation. The differential deformation velocity and differential residual topographic error are evaluated. Then the estimated linear model is subtracted and the APS is estimated from the phase residuals by graph inversion. It is also necessary to fix the velocity amount of at least one pixel (reference point). The APS estimation quality was assessed by an analysis of the temporal coherence of PSs after graph inversion and APS removal.
The next step was the selection of second order PS. The final processing for APS removal was conducted using the same parameters and the reference point similar to APS estimation. The final set of PS was then geocoded and displayed in Google Earth.
For the detailed study of Taftan volcano crater movement, two points were randomly selected and used through time series of ascending displacements.
SNAP and Stamps software was used for processing of 102 images of Sentinel 1-A, Sentinel 1-B in Taftan area acquired on the descending orbit all from the time interval in between October 2014 to June 2020.
Complete processing steps in SNAP, Stamps are as followings:
- Selection of master image and co-registration of master and slave images,
- Enhanced spectral diversity correction,
- Production of interferogram,
- Removing of topography phase,
- Selection of initial PS based on amplitude dispersion threshold,
- Estimation of phase noise and PS selection,
- Weeding of PS,
- Estimation of DEM error and resampling,
- Phase unwrapping,
- Estimation of spatially corrected look angle error and atmospheric phase,
- correction of atmospheric phase,
and Creation of PS displacement time series.
In order to remove the effect of topography phase, SRTM DEM with 30 meters spatial resolution has been used.
a |
b |
20150103 20170901
Figure 2. Images of star graphs depicting the temporal and perpendicular baseline of timeseries in the Taftan peak (a), SARAVAN station (b). Red circle is master image
2.2.2. Validation
One method of measuring behaviour of volcanoes is the use of GPS observations. The network of permanent stations in Sistan and Baluchistan province including 5 GPS stations are shown in Figure 3.
Considering not to be GPS station near Taftan volcano (this is a lack in this area), it was not possible to check the precision of PSI method in Taftan area so, necessarily for the validation of PSI method, one of the nearby stations was used. Saravan station was chosen due to the fact that Saravan station is within the urban area and is suitable for PSI method. 41 images of Sentinel 1-A acquired in 2014 to 2018 Saravan area were used and processed in SARPROZ by PSI method and the results of this step have been compared with GPS’s observation data of Saravan station. The three- dimensional position of the scatterere points with a millimeter displacement along LOS of the satellite in SARAVAN GPS station has been shown in figure 4. To calculate displacement accurately P1, P2, P3 points close to SARAVAN station were used.
Figure 3. The network of permanent GPS stations in Sistan and Baluchistan province
Figure 4. Displacement velocity (in LOS) in SARAVAN city area that consist of GPS station
3. Results
3.1. PSI Results
PSI technique gives the best results in urban and barren mountains areas. The three- dimensional position of the scattering points with a millimeter displacement along the line of sight of the satellite for ascending orbit in Taftan volcano area, has been calculated by the SARPROZ software. It has been shown in figure 5.
Figure 5. Displacement velocity (in LOS) calculated by SARPROZ software in Taftan area and displayed on Google Earth
As it is indicated in figure 5, the total velocity of displacement between -3 to +6 mm yr-1along LOS in the whole region and uplift about 0 to 3 mm yr-1in the Taftan peak in figure 5.
For images of descending orbit processed in SNAP,StaMPS in Tatan volcano area displacement velocity is from -5 to +5 mm yr-1along LOS in the whole area and uplift about 0 to 1 mm yr-1in the Taftan peak ( figure 6).
Figure 6. Displacement velocity (in LOS) calculated by SNAP, StaMPS software in Taftan area
3.2. Validation Results
In P1, P2, P3 points around Saravan GPS station, it is observed that the rate of displacement velocity (in LOS) is 0.9mm,0.2mm,0.5mm yr-1 respectively. Average displacement velocity is 0.5mm yr-1 (in LOS) in SARAVAN GPS station (figure 7, 8, 9).
Figure 10 shows time series of displacement in SARAVAN GPS station on basis on GPS data for three directions (East-West, North-South, Vertical). Based on displacements in these 3 directions, displacement velocities in SARAVAN station are as follows:
VE = 0.5 mm
VN = 0.5 mm
VH = 2 mm
Figure 7. Time series and displacement velocity (in LOS) P1 point are located in near SARAVAN GPS station
Figure 8. Time series and displacement velocity (in LOS) P2 point are located in near SARAVAN GPS station
Figure 9. Time series and displacement velocity (in LOS) P3 point are located in near SARAVAN GPS station
Figure 10. Time series of displacement in direction of north-south, east-west, vertical between 2014 to 2018 SARAVAN GPS station
To compare GPS to a single geometry, it is necessary to convert the GPS measurements to LOS (Hanssen, 2001). The following matrix relationship should be followed:
VN
VLOS = [ sin θ sin α - sin θ cos α cos θ] VE (4)
VH
The VE, VN, VH are the East-West, North- South, Vertical components from GPS station, θ is the incidence angle of satellite and α is the heading angle(azimuth angle). By substituting the values in equation 4, displacement velocity (in LOS) of SARAVAN GPS station is calculated about 1 mm per year. It is concluded that the amount of displacement velocity calculated using PSI is in agreement with displacement velocity of SARAVAN GPS station.
4. Discussion
As it is clear from figures 5 and 6, the amount of displacement in Taftan volcano crater in two ascending and descending modes using PSI method are consistent. For a more detailed investigation in the Taftan volcano crater, the results of time series analysis on displacement indicates an uplift about 0.5 and 0.9 mm yr-1 at point (figure 11). It is also revealed a good time coherency of 0.96 and 0.94 in this point, which shows the reliability of this point.
Figure 11. Time series and displacement velocity (in LOS) ascending mode at selected one point in Taftan peak calculate in SARPROZ software
Figure 12 shows time series of displacement in descending mode obtained using SNAP and StaMPS. It is worth mentioning that this result is in agreement with ascending mode (shown in Figures 11, 12). Despite of the result obtained in this study (as uplift), the result of the study conducted by Saadatmand (Saadatmand, 2018), was different (as subsidence). This may be due to different in data acquisition time, less number of images used by Saadatmand, less spatial resolution of Envisat data in comparisons with Sentinel1 data.
Figure 12. Time series and displacement descending mode in Taftan volcano area calculated in SNAP, StaMPS software
5. Conclusions
Both results of ascending and descending orbital modes show uplift in Taftan volcano crater. So, it concludes uplift between 2014 to 2020.
On the other hand, the obtained results using PSI technique and GPS data in SARAVAN GPS station had similar results. It results that the radar interferometric method can be used for calculating of displacement as an accurate method in areas where there is a lack of GPS stations.
It is suggested as a future research and for a more detailed investigation of the volcanic activity, in addition to using PSI technique, the surface temperature of Taftan volcano is checked using satellite images, and seismicity data are explored. finally, using composing these parameters, a more accurate analysis of volcanic activity will obtain.
References
Bathke, H., Shirzaei, M., Walter, T., R., 2011. Inflation and deflation at the steep‐sided Llaima stratovolcano (Chile) detected by using InSAR. Geophysical Research Letters, 38(10), doi.10.1029/2011GL047168.
Gorzi Zanjani, S., Mousavi, Z., Rezaeian, M., 2019. Survey of deformation around the Sabalan volcano using Interferometric Synthetic Aperture Radar. 18th Iran Geophysical conference, 1008-1011.
Hanssen, R., F., 2001. Radar interferometry: data interpretation and error analysis. Springer publishers, 308.
Manconi, A., Walter, T., R., Manzo, M., Zeni, G., Tizzani, P., Sansosti, E., Lanari, R., 2010. On the effects of 3‐D mechanical heterogeneities at Campi Flegrei caldera, southern Italy. Journal of Geophysical Research: Solid Earth, v. 115, no. B8.
Mao, W., Liu, G., Wang, X., Xiang, W., Wu, S., Zhang, B., Bao, J., Cai, J., Zhang, R., Pirasteh, S., 2021. An InSAR Ionospheric Correction Method Based on Variance Component Estimation with Integration of MAI and RSS Measurements. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing. 14: 1423-1433. doi.10.1109/JSTARS.2020.3045267.
Mao, W., Wang, X., Liu, G., Pirasteh, S., Ma, Z., Zhang, R., Lin, H., Xie, Y., Xiang, W., Wu, S., Ma, Z., Ma, P., 2023. Time Series InSAR Ionospheric Delay Estimation, Correction, and Ground Deformation Monitoring with Reformulating Range Split-Spectrum Interferometry. IEEE Transactions on Geoscience and Remote Sensing. 61 (99): 1-18 doi. 10.1109/TGRS.2023.3298919.
Hooper, A., J., 2006. Persistent scatter radar interferometry for crustal deformation studies and modeling of volcanic deformation. PHD thesis; Stanford University, 144
Massonnet, D., Briole, P., Arnaud, A., 1995. Deflation of Mount Etna monitored by spaceborne radar interferometry. Nature, v. 375, no. 6532, 567.
Saadatmand, Y., 2018. Surface deformation monitoring of volcano using InSAR,case study: Makran. Master thesis; University of Zanjan, 132.
Shirzaei, M., Walter, T., Nankali, H., Holohan, E., 2011. Gravity-driven deformation of Damavand volcano, Iran, detected through InSAR time series. Geology, v. 39, no. 3, 251- 254. doi.10.1130/G31779.
[1] Corresponding author. Tel:.+98-9924625843
E-mail address: mshirmohammadi.khu@gmail.com