# Through-Wall Image Enhancement Based on Singular Value Decomposition

Authors: Muhammad Mohsin Riaz; Abdul Ghafoor
Journal: International Journal of Antennas and Propagation
 (2012)

---

## Abstract

Singular value decomposition and information theoretic criterion-based image enhancement is proposed
for through-wall imaging. The scheme is capable of discriminating target, clutter, and noise
subspaces. Information theoretic criterion is used with conventional singular value decomposition to
find number of target singular values. Furthermore, wavelet transform-based denoising is performed (to
further suppress noise signals) by estimating noise variance. Proposed scheme works also for extracting
multiple targets in heavy cluttered through-wall images. Simulation results are compared on the basis of
mean square error, peak signal to noise ratio, and visual inspection.


---

## Body

## 1. Introduction

Through-wall imaging (TWI) is an active research area [1–11] due to its wide range of applications especially in rescue, military, surveillance, and remote sensing. TWI (seeing through opaque materials) gives the ability to examine structure layout of building, detection, and localization of target(s). As compared to other remote sensing techniques (ground penetrating radar and medical imaging), TWI has to deal with variety of challenges (like propagation environment, sensor positioning, and operational requirements). Moreover, propagation medium (which is often composed of multiple unknown and nonhomogenous walls) leads to multipaths and strong clutters which makes TWI a complex and challenging problem [1].TWI system works on RADAR principle [11]. Electromagnetic pulses of certain frequency are transmitted to get reflections from the target. The reflected electromagnetic pulses are received with some attenuation. Lower operating frequencies provide good penetration through walls (as compared to higher operating frequencies) but result in poor resolution and large antenna size. On the other hand, if we address detection through concrete walls, the upper frequency range is limited to 4 GHz [11].Image enhancement in TWI, has enjoyed an increasing interest over last few years [11–18]. Clutter and noise (due to antenna attenuation, cross talk, false targets, and wall reflections) result in degradation of image quality and appearance of false targets.Techniques for image enhancement in TWI includes, background subtraction [12], spatial filtering [13], wall parameter estimation/modeling based [14, 19], doppler domain filtering [17], image fusion [18, 20], and statistical methods [21–25].Main drawback of background subtraction technique is that it requires a surveillance mode of operation in which there is an access to the background (image scene that is free from targets) or reference [12, 25]. Spatial filtering relies on invariance of wall parameters (wall return remains same with changing antenna location). Moreover, this scheme works only for homogeneous (or near-homogeneous) walls at low operating frequencies [13, 25]. Limitation of wall parameter estimation/modeling-based approach is that it requires an accurate wall modeling and parameter estimation [14, 25]. Doppler domain filtering assumes that background is stationary and targets are moving [17]. Image fusion methods require multiple images of the same scene from different locations [18, 20].Statistical methods for TWI enhancement include: singular value decomposition (SVD), factor analysis (FA), principal component analysis (PCA) and independent Component Analysis (ICA) [21, 22]. Statistical methods (having less computational complexity) provide comparable results to other image enhancement methods. However these methods (SVD, PCA, FA, and ICA) require a priori knowledge of total number of targets and subjective (user defined) threshold value (for multiple targets).Information theoretic criterion (ITC) is a scheme used in array signal processing for determining number of target eigenvalues [26, 27]. It does not require subjective threshold setting and is never explored for TWI (to the best of author’s knowledge). In this paper, a TWI image enhancement method is proposed by modifying conventional SVD scheme. ITC is used to overcome subjective threshold requirement. SVD is chosen for its low complexity and simplicity over other methods (PCA, FA, and ICA). Proposed method successfully estimates total number of target singular values and results in improved image quality. Wavelet transform (WT) is used to further enhance image quality by estimating noise variance and thresholding WT coefficients. Visual inspection, mean square error (MSE) and peak signal to noise ratio (PSNR) comparison with conventional SVD is also provided.


## 2. Image Enhancement

### 2.1. Data Acquisition and Image Reconstruction

