Automated procedure for InSAR data inversion using

Download Report

Transcript Automated procedure for InSAR data inversion using

Automated procedure for InSAR data inversion using Finite Element Method

Currenti G.

1

, Del Negro C.

1

, Scandura D.

1

, Williams C.A.

2

1) Istituto Nazionale di Geofisica e Vulcanologia – Sezione Catania 2) GNS Science – New Zealand USEReST, Naples, 11-13 November 2008

Deformation model: Fault slip distribution

Dislocation source parameters determined from multiple types of static deformation data such as GPS displacements, InSAR imagery, and surface offset measurements suggest that slip along a fault is not uniform and is best described as a distribution of dislocation sources.

d i

 

j G ij s j

displacement Green’s function d

Gs slip

The procedure requires:

1)

subdividing the faults in a finite number of patches;

2)

computing the Green's Function for all the patches and the measurement

3)

points; solving a linear inversion problem to determine the slip distribution.

FEM generated Green’s function Medium Heterogeneity

At each node of the mesh different values of Young modulus and Poisson ratio are assigned on the basis of seismic tomography investigations (Patanè et al., 2006).

Young modulus

E

 5 / 6 

V p

2   [(

V p

/

V s

) 2  2 ] /[ 2 (

V p

/

V s

) 2  2 ]

11.5-133 GPa

Poisson ratio

0.12-0.32

Real Topography

LaGrit

Mesh refinement and smoothing are available to provide a mesh with more resolution in areas of interest.

www.meshing.lanl.gov

PyLith

PyLith is a finite element parallel code for the solution of

dynamic

and

quasi-static

deformation problems.

www.geodynamcis.org

FEM Inversion System of linear equations

d d

2 1  

G

11

s

1

G

12

s

1  

G

1

M

 

G

2

M s M s M d N

   

G N

1

s

1  

G NM s M

n x m patches 3 x n x m unknown variables!

(dip-slip, strike-slip, tensile)

j ≥ 3 x n x m

Non-uniquiness of the solution can arise from limited data and poor knowledge of the internal structure. This issue can be faced using: (i) as much data as possible and (ii) numerical technique to reduce the algebraic ambiguity.

Geodetic Data Integration

                       

U U d d U

 1

S N

   

U U

 

U x

N N y z y z N

1

x

1 1

S G G G

           

G G G

   

G G G N N N S S S G N N N S

 3    

S S S

3 1 , 1  

N N

3 ,

M

 1 , 1 2

G G N

, 1 3 , 1 

G

2  1 , 1 , 1 , 1           EDM, GPS, leveling data assure high measurements accuracy but the

G G G N G N G G N S G N N S N S S S

1 ,  

S M

,

M

 1 ,

M

 2 ,

M

 3 ,

M G N S

 3

N G

 2 ,

M

 3

N G

 1 ,

M

 3

N G

,

M

                  

s

 1

s M

    coverage area is usually limited. InSAR data, despite the lower accuracy, can provide a better deformation pattern thanks to the wide coverage.

 

d SAR d GPS

    

G G

overview

SAR GPS

 

s

A

W d

covariance matrix is introduced to weight the data depending on measurements uncertainties.

of the Linear least-squares methods for this problem require to incorporate regularization techniques in order to stabilize the problem and to reduce the set of likely solution.

min 

d

 min 1 2  (

G s

d

)

T

W

d T

W

d

(

G s

d

)     min   min  

d

 

r

FEM Inversion & Regularization

The inverse problem can be re-formulated as an optimization problem aimed at finding the unknown slip values s that minimize a data misfit and a smoothing functional: min   min    1 2  (

G s

d

)

T

W

d T

W

d

(

G s

d

)  

s

T

W

T

Ws

    ,

L

s

U

As smoothing functional, the Laplacian Operator was used to avoid large variations between neighboring dislocations.

The minimization of the quadratic functional φ subjected to bound constraints can be solved by using a Quadratic Programming algorithm based on an active set strategy (Gill et al., 1991): min   min 1 2 

s

T

Q s

b

T

s

  ,

L

s

U

where:

Q

G

T

W

T d

W

d

G

 

W

T

W b

G

T

W

d T

W

d

d

L-curve criterion

The L-curves of the data misfit versus the model norm as a function of the regularization parameter. The best value of regularization parameter lies on the corner of the L-curve.

Input Fault Slip Distribution φ d φ w

Procedure Scheme Ground deformation data (SAR, GPS, EDM) 1 – Domain and Fault discretization: LaGrit 2 – FEM generated Green’s function: PyLith 3 – Solving a linear inversion problem: QP algorithm

Analytical vs Numerical Solution

We considered a 3D FEM model reproducing a rectangular dislocation source in a homogeneous and isotropic half-space and compare it with the Okada’s solution.

A misfit index was used to quantify the discrepancy between the numerical and the analytical solution: 

i

(

U

) 

i N

  1

U i FE

U i AN

/

i N

  1

U i AN

Simulat.

Nodes Elements Quality # 1

2228 9858 0.3544

# 2

2513 11236 0.3083

# 3

2570 11505 0.3441

# 4

3195 14663 0.3109

# 5

3977 18888 0.3417

# 6

5101 25141 0.2960

# 7

8591 45021 0.3050

# 8

36196 205700 0.2454

Currenti et al. PEPI 2008

Tests on the Accuracy

