Local Statistical Filtering via Domain Dissection for Medical - - PowerPoint PPT Presentation

Local Statistical Filtering via Domain Dissection for Medical Imaging

GTC 2016 – San Jose, CA, USA

Alexandros-Stavros Iliopoulos1 Dimitris Floros2 Nikos Pitsianis2,1 Xiaobai Sun1 Fang-Fang Yin3 Lei Ren3

1Department of Computer Science, Duke University 2Department of Electrical and Computer Engineering, Aristotle University of Thessaloniki 3Department of Radiation Oncology, Duke University School of Medicine

April 6, 2016

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 1 / 33

1 Spatially variant signal-noise analysis: motivation needs & challenges LA-SAS contribution 2 Locally adaptive signal-noise analysis formulation example filters analytic advance & technical challenges 3 LA-SAS: design & development design principle: multi-layer configuration domain dissection: local adaptivity & global concurrency CUDA LA-SAS experimental results 4 Recap & discussion 5 References

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 2 / 33

1 Spatially variant signal-noise analysis: motivation needs & challenges LA-SAS contribution 2 Locally adaptive signal-noise analysis formulation example filters analytic advance & technical challenges 3 LA-SAS: design & development design principle: multi-layer configuration domain dissection: local adaptivity & global concurrency CUDA LA-SAS experimental results 4 Recap & discussion 5 References

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 3 / 33

Spatially variant signal-noise analysis: needs & challenges

• Noise is prevalent in medical images

– multiple sources (acquisition, processing, ...) – multiple types (Gaussian, Poisson, scatter, ...)

• Noise study: characterization & suppression

– critical to high-fidelity analysis (noise propagation in processing pipeline e.g. gradient calculation) – need effective tools for systematic investigation – speed important for on-board imaging applications

• Challenging conditions

– valuable low-contrast content (especially in CT) – acquisition constraints (resolution, imaging dose) – motion: nonlinear intensity-deformation relationship ⋆ spatial variance (w.r.t. material, density, acquisition set-up)

pelvis cone-beam OBI (125 kV, coronal projection) with spatially variant scattering

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 4 / 33

Spatially variant signal-noise analysis: contribution

• LA-SAS

(locally adaptive signal-noise analysis system)

– revealing local noise statistics and signal structure – filtering in adaptation to local structures – enabling effective noise suppression

• LA-SAS design and development

– basic operations – versatile filter composition – CUDA LA-SAS (efficiency)

Range bins 0.5 1 1.5 2 2.5 3 3.5 Frequency 1000 2000 3000 4000 5000 6000 7000 8000 9000

Range bins 1.6 1.8 2 2.2 2.4 2.6 2.8 3 Frequency 500 1000 1500 2000 2500 3000 3500 Range bins 1.9 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 Frequency 50 100 150 200 250 300 350 400 450 500

range histograms: global region (top)

• vs. nested sub-regions (bottom)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 5 / 33

1 Spatially variant signal-noise analysis: motivation needs & challenges LA-SAS contribution 2 Locally adaptive signal-noise analysis formulation example filters analytic advance & technical challenges 3 LA-SAS: design & development design principle: multi-layer configuration domain dissection: local adaptivity & global concurrency CUDA LA-SAS experimental results 4 Recap & discussion 5 References

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 6 / 33

Locally adaptive signal-noise analysis: problem description

• Dual task of an analysis/filtering mechanism (F)

I(x) = ˆ I(x) + η(x), x ∈ Ω ⊂ RD

– detect/reconstruct unknown signal, ˆ I – estimate/suppress unknown noise, η

• Adaptation to local variation

ˆ I(x) :=

∑︂

x′∈𝒪(x)

F

⎞

x′, I(x′); p𝒪(x)

)︂ – based on local statistics, pN(x), over spatial neighborhood, N(x) (mean, median, deviation, range distribution, etc) – preserving signal structure (smooth subregions, discontinuities at region boundaries, etc)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 7 / 33

Locally adaptive filtering example: median

ˆ I(x) = p𝒪(x) = median

𝒪(x) {I(x)}