TWI setup is shown in Figure1 while geometrical representation of TWI is shown in Figure 2. Let H transceivers be placed (parallel to the x-axis) in the x-y plane. Image region (located beyond the wall along the positive y-axis) is divided into grid of M×N pixels (m=1,2,3…,M and n=1,2,3…,N). Let θ(t) be a wideband transmitted signal then pixel value at location mn can be computed by weighted sum and delay beamforming [1]. Beamforming algorithm align different input signals (in time and spatial domain) by applying time delays and weights. Output ζmn(t) for target located in x-y plane at pixel location mn is given by the following:
(1)ζmn(t)=∑p,qHξ(p,q)ϑ(t+τ^mn(p,q)),
where ξ(p,q) are weights (normally based on the Kaiser or Hamming window) used to control side lobes and τ^mn(p,q), are applied focusing delays and can be calculated by various methods depending on the available wall information [1]. Received signal ϑ(t) is delayed version of transmitted signal θ(t) with some attenuation αmn(p,q) that is, ϑmn(t)=αmn(p,q)θ(t-τmn(p,q)), where τmn(p,q) are time delays. Let θ^(t)=θ(-t) be a filter matched to transmitted signal then the deconvolved output for pixel mn, xmn is given as follows:
(2)xmn=(ζmn(t)*θ^(t))|t=0=(∑p,qHαmn(p,q)ξ(p,q)θ(t-τmn(p,q)+τ^mn(p,q))*θ^(t)∑p,qHαmn(p,q)ξ(p,q)θ(t-τmn(p,q)+τ^mn(p,q)))|t=0.
Above process is repeated for each pixel location mn to obtain B-scan image as follows:
(3)X=[x11x12⋯x1Nx22x22⋯x2N⋮⋮⋱⋮xM1xM2⋯xMN].Figure 1
Geometrical representation of TWI.Figure 2
2D wavelet transform.


### 2.2. SVD- and ITC-Based Image Enhancement

Image enhancement in TWI can be performed by decomposingB-scan X into different spectral components using singular value decomposition, that is,
(4)X=USVT=s1[⋮u1⋮][⋯v1T⋯]+⋯+sM[⋮uM⋮][⋯vMT⋯],
where (for simplicity M≤N), U=[u1u2⋯uM] and V=[v1v2⋯vN] having dimensions M×M and N×N are called unitary matrices and computed as left XXT and right XTX eigenvectors, respectively. Let S=diag(s1,s2,…,sM) with s1≥s2≥⋯≥sM≥0, are singular values of X.ImageX can be decomposed into three spectral images (clutter, target, and noise). The spectral images Xm of X are,
(5)X=∑m=1MXm=∑m=1MsmumvmT=∑m=1k1smumvmT+∑m=k1+1k2smumvmT+∑m=k2+1MsmumvmT,
where the first k1 singular values that belong to wall clutters followed by k2-k1 singular values belong to target(s) and rest singular values represent noise. Verma et al. in [21] state that k1=1 for wall clutter and k2=2 for target subspaces (i.e., XSVD=s2u2v2T), and rest subspaces represents noise. However, we note that this statement (k2=2) is not true in case of multiple targets. In fact, the target subspace can be more than one dimensional even when only a single target is present in the scene. The problem is reformulated as follows:
(6)X^=Y+Z=∑m=2MsmumvmT,
where X^ is clutter reduced image, Y is target image and Z is noise image. Since rank[Y]=k2<M, SVD of Y is [28, 29]:
(7)Y=[UY1UY2][SY1000][VY1TVY1T],
where, UY1, UY2, VY1, VY1 are unitary matrix (i.e., VY1VY1T=I) containing left and right singular vectors and SY1 is diagonal matrix containing singular value of Y. Therefore,
(8)X^=UY1SY1VY1T+Z=[UY1SY1+ZVY1ZVY2][VY1TVY2T]=[U^1U^2][S^100S^2][V^1TV^2T],
where, σZ is noise variance, U^1=(UY1SY1+ZVY1)(S^Y12+σZIk2)-1/2, U^2=ZVY2, S^1=SY12+σZIk2, S^2=σZIM-k2, V^1T=VY1T and V^2T=VY2T The original image X in terms of wall clutter, target, and noise subspaces is:
(9)X=[u1U^1U^2][s1S^1S^2][v1TV^1TV^2T],
where S^1 and S^2 are diagonal matrices containing singular values of target (s2,s3,…sk2) and noise (sk2+1,sk2+2,…sM), respectively. Note that sk2>sk2+1 and sk2+1≃sk2+2≃⋯sM≃σZ. Therefore, some statistical analysis needs to be performed in order to determine value for k2. In this regard some schemes in literature include difference of singular values (sm-sm+1), ratio of singular values (sm/sm+1) and percent of total power in an singular value (sm/tr[X]) [26]. However, these schemes do not always provide satisfactory results (and sometimes do not overcome the requirement of user defined threshold). It is observed that difference between noise singular values are relatively smaller than target singular values. To accurately determine the number of target singular values ITC methods are explored. These methods (ITC) do not require knowledge of a threshold value. Two well-known ITC methods are akaike information Criterion (AIC) and minimum description length (MDL). AIC and MDL utilize measures of the relative cross entropy between target and noise singular values [27]. AIC criterion is given as follows:
(10)AIC(k2)=Nln⁡[[1/(M-(k2+1))∑m=k2+1Msm]M-(k2+1)∏m=k2+1Msm]+k2(2M-k2).
Number of target singular values is determined by value of k2 for which AIC is minimized, that is, k2AIC=argmink2{AIC(k2)} as follows:
(11)XtarAIC=∑m=2k2AICsmumvmT.
Conceptually, MDL criterion is same as AIC criterion and simply adds a correction term in the minimization function as follows:
(12)MDL(k2)=Nln⁡[[1/(M-(k2+1))∑m=k2+1Msm]M-(k2+1)∏m=k2+1Msm]+12k2(2M-k2)ln⁡N.
Number of target singular values is determined by value of k2 for which MDL criterion is minimized that is, k2MDL=argmink2{MDL(k2)} as follows:
(13)XtarMDL=∑m=2k2MDLsmumvmT.
AIC and MDL works on the inequality of geometric and arithmetic mean of singular values as follows:
(14)1M-(k2+1)∑m=k2+1Msm≥(∏m=k2+1M‍sm)1/(M-(k2+1)).
In (14), equality holds only if sk2+1=sk2+2=⋯=sM. As a consequence AIC and MDL is minimized. Note here that for noise only singular values sk2+1≃sk2+2≃⋯≃sM≃σZ.