We compare the ground deformation achieved as the sum of the displacements generated by each patch and those obtained by the overall fault. In these computations a uniform slip distribution on the patch is assumed. The error has a RMSE (root mean square error) of the order of 10 -7 m, which is well below the measurement error.

The accuracy of the solution is dependent on the number of the observation points.

A uniform distribution of observation points within a regular grid was assumed, centered on the ground projection of the sources. The RMSE between the assumed and the inverted slip is computed as the number measurement points increases.

of

A case study: Northeastern flank movement at Etna volcano in 2002

On 22 September 2002, 1 month before the beginning of the flank eruption on the NE Rift, an M-3.7 earthquake struck the northeastern part of Mt. Etna, on the westernmost part of the Pernicana fault.

The GPS surveys carried out in September and July 2002 shows a ground deformation pattern that affects the northeastern flank clearly shaped by the Pernicana fault.

Differential interferogram for ascending scene pair 31 July 2002 to 09 October 2002. The scale indicates the phase variation along the LOS (negative values correspond to the approaching of the surface to the sensor)

Bonforte et al., Bull Volc. 2007

Source from GPS data Inversion

A numerical model was set up using the dislocation sources from GPS data inversion.

Topography and medium heterogeneity are also taken into account.

Bonforte et al., Bull Volc. 2007

3D FEM of Mt. Etna: slip on Northeastern Flank Mesh Domain

:100 x 100 x 50 km

Nodes

: 129253

Tetrahedra Elements

: 744886

Computing Time

: 15 minutes

Sar Data and Source Discretization

Line of sight change map calculated by the unwrapped interferograms, processed using ROI_PAC software, developed by Jet Propulsion Laboratory & Caltech; Courtesy of F. Guglielmino ed A. Spata.

A deformation pattern is clearly seen around the Pernicana subdivided in patches to obtain the slip distribution from the inversion procedure.

Green Functions Computations

248 fault patches

FEM simulations

1 2 87

248 x 3 = 744 FEM simulations

744 Dip-slip, strike-slip and tensile dislocation were assigned to the nodes lying along the fault surface. For each patches, in which the fault is discretized, Green’s functions were computed using PyLith. For each each accuracy simulations, was the warranted checking the convergency of the FEM solution. The iteration of GMRES solver is stopped when a threshold of 10 -9 is reached or the number of iterations is higher than 200.

By a linear speedup on a cluster of 20 nodes the computing reduces from 10 days to time 10 hours

.

Inversion Results

19 18.5

18 17.5

17 16.5

16 0 1 2 3 4 Smoothing Functional 5 6 1000 7 x 10 -9

SAR data

The kinematic of the faults was constrained to those derived from GPS data inversion. The solution in correspondence of the corner of the L-curve provides a deformation pattern which resembles the SAR data and also matches the GPS obtained from the analytical solution.

Inverted Solution

Slip Distribution

Strike-slip, dip-slip and tensile dislocation distribution.

Model Comparison SAR data

The GPS data inversion provides a deformation pattern which misses the clear anomaly around the Pernicana fault system.

A heterogeneous distribution of the slip along the structure is able to better justify the SAR data. The FEM model based on the numerical inversion provides a more complex pattern which is likely to be expected.

Analytical Model Numerical Model

Conclusions

An automated procedure for geodetic data inversion is proposed to estimate slip distribution along fault interfaces. 3D finite element models (FEMs) are implemented to compute synthetic Green’s functions for static displacement.

FEM-generated Green’s functions computed using PyLith, a source-free parallel finite element code, are Programming combined algorithm with to a invert Quadratic ground deformation data and obtain an estimate of slip distributions along seismogenic structures.

The procedure was applied to study the ground deformation preceding the 2002-03 Etna eruption.

SAR images showed a significant deformation pattern of the north-eastern flank of the volcano involving the main local volcanic edifice features on the northeastern flank. The numerical model highlighted a heterogeneity slip distribution along the Pernicana fault with a predominant strike-slip mechanism associated with a dip-slip movement in the western part.

The time consuming computation of the procedure is related to the generation of the Green’s functions. For the main volcano edifice features and the seismogenic structures they has to be computed only once and stored in a database.

Numerical Solution

Several computations were performed to better understand the effects induced by the topography and the heterogeneity.

medium The topography significantly alters the general pattern of the ground deformation especially near the volcano summit. On the contrary, the medium heterogeneity does not strongly affect the expected deformation.

Currenti et al. PEPI 2008

Smoothing Functional

min 1 2 

G

 

s

1 2   (

d G

T

s

d

T

)

T

W

d T

(

G

d

s

d

(

d

)

s

  (

s

)

s

 0

s

)

L L T T L L

s

(

s

  

s

0 )   The Laplacian Operator is introduced to avoid large variations between neighboring dislocations  2

s

 0 

s i

 1 ,

j

 2

s i

,

j

( 

x

) 2 

s i

 1 ,

j

s i

,

j

 1  2

s i

,

j

( 

y

) 2 

s i

,

j

 1 The matrix the

n-th

W L

row is constructed so that in

W L

contains coefficients in the equation above for columns corresponding neighboring source.

to the

W L

      

W W W W W

( ( ( (

j j j j

 1 )  1 )  1 )  2

jnc

i

)

nc nc nc nc

     ( (

i i i i

  1 )  1 )  

y

2     

   

2 

y x

   

 2   2 2 

y

 2 