BLOCK COMPRESSED SENSING OF IMAGES USING DIRECTIONAL TRANSFORMS
Sungkwang Mun and James E.Fowler
Department of Electrical and Computer Engineering,Geosystems Research Institute,
Mississippi State University,USA
ABSTRACT
Block-based random image sampling is coupled with a projection-driven compressed-sensing recovery that encourages sparsity in the domain of directional transforms simultaneously with a smooth reconstructed image.Both contourlets as well as complex-valued dual-tree wavelets are considered for their highly directional rep-resentation,while bivariate shrinkage is adapted to their multi-scale decomposition structure to provide the requisite sparsity con-straint.Smoothing is achieved via a Wienerfilter incorporated into iterative projected Landweber compressed-sensing recovery, yielding fast reconstruction.The proposed approach yields im-ages with quality that matches or exceeds that produced by a pop-ular,yet computationally expensive,technique which minimizes to-tal variation.Additionally,reconstruction quality is substantially superior to that from several prominent pursuits-based algorithms that do not include any smoothing.
Index Terms—Compressed sensing,contourlets,dual-tree dis-crete wavelet transform,bivariate shrinkage
1.INTRODUCTION
Recent years have seen significant interest in the paradigm of com-pressed sensing(CS)[1–3]which permits,under certain condi-tions,signals to be sampled at sub-Nyquist rates via linear pro-jection onto a random basis while still enabling exact reconstruc-tion of the original signal.As applied to2D images,however, CS faces several challenges including a computationally expen-sive reconstruction process and huge memory required to store the random sampling operator.Recently,several fast , [4–6])have been developed for CS reconstruction,while the lat-ter challenge was addressed in[7]using a block-based sampling operation.Additionally in[7],projection-based Landweber iter-ations were proposed to accomplish fast CS reconstruction while simultaneously imposing smoothing with the goal of improving the reconstructed-image quality by eliminating blocking artifacts.
In this paper,we adopt this same basic framework of block-based CS sampling of images coupled with iterative projection-based reconstruction with smoothing.Our contribution lies in that we cast the
reconstruction in the domain of recent transforms that feature a highly directional decomposition.These transforms—specifically,contourlets[8]and complex-valued dual-tree wavelets [9]—have shown promise to overcome deficiencies of widely-used wavelet transforms in several application areas.In their applica-tion to iterative projection-based CS recovery,we adapt bivariate shrinkage[10]to their directional decomposition structure to pro-vide sparsity-enforcing thresholding,while a Wiener-filter step en-courages smoothness of the result.In experimental simulations, wefind that the proposed CS reconstruction based on directional transforms outperforms equivalent reconstruction using common wavelet and cosine transforms.Additionally,the proposed tech-nique usually matches or exceeds the quality of total-variation(TV) reconstruction[11],a popular approach to CS recovery for images whose gradient-based operation also promotes smoothing but runs several orders of magnitude slower than our proposed algorithm.
2.BACKGROUND
Suppose that we want to recover real-valued signal x with length N from M samples,M≪,we want to recover x from y=Φx,where y has length M,andΦis an M×N measure-ment matrix.Because the number of unknowns is much larger than observations,recovering every x∈ℜN from its corresponding y is impossible in general;however,if x is sufficiently sparse,exact recovery is possible—
this is the fundamental tenant of CS theory; ,[3],for a more complete overview.The usual choice for the measurement basisΦis a random matrix;here,we further assume thatΦis orthonormal such thatΦΦT=I.
Quite often,the requisite sparsity will exist with respect to some transformΨ.In this case,the key to CS recovery is the production of a sparse set of significant transform coefficients,ˇx=Ψx,and the ideal recovery procedure searches for theˇx with the smallestℓ0norm consistent with the observed y.How-ever,thisℓ0optimization being NP-complete,several alternative solution procedures have been proposed.Perhaps the most promi-nent of these is basis pursuit(BP)[12]which applies a convex relaxation to theℓ0problem resulting in anℓ1optimization,ˇx=arg min
ˇx ˇx 1,such that y=ΦΨ
−1ˇx,(1)
whereΨ−1is the inverse transform.Although BP can be im-plemented effectively with linear programming,its computational complexity is often high,leading to recent interest in reduced-complexity ,gradient projection for sparse recon-struction(GPSR)[4])as well as in greedy BP variants,includ-ing matching pursuits,orthogonal matching pursuits,and,recently, sparsity
adaptive matching pursuits(SAMP)[5].Such algorithms significantly reduce computational complexity at the cost of lower reconstruction quality.
As an alternative to the pursuits class of CS reconstruction, techniques based on projections have been proposed , [6]).Algorithms of this class formˇx by successively projecting and thresholding;for example,the reconstruction in[6]starts from some initial approximationˇx(0)and forms the approximation at iteration i+1astransform英文
ˇx(i)=ˇx(i)+1
γ
ΨΦT“y−ΦΨ−1ˇx(i)”,(2)ˇx(i+1)=(ˇx(i),˛˛˛ˇx(i)˛˛˛≥τ(i),
0else.
(3)
Here,γis a scaling factor([6]uses the largest eigenvalue ofΦTΦ), whileτ(i)is a threshold set appropriately at each iteration.It is straightforward to see that this procedure is a specific instance of
a projected Landweber(PL)algorithm[13].Like the greedy algo-rithms of the pursuits class,PL-based CS reconstruction also pro-vides reduced computational complexity.Additionally,and per-haps more importantly,the PL formulation offers the possibility of easily incorporating additional optimization criteria.For exam-ple,the technique that we overview in the next section incorporates Wienerfiltering into the PL iteration to search for a CS reconstruc-tion simultaneously achieving sparsity and smoothness.
3.BLOCK-BASED CS WITH SMOOTHED PL
RECONSTRUCTION
In[7],a paradigm for the CS of2D images was proposed.In this technique,the sampling of an image is driven by random matri-ces applied on a block-by-block basis,while the reconstruction is a variant of the PL reconstruction of(2)–(3)that incorporates a smoothing operation.Due to its combination of block-based CS (BCS)sampling and smoothed-PL(SPL)reconstruction,we refer to the overall technique as BCS-SPL.We overview its constituent components below.
3.1.BCS—Block-Based CS Sampling
In BCS,an image is divided into B×B blocks and sampled using an appropriately-sized measurement matrix.That is,suppose that x j is a vector representing,in raster-scan fashion,block j of input image x.The corresponding y j is then y j=ΦB x j,whereΦB is an M B×B2orthonormal measurement matrix with M B=⌊M N B2⌋.
Using BCS rather than random sampling applied to the entire image x has several merits[7].First,the measurement operator ΦB is conveniently stored and employed because of its compact size.Second,the encoder does not need to wait until the entire image is measured,but may send each block after its linear pro-jection.Last,an initial approximation x(0)with minimum mean squared error can be feasibly calculated due to the small size of ΦB[7].As done in[7],we employ blocks of size B=32.
3.2.SPL—A Smoothed PL Variant
In[7],Wienerfiltering was incorporated into the basic PL frame-work described in Sec.2in order to remove blocking artifacts.In essence,this operation imposes smoothness in addition to the spar-sity inherent to PL.Specifically,in[7],a Wiener-filtering step was interleaved with the PL projection of(2)–(3);thus,the approxima-tion to the image at iteration i+1,x(i+1),is produced from x(i) as:
function x(i+1)=SPL(x(i),y,ΦB,Ψ,λ)
ˆx(i)=Wiener(x(i))
for each block j
ˆˆx(i) j =ˆx(i)j+ΦT B(y−ΦBˆx(i)j)
ˇx(i)=Ψˆˆx(i)
ˇx(i)=Threshold(ˇx(i),λ)
¯x(i)=Ψ−1ˇx(i)
for each block j
x(i+1)
j
=¯x(i)j+ΦT B(y−ΦB¯x(i)j)
Here,Wiener(·)is pixelwise adaptive Wienerfiltering using a neighborhood of3×3,while Threshold(·)is a thresholding pro-cess as discussed below.In our use of SPL,we initialize with x(0)=ΦT y and terminate when|D(i+1)−D(i)|<10−4,where D(i)=1√
N x(i)−
ˆˆx(i−1) 2.
4.DIRECTIONAL TRANSFORMS AND BCS-SPL
4.1.Transforms
In[7],several iterations of the SPL(·)procedure described above are used as an initial step in a dual-stage algorithm for CS recon-struction.The stages employ PL iterations in the form of(2)–(3) using several different transforms,including a block-based lapped cosine transform as well as a redundant wavelet transform.For reasons of simplicity,we now depart from this methodology and instead focus on a single stage of SPL(·)iterations.This allows us to incorporate several prominent directional transforms into the basic SPL formulation to evaluate their relative efficacy at CS re-construction.Although we do not pursue it here,multiple SPL stages in the style of[7]could be employed along with these di-rectional transforms to potentially refine performance.
Although the discrete wavelet transform(DWT)is widely used for image compression,DWTs in their traditional critically sam-pled form are known to be somewhat deficient in several charac-teristics,lacking such properties as shift invariance and significant directional selectivity.As a result,there have been several recent proposals made for transforms that feature a much higher degree of directional representation than is obtained with traditional DWTs. Two prominent families of such directional transforms are con-tourlets and complex-valued DWTs.
The contourlet transform(CT)[8]preserves interesting fea-tures of the traditional DWT,namely multiresolution and local characteristics of the signal,and,at the expense of a spatial re-dundancy,it better represents the directional features of the image. The CT couples a Laplacian-pyramid decomposition with direc-tionalfilterbanks,inheriting the redundancy of the Laplacian ,4/3).
Alternatively,complex-valued wavelet transforms have been proposed to improve upon DWT deficiencies,with the dual-tree DWT(DDWT)[9]becoming a preferred approach due to the ease of its implementation.In the DDWT,real-valued waveletfilters produce the real and imaginary parts of the transform in parallel decomposition trees.DDWT yields a decomposition with a much higher degree of directionality than that possessed by the tradi-tional DWT;however,since both trees of the DDWT a
re them-selves orthonormal or biorthogonal decompositions,the DDWT taken as a whole is a redundant tight frame.
Albeit redundant,both the CT and DDWT have been effec-tively used for image ,[14–16]).The experi-mental results below explore the efficacy of these directional trans-forms in the SPL-based CS reconstruction of Sec.3.
4.2.Thresholding
As originally described in[7],SPL(·)used hard thresholding in the form of(3).To set a properτfor hard thresholding,we employ the universal threshold method of[17].Specifically,in(3),
τ(i)=λσ(i)p2log K,(4) whereλis a constant control factor to manage convergence,and K is the number of the transform coefficients.As in[17],σ(i)is estimated using a robust median estimator,
σ(i)=
median(|ˇx(i)|)
0.6745
.(5)
Hard thresholding inherently assumes independence between coefficients.However,bivariate shrinkage [10]is better suited to directional transforms in that it exploits statistical dependency between transform coefficients and their respective parent coeffi-cients,yielding performance superior to that of hard thresholding.In [10],a non-Gaussian bivariate distribution was proposed for the current coefficient and its lower-resolution parent coefficient based on an empirical joint histogram of DWT coefficients.However,it is straightforward to apply this process to any transform having a multiple-level decomposition,such as the directional transforms we consider here.Specifically,given a specific transform coeffi-cient ξand its parent coefficient ξp in the next coarser scale,the Threshold(·)operator in SPL is the MAP estimator of ξ,
Threshold(ξ,λ)=…p ξ2+ξ2p
−λ√3σ(i )σξ
«+
p ξ2+ξ2p
·
ξ,(6)
where (g )+=0for g <0,(g )+=g else;σ(i )is the median estimator of (5)applied to only the finest-scale transform coeffi-cients;and,again,λis a convergence-control factor.Here,σ2ξ
is the marginal variance of coefficient ξestimated in a local 3×3neighborhood surrounding ξas in [10].
5.EXPERIMENTAL RESULTS
To evaluate directional transforms for CS reconstruction,we de-ploy both the CT and DDWT within the BCS-SPL framework de-scribed in Sec.3.We refer to the resulting implementations as BCS-SPL-CT and BCS-SPL-DDWT,respectively.To evaluate the effectiveness of the increased directionality of the CT and DDWT,we compare to BCS-SPL-DWT,an equivalent approach using the ubiquitous biorthogonal 9-7DWT.We also compare to a BCS-SPL variant using a block-based DCT for SPL reconstruction;the resulting algorithm (BCS-SPL-DCT)is similar to that initially pro-posed as the first-stage reconstruction in [7],although a lapped transform was used there.In SPL,we use bivariate shrinkage (6)with λ=10,25,and 25,respectively,for BCS-SPL-CT,BCS-SPL-DDWT,and BCS-SPL-DWT.Lacking parent-child relations,BCS-SPL-DCT uses hard thresholding (4)with λ=6.
We compare also to BCS-TV wherein the block-based BCS is still used for image sampling while the SPL reconstruction is replaced with the minimum TV optimization of [11].Like SPL,such TV-based reconstruction also imposes sparsity and smooth-ness constraints;unlike the explicit smoothing of SPL’s Wiener filtering,however,smoothness in TV is implicit in that the solution is sparse in a gradient space.We use ℓ1-MAGIC 1in the BCS-TV implementation.Finally,as representative of the pursuits class of CS reconstruction,we compare to GPSR 2[4]as well as SAMP 3[5];for these,we use implementations provided by their respective authors.
Table 1compares PSNR for several 512×512images at sev-eral measurement ratios,M/N .We note that,since the quality of reconstruction can vary due to the randomness of the measure-ment matrix Φ(or ΦB ),all PSNR figures are averaged over 5independent trials.The results indicate that BCS-SPL with the directional transforms achieves the best performance at low mea-surement rates.At higher measurement rates,performance is more
1www.acm.caltech.edu/l1magic/2www.lx.it.pt/
~mtf/GPSR/
3lepages/samp_intro/
varied—BCS-TV is more competitive;however,the directional
BCS-SPL techniques usually produce PSNR close to that of the TV-based algorithm.However,the ℓ1-MAGIC implementation of TV is quite slow—BCS-TV takes 3–4hours for each trial,whereas the BCS-SPL implementations run for only 1–5minutes depend-ing on the complexity of the transform used.These later times are in line with GPSR (less than 60seconds)and SAMP (several minutes).Times are for a 3.2-GHz dual-core processor.
Fig.1illustrates example visual results.We note that,despite the smoothing inherent to TV reconstruction,blocking artifacts are apparent.On the other hand,the smoothing of the BCS-SPL-DDWT reconstruction (BCS-SPL-CT yields similar visual quality)eliminates blocking while the enhanced directionality of the trans-form provides better quality than the DWT-and DCT-based tech-niques.Finally,we note that the pursuits-based algorithms which do not factor in any sort of smoothing—GPSR and SAMP—yield images of noticeably deficient visual quality.
6.CONCLUSIONS
In this paper,we examined the use of recently proposed direc-tional transforms in the CS reconstruction of images.We adopted the general paradigm of block-based random image sampling
cou-pled with a projection-based reconstruction promoting not only sparsity but also smoothness of the reconstruction.This frame-work facilitates the incorporation into the CS-recovery process of directional transforms based on contourlets and complex-valued dual-tree wavelets.The resulting algorithms inherit the fast exe-cution speed of the projection-based CS reconstruction while the enhanced directionality coupled with a smoothing step encourages superior image quality,particularly at low sampling rates.
7.REFERENCES
[1]  E.Cand`e s and T.Tao,“Near-optimal signal recovery from
random projections:Universal encoding strategies?”IEEE Transactions on Information Theory ,vol.52,no.12,pp.5406–5425,December 2006.[2]  D.L.Donoho,“Compressed sensing,”IEEE Transactions
on Information Theory ,vol.52,no.4,pp.1289–1306,April 2006.[3]  E.J.Cand`e s and M.B.Wakin,“An introduction to compres-sive sampling,”IEEE Signal Processing Magazine ,vol.25,
no.2,pp.21–30,March 2008.[4]M.A.T.Figueiredo,R.D.Nowak,and S.J.Wright,“Gradi-ent projection for s
parse reconstruction:Application to com-pressed sensing and other inverse problems,”IEEE Journal on Selected Areas in Communications ,vol.1,no.4,pp.586–597,December 2007.[5]T.T.Do,L.Gan,N.Nguyen,and T.D.Tran,“Sparsity adap-tive matching pursuit algorithm for practical compressed sensing,”in Proceedings of the 42th Asilomar Conference on Signals,Systems,and Computers ,Pacific Grove,California,October 2008,pp.581–587.[6]J.Haupt and R.Nowak,“Signal reconstruction from noisy
random projections,”IEEE Transactions on Information Theory ,vol.52,no.9,pp.4036–4048,September 2006.
Figure 1:Lenna for M/N =20%.Left to right:BCS-SPL-DDWT,PSNR =31.37dB;BCS-SPL-DCT,PSNR =30.45dB;BCS-TV ,PSNR =30.60dB;SAMP,PSNR =28.54dB.
[7]L.Gan,“Block compressed sensing of natural images,”in
Proceedings of the International Conference on Digital Sig-nal Processing ,Cardiff,UK,July 2007,pp.403–406.[8]M.N.Do and M.Vetterli,“The contourlet transform:An
efficient directional multiresolution image representation,”IEEE Transactions on Image Processing ,vol.14,no.12,pp.2091–2106,December 2005.[9]N.G.Kingsbury,“Complex wavelets for shift invariant anal-ysis and filtering of signals,”Journal of Applied Computa-tional Harmonic Analysis ,vol.10,pp.234–253,May 2001.[10]L.S ¸endur and I.W.Selesnick,“Bivariate shrinkage func-tions for wavelet-based denoising exploiting interscale de-pendency,”IEEE Transactions on Signal Processing ,vol.50,
no.11,pp.2744–2756,November 2002.[11]  E.Cand`e s,J.Romberg,and T.Tao,“Stable signal recovery
from incomplete and inaccurate measurements,”Communi-cations on Pure and Applied Mathematics ,vol.59,no.8,pp.1207–1223,August 2006.[12]S.S.Chen,D.L.Donoho,and M.A.Saunders,“Atomic
decomposition by basis pursuit,”SIAM Journal on Scientific Computing ,vol.20,no.1,pp.33–61,August 1998.[13]M.Bertero and P.Boccacci,Introduction to Inverse Problems
in Imaging .Bristol,UK:Institute of Physics Publishing,1998.[14]M.Trocan,B.Pesquet-Popescu,and J.E.Fowler,“Graph-cut rate distortion algorithm for contourlet-based image com-pression,”in Proceedings of the International Conference on Image Processing ,vol.3,San Antonio,TX,September 2007,pp.169–172.[15]J.E.Fowler,J.B.Boettcher,and B.Pesquet-Popescu,“Im-age coding using a complex dual-tree wavelet transform,”in Proceedings of the European Signal Processing Conference ,Pozna´n ,Poland,September 2007.[16]J.B.Boettcher and J.E.Fowler,“Video coding using a com-plex wavelet transform and set partitioning,”IEEE Signal Processing Letters ,vol.14,no.9,pp.633–636,September 2007.[17]  D.L.Donoho,“De-noising by soft-thresholding,”IEEE
Transactions on Information Theory ,vol.41,no.3,pp.613–627,May 1995.
Table 1:PSNR Performance in dB Measurement Rate (M/N )Algorithm 0.10.20.30.4
0.5Lenna
BCS-SPL-CT 28.1731.0232.9934.6836.25BCS-SPL-DDWT 28.3131.3733.5035.2036.78BCS-SPL-DWT 27.8130.8932.9434.6136.15BCS-SPL-DCT 27.7030.4532.4634.1935.77BCS-TV 27.8630.6032.5634.2535.89SAMP 25.9428.5432.0433.9335.37GPSR 24.6928.5431.5333.69
35.82Barbara
BCS-SPL-CT 22.7524.3325.9027.5429.38BCS-SPL-DDWT 22.8524.2925.9227.5029.12BCS-SPL-DWT 22.6223.9425.2026.5628.05BCS-SPL-DCT 22.7624.3825.9127.4229.05BCS-TV 22.4523.6024.5725.5726.73SAMP 20.9722.8325.0427.6830.08GPSR 20.2322.6624.9927.42
30.15Peppers
BCS-SPL-CT 28.5631.0432.5733.7734.88BCS-SPL-DDWT 28.8831.4432.8934.0635.18BCS-SPL-DWT 28.6931.0432.4833.6334.74BCS-SPL-DCT 27.8830.4131.9033.0834.20BCS-TV 28.5231.2132.7433.9635.17SAMP 25.9428.6130.6931.7132.42GPSR 24.5828.1930.1931.76
33.21Mandrill
BCS-SPL-CT 22.8724.9726.9528.9030.93BCS-SPL-DDWT 22.9424.8726.6928.4230.28BCS-SPL-DWT 22.5424.3326.0327.6829.38BCS-SPL-DCT 22.3124.1525.9427.7729.68BCS-TV 22.3124.3426.0827.7729.45SAMP 20.2021.5624.0027.1030.29GPSR 20.1122.2324.3726.99
30.06Goldhill
BCS-SPL-CT 26.8528.9530.4831.9233.28BCS-SPL-DDWT 26.9628.9330.4531.7933.11BCS-SPL-DWT 26.7128.6830.1331.5332.85BCS-SPL-DCT 26.1028.3229.6330.9832.57BCS-TV 26.5328.8530.
5632.0933.61SAMP 24.3126.3028.0729.4530.86GPSR
23.6326.1428.0930.02
31.72

版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。