TBMI 02 ¨ Linkoping University Mats Andersson Medical image analysis TBMI02 Mats Andersson Mini project: Adaptive filtering Mini project Mats Andersson teacher for: TBMI 02 ¨ Linkoping University - Filter optimization(2D). ([email protected]), - Tensor representation of local structure. your - Adaptive image enhancement - Mini project - MRI lab 2 classes 3 laborations Lab-PM So what is this mini project about? Biomedical Engineering November 14, 2014 Page 1 TBMI 02 ¨ Linkoping University Mats Andersson November 14, 2014 ¨ Linkoping University 2D Enhancement of MR-image Biomedical Engineering TBMI 02 Page 3 Mats Andersson Mini project: Examination Result from last years students 1. Demonstrate working code during lab. 2. Written report. - Answers to the questions in the lab-PM. - Lots of images of filters and intermediate results. - Motivate you settings and parameters. Last date to hand in the report Dec 1 Test a lot at the laborations - plenty of time. Work in groups of one or two people. Enhanced image, Original image. November 14, 2014 Biomedical Engineering Page 2 November 14, 2014 Biomedical Engineering Page 4 TBMI 02 ¨ Linkoping University Mats Andersson TBMI 02 ¨ Linkoping University Mini project: Literature Mats Andersson Mini project: Filter optimization Granlund, Knutsson: Signal Processing for Computer Vision Filter optimization as a LS-problem using ideal and weight functions. ch 4 the Fourier transform. • Fi (u) – frequency ideal function. ch 5 Filter optimization + Advanced filter opt paper, sec 1-5. • Fw (u) – frequency weight function. ch 6 Orientation and tensor parts. • fi (x) – spatial ideal function. ch 10 Adaptive filtering. • fw (x) – spatial weight function. Note lower case letters for SD and upper case letters for FD. 2 = Fw (Fi − B f ) 2 + α fw (fi − f ) 2 ∂ 2 =0 ∂f Use matlab script krnopt Biomedical Engineering November 14, 2014 TBMI 02 ¨ Linkoping University Page 5 Mats Andersson Biomedical Engineering November 14, 2014 ¨ Linkoping University Mini project: Filter optimization TBMI 02 Page 7 Mats Andersson Mini project: Filter optimization Filter optimization wish list: Fi (u) – ideal function in the FD. Use e.g the matlab scripts quadrature, enhlp, enhbp, . . . to generate different Fi (u). - Minimal deviation from an ideal frequency function (e.g. lognormal quadrature). - Small spatial support. - Adapted to the image spectrum. Fw (u) – weight function in FD. The spectrum of the image defines the importance of a close fit to Fi (u). Use F weightgensnr to generate Fw (u). - Adapted to the noise spectrum. Conflicting demands that relates back to the uncertainty principle. Some demands are easy to formulate in the Fourier domain (FD). Others are best formulated in the spatial domain (SD). November 14, 2014 Biomedical Engineering Page 6 November 14, 2014 Biomedical Engineering Page 8 TBMI 02 ¨ Linkoping University Mats Andersson Mini project: Filter optimization Mats Andersson Mini project: Filter optimization Coordinates and sampling in FD fi (x) – ideal function SD. Not F −1 (Fi (u))! 1 fi (x) = 0 TBMI 02 ¨ Linkoping University n−1 2 x=0 x=0 f (x) e−i u x F (u) x=− n−1 2 The most compact filter there is. f (x) is a discrete function. F (u) is continuous and periodic (2 π ). fw (x) – weight function in SD. Preserve spatial locality. fw (x) ≈ x2 . How to sample used? Use goodsw to generate fw (x). F (u) and how many samples, N , should be KRNOPT VIDEO Biomedical Engineering November 14, 2014 Page 9 TBMI 02 ¨ Linkoping University Mats Andersson Always use an odd number n • bandwidth β = Sample F (u) using FFT coordinates (N is even) = 7, 9, 11, . . . u∈ − √ ln(uh /ul ) [Octaves] ln 2 uo corresponds to the wave length λ0 = ≥ 1.5 − 2 λ0 pixels. ex :N =4 Sample F (u) using KRNOPT coordinates (N is odd) 2π u0 u0 π/2 π/4 N N 2π : −1 2 2 N • One sample in the origin is good! • Impossible to have different values close to π and −π uh ul • where uh and ul are 6dB cut off frequencies. Let n Mats Andersson Mini project: Filter optimization × n pixels Ex. lognormal quadrature filter uo , β • center frequency uo = Page 11 TBMI 02 ¨ Linkoping University Mini project: Filter optimization Decide spatial support size n Biomedical Engineering November 14, 2014 u∈ − λ0 4 8 n 7, 9 15 N − 1 N − 1 2π : 2 2 N ex :N =5 • One sample in the origin is good! • symmetric for values close to π and −π LIFT: Local Interactive FT Demo November 14, 2014 Biomedical Engineering Page 10 November 14, 2014 Biomedical Engineering Page 12 TBMI 02 ¨ Linkoping University TBMI 02 ¨ Linkoping University Mats Andersson Mini project: Filter optimization Mats Andersson Mini project: Filter optimization Filter optimization hints How many samples N to use in the FD? • Use no spatial weight first time We sample F (u) and measure the error F (u) − Fi (u) using N sample points. For small N (N = n) the error will be small at the sample points but may be very large in between! Use as large N as possible • Compute error using filterr F (F −F ) w i Most important FD-error δ2 = Fw Fi DC-error (weighted and unweighted) and • plot filters in both domains, use wamesh N ≥ 2 n + 1 (in each dimension) - No ripple in the borders of SD kernel - Note shape and magnitude in the FD SD: 15 × 15 pixels → FD: 31 × 31 samples in ] − π, π[ • Increase spatial weight gradually and observe the error and the filter shapes. Biomedical Engineering November 14, 2014 Page 13 TBMI 02 ¨ Linkoping University Biomedical Engineering November 14, 2014 TBMI 02 ¨ Linkoping University Mats Andersson Mini project: Filter optimization Mini project: Local structure krnopt matlab script Frequency ideal filter Page 15 Mats Andersson [f, F ] = krnopt(Fi , Fw , fm , fi , fw ) FD SD 1 All in and outputs are wa (weight arrays) wa = 0.8 Data: multidim array FD/SD-flag coordinate system (SD: pixels – Defined in polar separable coordinates. 0.6 – How many filters do we need? 0.4 FD: ] − π, π[ ) – How should they be oriented? 0.2 0 π – Phase invariant magnitude? π π/2 Most element wise operations work for wa use .∗ and ./ π/2 0 0 −π/2 y−axis Many functions for generation of weight and ideal functions exist in the kerngen library e.g. quadrature, F weightgensnr, goodsw, see labPM. November 14, 2014 Biomedical Engineering Page 14 −π/2 −π −π Quadrature filter u0 November 14, 2014 x−axis = π/3, β = 2 octaves. Biomedical Engineering Page 16 TBMI 02 ¨ Linkoping University Mats Andersson TBMI 02 ¨ Linkoping University Mini project: Local structure Mats Andersson Mini project: Local structure Ideal enhancement filters 1 1 0.8 0.8 0.6 0.6 0.4 0.4 0.2 0.2 0 π 0 π π π/2 −π/2 −π 0 −π/2 −π/2 −π π/2 0 0 y−axis π π/2 π/2 0 y−axis x−axis −π/2 −π −π x−axis The LP and one of the broadband BP-filters in the enhancement filterset. Testimage ploop Use the same number an orientations for the BP-filters as for the quadrature filters. Biomedical Engineering November 14, 2014 ¨ Linkoping University TBMI 02 Page 17 Mats Andersson Biomedical Engineering November 14, 2014 TBMI 02 ¨ Linkoping University Mini project: Local structure Page 19 Mats Andersson Mini project: Local structure Filter response from a quadrature filter Compute quadrature filter responses using convolution with edge extraction to avoid large artifacts at the image border. q = conv ext2(myimage, f ) q = imfilter(myimage, getdata(f ), ’replicate’) Real (even) part. November 14, 2014 Biomedical Engineering Page 18 November 14, 2014 Odd (imaginary) part. Biomedical Engineering Page 20 TBMI 02 ¨ Linkoping University Mats Andersson Mini project: Local structure t1 t2 Mats Andersson t2 t3 k=1 t3 red – negative Biomedical Engineering November 14, 2014 Page 23 TBMI 02 ¨ Linkoping University T = tn – images M k – 2 × 2 matrix (2D) qk M k k 4 1 n ˆkn ˆ Tk − I Nf Nf Mats Andersson A 2nd order tensor implies an ellipse model of the local energy contribution in an image Nf = Mk = Mini project: Local structure Local tensor T : symmetric, positive definite matrix. t1 t2 qk M k t2 Mini project: Local structure T = = green – positive Page 21 TBMI 02 ¨ Linkoping University t2 t3 t1 q Biomedical Engineering November 14, 2014 Mats Andersson Mini project: Local structure T = Magnitude response from a quadrature filter TBMI 02 ¨ Linkoping University t1 t2 t2 t3 = λ1 e ˆ1 e ˆT1 + λ2 e ˆ2 e ˆT2 eigen values λ1 ≥ λ2 ≥ 0 eigen vectors e ˆ1 ⊥ e ˆ2 • The tensor has orientation but no direction, sign of e ˆ1 and e ˆ2 is arbitrary. Projection tensors M k 4 1 Mk = n ˆkn ˆ Tk − I Nf Nf November 14, 2014 • Noise makes the tensor more isotropic. where nk I Nf – filter orienting vectors – isotropic term – number of quadrature filters Nf Biomedical Engineering • For an isotropic tensor the eigen vectors can be chosen arbitrary. ≥3 Page 22 November 14, 2014 Biomedical Engineering Page 24 TBMI 02 ¨ Linkoping University Mats Andersson TBMI 02 ¨ Linkoping University Mini project: Local structure Mats Andersson Visualization of local structure T = Why are tensors superior to vectors for representing orientation? t1 t2 t2 t3 = 1 2 1 − cos(2ϕ) sin(2ϕ) sin(2ϕ) 1 + cos(2ϕ) Compute G(ϕ) = t3 − t1 + i 2 t2 = cos(2ϕ) + i sin(2ϕ) A computer can look at three tensor components simultanously but it is hard for us humans. Compute one complex valued image from T = t1 t2 t2 t3 G(ϕ) = G(ϕ + π) continuous representation! Note that G(ϕ) contains less information compared to T . Biomedical Engineering November 14, 2014 TBMI 02 ¨ Linkoping University Page 25 Mats Andersson Visualization of local structure Biomedical Engineering November 14, 2014 ¨ Linkoping University TBMI 02 Page 27 Mats Andersson Orientation estimation from tensor components Compute one orientation image from the 3 tensor components. Simple case λ1 =1 λ2 = 0 e1 = T = λ1 e1 e1 T = T = 1 2 sin(ϕ) cos(ϕ) sin2 (ϕ) sin(ϕ) cos(ϕ) sin(ϕ) cos(ϕ) cos2 (ϕ) 1 − cos(2ϕ) sin(2ϕ) sin(2ϕ) 1 + cos(2ϕ) G(ϕ) = t3 − t1 + i 2 t2 = cos(2ϕ) + i sin(2ϕ) November 14, 2014 Biomedical Engineering Page 26 November 14, 2014 Biomedical Engineering Page 28 TBMI 02 ¨ Linkoping University Mats Andersson Mini project: Adaptive filtering t1 t2 t2 t3 C= c1 c2 c2 c3 Mats Andersson Mini project: Adaptive filtering Compute the Control tensor C T = TBMI 02 ¨ Linkoping University T T We need to compute λ1 , λ2 and e ˆ1 e ˆ1 , e ˆ2 e ˆ2 Solve the secular equation: = λ1 e ˆ1 e ˆT1 + λ2 e ˆ2 e ˆT2 λ1 = 1 2 t1 + t3 + (t1 − t3 )2 + 4 t22 λ2 = 1 2 t1 + t3 − (t1 − t3 )2 + 4 t22 = γ1 e ˆ1 e ˆT1 + γ2 e ˆ2 e ˆT2 • The Local tensor and the control tensor share the same eigensystem. Use λ2 • γ1 and γ2 are mappings of λ1 and λ2 . Parallel computation, never loop over pixels in matlab! Biomedical Engineering November 14, 2014 ¨ Linkoping University Page 29 TBMI 02 Mats Andersson Biomedical Engineering November 14, 2014 Mats Andersson Mini project: Adaptive filtering C = γ1 e ˆ1 e ˆT1 + γ2 e ˆ2 e ˆT2 T T Compute e ˆ1 e ˆ1 and e ˆ2 e ˆ2 T = λ1 e ˆ1 e ˆT1 + λ2 e ˆ2 e ˆT2 γ1 = 1 γ2 = 0 Page 31 TBMI 02 ¨ Linkoping University 2D adaptive filters in Fourier domain Control tensor = |λ2 | to ensure a positive definite tensor. I=e ˆ1 e ˆT1 + e ˆ2 e ˆT2 T = (λ1 − λ2 ) e ˆ1 e ˆT1 + λ2 I γ1 = 1 γ2 = 1 (λ1 − λ2 ) e ˆ1 e ˆT1 = T − λ2 I Normalize, do not divide! γ1 = γ2 = 0 T Finally e ˆ2 e ˆ2 γ1 = 0.5 γ2 = 0 November 14, 2014 Biomedical Engineering γ1 = 0.5 γ2 = 0.5 Page 30 November 14, 2014 =I −e ˆ1 e ˆT1 Biomedical Engineering Page 32 TBMI 02 ¨ Linkoping University Mats Andersson The control tensor C , mapping of local energy T = λ21 + λ22 x = x = x/ max(x) xβ xβ+α + σ β γ1 = m(x, σ, α, β, j) = TBMI 02 ¨ Linkoping University Mats Andersson The control tensor C , anisotropy 1/j λ2 x= λ1 µ2 (x, α, β, j) = (x (1 − α))β (x (1 − α))β + (α (1 − x))β m−func s=0.1 a=0.3 b=3.0 j=1.0 1/j γ2 = µ2 γ1 mu−func a=0.4 b=3.0 j=1.0 1.4 1 1.2 0.8 1 0.8 0.6 0.6 0.4 0.4 0.2 0.2 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 0 1 Biomedical Engineering November 14, 2014 ¨ Linkoping University TBMI 02 0 0.1 anisotropic (λ2 m func example Page 33 Mats Andersson 0.2 λ1 ) 0.3 0.4 0.5 0.6 mu func example 0.7 0.8 0.9 isotropic (λ2 1 = λ1 ) Biomedical Engineering November 14, 2014 ¨ Linkoping University TBMI 02 Page 35 Mats Andersson The control tensor C , mapping of γ1 The control tensor C , mapping of µ2 Tune the m func parameters such that γ1 is large where the image has a well defined orientation and dark in noisy areas without structure. Tune the mu func parameters such that µ2 is dark when the orientation is well defined and bright at the isotropic areas. Note that µ2 only depends on the shape of the tensor (and not the energy). November 14, 2014 Biomedical Engineering Page 34 November 14, 2014 Biomedical Engineering Page 36 TBMI 02 ¨ Linkoping University Mats Andersson C= c1 c2 4 1 n ˆkn ˆ Tk − I= Mk = Nf Nf c2 c3 mk1 mk2 Mats Andersson Iterative adaptive filtering Computation of filter weights for enhancement filters The control tensor TBMI 02 ¨ Linkoping University • Apply the enhancement filters multiple times • Use the original control tensor do not recompute C mk2 mk3 • Let j in the m-map and µ-map correspond to the number of iterations. Weight for enhancement filter k , wk m= wk = C • M k = c1 mk1 + 2 c2 mk2 + c3 mk3 xβ β+α x + σβ 1/j µ2 = (x (1 − α))β (x (1 − α))β + (α (1 − x))β 1/j Nf enhanced image: enh = enhlp + – Bring your own images or use tha camera in the lab wk enhbp (k) k=1 Biomedical Engineering November 14, 2014 Page 37 TBMI 02 ¨ Linkoping University Mats Andersson Signal to noise relation (SNR) for images Compute standard deviation for the image and the noise im noise σim σnoise (σploop = 0.35) A testimage with a desired SNR is computed as: noiseim = α (im − im) + (1 − α) (noise − noise) + im where: SNR = 20 log10 α σim (1 − α) σnoise [dB] Compute α that corresponds to 5dB and 20dB. November 14, 2014 Biomedical Engineering Page 38 November 14, 2014 Biomedical Engineering Page 39
© Copyright 2024 ExpyDoc