• basic denoising &

processing sub-module

(regional dynamic range) median Ąlter output (5 × 5) residual image

Chung et al. NSS/MIC, 2010 Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 8 / 33

Locally adaptive filtering example: entropy

p𝒪(x) = Pr𝒪(x){I(x)} H(x) = −

∑︂

p𝒪(x) log(p𝒪(x))

• multimodal registration
• basic step for other

processing modules (e.g. segmentation, histogram equalization)

(regional dynamic range) local entropy map (9 × 9)

Zhang et al. ICBBE, 2008 Pluim et al. IEEE TMI (22), 2003 Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 9 / 33

Locally adaptive filtering example: histogram equalization (HE)

p𝒪(x) = hist[r, I(N(x))] where hist: local histogram r: quantized ranges

• local contrast

enhancement

• local + global

distribution information

(global dynamic range) global HE local HE (adapthisteq) to be replaced with

• verlapping LHE

Zhu et al. CVIA (73), 1999 Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 10 / 33

Locally adaptive filtering example: bilateral filter (BF)

p𝒪(x) = σr(x) ks(x, x′) = e

⊗ ‖x⊗x′‖2

2 σ2 s

kr(I(x), I(x′)) = e

⊗ ‖I(x)⊗I(x′)‖2

2 σ2 r (x)

(space- and range-kernels)

• boundary-preserving

denoising

• local adaptation

to boundary “jumps”

(regional dynamic range) global BF σs = 1.5, σr = 0.157 locally adaptive BF σs = 1.5, σr(x) = mad(𝒪9×9(x))

Tomasi & Manduchi. ICCV, 1998 Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 11 / 33

Locally adaptive filtering example: bilateral filter (BF)

p𝒪(x) = σr(x) ks(x, x′) = e

⊗ ‖x⊗x′‖2

2 σ2 s

kr(I(x), I(x′)) = e

⊗ ‖I(x)⊗I(x′)‖2

2 σ2 r (x)

(space- and range-kernels)

• boundary-preserving

denoising

• local adaptation

to boundary “jumps”

(regional dynamic range) global BF σs = 1.5, σr = 0.157 (residual image) locally adaptive BF σs = 1.5, σr(x) = mad(𝒪9×9(x)) (residual image)

Tomasi & Manduchi. ICCV, 1998 Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 11 / 33

Locally adaptive filtering: analytic advance & technical challenges

• Reveal and preserve spatially variant signal

structure (same as in conventional methods with spatial adaptivity)

• Permit spatially inhomogeneous noise behavior

(often observed in medical imaging)

• Depart from filtering algorithms with

predetermined, global parameters (including histograms and some bilateral filters)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 12 / 33

Locally adaptive filtering: analytic advance & technical challenges

• Reveal and preserve spatially variant signal

structure (same as in conventional methods with spatial adaptivity)

• Permit spatially inhomogeneous noise behavior

(often observed in medical imaging)

• Depart from filtering algorithms with

predetermined, global parameters (including histograms and some bilateral filters)

• Challenge traditional parallel primitives in

multiple aspects (algorithmic complexity, concurrency, numerical behavior)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 12 / 33

Locally adaptive filtering: analytic advance & technical challenges

• Reveal and preserve spatially variant signal

structure (same as in conventional methods with spatial adaptivity)

• Permit spatially inhomogeneous noise behavior

(often observed in medical imaging)

• Depart from filtering algorithms with

predetermined, global parameters (including histograms and some bilateral filters)

• Challenge traditional parallel primitives in

multiple aspects (algorithmic complexity, concurrency, numerical behavior)

• Fully parallel but highly redundant
• Work-efficient may become

– highly sequential – highly divergent – numerically unstable

• Compromises:

increase redundancy to stabilize numerics and expose parallelism

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 12 / 33

1 Spatially variant signal-noise analysis: motivation needs & challenges LA-SAS contribution 2 Locally adaptive signal-noise analysis formulation example filters analytic advance & technical challenges 3 LA-SAS: design & development design principle: multi-layer configuration domain dissection: local adaptivity & global concurrency CUDA LA-SAS experimental results 4 Recap & discussion 5 References

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 13 / 33

