Mapping Volcano Deformation using InSAR
- Arpit Shah
- Feb 28, 2022
- 12 min read
Updated: 10 minutes ago
1. Introduction
I tend to recall this comic strip whenever I encounter the word Volcano-
Sabu is a giant alien from Jupiter who resides on Earth with a clever old man Chacha Chaudhary - also the name of this highly popular Indian comic book series. Together, they form a perfect pair of brawn and wit who help get rid of unscrupulous elements from the society. While being proficient at mauling thugs, Sabu does tend to struggle against a devious or occasionally, an even stronger adversary. When the lives of Chacha and his family are threatened in particular, Sabu enters into a Hulk-esque rage which increases his strength manifold and is a forewarning about the destruction that is to follow which Pran, the creator of the comic series, signals to the readers using the erupting volcano visualization.
Post Credits: RUS Copernicus & EO College
HYPERLINKS TO SECTIONS
2. Setting the Context
Humans have feared and revered volcanoes through time. While India has only one active volcano currently, and that too far away from the mainland, we do have a striking name for this phenomenon since ancient times - Jwalamukhi (Sanskrit) which translates to the mouth that spurts fire.
Much like a human mouth, a volcano too breathes in and out, albeit at a longer timescale. Its size, shape and volume expands and contracts due to the changes in pressure exerted by the hot magma within, which the surrounding rocks accommodate by deforming. The trend of deformation is what volcanologists are keen to observe as it is a leading indicator of eruption potential - the more restless the magmatic system becomes, it causes an increase in both, the frequency as well as the intensity of ground displacement through uplift (ground level swelling) and/or subsidence (ground level sinking).
Not discernible to the naked eye, researchers use geodetic monitoring techniques such as Leveling, GNSS and InSAR to detect and measure deformation. In this post, I will demonstrate the method of mapping Ground Deformation and deriving the rate of Displacement using the latter - Interferometry (In) applied on Synthetic Aperture Radar (SAR) Satellite Imagery at three recent eruption sites-
Monitoring Deformation at a subsurface level closer to the source of volcanic activity - the Magma Chamber located 1-10 kilometers below the Earth's surface - would be more pertinent than ground level-studies near the eruption outlet, however, the technology at our disposal is able to permeate only a few meters of rocks that lie between the surface and the underground chamber😓.
If there is a layer of water in between, then even doing that is nigh impossible. Our inability to monitor deformation in submarine conditions is worrisome as this is where most of our planet's productive volcanic systems lie - we had no inkling that Hunga Tonga–Hunga Haʻapai was to erupt in a violent expression of fury on 15th January 2022.
3. The Processing Chain and Video Walkthrough
As indicated earlier, I will apply Interferometry on Sentinel-1 SAR Imagery using SNAP tool to map Ground/Surface-level Deformation and to derive the rates of Displacement (Uplift and Subsidence).
The processing chain of doing so on the software is extensive, here is the graphical representation-
I have prepared a step-by-step video walkthrough of this workflow below, the area of study being Cumbre Vieja Volcano in Canary Islands, Africa which erupted in 2021. Those who prefer to watch content over reading it will stand to benefit from this visual guide. However, I'll recommend that you read this post as it is more detailed and subsequently choose to view this video demonstration.
4. Step-by-Step Guide on InSAR Interferometry to map Deformation
Let me begin to explain the InSAR processing chain as depicted in Figure 3 - you may download the image to your desktop in order to refer to it in a more convenient manner (left-click to expand the image and then right-click it to save it).
PS: In the visual that accompanies each processing step note, I often alternate between the three eruption events I have chosen to study in this post - the image description specifies the study area
Read & Read(2) - Two sets of SAR imagery are used as the input datasets in the workflow - one acquired prior to the eruption and the other after the eruption event. Essentially, I am asking SNAP to read these input datasets which I will subsequently compare/perform change detection on, in order to map deformation and measure the rate of displacement.
TOPSAR-Split - Satellite Imagery datasets tend to be large - each Sentinel-1 SAR dataset that I am using is around 1 GB in size and has a coverage (swath width) of 250 kilometers and spatial resolution is 5 x 20 m per pixel. The TOP in TOPSAR refers to Terrain Observation with Progressive scans - a method of illuminating targets which Sentinel-1 satellite uses, and which is particularly useful for Interferometry and Time Series analysis. In order to reduce computational load and save time during subsequent processing, I would like to restrict the spatial extent of the Imagery dataset to just the the region of my interest i.e. the volcanic eruption site, as well as restrict the number of spectral bands containing reflectance information (backscatter) to just the one which I need for this demonstration - VV polarization. The TOPSAR-Split operator helps me to perform these subsetting operations.
Apply Orbit Files - A Satellite captures data as it navigates in its Orbit. The Orbit Files within the SAR Imagery dataset package contain the satellite's velocity and position metadata at the time of Imagery acquisition - this radar geometry information helps in the precise assignment of acquired backscatter values to the individual pixels in the Imagery dataset. It is necessary to perform this step is because refined Orbit Files take up to a few weeks to be prepared by the satellite operator - the original ones that come bundled with the package may not be entirely accurate, particularly if the Imagery has been acquired very recently and using these would adversely impact the change detection analysis and output interpretation.
Consider the visual output of this processing step the same as Figure 5 as there is no discernible change in it after running the Apply Orbit Files operator.
Back Geocoding - Now that I have updated the Orbit files, the Back Geocoding operator will help me to assign the backscatter readings acquired by the satellite's signal receiver instrument to coordinates on the Earth's surface. This conversion i.e reprojection of data from radar coordinates to geographic coordinates necessitates the use of a Digital Elevation Model. A DEM dataset contains the elevation values of the Earth's bare surface - it excludes the elevation values of natural and built features above the bare surface such as vegetation and buildings respectively. Besides assigning surface reflectance information to pixels (aided by interpolation), this operator also bundles both the reprojected Imagery datasets in a single Product Stack - as can be observed in the Product Explorer box in Figure 6 below. This bundling step is individually known as Coregistration and has a separate operator within SNAP should one want to perform just this operation alone. Coregistering both the processed SAR Imagery datasets is absolutely necessary - only upon doing so will SNAP be able to compare them, and this is vital for our deformation mapping workflow which entails change detection.
Figure 6: Back-Geocoded Imagery Stack factoring in elevation values from a DEM dataset. This stack contains both pre (19 May 2021) and post (12 June 2021) Imagery datasets over Mt. Nyiragongo Volcano eruption site in Congo.
Enhanced Spectral Diversity - This operator refines the Back-Geocoding output by enhancing the geometric precision of the backscatter assignment to pixels using an optimization algorithm. This is the last step in the Coregistration portion of the workflow. The output - a Coregistered ESD Stack - contains precise and comparable sets of Imagery datasets pre and post eruption - ripe for applying Interferometry to extract deformation and displacement information.
Interferogram - This operator is the first step in my quest to detect Deformation and measure the Ground Displacement rate between the pre and post eruption SAR Imagery datasets, with centimeter-level accuracy. The three interferometric outputs generated by this operator are the Intensity, Phase and Coherence bands-
Intensity - A pixel in the Intensity band of a single, raw Imagery dataset (the 'i_' band seen in Figure 4) contains a data value that is the square of the Amplitude. Amplitude is the value of the acquired reflected microwave energy from the Earth's surface, originally transmitted from the satellite itself. However, the Intensity band generated in the interferometric output is formed in a different way - through the multiplication of the Amplitude value of the same pixels across both the Imagery datasets i.e. Intensity of Pixel 1 in interferometric output of Coregistered Stack = Amplitude of Pixel 1 in Imagery 1 x Amplitude of Pixel 1 in Imagery 2.
Phase - Likewise, a pixel in the Phase band of a single, raw Imagery dataset (the 'q_' band seen in Figure 4), contains a data value that is a modified representation of the distance between the satellite's antenna and the ground target. As I am seeking to measure Displacement (rate of Uplift or Subsidence), it the change in Phase value across both the datasets that would matter to me. The Phase band generated in the interferometric output is exactly this - it captures the Phase difference i.e. Phase of Pixel 1 in interferometric output = Phase of Pixel 1 in Imagery 1 − Phase of Pixel 1 in Imagery 2. This interferometric Phase output is what is known as an Interferogram.
Coherence - Coherence is a measure of interferometric quality and is a normalized value ranging from 0 to 1 derived by combining the default Amplitude and Phase values across both the SAR Imagery datasets. If the Amplitude values of the same pixel across both the Imagery datasets are similar and the Phase values of the same pixel are similar too, then the computed Coherence takes a value closer to 1 i.e. very high Coherence. If the Amplitude and/or Phase values across both the datasets are dissimilar, then the Coherence value diminishes towards 0 i.e. very low Coherence. What a deformation researcher is keen to observe in the interferometric output is cluster(s) of pixels with low Coherence, a telltale sign of Deformation.
TOPSAR Deburst - TOPSAR mode of Imaging used by Sentinel-1 satellite involves the transmission of microwave energy in short bursts which illuminates sections of the ground target. The data discontinuity between two bursts result in the formation of horizontal black lines (no data) between the Imagery sections/sub-swaths (as can be seen in Figure 4). The TOPSAR Deburst operator fills in this missing data by estimating the values (technique described here), thereby making the dataset seamless.
Topographic Phase Removal - The difference in Phase values between both the pre and post Imagery datasets can occur due to two reasons-
change in Topographic features (eg. Flooding, Harvesting, Wind reorienting crops etc.)
movement (Uplift/Subsidence) of ground surface i.e. Displacement
In order to visualize Deformation, it becomes important to remove the former - the Topographic Phase, so that what is left is just the Phase (difference) due to Displacement.
Multilooking - The pixel sizes within the Interferogram output are not equal, due to the differing angles of the incoming energy signal (backscatter) - this creates speckle noise which contributes to a grainy image. The Multilooking operator utilizes an algorithm which standardizes all the pixel sizes into square shapes, thereby reducing the speckle and improving the radiometric resolution of the Interferogram albeit at the expense of degrading its spatial resolution. As evident in Figure 11 below, the clarity of the Interferogram improves after Multilooking, allowing for easier loss-of-Phase (Deformation) detection.
Goldstein Phase Filtering - This and the next step is where it gets too technical for me to both, understand as well as explain what is happening in a simple paragraph or two. To summarize from what I have gathered - like Multilooking, Phase Filtering is also used to reduce the noise in the Interferogram output (compare Figure 12 with Figure 11). Phase Filtering also helps the next processing step involving Phase Unwrapping to generate a more accurate output.
Snaphu Export - Snaphu (Phase Unwrapping) - Snaphu Import - Phase to Displacement - SNAPHU (Statistical-cost, Network-flow Algorithm for Phase Unwrapping) is a critical and a complex step in the interferometric processing chain. While Multilooking and Goldstein Phase Filtering helped me to ascertain the presence of Deformation (by observing the loss of Phase/presence of Phase fringes), Phase Unwrapping using SNAPHU and the Phase to Displacement operators help derive the precise rate of Displacement between the pre and post eruption Imagery datasets.
Range-Doppler Terrain Correction - Besides using a Digital Elevation Model to remove water bodies (and their misrepresentative displacement rates) from the geographic extent, the Range-Doppler Terrain Correction operator's primary objective is to make the Interferogram correspond geometrically to the Earth's globular shape, technically known as orthorectification. This is because the tilt of the satellite and of its sensor, the curvature of earth, among other factors have not been accounted for until now. So while I had converted the Imagery Dataset from Radar Coordinates to Geographic Coordinates (as they appear on a 2D Map of Earth) using the Back-Geocoding operator earlier, the Range-Doppler Terrain Correction operator adjusts the spatial relationship between the pixels so that the output appears as it would on a 3D Globe of Earth.
Write (& Export) - This operator enables me to download the output on my system.
Visualizing the Interferometric Phase and Displacement Output on Google Earth
Interferometric Phase Visualization on Google Earth
Q. What could be the reason behind the Phase fringes found far away from the Volcano cone?
A. These Phase fringes correspond to persistent Seismic activity (earthquakes) following the eruption. Deformation Mapping through Remote Sensing helps to ascertain the occurrence of Earthquakes, Landslides or any phenomenon which involves the regional-scale alteration of surface features.
Displacement Visualization on Google Earth
Thank you for reading this post and I hope that you found it as interesting as a blazing Volcano!🤯😁. Feel free to share your feedback. Do watch the Video walkthrough in case you'd like to see the usage of SNAP software to execute the InSAR processing chain.
Let me leave you with a spooky possibility so that you do not crave for volcanoes anymore!
ABOUT US
Intelloc Mapping Services, Kolkata | Mapmyops.com offers Mapping services that can be integrated with Operations Planning, Design and Audit workflows. These include but are not limited to Drone Services, Subsurface Mapping Services, Location Analytics & App Development, Supply Chain Services, Remote Sensing Services and Wastewater Treatment. The services can be rendered pan-India and will aid your organization to meet its stated objectives pertaining to Operational Excellence, Sustainability and Growth.
Broadly, the firm's area of expertise can be split into two categories - Geographic Mapping and Operations Mapping. The Infographic below highlights our capabilities-
Our Mapping for Operations-themed workflow demonstrations can be accessed from the firm's Website / YouTube Channel and an overview can be obtained from this brochure. Happy to address queries and respond to documented requirements. Custom Demonstration, Training & Trials are facilitated only on a paid-basis. Looking forward to being of service.
Regards,