Cascadic Multilevel Methods for Cascadic Multilevel Methods for Large-Scale Ill-Posed Problems Large-Scale Ill-Posed Problems Fiorella Sgallari University of Bologna , Italy Faculty of Engineering Department of Mathematics - CIRAM Joint work with Serena Morigi, Lothar Reichel, Andriy Shyshkov GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 GAMM 08 nd 2008 1
Outline 1. Inverse and ill-posed problems in image processing 2. Linear and nonlinear regularization methods 3. Cascadic multilevel method for image deblurring 4. Examples CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 2
Inverse problem in image processing Inverse problems arise when one seeks to determine the cause of an observed effect Image restoration: Determine the unavailable exact image from an available contaminated version Ill-posed problems (in the sense of Hadamard) CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 3
Ill-posed problem in image processing U is the representation of the real world into the image plane F is the observed corrupted image Let ϕ be the degradation model function F = ϕ (U) ∀ϕ (U) = U + noise ∀ϕ (U) = K * U + noise •…………………….. U = ϕ -1 (F) CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 4
Different type of noise Original salt & pepper speckle Gauss-add Gauss-molt CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008
Degradation model Continuous degradation model: ∫∫ = + η ∈ Ω ˆ f x ( ) h x y u y dy ( , ) ( ) ( ) x x Ω Data noise Perturbed observed image Blur and noise-free image Point Spread Function Integral equation can be expressed as = ∗ + η ˆ f h u Discretization yields b = Au with matrix A block Toeplitz with Toeplitz blocks CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 6
Linear discrete ill-conditioned problems Au = b arise from the discretization of linear ill-posed problems • The matrix A is of ill-determined rank, possibly singular. A may be indefinite. • The right-hand side b represents available data that generally is contaminated by an error. CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 7
Linear discrete ill-conditioned problems Available contaminated, possibly inconsistents, linear system = Au b (1) Unavailable associated consistent linear system with error-free right-hand side ˆ = Au b (2) ˆ u Let denote the desired solution of (2), e.g., the minimal-norm solution. Task: Determine an approximate solution of (1) that is a good ˆ u approximation of CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 8
u=A -1 b Shaw.m Solution Au=b : add 0.1% noise to rhs ˆ = + b b e ˆ ˆ = − + = − + − = + − 1 1 1 1 ˆ u A ( b e ) A b A e u A e u=A -1 b CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 9
Popular regularization methods 1. Tikhonov regularization 2. Iterative regularization methods: regularization by truncated iteration CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 10
Tikhonov regularization Solve the Tikhonov minimization problem { } 2 2 + µ min Au-b Lu , u where ∈ mxn A R ; ∈ ≤ pxn L R , p n, is the regularization operator. Common choices: L = I or a finite difference operator ; μ > 0 is the regularization parameter to be determined. It is important to determine a suitable value; see Engl, Hanke, Neubauer. CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 11
Regularization by truncated iteration = A Au * A b * CGLS: CG applied to Define the Krylov subspace { K ( A * A A , * b ) span A * ,( b A * A A ) * ,.., b = k } ( A * A ) A * ,( b A * A ) A * b k 2 k 1 − − Then u K ( A * A A , * b ) and ∈ k k Minimum-Residual Method Au b min Au b − = − k u K ( A A A b * , * ) ∈ k k Therefore discrepancy d b Au satisfies = − j j b d ... d ≥ ≥ ≥ 1 k CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 12
Stopping Criterion Discrepancy principle ˆ = − = δ e b b Let α > 1 be fixed, The iterate u k satisfies the discrepancy principle if − ≤ α δ Au b k (Stopping rule) Terminate the iterations as soon as iterate u k satisfies Au b − ≤ α δ k Au b − > α δ k 1 − Denote the termination index by k δ CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 13
----- exact solution Hansen, 2005 CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 14
An iterative method is a regularization method if − = ˆ limsup u u 0 k δ δ → 0 ≤ δ e Minimum-residual methods and their Krylov subspaces Matrix Algorithm Krylov subspace κ ( A,b ) Symmetric MINRES k κ MR-II ( A, Ab ) k κ ( A,b ) Nonsymmetric GMRES k κ and square RRGMRES ( A, Ab ) k Any CGLS, LSQR κ T T ( A A, A b ) k [ see Nemirovskii,Hanke,Reichel] CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 15
CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 16
CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 17
Nonlinear denoising-deblurring methods Return to Tikhonov regularization: ∫ 2 − + µ min ( h* u f ) R( u )dx u Ω Euler-Lagrange equations with gradient descent: ∂ u = − − + µ = h* ( h* u f ) D( u ), u f 0 ∂ t 2 = ∇ ⇒ = ∆ R( u ) u D( u ) u Tikhonov ∇ u = ∇ ⇒ = R( u ) u D( u ) div( ) TV ∇ u CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 18
Edge-enhanched regularization Tikhonov regularization TV regularization CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 19
Edge-enhanched regularization Tikhonov regularization TV regularization CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 20
Perona-Malik diffusivity 2 = ∇ ∇ = + κ ρ > D( u ) div( g( u ) u ) g( s ) 1 /( 1 s / ), 0 Perona-Malik denoising model g(s) ∂ u 2 − ∇ • ∇ ∇ = ( g( u ) u ) 0 ρ ∂ t s [ Perona Malik ‘87, Cattè Lions Morel ‘92] CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 21
Computational effort ∂ u 2 = − − + µ ∇ ⋅ ∇ ∇ = h* ( h* u f ) ( g( u ) u ), u f 0 ∂ t = u 0 1 = g( s ) + s / κ 1 Space discretization Explicit time schemes: •Easy implementation •CFL stability conditions (small step, many iterations) ⇒ method expensive Semi-implicit time schemes: •Solve a linear system at each time step •Stability conditions due to the blur term ⇒ method expensive CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 22
Combine: Computational efficiency of truncated iteration for linear systems Edge-preserving property of nonlinear models = Au b CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 23
Multilevel methods n nested linear subspaces ⊂ ⊂ ⊂ ⊂ ℜ W W ... W 1 2 k W 1 coarsest level W 1 finest level = ℜ n W k = dimW n j j = n 4 n + j 1 j W k CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 24
Multilevel methods : Res triction Operator n ℜ → R W i i → M :W W W 1 i+1 i + i 1 = = b M b R M M L M + + + + i i 1 i 1 i i 1 i 2 k * = ≤ < A R AR 1 i k i i i • averaging (replaces groups of four adjacent by one pixel, whose value is the average of the values of the four pixels it replaces) • local weighted least-squares W k approximation CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 25
Multilevel methods Prolongation Oper ato r P i 1 <i k → ≤ P :W W W 1 i-1 i i • Piecewise Linear Interpolation • Nonlinear prolongation W k CIRAM GAMM 08 Hamburg- September 12 Hamburg- September 12 nd 2008 University of Bologna GAMM 08 nd 2008 26
Recommend
More recommend