LA-SAS design principles

• Multi-layer configuration between

– filtering at the application level (top) – parallel computing at the architecture level (bottom)

(not one-to-all cross-bar configuration)

• Abstraction of basic LA-SAS operations &

versatile filter composition

• Enabling systematic investigation &

efficient implementation

• Domain dissection
• Locality adaptation

– Local range per tile

• Redundancy reduction

– Single-pass precomputation – Common sub-regions

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 14 / 33

LA-SAS MATLAB Syntax and Pyramid Data Structure

locStat = local

∏︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⨄︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋃︂

Mean Std Min Max Median Mad Hist Entropy

∫︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋀︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋁︂ ⋂︂ ⋃︁ ⋁︁ ⋁︁ ⨄︁

(Image, b) Tile(Image, b, t) Tile1x1(Image, b)

⋂︁ ⎥ ⎥ ⋀︁

localStat image integral integral2 minLocal maxLocal integralHistogram min max tile{1} size integral integral2 minLocal maxLocal integralHistogram min max tile{2} size etc Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 15 / 33

LA-SAS Dependency Graph and Arithmetic Complexity

Entropy : 2b Median : 14w Min Max : 3 Local Histogram : 3b Integral x2 : 3 Mean : 4 Integral : 2 StD : 3 MAD : 14w √︂ x2 : 3 Integral Histogram : 2b

• perations per image pixel

w : number of pixels in a local neighborhood b : number of histogram bins

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 16 / 33

LA-SAS primitives: local grid for local histograms in constant time

(D + 1)-dimensional grid embedding, with local dynamic range bins Histograms in constant time from integral histogram representation (Poostchi et al, ACCVW, 2012))

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 17 / 33

Domain dissection: local adaptivity & global concurrency

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 18 / 33

Domain dissection: local adaptivity & global concurrency

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 18 / 33

Domain dissection: local adaptivity & global concurrency

Sliding window operation

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 19 / 33

M + 2bM

LA-SAS Implementation

• Utilize MATLAB built-ins
• Data in gpuArray
• CUDA kernels for the rest
• Compute all you can in each pass
• Cache computed results

in cell array of struct to reuse (memoization)

• User responsible for caching

and data motion to/from GPU

• Three CUDA kernels:

– localStats – histPrefixSum – localMinMax

• Tiles mapped on grid
• Thread blocks working locally
• Common neighborhood shared

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 20 / 33

Catphan Phantom

Catphan 504 phantom by Varian Medical Systems CTP515 low contrast module with supra-slice and subslice contrast targetss CTP528 High resolution module with 21 line pair per cm gauge and point source Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 21 / 33

Local statistics exploration

Catphan cone-beam 125kV X-ray (θ = 21o)

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Normalized intensity range bins (r)

0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5

• p r log 2p r

global and ROI-specific 32-bin Shannon distributions

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

Catphan cone-beam projection example histogram windows (7 × 7)

0.53 0.54 0.55 0.56 0.57 0.58 5 0.94 0.95 0.96 0.97 10 20 0.121 0.122 0.123 0.124 0.125 0.126 5 10

(7 × 7)

0.5 0.52 0.54 0.56 0.58 0.6 10 20 0.92 0.94 0.96 0.98 1 20 40 0.12 0.122 0.124 0.126 20 40

(15 × 15)

example 32-bin local histograms

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 mean 15 × 15 mean

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 standard deviation 15 × 15 standard deviation

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 median 15 × 15 median

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 median absolute deviation (MAD) 15 × 15 median absolute deviation (MAD)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 dynamic range (log-scale) 15 × 15 dynamic range (log-scale)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 entropy via locally adaptive 32-bin histograms 15 × 15 entropy via locally adaptive 32-bin histograms

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

7 × 7 entropy via uniform 256-bin histogram 15 × 15 entropy via uniform 256-bin histogram

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration

global SVD filtering (50 out of 384 components) global SVD residual (50 out of 384 components)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 22 / 33

