Upgrade to Pro — share decks privately, control downloads, hide ads and more …

Inversão robusta do gradiente da gravidade através da plantação de anomalias de densidade

Leonardo Uieda
September 21, 2011

Inversão robusta do gradiente da gravidade através da plantação de anomalias de densidade

Seminário anual (2011) da pós-graduação em geofísica do Observatório Nacional.

Leonardo Uieda

September 21, 2011
Tweet

More Decks by Leonardo Uieda

Other Decks in Science

Transcript

  1. Robust 3D gravity gradient inversion
    by planting anomalous densities
    Leonardo Uieda
    Valéria C. F. Barbosa
    2011
    Observatório Nacional

    View full-size slide

  2. Forward Problem
    Outline

    View full-size slide

  3. Inverse Problem
    Forward Problem
    Outline

    View full-size slide

  4. Inverse Problem Planting Algorithm
    Forward Problem
    Inspired by René (1986)
    Outline

    View full-size slide

  5. Inverse Problem Planting Algorithm
    Synthetic Data
    Forward Problem
    Inspired by René (1986)
    Outline

    View full-size slide

  6. Inverse Problem Planting Algorithm
    Synthetic Data Real Data
    Forward Problem
    Inspired by René (1986)
    Outline
    Quadrilátero Ferrífero

    View full-size slide

  7. Forward problem

    View full-size slide

  8. Observed g
    αβ
    gαβ

    View full-size slide

  9. Observed g
    αβ
    gαβ

    View full-size slide

  10. Observed g
    αβ
    gαβ
    Anomalous density

    View full-size slide

  11. Observed g
    αβ
    gαβ
    Anomalous density
    Want to model this

    View full-size slide

  12. Interpretative model

    View full-size slide

  13. Interpretative model
    Right rectangular prism
    Δρ=p
    j

    View full-size slide

  14. Prisms with
    not shown
    p
    j
    =0

    View full-size slide

  15. p=
    [p
    1
    p
    2

    p
    M
    ]
    Prisms with
    not shown
    p
    j
    =0

    View full-size slide

  16. Predicted g
    αβ
    dαβ
    p=
    [p
    1
    p
    2

    p
    M
    ]
    Prisms with
    not shown
    p
    j
    =0

    View full-size slide

  17. Predicted g
    αβ
    dαβ
    p=
    [p
    1
    p
    2

    p
    M
    ]
    dαβ=∑
    j=1
    M
    p
    j
    a
    j
    αβ
    Prisms with
    not shown
    p
    j
    =0

    View full-size slide

  18. Predicted g
    αβ
    dαβ
    p=
    [p
    1
    p
    2

    p
    M
    ]
    dαβ=∑
    j=1
    M
    p
    j
    a
    j
    αβ
    Contribution of jth prism
    Prisms with
    not shown
    p
    j
    =0

    View full-size slide

  19. dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    More components:

    View full-size slide

  20. d
    dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    More components:

    View full-size slide

  21. d
    dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    =∑
    j=1
    M
    p
    j
    a
    j
    More components:

    View full-size slide

  22. d
    dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    =∑
    j=1
    M
    p
    j
    a
    j
    = A p
    More components:

    View full-size slide

  23. d
    dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    =∑
    j=1
    M
    p
    j
    a
    j
    = A p
    Jacobian (sensitivity) matrix
    More components:

    View full-size slide

  24. d
    dxx
    dxy
    dxz
    dyy
    dyz
    dzz
    =∑
    j=1
    M
    p
    j
    a
    j
    = A p
    Column vector of
    More components:
    A

    View full-size slide

  25. Forward problem:
    p d
    d=∑
    j=1
    M
    p
    j
    a
    j

    View full-size slide

  26. ̂
    p g
    ?
    Inverse problem:

    View full-size slide

  27. Inverse problem

    View full-size slide

  28. Minimize difference between and
    g d
    r=g−d
    Residual vector

    View full-size slide

  29. Minimize difference between and
    g d
    r=g−d
    Residual vector
    Data­misfit function:
    ϕ( p)=∥r∥2
    =
    (∑
    i=1
    N
    (g
    i
    −d
    i
    )2
    )1
    2

    View full-size slide

  30. Minimize difference between and
    g d
    r=g−d
    Residual vector
    Data­misfit function:
    ϕ( p)=∥r∥2
    =
    (∑
    i=1
    N
    (g
    i
    −d
    i
    )2
    )1
    2
    ℓ2­norm of r

    View full-size slide

  31. Minimize difference between and
    g d
    r=g−d
    Residual vector
    Data­misfit function:
    ϕ( p)=∥r∥2
    =
    (∑
    i=1
    N
    (g
    i
    −d
    i
    )2
    )1
    2
    ℓ2­norm of r
    Least­squares fit

    View full-size slide

  32. Minimize difference between and
    g d
    r=g−d
    Residual vector
    Data­misfit function:
    ϕ( p)=∥r∥2
    =
    (∑
    i=1
    N
    (g
    i
    −d
    i
    )2
    )1
    2
    ℓ2­norm of r
    Least­squares fit
    ϕ( p)=∥r∥1
    =∑
    i=1
    N
    ∣g
    i
    −d
    i

    ℓ1­norm of r

    View full-size slide

  33. Minimize difference between and
    g d
    r=g−d
    Residual vector
    Data­misfit function:
    ϕ( p)=∥r∥2
    =
    (∑
    i=1
    N
    (g
    i
    −d
    i
    )2
    )1
    2
    ℓ2­norm of r
    Least­squares fit
    ϕ( p)=∥r∥1
    =∑
    i=1
    N
    ∣g
    i
    −d
    i

    ℓ1­norm of r
    Robust fit

    View full-size slide

  34. ill­posed problem
    non­existent
    non­unique
    non­stable

    View full-size slide

  35. ill­posed problem
    non­existent
    non­unique
    non­stable
    constraints

    View full-size slide

  36. ill­posed problem
    non­existent
    non­unique
    non­stable
    well­posed problem
    exist
    unique
    stable
    constraints

    View full-size slide

  37. Constraints:
    1. Compact

    View full-size slide

  38. Constraints:
    1. Compact no holes inside

    View full-size slide

  39. Constraints:
    1. Compact no holes inside
    2. Concentrated around “seeds”

    View full-size slide

  40. Constraints:
    1. Compact no holes inside
    2. Concentrated around “seeds”

    User­specified prisms

    Given density contrasts

    Any # of ≠ density contrasts
    ρs

    View full-size slide

  41. Constraints:
    1. Compact no holes inside
    2. Concentrated around “seeds”

    User­specified prisms

    Given density contrasts
    3. Only

    Any # of ≠ density contrasts
    or
    p
    j
    =0 p
    j
    =ρs
    ρs

    View full-size slide

  42. Constraints:
    1. Compact no holes inside
    2. Concentrated around “seeds”

    User­specified prisms

    Given density contrasts
    3. Only

    Any # of ≠ density contrasts
    or
    p
    j
    =0 p
    j
    =ρs
    ρs
    4. of closest seed
    p
    j
    =ρs

    View full-size slide

  43. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)

    View full-size slide

  44. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    Data­misfit function

    View full-size slide

  45. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    (Tradeoff between fit and regularization)
    Regularizing parameter

    View full-size slide

  46. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    Regularizing function
    θ( p)=∑
    j=1
    M p
    j
    p
    j

    l
    j
    β

    View full-size slide

  47. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    Regularizing function
    θ( p)=∑
    j=1
    M p
    j
    p
    j

    l
    j
    β
    Similar to
    Silva Dias et al. (2009)

    View full-size slide

  48. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    Regularizing function
    θ( p)=∑
    j=1
    M p
    j
    p
    j

    l
    j
    β
    Similar to
    Silva Dias et al. (2009)
    Distance between
    jth prism and seed

    View full-size slide

  49. Well­posed problem: Minimize goal function
    Γ( p)=ϕ( p)+μθ( p)
    Regularizing function
    θ( p)=∑
    j=1
    M p
    j
    p
    j

    l
    j
    β
    Similar to
    Silva Dias et al. (2009)
    Distance between
    jth prism and seed
    Imposes:

    Compactness ●
    Concentration around seeds

    View full-size slide

  50. Constraints:
    1. Compact
    2. Concentrated around “seeds”
    Regularization
    3. Only or
    p
    j
    =0 p
    j
    =ρs
    4. of closest seed
    p
    j
    =ρs

    View full-size slide

  51. Constraints:
    1. Compact
    2. Concentrated around “seeds”
    Regularization
    Algorithm
    3. Only or
    p
    j
    =0 p
    j
    =ρs
    4. of closest seed
    p
    j
    =ρs Based on René (1986)

    View full-size slide

  52. Planting Algorithm

    View full-size slide

  53. Setup: g = observed data

    View full-size slide

  54. Setup:
    Define interpretative model
    Interpretative model
    g = observed data

    View full-size slide

  55. Setup:
    Define interpretative model
    All parameters zero
    Interpretative model
    g = observed data

    View full-size slide

  56. Setup:
    seeds
    N
    S
    Define interpretative model
    All parameters zero
    Interpretative model
    g = observed data

    View full-size slide

  57. Setup:
    seeds
    N
    S
    Define interpretative model
    All parameters zero
    Include seeds
    Prisms with
    not shown
    p
    j
    =0
    g = observed data

    View full-size slide

  58. Setup:
    seeds
    N
    S
    Define interpretative model
    All parameters zero
    Include seeds
    Compute initial residuals
    r(0)=g−d(0)
    Prisms with
    not shown
    p
    j
    =0
    g = observed data
    d = predicted data

    View full-size slide

  59. Setup:
    seeds
    N
    S
    Define interpretative model
    All parameters zero
    Include seeds
    Compute initial residuals
    r(0)=g−d(0)
    Prisms with
    not shown
    p
    j
    =0
    g = observed data
    Predicted by seeds
    d = predicted data

    View full-size slide

  60. Setup:
    seeds
    N
    S
    Define interpretative model
    All parameters zero
    Include seeds
    Compute initial residuals
    Prisms with
    not shown
    p
    j
    =0
    g = observed data
    d = predicted data
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    )

    View full-size slide

  61. seeds
    N
    S
    Define interpretative model
    All parameters zero
    Include seeds
    Compute initial residuals
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    )
    Prisms with
    not shown
    g = observed data
    d = predicted data
    Neighbors
    Find neighbors of seeds
    p
    j
    =0
    Setup:

    View full-size slide

  62. Prisms with
    not shown
    Try accretion to sth seed:
    p
    j
    =0
    Growth:

    View full-size slide

  63. Prisms with
    not shown
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    Choose neighbor:
    p
    j
    =0
    Growth:

    View full-size slide

  64. Prisms with
    not shown
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Choose neighbor:
    p
    j
    =0
    Growth:
    (New elements)
    j

    View full-size slide

  65. Prisms with
    not shown
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Choose neighbor:
    p
    j
    =0
    Growth:
    (New elements)
    j
    Update residuals
    r(new)=r(old )− p
    j
    a
    j

    View full-size slide

  66. Prisms with
    not shown
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Choose neighbor:
    p
    j
    =0
    Growth:
    (New elements)
    j
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Contribution of j

    View full-size slide

  67. Prisms with
    not shown
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Choose neighbor:
    p
    j
    =0
    Growth:
    (New elements)
    j
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Variable sizes
    None found = no accretion

    View full-size slide

  68. Prisms with
    not shown
    None found = no accretion
    N
    S
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Choose neighbor:
    p
    j
    =0
    Growth:
    (New elements)

    View full-size slide

  69. Prisms with
    not shown
    None found = no accretion
    N
    S
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Choose neighbor:
    p
    j
    =0
    Growth:
    j
    (New elements)

    View full-size slide

  70. Prisms with
    not shown
    None found = no accretion
    N
    S
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Choose neighbor:
    At least one seed grow?
    p
    j
    =0
    Growth:
    j
    (New elements)

    View full-size slide

  71. Prisms with
    not shown
    None found = no accretion
    N
    S
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Choose neighbor:
    At least one seed grow?
    Yes
    p
    j
    =0
    Growth:
    j
    (New elements)

    View full-size slide

  72. Prisms with
    not shown
    None found = no accretion
    N
    S
    Try accretion to sth seed:
    1. Reduce data misfit
    2. Smallest goal function
    p
    j

    s
    j = chosen
    Update residuals
    r(new)=r(old )− p
    j
    a
    j
    Choose neighbor:
    At least one seed grow?
    Yes No
    p
    j
    =0
    Growth:
    Done!
    j
    (New elements)

    View full-size slide

  73. Advantages:
    Compact & non­smooth
    Any number of sources
    Any number of different density contrasts
    No large equation system
    Search limited to neighbors

    View full-size slide

  74. Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρs
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector

    View full-size slide

  75. Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector
    No matrix multiplication (only vector +)

    View full-size slide

  76. No matrix multiplication (only vector +)
    Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector
    Only need some columns of A

    View full-size slide

  77. No matrix multiplication (only vector +)
    Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector
    Only need some columns of A
    Calculate only when needed

    View full-size slide

  78. No matrix multiplication (only vector +)
    Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector
    Only need some columns of A
    Calculate only when needed & delete after update

    View full-size slide

  79. No matrix multiplication (only vector +)
    Remember equations:
    r(0)=g−
    (∑
    s=1
    N
    S
    ρ
    s
    a
    j
    S
    ) r(new)=r(old)− p
    j
    a
    j
    Initial residual Update residual vector
    Only need some columns of A
    Calculate only when needed
    Lazy evaluation
    & delete after update

    View full-size slide

  80. Advantages:
    Compact & non­smooth
    Any number of sources
    Any number of different density contrasts
    No large equation system
    Search limited to neighbors

    View full-size slide

  81. Advantages:
    Compact & non­smooth
    Any number of sources
    Any number of different density contrasts
    No large equation system
    Search limited to neighbors
    No matrix multiplication (only vector +)
    Lazy evaluation of Jacobian

    View full-size slide

  82. Advantages:
    Compact & non­smooth
    Any number of sources
    Any number of different density contrasts
    No large equation system
    Search limited to neighbors
    No matrix multiplication (only vector +)
    Lazy evaluation of Jacobian
    Fast inversion + low memory usage

    View full-size slide

  83. Synthetic Data

    View full-size slide

  84. Data set:

    3 components

    51 x 51 points

    2601 points/component

    7803 measurements

    5 Eötvös noise

    View full-size slide

  85. Model: ●
    11 prisms

    View full-size slide

  86. Model: ●
    11 prisms ●
    4 outcropping

    View full-size slide

  87. Model: ●
    11 prisms ●
    4 outcropping

    View full-size slide

  88. Model: ●
    11 prisms ●
    4 outcropping

    View full-size slide


  89. Strongly interfering effects

    View full-size slide


  90. Strongly interfering effects

    What if only interested in these?

    View full-size slide


  91. Common scenario

    View full-size slide


  92. Common scenario

    May not have prior information

    Density contrast

    Approximate depth

    View full-size slide


  93. Common scenario

    May not have prior information

    Density contrast

    Approximate depth

    No way to provide seeds

    View full-size slide


  94. Common scenario

    May not have prior information

    Density contrast

    Approximate depth

    No way to provide seeds

    Difficult to isolate effect of targets

    View full-size slide

  95. Robust procedure:

    View full-size slide

  96. Robust procedure:

    Seeds only for targets

    View full-size slide

  97. Robust procedure:

    Seeds only for targets

    View full-size slide

  98. Robust procedure:

    Seeds only for targets

    ℓ1
    ­norm to “ignore” non­targeted

    View full-size slide

  99. Robust procedure:

    Seeds only for targets

    ℓ1
    ­norm to “ignore” non­targeted

    View full-size slide

  100. Inversion: ●
    13 seeds ●
    7,803 data

    View full-size slide

  101. Inversion: ●
    13 seeds ●
    7,803 data

    View full-size slide

  102. Inversion: ●
    13 seeds ●
    7,803 data

    View full-size slide

  103. Inversion: ●
    13 seeds ●
    7,803 data

    View full-size slide

  104. Inversion: ●
    13 seeds ●
    7,803 data

    View full-size slide

  105. Inversion: ●
    13 seeds ●
    7,803 data ●
    37,500 prisms

    View full-size slide

  106. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Only prisms with zero
    density contrast not shown

    View full-size slide

  107. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Only prisms with zero
    density contrast not shown

    View full-size slide

  108. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Only prisms with zero
    density contrast not shown

    View full-size slide

  109. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Only prisms with zero
    density contrast not shown

    View full-size slide

  110. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Only prisms with zero
    density contrast not shown

    View full-size slide

  111. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds

    Recover shape of targets
    Only prisms with zero
    density contrast not shown

    View full-size slide

  112. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds

    Recover shape of targets

    Total time = 2.2 minutes (on laptop) Only prisms with zero
    density contrast not shown

    View full-size slide

  113. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Predicted data in contours

    View full-size slide

  114. Inversion: ●
    7,803 data ●
    37,500 prisms

    13 seeds
    Predicted data in contours
    Effect of true targeted sources

    View full-size slide

  115. Data:

    3 components

    FTG survey

    Quadrilátero Ferrífero, Brazil

    View full-size slide

  116. Data:

    3 components

    FTG survey

    Quadrilátero Ferrífero, Brazil
    Targets:

    Iron ore bodies

    BIFs of Cauê Formation

    View full-size slide

  117. Data:

    3 components

    FTG survey

    Quadrilátero Ferrífero, Brazil
    Targets:

    Iron ore bodies

    BIFs of Cauê Formation

    View full-size slide

  118. Data:
    Seeds for iron ore:

    Density contrast 1.0 g/cm3

    Depth 200 m

    View full-size slide

  119. Inversion:
    Observed
    Predicted

    46 seeds ●
    13,746 data

    View full-size slide

  120. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  121. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  122. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  123. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  124. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  125. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  126. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms
    Only prisms with zero
    density contrast not shown

    View full-size slide

  127. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms
    Only prisms with zero
    density contrast not shown

    View full-size slide

  128. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms
    Only prisms with zero
    density contrast not shown

    View full-size slide

  129. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  130. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    View full-size slide

  131. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    Agree with previous interpretations
    (Martinez et al., 2010)

    View full-size slide

  132. Inversion: ●
    46 seeds ●
    13,746 data ●
    164,892 prisms

    Agree with previous interpretations
    (Martinez et al., 2010)

    Total time = 14 minutes
    (on laptop)

    View full-size slide


  133. New 3D gravity gradient inversion

    Multiple sources

    Interfering gravitational effects

    Non­targeted sources

    No matrix multiplications

    No linear systems

    Lazy evaluation of Jacobian matrix
    Conclusions

    View full-size slide


  134. Estimates geometry

    Given density contrasts

    Ideal for:

    Sharp contacts

    Well­constrained physical properties
    – Ore bodies
    – Intrusive rocks
    – Salt domes
    Conclusions

    View full-size slide

  135. Mestrado

    2010: Cumprir disciplinas

    2010­2011: 7 trabalhos em congresso (5 primeiro autor)

    10/2011: Submeter artigo para Geophysics

    11/2011: Defesa da dissertação de mestrado
    Continuação

    Adaptar para gravimetria e magnetometria

    Disponibilizar software Open Source
    Cronograma

    View full-size slide