### 2.3. SVD- and WT-Based Denoising

WT (multiresolution analysis) localize image (in both space and scale) using scaled and translated copies of a finite-length waveform (mother wavelet). Image on larger scales provide gross features while small scales provide detail features. WT have advantages over other (spatial and fourier) transforms for: accurate representation of functions (having discontinuities and sharp peaks); data compression; noise reduction; probability density function estimation [30, 31]. In WT image is represented by linear combination of the wavelet functions and coefficients. Different operations (such as truncation of wavelet coefficients below a certain threshold are used for data compression and noise reduction) may be performed on the corresponding wavelet coefficients [30–33]. For overlapping boundaries of noise and target signals WT denoising is proposed after ITC for more enhancement of target(s). WT performs better denoising (since it performs denoising without smoothing out the edges) as compared to other filtering techniques (spatial filtering, fourier filtering, etc.) [30].WT for imageXMDL is given as follows:
(15)Wφ(j0,a,b)=1MN∑m=0M-1‍∑n=0N-1XMDL(m,n)φj0,a,b(m,n),Wψi(j,a,b)=1MN∑m=0M-1‍∑n=0N-1XMDL(m,n)ψj,a,bi(m,n),
where i=H,V,D and Wφ, WψH, WψV, WψD, are approximation, horizontal details, vertical details, and diagonal detail coefficients, respectively. φ and ψ are wavelet and scaling coefficients, respectively, and are determined by type of mother wavelet. Generally, j≥j0, where j0 is minimum resolution level and j=0,1,2,…,J-1, J=⌊log2(min(M,N))⌋, a,b=0,1,2,…,2j-1.Figure2 shows 2D iterative WT while Figure 3 shows wavelet decomposition structure. Noise reduction using WT is perform by thresholding detail coefficients at certain decomposition level.Figure 3
2D wavelet decomposition.


#### 2.3.1. Wavelet Parameters Selection

Various mother wavelet functions (Daubechies, Haar, Maxican Hat, Symlets, Morlet, etc.) may be used to calculate WT. These wavelets are different due to their complexity, accuracy, and time frequency analysis. Wavelet and scaling coefficients of Daubechies wavelet is shown in Figure4. Mother wavelet function, filter order, and decomposition level are important parameters for achieving better wavelet denoising results. Filter order for better performance (smoothness, localization of targets, and robustness against noise) is high for high and low for low resolution images. Generally, images in TWI are low resolution so a low filter order is suitable to avoid target delocalization and blurring effects. Decomposition level is also selected by same ideology as that of filter order [30, 32, 33].Daubechies wavelet and scaling coefficients.
(a)
Wavelet coefficients
(b)
Scaling coefficients


#### 2.3.2. Threshold Selection