Local statistics exploration: ROI detail

Catphan cone-beam 125kV X-ray (θ = 21o)

0.42 0.44 0.46 0.48 0.5 0.52 0.54 0.56 0.58 0.6

Normalized intensity range bins

200 400 600 800 1000 1200 1400 1600

Pixel frequency

global and ROI-specific 32-bin histograms

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 mean 15 × 15 mean

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 standard deviation 15 × 15 standard deviation

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 median 15 × 15 median

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 median absolute deviation (MAD) 15 × 15 median absolute deviation (MAD)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 dynamic range (log-scale) 15 × 15 dynamic range (log-scale)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 entropy via locally adaptive 32-bin histograms 15 × 15 entropy via locally adaptive 32-bin histograms

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

7 × 7 entropy via uniform 256-bin histogram 15 × 15 entropy via uniform 256-bin histogram

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Local statistics exploration: ROI detail

regional SVD filtering (20 out of 80 components) regional SVD residual (20 out of 80 components)

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 23 / 33

Accuracy Results

CPU: 4 x AMD Opteron(tm) Processor 6376 @ 2.3 GHz (4 x 16 cores), 128 GB DDR3 GPU: NVIDIA Tesla K20c, 13 SMs @ 0.7 GHz (13 x 192 = 2496 cores), 5 GB GDDR5

Single precision error in Mean and STD for 1000x1000 image

Tile Size

10 1 10 2 10 3

Error

10 -3 10 -2 10 -1 10 0 Mean STD

Histogram difgerence for same ROI

Frequency

0.05 0.1 0.15 0.2 0.25 0.3

Range values

0.34 0.36 0.38 0.4 0.42 0.44 0.46 0.48 0.5

Window size 7 × 7

Local Histogram Tile 8x8 Tile 16x16 Whole Image

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 24 / 33

Timing Results

LA-SAS kernels precomputations with different image-tile sizes

Timing results for difgerent tile sizes

Image Size

500 1000 1500 2000 2500 3000

Time (sec)

0.05 0.1 0.15 0.2 0.25 0.3 10 20 25 50

Timing results for difgerent image sizes

Tile Size

10 15 20 25 30 35 40 45 50

Time (sec)

0.05 0.1 0.15 0.2 0.25 0.3 500 1000 2000 3000

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 25 / 33

Timing Results

Throughput results for difgerent image sizes

Image Size

400 600 800 1000 1200 1400 1600 1800 2000

Throughput (MPix/s)

30 40 50 60 70 80 90 100 110 120 T = 10 × 10 T = 20 × 20 T = 25 × 25 Whole Image

Timing results for difgerent image sizes

Image Size

400 600 800 1000 1200 1400 1600 1800 2000

Time (msec)

20 40 60 80 100 120 140 T = 10 × 10 T = 20 × 20 T = 25 × 25 T = 50 × 50 Whole Image

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 26 / 33

localStat: vs MATLAB GPU

localStats Kernel over image size

Image Size

500 1000 1500 2000 2500 3000

Time (msec)

20 40 60 80 100 120 140 160 T = 8 × 8 localStat T = 8 × 8 MATLAB

localStats Kernel over image size

Image Size

500 1000 1500 2000 2500 3000

Time (msec)

10 20 30 40 50 60 70 T = 16 × 16 localStat T = 16 × 16 MATLAB

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 27 / 33

Timing Results

localStats Kernel over image size

Image Size

500 1000 1500 2000 2500 3000

Time (msec)

5 10 15 20 25 30 35 40 45 50 55 T = 20 × 20 localStat T = 20 × 20 MATLAB

localStats Kernel over image size

Image Size

500 1000 1500 2000 2500 3000

Time (msec)

5 10 15 20 25 30 35 40 45 T = 30 × 30 localStat T = 30 × 30 MATLAB

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 28 / 33

Timing Results

Timing ratio between whole image and tile for localStats

Tile Size

10 15 20 25 30 35 40 45 50 55 60

Ratio (%)

