Document Sample

Structured Pursuits for Geometric Super-Resolution e St´ phane M ALLAT and Guoshen Y U CMAP, Ecole Polytechnique, 91128 Palaiseau Cedex, France Abstract—Super-resolution image zooming is possible when the II. D IRECTIONAL I NTERPOLATIONS image has some geometric regularity. We introduce a general class of non-linear inverse estimators, which combines linear Let us consider a subsampling operator U by a factor 2. estimators with mixing weights in a frame providing a sparse The measured low-resolution image y is representation. Mixing weights are computed with a block decom- position, which minimizes a Tikhonov energy penalized by an l1 y[n1 , n2 ] = f [2n1 , 2n2 ] + w[n1 , n2 ] (2) norm of the mixing weights. A fast orthogonal matching pursuit and if we neglect the noise w its Fourier transform is algorithm computes the mixing weights. Adaptive directional image interpolations are calculated with mixing weights in a 1 ˆ ˆ wavelet frame. y (ω1 , ω2 ) = f (ω1 , ω2 ) + f (ω1 + π, ω2 ) ˆ 4 I. I NTRODUCTION ˆ ˆ + f (ω1 , ω2 + π) + f (ω1 + π, ω2 + π) . (3) We consider linear inverse problems where the measure- Fig. 1(b) shows such an aliased Fourier transform. ments y of a signal f are obtained through a non-invertible operator U to which is added a noise w: y = U f + w. Image interpolation is an example of such ill-posed inverse −pi −pi problems, where U is a subsampling operator. This paper studies non-linear inverse estimators, obtained as a mixture of linear inverse estimators, for image interpolation. 0 0 Increasing the size of images is needed for many image or video screens. Most high deﬁnition televisions perform pi −pi 0 pi pi −pi 0 pi a scaling of input standard deﬁnition images, with a spatial a b c d interpolation. When images are aliased, linear interpolations introduce artifacts such as Gibbs oscillations or zigzag along Fig. 1. (a): a high-resolution edge image f [n] and its Fourier transform. (b): a low-resolution image y[n] = f [2n] and its aliased edges. Better estimations can be obtained with non-linear Fourier transform. (c): 2D cubic-spline interpolation (50.45 dB). (d): estimators which take into account local image properties. super-resolution estimation with directional interpolation along the Suppose that we are given a family of linear estimators contour direction (72.72 dB). The ellipse in (b) gives the frequency ˜ domain preserved by the directional interpolator. fθ = Iθ y that are optimized with speciﬁc priors on f . For ˆ To recover f (ω1 , ω2 ), a ﬁlter h is adjusted to reduce the image interpolation, Iθ is a directional interpolation along impact of the three aliased components: a direction θ, which is precise if the image has regular ˜ f (n) = Iy(n) = y(q) h(n − 2q) for n = (n1 , n2 ). variations along θ. Mixing linear estimators is very effective q for noise removal [3]. However, these approaches do not apply to ill-posed inverse problems because they rely on a Stein Interpolation ﬁlters satisfy h(0) = 1 and h(2n) = 0 for n = 0, ˜ so that f (2q) = f (2q) if w = 0. If f is uniformly regular unbiased empirical estimator of the risk, and there is no such estimator for ill-posed inverse problems. In the context of then its Fourier transform has an energy mostly concentrated image interpolation, this paper introduces a new approach to in [−π/2, π/2]2 and h is chosen to be a low-pass ﬁlter. deﬁne an adaptive inverse estimator which combines linear Among all linear ﬁlters, a cubic spline interpolation ﬁlter estimators. The measurement y is decomposed in a basis or h = hs deﬁnes an interpolator Is which yield nearly the a frame {φp }p∈G providing a sparse representation, and a best PSNR for most images [5]. When the frequency support ˆ of f exceeds [−π/2, π/2]2 then a low frequency interpolator family of linear estimators {Iγ }γ∈Γ are combined with mixing weights a(γ, p): produces Gibbs oscillations, zigzag jaggies along contours. Fig. 1(c) gives an example with a cubic spline ﬁlter. The optimization of h can be adapted to prior information ˜ f= Iγ ˜ a(γ, p) y, φp φp . (1) on f with a quadratic Tikhonov regularization. Suppose that f γ∈Γ p∈G has some form of regularity, which is measured by a regularity Mixing weights are computed for directional image interpo- norm Rf , where Rf is typically a high-pass ﬁltering. lations with an adaptive block transform in a wavelet frame. Sobolev norms are particular instances of such regularity Directional blocks are adapted to the directional regularity of norms. Let σ 2 be the variance of the noise w. A Tikhonov wavelet coefﬁcients, with a block matching pursuit algorithm. ˜ estimator computes f = Iy by 2 2 minimizing RIy subject to U Iy − y≤ σ 2 . (4) ˜ We shall neglect the noise for image interpolation, so U f = y. Minimizing Rθ f ˜ ˜ where Rθ computes variations of f in the direction θ is a prior regularity condition which corresponds to ˆ images where f (ω1 , ω2 ) has a frequency support concentrated in a narrow band around a line of angle θ+π/2, as illustrated in Fig. 1(a). The resulting Tikhonov interpolator Iθ in Fig. 1(d), Fig. 2. A block Bθ,q of orientation θ regroups wavelet coefﬁcients that are aligned with the directional image regularity in a neighbor- performs an interpolation in the direction θ with a ﬁlter hood of a point q, in the 3 wavelet coefﬁcient images corresponding hθ which eliminates the aliased frequency components while to the 3 wavelet directions d = 1, 2, 3. Three blocks are shown in ˆ preserving the support of f (ω1 , ω2 ). this ﬁgure. The direction of image regularity changes in space and it is therefore necessary to adapt locally the interpolation angle θ. function of a block, which does not depend on d = 1, 2, 3. Many researchers and engineers in industry have studied such The block position q is sampled along a grid Gθ ⊂ G so adaptive interpolators. The interpolation direction is adapted that G = ∪q∈Gθ Bθ,q and each block overlaps with only two along edges, using ad-hoc edge detectors [2] or gradient based other blocks of same orientation. Among all Bθ,q , the blocks orientation estimators [4], [6]. These algorithms can provide selected by interpolator Iθ are located where the image is better visual quality then a linear interpolation but the PSNR regular in the direction θ, as illustrated in Fig. 2. This selection is hardly improved. is performed by coefﬁcients (θ, q). They should be close to 1 when the image is regular in the direction θ and otherwise III. B LOCK A DAPTIVE T IKHONOV E STIMATIONS close to 0. Selection coefﬁcients decompose yh according to Adaptive interpolation amounts to combine a family of its directional regularity plus a residue: directional interpolators Iθ and the cubic spline interpolator Is . This section introduces an algorithm which computes yh = (θ, q) yθ,q + yr . (6) θ∈Θ q∈Gθ mixing weights in a wavelet frame, according to the general formulation (1), with an adaptive block decomposition. The residue yr has no directional regularity and is therefore The subsampled image y[n] for n ∈ G is decomposed in a interpolated with a cubic spline interpolator Is . The resulting d estimator is: translation invariant wavelet tight frame {ψp }0≤d≤3,p∈G on 1 scale. It separates low frequencies yl over the low frequency 0 ˜ scaling ﬁlters {ψp }p∈G and high frequencies yh over ﬁne scale f = Is (yl ) + Iθ (θ, q) yθ,q + Is (yr ). (7) d wavelets {ψp }1≤d≤3,p∈G : θ∈Θ q∈Gθ 3 This interpolator is a mixture of the interpolators {Iθ , Is }θ∈Θ 0 d yl = c(0, p) ψp and yh = c(d, p) ψp , (5) where the image low frequency and the high frequency residue p∈G d=1 p∈G are interpolated with a cubic spline: ˜d where the wavelet coefﬁcients c(d, p) = y, ψp are computed 3 with the dual wavelet frame. The low frequency image yl has ˜ f= Iθ a(θ, p) c(d, p) ψp d (8) little aliasing and can thus be precisely interpolated with a θ∈Θ p∈G d=1 cubic spline interpolator Is . The high frequency image yh is 3 interpolated by combining directional interpolators and a cubic + Is 0 c(0, p) ψp + a(s, p) c(d, p) ψp . d spline interpolator. p∈G p∈G d=1 Let Θ be a set of angles θ uniformly discretized between 0 and π. An adaptive interpolator is implemented over ﬁne The mixing weights are given by d scale wavelets ψp , by mixing the interpolators {Iθ , Is }θ∈Θ a(θ, p) = (θ, q)1Bθ,q (p) and a(s, p) = 1− a(θ, p) . with weights which depend on the image regularity. The ap- q∈Gθ θ∈Θ pendix speciﬁes the directional interpolators used in numerical (9) experiments. Selection coefﬁcients (θ, q) are computed by minimizing To construct a mixing operator, wavelets are regrouped in a criterion which multiplies these coefﬁcients by a Tikhonov blocks for which a speciﬁc interpolator is chosen. Let us regularity term in the corresponding blocks, while minimizing denote the norm of the residue yr . The norm of the residue is d yθ,q = c(d, p) ψp controlled by their wavelet coefﬁcient energy. The regularity (d,p)∈Bθ,q in the direction θ is calculated with an operator Rθ applied to wavelet coefﬁcients, which computes their variation in the the image projection over wavelets in a spatial block Bθ,q . direction θ: This block is a rectangle elongated in the direction θ, and 2 located in the spatial neighborhood of q for the three directions Rθ c Bθ,q = |c(d, p) − Aθ c(d, p)|2 , d = 1, 2, 3. We write 1Bθ,q (d, p) = 1Bθ,p (p) the indicator (d,p)∈Bθ,q where Aθ c(d, p) is constant over each line of angle θ included The overall algorithm thus proceeds as follow: in Bθ,q and equal to the average of c(d, p) over each of these • Compute the wavelet transform of y. lines, for each d = 1, 2, 3. Selection coefﬁcients are thus • Compute selection coefﬁcients with the orthogonal computed by minimizing pursuit. 3 • Apply the reconstruction (8) with the mixing weights (9). 1 2 L( ) = |c(d, p)|2 1 − (θ, q) 1Bθ,q (p) Fig. 3 compares a separable cubic spline interpolation with 2 d=1 p∈G θ∈Θ q∈Gθ this adaptive interpolation computed by mixing directional +λ (θ, q) Rθ c 2 Bθ,q , (10) interpolations. The mixing algorithm improves the PSNR and θ∈Θ q∈Gθ the visual image quality where the image is geometrically regular. This clearly appears in the straws, the hat border, the with 2 cloth and the hairs of various directions. c Bθ,q = |c(d, p)|2 . (d,p)∈Bθ,q A PPENDIX The energy L( ) has a quadratic term in , which is penalized We describe a directional interpolator Iθ for upscaling by a weighted l1 norm. It is minimized by standard numerical images by a factor 2, which is illustrated in Fig. 4. Origi- algorithms such as an iterated thresholding [1]. nal samples are shown as crosses. Step 1 computes a one- dimensional interpolations in the direction θ. We consider all IV. O RTHOGONAL B LOCK M ATCHING P URSUIT lines of angle θ that intersect original image samples (crosses An l1 minimization is computationally too expensive for in Fig. 4) and we compute mid-points (circles) with a cubic real time image processing applications. In the following, we spline interpolation. It oversamples by a factor two either develop a fast orthogonal matching pursuit algorithm, which the image rows, or the image columns or the diagonals of computes an approximate solution. angle ±π/4. The missing coefﬁcients are shown as squares An orthogonal block matching pursuit is a greedy algorithm in Fig. 4. Step 2 computes new samples (dots) with a cubic which selects progressively non-overlapping blocks which spline interpolation along these oversampled rows, columns or minimize the energy (10). If a single block Bθ,q is selected diagonals. The position of these new samples are chosen so then one can verify that (10) is minimized for that any missing coefﬁcient (square) is a mid-point between Rθ c 2 two dots on a line of angle θ. Step 3 computes missing samples Bθ,q (θ, q) = max 1 − λ 2 ,0 . (11) (squares) with a cubic spline linear interpolation along the c Bθ,q direction θ from the previously calculated new samples (dots). Setting (θ, q) to this optimal value reduces L( ) by 2 2 c Bθ,q Rθ c Bθ,q 2 e(θ, q) = max 1 − λ 2 ,0 . (12) 2 c Bθ,q The best block Bθ,q is the one which minimizes L( ) and hence which maximizes e(θ, q). An orthogonal matching pursuit algorithm selects blocks Fig. 4. Directional interpolation in three steps: a one-dimensional that do not intersect and which reduce as much as possible interpolations along the angle θ, a vertical interpolation and another L( ). At each iteration, among all blocks that do not intersect interpolation along θ. with the previously selected blocks, it ﬁnds the block which R EFERENCES maximizes e(θ, q). It proceeds as follows. 1) Initialize Ω0 = {(θ, q) : θ ∈ Θ , q ∈ Gθ }, (θ, q) = 0 [1] I. Daubechies. and M. Defrise, and C. De Mol, “An Iterative Thresholding Algorithm for Linear Inverse Problems with a for all (θ, q) ∈ Ω0 and l = 0. Sparsity Constraint”, Comm. on Pure Appl. Math., 57:1413– 2) Find (θl , ql ) ∈ Ωl which maximizes e(θ, q) in (12). and 1457, 2004. set (θl , ql ) with (11). [2] K. Jensen and D. Anastassiou, “Subpixel edge localization 3) Set Ωl+1 = {(θ, q) ∈ Ωl : Bθ,q ∩ Bθl ,ql = 0} / and the interpolation of still images”, IEEE Trans. on Image 4) If e(θl , ql ) > 0 set l = l + 1 and go back to (2). Processing, 4(3):285–295, 1995. [3] G. Leung and Barron, A.R., “Information theory and mixing Otherwise stop. least-squares regressions”, IEEE Trans. on Information Theory, This algorithm computes selection coefﬁcients (θ, q) for 52(8):3396–3410, 2006. all θ ∈ Θ and q ∈ Gθ . A translation invariant estimator [4] X. Li and M. T. Orchard, “New edge-directed interpolation”, computes selection coefﬁcients (θ, q) for all q ∈ G. This IEEE Trans. on Image Processing, 10(10):1521–1527, 2001. [5] M. Unser, “Splines: A perfect ﬁt for signal and image process- is done by covering the original image grid G by a union ing”, IEEE Signal Processing Magazine, 16(6):22–38, 1999. of translated grids Gθ,τ = Gθ + τ . The translation invariant [6] Q. Wang and R.K. Ward, “A New Orientation-Adaptive Interpo- selection coefﬁcients (θ, q) are the average of the selection lation Method” IEEE Trans. on Image Processing, 16(4):889– coefﬁcients calculated over each translated grid Gθ,τ , with the 900, 2007. orthogonal matching pursuit algorithm. High-resolution image Cubic spline 21.37 dB Adaptive interpolation 23.82 dB High-resolution image Cubic spline 34.02 dB Adaptive interpolation 35.84 dB High-resolution image Cubic spline 29.89 dB Adaptive interpolation 30.49 dB High-resolution image Cubic spline 28.66 dB Adaptive interpolation 29.27 dB Fig. 3. From left to right: high-resolution image, cubic-spline interpolation, adaptive directional interpolation with mixing weights.

DOCUMENT INFO

Shared By:

Categories:

Tags:
signal processing, video coding, image processing, compressed sensing, compressive sensing, united kingdom, motion estimation, ieee trans, matching pursuit, ieee international conference, stéphane mallat, richard baraniuk, technical report, wavelet coefﬁcients, image segmentation

Stats:

views: | 9 |

posted: | 8/19/2010 |

language: | English |

pages: | 4 |

OTHER DOCS BY nne25858

How are you planning on using Docstoc?
BUSINESS
PERSONAL

By registering with docstoc.com you agree to our
privacy policy and
terms of service, and to receive content and offer notifications.

Docstoc is the premier online destination to start and grow small businesses. It hosts the best quality and widest selection of professional documents (over 20 million) and resources including expert videos, articles and productivity tools to make every small business better.

Search or Browse for any specific document or resource you need for your business. Or explore our curated resources for Starting a Business, Growing a Business or for Professional Development.

Feel free to Contact Us with any questions you might have.