Thresholding techniques are broadly classified into soft and hard [32]. Hard thresholding preserve edges but is less accurate for image denoising [32]. On the other hand soft thresholding reduces noise significantly but is not good in preserving edges [32]. Hard thresholding is
(16)W^ψi(j)={Wψi(j)if‖Wψi‖>ρ0otherwise,
and soft thresholding is
(17)W^ψi={sign⁡(Wψi)(‖Wψi(j)‖-ρ)if‖Wψi‖>ρ0otherwise,
where W^ψi thresholded coefficients, ρ is threshold value and can be determined by various methods including universal or global thresholding, Visu Shrink, Sure Shrink, Bayes Shrink, and so forth [32, 33]. We have used universal thresholding method for calculating value of ρ. 										(18)ρ=2σ^Zlog⁡10(MN).
Noise variance σ^Z is estimated using k2MDL obtained from MDL, that is,
(19)σ^Z=1M1∑m=k2MDL+1Msm,
where M1=M-k2MDL-1 represents number of singular values belonging to noise. Enhanced image XWT,MDL is obtained by performing inverse WT on the truncated (thresholded) coefficients as follows:
(20)XWT,MDL(m,n)=1MN∑a,b‍Wφ(j0,a,b)φj0,a,b(m,n)+1MN∑i=H,V,D∑j=j0∞∑a,bW^ψ(j,a,b)ψj,a,bi(m,n).


### 2.4. Performance of SVD and ITC

ITC (AIC and MDL) is applied on singular values ofB-scan image. SVD will change with change in the input image, consequently it changes ITC (AIC and MDL) graphs. It is important to note that, its not only the number of targets that determines the value of k2, rather it is the location, shape, size, and reflectiveness of targets also. This opens a new research direction. In some cases a single target may depend on more than one singular value. On the other hand it is also possible that one singular value represents more than one target. Distribution of singular value changes (due to change in target location, size and shape) AIC and MDL plots and consequently the value of k2. The image formed by using k2 (obtained using AIC or MDL) exactly extracts all targets.


## 2.1. Data Acquisition and Image Reconstruction

TWI setup is shown in Figure1 while geometrical representation of TWI is shown in Figure 2. Let H transceivers be placed (parallel to the x-axis) in the x-y plane. Image region (located beyond the wall along the positive y-axis) is divided into grid of M×N pixels (m=1,2,3…,M and n=1,2,3…,N). Let θ(t) be a wideband transmitted signal then pixel value at location mn can be computed by weighted sum and delay beamforming [1]. Beamforming algorithm align different input signals (in time and spatial domain) by applying time delays and weights. Output ζmn(t) for target located in x-y plane at pixel location mn is given by the following:
(1)ζmn(t)=∑p,qHξ(p,q)ϑ(t+τ^mn(p,q)),
where ξ(p,q) are weights (normally based on the Kaiser or Hamming window) used to control side lobes and τ^mn(p,q), are applied focusing delays and can be calculated by various methods depending on the available wall information [1]. Received signal ϑ(t) is delayed version of transmitted signal θ(t) with some attenuation αmn(p,q) that is, ϑmn(t)=αmn(p,q)θ(t-τmn(p,q)), where τmn(p,q) are time delays. Let θ^(t)=θ(-t) be a filter matched to transmitted signal then the deconvolved output for pixel mn, xmn is given as follows:
(2)xmn=(ζmn(t)*θ^(t))|t=0=(∑p,qHαmn(p,q)ξ(p,q)θ(t-τmn(p,q)+τ^mn(p,q))*θ^(t)∑p,qHαmn(p,q)ξ(p,q)θ(t-τmn(p,q)+τ^mn(p,q)))|t=0.
Above process is repeated for each pixel location mn to obtain B-scan image as follows:
(3)X=[x11x12⋯x1Nx22x22⋯x2N⋮⋮⋱⋮xM1xM2⋯xMN].Figure 1
Geometrical representation of TWI.Figure 2
2D wavelet transform.


## 2.2. SVD- and ITC-Based Image Enhancement

Image enhancement in TWI can be performed by decomposingB-scan X into different spectral components using singular value decomposition, that is,
(4)X=USVT=s1[⋮u1⋮][⋯v1T⋯]+⋯+sM[⋮uM⋮][⋯vMT⋯],
where (for simplicity M≤N), U=[u1u2⋯uM] and V=[v1v2⋯vN] having dimensions M×M and N×N are called unitary matrices and computed as left XXT and right XTX eigenvectors, respectively. Let S=diag(s1,s2,…,sM) with s1≥s2≥⋯≥sM≥0, are singular values of X.ImageX can be decomposed into three spectral images (clutter, target, and noise). The spectral images Xm of X are,
(5)X=∑m=1MXm=∑m=1MsmumvmT=∑m=1k1smumvmT+∑m=k1+1k2smumvmT+∑m=k2+1MsmumvmT,
where the first k1 singular values that belong to wall clutters followed by k2-k1 singular values belong to target(s) and rest singular values represent noise. Verma et al. in [21] state that k1=1 for wall clutter and k2=2 for target subspaces (i.e., XSVD=s2u2v2T), and rest subspaces represents noise. However, we note that this statement (k2=2) is not true in case of multiple targets. In fact, the target subspace can be more than one dimensional even when only a single target is present in the scene. The problem is reformulated as follows:
(6)X^=Y+Z=∑m=2MsmumvmT,
where X^ is clutter reduced image, Y is target image and Z is noise image. Since rank[Y]=k2<M, SVD of Y is [28, 29]:
(7)Y=[UY1UY2][SY1000][VY1TVY1T],
where, UY1, UY2, VY1, VY1 are unitary matrix (i.e., VY1VY1T=I) containing left and right singular vectors and SY1 is diagonal matrix containing singular value of Y. Therefore,
(8)X^=UY1SY1VY1T+Z=[UY1SY1+ZVY1ZVY2][VY1TVY2T]=[U^1U^2][S^100S^2][V^1TV^2T],
where, σZ is noise variance, U^1=(UY1SY1+ZVY1)(S^Y12+σZIk2)-1/2, U^2=ZVY2, S^1=SY12+σZIk2, S^2=σZIM-k2, V^1T=VY1T and V^2T=VY2T The original image X in terms of wall clutter, target, and noise subspaces is:
(9)X=[u1U^1U^2][s1S^1S^2][v1TV^1TV^2T],
where S^1 and S^2 are diagonal matrices containing singular values of target (s2,s3,…sk2) and noise (sk2+1,sk2+2,…sM), respectively. Note that sk2>sk2+1 and sk2+1≃sk2+2≃⋯sM≃σZ. Therefore, some statistical analysis needs to be performed in order to determine value for k2. In this regard some schemes in literature include difference of singular values (sm-sm+1), ratio of singular values (sm/sm+1) and percent of total power in an singular value (sm/tr[X]) [26]. However, these schemes do not always provide satisfactory results (and sometimes do not overcome the requirement of user defined threshold). It is observed that difference between noise singular values are relatively smaller than target singular values. To accurately determine the number of target singular values ITC methods are explored. These methods (ITC) do not require knowledge of a threshold value. Two well-known ITC methods are akaike information Criterion (AIC) and minimum description length (MDL). AIC and MDL utilize measures of the relative cross entropy between target and noise singular values [27]. AIC criterion is given as follows:
(10)AIC(k2)=Nln⁡[[1/(M-(k2+1))∑m=k2+1Msm]M-(k2+1)∏m=k2+1Msm]+k2(2M-k2).
Number of target singular values is determined by value of k2 for which AIC is minimized, that is, k2AIC=argmink2{AIC(k2)} as follows:
(11)XtarAIC=∑m=2k2AICsmumvmT.
Conceptually, MDL criterion is same as AIC criterion and simply adds a correction term in the minimization function as follows:
(12)MDL(k2)=Nln⁡[[1/(M-(k2+1))∑m=k2+1Msm]M-(k2+1)∏m=k2+1Msm]+12k2(2M-k2)ln⁡N.
Number of target singular values is determined by value of k2 for which MDL criterion is minimized that is, k2MDL=argmink2{MDL(k2)} as follows:
(13)XtarMDL=∑m=2k2MDLsmumvmT.
AIC and MDL works on the inequality of geometric and arithmetic mean of singular values as follows:
(14)1M-(k2+1)∑m=k2+1Msm≥(∏m=k2+1M‍sm)1/(M-(k2+1)).
In (14), equality holds only if sk2+1=sk2+2=⋯=sM. As a consequence AIC and MDL is minimized. Note here that for noise only singular values sk2+1≃sk2+2≃⋯≃sM≃σZ.


## 2.3. SVD- and WT-Based Denoising

WT (multiresolution analysis) localize image (in both space and scale) using scaled and translated copies of a finite-length waveform (mother wavelet). Image on larger scales provide gross features while small scales provide detail features. WT have advantages over other (spatial and fourier) transforms for: accurate representation of functions (having discontinuities and sharp peaks); data compression; noise reduction; probability density function estimation [30, 31]. In WT image is represented by linear combination of the wavelet functions and coefficients. Different operations (such as truncation of wavelet coefficients below a certain threshold are used for data compression and noise reduction) may be performed on the corresponding wavelet coefficients [30–33]. For overlapping boundaries of noise and target signals WT denoising is proposed after ITC for more enhancement of target(s). WT performs better denoising (since it performs denoising without smoothing out the edges) as compared to other filtering techniques (spatial filtering, fourier filtering, etc.) [30].WT for imageXMDL is given as follows:
(15)Wφ(j0,a,b)=1MN∑m=0M-1‍∑n=0N-1XMDL(m,n)φj0,a,b(m,n),Wψi(j,a,b)=1MN∑m=0M-1‍∑n=0N-1XMDL(m,n)ψj,a,bi(m,n),
where i=H,V,D and Wφ, WψH, WψV, WψD, are approximation, horizontal details, vertical details, and diagonal detail coefficients, respectively. φ and ψ are wavelet and scaling coefficients, respectively, and are determined by type of mother wavelet. Generally, j≥j0, where j0 is minimum resolution level and j=0,1,2,…,J-1, J=⌊log2(min(M,N))⌋, a,b=0,1,2,…,2j-1.Figure2 shows 2D iterative WT while Figure 3 shows wavelet decomposition structure. Noise reduction using WT is perform by thresholding detail coefficients at certain decomposition level.Figure 3
2D wavelet decomposition.


### 2.3.1. Wavelet Parameters Selection

Various mother wavelet functions (Daubechies, Haar, Maxican Hat, Symlets, Morlet, etc.) may be used to calculate WT. These wavelets are different due to their complexity, accuracy, and time frequency analysis. Wavelet and scaling coefficients of Daubechies wavelet is shown in Figure4. Mother wavelet function, filter order, and decomposition level are important parameters for achieving better wavelet denoising results. Filter order for better performance (smoothness, localization of targets, and robustness against noise) is high for high and low for low resolution images. Generally, images in TWI are low resolution so a low filter order is suitable to avoid target delocalization and blurring effects. Decomposition level is also selected by same ideology as that of filter order [30, 32, 33].Daubechies wavelet and scaling coefficients.
(a)
Wavelet coefficients
(b)
Scaling coefficients


### 2.3.2. Threshold Selection

Thresholding techniques are broadly classified into soft and hard [32]. Hard thresholding preserve edges but is less accurate for image denoising [32]. On the other hand soft thresholding reduces noise significantly but is not good in preserving edges [32]. Hard thresholding is
(16)W^ψi(j)={Wψi(j)if‖Wψi‖>ρ0otherwise,
and soft thresholding is
(17)W^ψi={sign⁡(Wψi)(‖Wψi(j)‖-ρ)if‖Wψi‖>ρ0otherwise,
where W^ψi thresholded coefficients, ρ is threshold value and can be determined by various methods including universal or global thresholding, Visu Shrink, Sure Shrink, Bayes Shrink, and so forth [32, 33]. We have used universal thresholding method for calculating value of ρ. 										(18)ρ=2σ^Zlog⁡10(MN).
Noise variance σ^Z is estimated using k2MDL obtained from MDL, that is,
(19)σ^Z=1M1∑m=k2MDL+1Msm,
where M1=M-k2MDL-1 represents number of singular values belonging to noise. Enhanced image XWT,MDL is obtained by performing inverse WT on the truncated (thresholded) coefficients as follows:
(20)XWT,MDL(m,n)=1MN∑a,b‍Wφ(j0,a,b)φj0,a,b(m,n)+1MN∑i=H,V,D∑j=j0∞∑a,bW^ψ(j,a,b)ψj,a,bi(m,n).


## 2.3.1. Wavelet Parameters Selection

Various mother wavelet functions (Daubechies, Haar, Maxican Hat, Symlets, Morlet, etc.) may be used to calculate WT. These wavelets are different due to their complexity, accuracy, and time frequency analysis. Wavelet and scaling coefficients of Daubechies wavelet is shown in Figure4. Mother wavelet function, filter order, and decomposition level are important parameters for achieving better wavelet denoising results. Filter order for better performance (smoothness, localization of targets, and robustness against noise) is high for high and low for low resolution images. Generally, images in TWI are low resolution so a low filter order is suitable to avoid target delocalization and blurring effects. Decomposition level is also selected by same ideology as that of filter order [30, 32, 33].Daubechies wavelet and scaling coefficients.
(a)
Wavelet coefficients
(b)
Scaling coefficients


## 2.3.2. Threshold Selection

Thresholding techniques are broadly classified into soft and hard [32]. Hard thresholding preserve edges but is less accurate for image denoising [32]. On the other hand soft thresholding reduces noise significantly but is not good in preserving edges [32]. Hard thresholding is
(16)W^ψi(j)={Wψi(j)if‖Wψi‖>ρ0otherwise,
and soft thresholding is
(17)W^ψi={sign⁡(Wψi)(‖Wψi(j)‖-ρ)if‖Wψi‖>ρ0otherwise,
where W^ψi thresholded coefficients, ρ is threshold value and can be determined by various methods including universal or global thresholding, Visu Shrink, Sure Shrink, Bayes Shrink, and so forth [32, 33]. We have used universal thresholding method for calculating value of ρ. 										(18)ρ=2σ^Zlog⁡10(MN).
Noise variance σ^Z is estimated using k2MDL obtained from MDL, that is,
(19)σ^Z=1M1∑m=k2MDL+1Msm,
where M1=M-k2MDL-1 represents number of singular values belonging to noise. Enhanced image XWT,MDL is obtained by performing inverse WT on the truncated (thresholded) coefficients as follows:
(20)XWT,MDL(m,n)=1MN∑a,b‍Wφ(j0,a,b)φj0,a,b(m,n)+1MN∑i=H,V,D∑j=j0∞∑a,bW^ψ(j,a,b)ψj,a,bi(m,n).


## 2.4. Performance of SVD and ITC

ITC (AIC and MDL) is applied on singular values ofB-scan image. SVD will change with change in the input image, consequently it changes ITC (AIC and MDL) graphs. It is important to note that, its not only the number of targets that determines the value of k2, rather it is the location, shape, size, and reflectiveness of targets also. This opens a new research direction. In some cases a single target may depend on more than one singular value. On the other hand it is also possible that one singular value represents more than one target. Distribution of singular value changes (due to change in target location, size and shape) AIC and MDL plots and consequently the value of k2. The image formed by using k2 (obtained using AIC or MDL) exactly extracts all targets.


## 3. Simulation and Results

Experimental setup (constructed using [4]) for TWI shown in Figure 5 (physical elements of experimental setup are shown in Figure 6). Agilents’s vector network analyzer (VNA) in the range of 300KHz to 3GHz (1GHz band width (BW)) is used to generates a stepped frequency 2 GHz-3 GHz waveform having step size Δf=5 MHz and Nf=201. Maximum range Rmax is calculated as follows:
(21)Rmax⁡=c(Nf-1)2BW=30m.
The range resolution ΔR is
(22)ΔR=c2NfΔf=0.37m.
Directional and broadband horn antenna with 12dB gain is used in monostatic mode (for transmitting and receiving signals). Antenna is mounted on 2D-scanning frame (having dimensions width 2.4m and height 3m) which can slide along cross range and height. Rear and side walls are covered with pyramidal radar absorbable modules. Scanning is controlled by microcontroller and at each point scattering parameters (magnitude and phase) are recorded by VNA and transferred to local computer. Wood wall is constructed having thickness 5cm, relative permittivity (approximately) equals to 2.3 and relative permeability (approximately) equals to 1. The antenna is positioned 0.03m from the wall. Received data is converted from frequency domain to time domain using inverse fourier transform. Time delays and weights are fed into beamforming algorithm for image reconstruction.Figure 5
TWI setup.Physical elements of experimental setup at Microwave Engineering Laboratory, College of Signals, NUST.
(a)
(b)Image enhancement algorithms based on conventional SVD and proposed schemes are simulated in MATLAB. For wavelet denoising we have used Daubechies wavelet, fourth order filter, third level decomposition and soft thresholding technique. Background subtracted imageXbs is constructed using the difference of two images (i.e., image with target and image without target) [4, 12]. This background image is used as a comparison measure for proposed and existing algorithm. Simulation results are compared on the basis of MSE, PSNR, and visual inspection as follows:
(23)MSE=1M×N∑m=1M‍∑n=1N(Xbs(m,n)-Xtar(m,n))2,PSNR(dB)=10log⁡101MSE,
where, Xtar∈{XSVD,XAIC,XMDL,XWT,MDL}.Example 1.
A single metallic target is placed approximately0.8m away from wall. Figure 7 shows different analysis schemes versus number of singular values. Figure 8 shows different spectral components (ranging from one to nine) for single target. It is observed that more than one spectral component contains target image. Note that predicted subspace dimension is one and k2AIC=k2MDL=3 from Figure 7 which leads to estimated target subspace dimension equals to two (k2AIC-k1=k2MDL-k1=2). An explanation for this value of k2 appears in Section 2.4. Figure 9 shows the performance of conventional SVD and proposed (AIC-based SVD, MDL-based SVD and WT MDL-based SVD) schemes for single target. Note that conventional and proposed schemes provide comparable results.Example1: Single target (a) amplitude of singular values (sm), (b) difference of singular values (sm-sm+1), (c) ratio of singular values (sm/sm+1), (d) power in singular values (sm/tr[X]), (e) AIC function of singular values (10), and (f) MDL function of singular values (12).
(a)
(b)
(c)
(d)
(e)
(f)Example1: Single target: different spectral components.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Example1: Single target (a) original image, (b) conventional SVD, (c) AIC-based SVD, (d) MDL-based SVD, (e) WT MDL-based SVD, and (f) background subtracted image.
(a)
(b)
(c)
(d)
(e)
(f)Example 2.
Three metallic targets are placed approximately0.15m, 1.25m, and 1.5m away from wall. Figure 10 shows different analysis schemes versus number of singular values. Figure 11 shows different spectral components (ranging from one to nine) for three targets. It is observed that more than one spectral component contains target image. Note that predicted subspace dimension is three for three targets and k2AIC=k2MDL=5 from Figure 11 which leads to estimated target subspace dimension equals to four (k2AIC-k1=k2MDL-k1=4). An explanation for this value of k2 appears in Section 2.4. Figure 12 shows the performance of conventional SVD and proposed (AIC-based SVD, MDL-based SVD, and WT MDL-based SVD) schemes for multiple (three) targets. Note that proposed schemes provide better results as compared to the conventional SVD scheme and detects all three targets, where as conventional SVD scheme fails to detects all targets.Example2: Multiple (three) targets (a) amplitude of singular values (sm), (b) difference of singular values (sm-sm+1), (c) ratio of singular values (sm/sm+1), (d) power in singular values (sm/tr[X]), (e) AIC function of singular values (10), and (f) MDL function of singular values (12).
(a)
(b)
(c)
(d)
(e)
(f)Example2: Multiple (three) targets: different spectral components.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Example2: Multiple (three) targets (a) original image, (b) conventional SVD, (c) AIC-based SVD, (d) MDL-based SVD, (e) WT MDL-based SVD, and (f) background subtracted image.
(a)
(b)
(c)
(d)
(e)
(f)Example 3.
Five metallic targets are placed approximately0.15m, 0.3m, 0.8m, 1.4m, and 1.7m away from wall. Figure 13 shows different analysis schemes versus number of singular values. Figure 14 shows different spectral components (ranging from one to nine) for five targets. It is observed that more than one spectral component contains target image. Note that predicted subspace dimension is five and k2AIC=6 and k2MDL=7 from Figure 13 which leads to estimated target subspace dimension equals to five by AIC (k2AIC-k1=5) and six by MDL (k2MDL-k1=6). An explanation for this value of k2 appears in Section 2.4. Figure 15 shows the performance of conventional SVD and proposed (AIC-based SVD, MDL-based SVD and WT MDL-based SVD) schemes for multiple (five) targets. Note that proposed schemes provide better results as compared to the conventional SVD scheme and detects all five targets, where as conventional SVD scheme fails to detects all targets.Example3: Multiple (five) targets, (a) amplitude of singular values (sm), (b) difference of singular values (sm-sm+1), (c) ratio of singular values (sm/sm+1), (d) power in singular values (sm/tr[X]), (e) AIC function of singular values (10), and (f) MDL function of singular values (12).
(a)
(b)
(c)
(d)
(e)
(f)Example3: Multiple (five) targets: different spectral components.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Example3: Multiple (five) targets (a) original image, (b) conventional SVD, (c) AIC-based SVD, (d) MDL-based SVD, (e) WT MDL-based SVD, and (f) background subtracted image.
(a)
(b)
(c)
(d)
(e)
(f)Example 4.
Five metallic targets are placed approximately0.15m, 0.15m, 0.15m, 0.8m, and 1.6m away from wall. Figure 16 shows different analysis schemes versus number of singular values. Figure 17 shows different spectral components (ranging from one to nine) for five targets. It is observed that more than one spectral component contains target image. Note that predicted subspace dimension is five for five targets and k2AIC=k2MDL=5 from Figure 16 which leads to estimated target subspace dimension equals to four (k2AIC-k1=k2MDL-k1=4). An explanation for this value of k2 appears in Section 2.4. Figure 18 shows the performance of conventional SVD and proposed (AIC based SVD, MDL based SVD and WT MDL based SVD) schemes for multiple (five) targets. Note that proposed schemes provide better results as compared to the conventional SVD scheme and detects all five targets, where as conventional SVD scheme fails to detects all targets.
Note that, proposed (AIC based SVD, MDL based SVD and WT MDL based SVD) schemes successfully detects all targets and cause some target spread in cross-range direction. This is due to the fact that we are only consideringk2-k1 subspaces (instead of N subspaces). This spread may easily be removed by applying some intensity thresholding technique.
Figure19 summarizes MSE and PSNR plots for different number of targets placed behind wall. For single target conventional SVD and proposed schemes have comparable results but as number of target increases proposed schemes (AIC based SVD, MDL based SVD and WT MDL based SVD) perform significantly better (as compared to conventional SVD scheme).Example4: Multiple (five) targets (a) amplitude of singular values (sm), (b) difference of singular values (sm-sm+1), (c) ratio of singular values (sm/sm+1), (d) power in singular values (sm/tr[X]), (e) AIC function of singular values (10), and (f) MDL function of singular values (12).
(a)
(b)
(c)
(d)
(e)
(f)Example4: Multiple (five) targets: different spectral components.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)Example4: Multiple (five) targets (a) original image, (b) conventional SVD, (c) AIC-based SVD, (d) MDL-based SVD, (e) WT MDL-based SVD, and (f) background subtracted image.
(a)
(b)
(c)
(d)
(e)
(f)MSE and PSNR comparison of conventional and proposed schemes against number of targets.
(a)
MSE Comparison
(b)
PSNR Comparison


## 4. Conclusion

SVD-based image enhancement for TWI using ITC and WT is proposed. The scheme is capable of discriminating between target and noise subspaces. The limitation of subjective threshold setting in conventional SVD is overcome using ITC scheme. Furthermore, noise is suppressed using WT-based denoising. Proposed method increases accuracy of conventional SVD-based TWI image enhancement scheme. Both AIC and MDL technique detect multiple targets (so one can use either of these). Proposed scheme can easily be modified for PCA, FA, and ICA methods to get better accuracy.



---

*Source: 961829-2012-07-30.xml*