100 150 200 250 300 350 I = 500 × 500 Calculated Redundancy I = 1000 × 1000 Calculated Redundancy I = 1500 × 1500 Calculated Redundancy I = 2000 × 2000 Calculated Redundancy I = 2500 × 2500 Calculated Redundancy I = 3000 × 3000 Calculated Redundancy Theoretical

Timing ratio between whole image and tile for localStats

Image Size

500 1000 1500 2000 2500 3000

Ratio (%)

100 150 200 250 300 350 T = 8 × 8 T = 10 × 10 T = 16 × 16 T = 20 × 20 T = 30 × 30 T = 32 × 32 T = 40 × 40 T = 50 × 50 T = 64 × 64

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 29 / 33

Timing Results

Multiple kernels all timings

Image Size

500 1000 1500 2000 2500 3000

Time (msec)

20 40 60 80 100 120 140 160 T = 8 × 8 localStat T = 10 × 10 localStat T = 16 × 16 localStat T = 20 × 20 localStat T = 30 × 30 localStat T = 32 × 32 localStat T = 40 × 40 localStat T = 50 × 50 localStat T = 64 × 64 localStat T = 8 × 8 MATLAB T = 10 × 10 MATLAB T = 16 × 16 MATLAB T = 20 × 20 MATLAB T = 30 × 30 MATLAB T = 32 × 32 MATLAB T = 40 × 40 MATLAB T = 50 × 50 MATLAB T = 64 × 64 MATLAB

Timing for histogram preĄx sum

Image Size

10 15 20 25 30 35 40 45 50 55 60

Time (msec)

50 100 150 200 250 300 350 400 450 500 550 I = 500 × 500 localStat I = 1000 × 1000 localStat I = 1500 × 1500 localStat I = 2000 × 2000 localStat I = 500 × 500 MATLAB I = 1000 × 1000 MATLAB I = 1500 × 1500 MATLAB I = 2000 × 2000 MATLAB

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 30 / 33

Acknowledgements

• ARO Equipment Grant #W911NF-13-l-0344
• NVIDIA Academic Research Equipment Support

Iliopoulos, Floros, Pitsianis, Sun, Yin, Ren (Duke|AUTh) Locally Adaptive Signal-Noise Analysis GTC-2016 Apr 6, 2016 31 / 33

References I

• A. Alekseychuk. Hierarchical recursive running median. In Proceedings - International Conference on Image Processing, ICIP, pages 109Ű112, 2012. doi:

10.1109/ICIP.2012.6466807.

• B. Bilgic, B. K. P. Horn, and I. Masaki. Effjcient integral image computation on the GPU. IEEE Intelligent Vehicles Symposium, Proceedings, pages 528Ű533,

2010.

• C. Chan, R. Fulton, D. D. Feng, and S. Meikle. Median non-local means Ąltering for low SNR image denoising: application to PET with anatomical knowledge. In

IEEE Nuclear Science Symposium Conference Record, pages 3613Ű3618, Knoxville, TN, USA, Oct. 2010. ISBN 978-1-4244-9106-3. doi: 10.1109/NSSMIC.2010.5874485. URL http://ieeexplore.ieee.org/lpdocs/epic03/wrapper.htm?arnumber=5874485.

• J. Chen, S. Paris, and F. Durand. Real-time Edge-Aware Image Processing with the Bilateral Grid. ACM Transactions on Graphics (SIGGRAPH), (3):1Ű9, 2007.
• C. C. Cheng, F. C. Cheng, P. H. Lin, and S. C. Huang. A cloud-computing local histogram construction algorithm for big image data. In Proceedings - 2014 IEEE

International Congress on Big Data, BigData Congress 2014, pages 200Ű203, 2014.

• F. Durand and J. Dorsey. Fast bilateral Ąltering for the display of high-dynamic-range images. ACM Transactions on Graphics, 21(3):257Ű266, 2002.
• M. Kachelrieß. Branchless vectorized median Ąltering. IEEE Nuclear Science Symposium Conference Record, pages 4099Ű4105, 2009.
• M. Kass and J. Solomon. Smoothed local histogram Ąlters. ACM Transactions on Graphics, 29(4):1, 2010.
• Z. Ma, K. He, Y. Wei, J. Sun, and E. Wu. Constant time weighted median Ąltering for stereo matching and beyond. Proceedings of the IEEE International

Conference on Computer Vision, pages 49Ű56, 2013.

• J. Malcolm. Gpu and accelerator software blog - case study in cuda optimization.

http://blog.accelereyes.com/blog/2010/02/20/a-case-study-in-cuda-optimization, 2010.

• U. Milic, I. Gelado, N. Puzovic, A. Ramirez, and M. Tomasevic. Parallelizing general histogram application for CUDA architectures. Proceedings - 2013

International Conference on Embedded Computer Systems: Architectures, Modeling and Simulation, IC-SAMOS 2013, pages 11Ű18, 2013.

• S. Perreault and P. Hebert. Median Filtering in Constant Time. IEEE Transactions on Image Processing, 16(9):2389Ű2394, 2007.
• G. Perrot, S. Domas, and R. Couturier. Fine-tuned high-speed implementation of a GPU-based median Ąlter. Journal of Signal Processing Systems, 75(3):

185Ű190, 2014.

References II

• J. P. Pluim, J. A. Maintz, and M. A. Viergever. Mutual-information-based registration of medical images: a survey. IEEE Transactions on Medical Imaging, 22(8):

986Ű1004, Aug. 2003. ISSN 0278-0062. doi: 10.1109/TMI.2003.815867. URL http://ieeexplore.ieee.org/lpdocs/epic03/wrapper.htm?arnumber=1216223.

• M. Poostchi, K. Palaniappan, F. Bunyak, M. Becchi, and G. Seetharaman. Effjcient GPU implementation of the integral histogram. Lecture Notes in Computer

Science (including subseries Lecture Notes in Artificial Intelligence and Lecture Notes in Bioinformatics), 7728 LNCS(PART 1):266Ű278, 2013.

• F. Porikli. Constant time O (1) bilateral Ąltering. Computer Vision and Pattern Recognition, 2008. CVPR 2008. IEEE Conference on, (1):1Ű8, 2008.
• R. M. Sánchez and P. A. Rodríguez. Highly parallelable bidimensional median Ąlter for modern parallel programming models. Journal of Signal Processing

Systems, 71(3):221Ű235, 2013.

• J. a. Stark. Adaptive image contrast enhancement using generalizations of histogram equalization. Image Processing, IEEE Transactions on, 9(5):889Ű896, 2000.
• C. Tomasi and R. Manduchi. Bilateral Filtering for Gray and Color Images. International Conference on Computer Vision, 1998.
• P. Viola and M. Jones. Rapid object detection using a boosted cascade of simple features. Computer Vision and Pattern Recognition (CVPR), 1, 2001.
• P. Viola and M. J. Jones. Robust real-time face detection. International Journal of Computer Vision, 57(2):137Ű154, 2004.
• Y. Wei and L. Tao. Effjcient histogram-based sliding window. Proceedings of the IEEE Computer Society Conference on Computer Vision and Pattern

Recognition, (February):3003Ű3010, 2010.

• Q. Yang, K. H. Tan, and N. Ahuja. Real-time O(1) bilateral Ąltering. 2009 IEEE Computer Society Conference on Computer Vision and Pattern Recognition

Workshops, CVPR Workshops 2009, (1):557Ű564, 2009.

• Q. Yang, N. Ahuja, and K. H. Tan. Constant Time Median and Bilateral Filtering. International Journal of Computer Vision, (October 2013), 2014.

X.-Y. Zhang, L. Ge, and T.-F. Wang. Entropy-Based Local Histogram Equalization for Medical Ultrasound Image Enhancement. 2008 2nd International Conference on Bioinformatics and Biomedical Engineering, (c):2427Ű2429, 2008.

• H. Zhu, F. H. Chan, and F. Lam. Image Contrast Enhancement by Constrained Local Histogram Equalization. Computer Vision and Image Understanding, 73(2):

281Ű290, 1999.