A Cross Correlation based Back-Projection Imaging Method for Through-Wall Imaging

A back-projection (BP) algorithm based on cross-correlation (CBP), which is both fast and good at suppressing artifacts, is adopted here for through wall imaging (TWI). Firstly, by “moving” the wall to the image point, the three layer TWI modal is equivalent to the two layer TWI modal. By doing so, only half of the reflection points are needed to be calculated; secondly, to further reduce the calculation time, on the one hand, the symmetric-invariant property as well as the shift-invariant property of the ray propagation is used to calculate the look-up table, on the other hand, the delay-multiply-sum equation in CBP is rewrite into a new form with less addition and multiply computation. Simulation results show the superiority of the proposed CBP algorithm both in operation speed and artifacts suppression.


Introduction
Through-wall imaging (TWI) radar is a very promising technique for detecting and identifying the unknown objects hidden behind the walls and it has become a hot research topic in recent years for it is urgently used in counterterrorism, urban warfare and calamity rescue and so on [1][2][3][4][5][6].
Many imaging algorithms have been proposed for through-wall imaging. In [1], a compressive sensing based imaging method is proposed for high-resolution ultra-wide band TWI radar. Though it can save the data collection time by using much less data than Nyquist sampling rule, a L1 norm optimization problem must be solved, which is time consuming. Inverse scattering method by differential evolution [2] can reconstruct the shape of the targets behind the wall, but the forward scattering problem must be solved in each cost function evaluation. Learning by example method [3] can image the target in real time when the support vector machine is properly trained, but the construction of a proper training data set needs much prior knowledge about the scene. The subspace optimization method [4] is quite tolerant with noise, but it is very difficult to extend to threedimension imaging for the huge amount of the layered green function calculation. The autofocusing method by time reverse mirror proposed in [5] is capable of imaging the targets without knowing the wall parameters, but still the layered green function needs to be calculated. The tomographic diffraction method [6] is suitable for on-site applications as FFT/IFFT is implemented, but it is unsuitable for strong dielectric scatters or metallic targets as it is based on Born approximation.
Because the refraction phenomenon is seriously considered and the distortion of the wave path is accurately compensated, the back projection (BP) algorithm [7-10] is one of the most practical imaging methods for its convenience and robustness, particularly when the imaging scene can be modeled as layered mediums. However, there are two main disadvantages of the BP algorithm. The first one is the high computation burden. In order to calculate the positions of refraction points, a quartic equation should be solved for every point behind the wall to be imaged. It takes lots of time to ISEITCE 2020 Journal of Physics: Conference Series 1607 (2020) 012020 IOP Publishing doi:10.1088/1742-6596/1607/1/012020 2 solve the equation and repeat it all over the imaging scene. The second disadvantage is that there are too much artifacts with high energy in the imaging result. The existence of artifacts decreases the contrast between the targets and the background, which will affect the following target identifying work. In paper [10], a fast back projection algorithm based on cross correlation (CBP) is proposed to deal with these two disadvantages for the ground penetrating radar imaging (GPR). To save the computation, on the one hand, an approved approximation method is proposed to calculate the intersection of the wave path; on the other hand, a lookup table is also used to reduce the redundancy in BP algorithm. What is more, an additional cross-correlation procedure in the original algorithm to suppress the artifacts.
In this paper, these three ideas in CBP for GPR are introduced to TWI. The work lies in two aspects: firstly, by "moving" the wall to the measurement line, the three layer TWI modal is equivalent to the two layer TWI modal which is quite similar to GPR imaging, and thus these ideas can be adopted in TWI; secondly, to further reduce the calculation time, the delay-multiply-sum equation in CBP is rewrite into a new form with less multiply computation. The simulation results show the work is valid.

Classic BP algorithm
The 2-D imaging configuration of the TWI system is shown in Fig.1(a). The scene is divided into three regions by the wall. The upper region and the lower region are the air with the permittivity and permeability , where ε and μ are the corresponding parameters in free space. The center region is homogeneous wall with the thickness of , permittivity and permeability , where is the relative permittivity of the wall.
(b). Fig.1. the geometry of the TWI: (a). Three-layer TWI modal; (b). Two-layer TWI modal which is equivalence to (a) when calculating the round-trip time delays. The antennas are in front of the wall and work in the bi-static way. They transmit and receive signals in each of the M synthetic aperture positions on the measurement line which parallels to the wall to collect the B-Scan data. At each synthetic aperture position, for a given point A with coordinates , in the imaging space behind the wall, the transmitting signal from the transmitter at , ℎ travels to , , with the first turning at the inflection point , 0 on the front wall surface and the second turning at the second inflection point , in the other wall surface, and returns along a different path with two tunings at the inflection point , and , 0 on each wall surface respectively to the receiver at , ℎ . These four inflection points can be exactly calculated by Snell's law. After finding them, the distances : can be calculated, and then the round-trip time delay from the imaging point A to the kth antenna position will be: where c is the wave propagation velocity in the air.
For a given point A in the imaging zone, the first step of BP algorithm is to calculate the round-trip time delays from A to each of the antenna pair positions . Suppose that at the Then, the value of point A in the final image can be formulated as: To get the final through wall image, all the points in the imaging space are considered and their image values are calculated by repeating aforementioned steps.

fast calculation of the round-trip time delay
From above description, it can be seen that the calculation of the round-trip time delay is of key importance for the BP algorithm. In this subsection, the fast calculation of the round-trip time delay is demonstrated as follow: Firstly, the equivalent two-layer TWI modal is introduced. For the three layered through wall modal, four inflection points must be calculated for each image point at each antenna position, which is difficult and time-consuming. To make it easier, we make an equivalent two-layer modal by "moving" the wall to such place where the back surface of the wall is on the imaging point (as shown in Fig.1(b)). Since the foreside of the wall is parallel to the backside of the wall, according to the ray theory view of electromagnetic wave propagation, this "moving" will not change the round-trip time delay. By doing so, only two inflection points , and , are needed to be calculate, which is half of the three-layer TWI modal.
Since the calculation of T is a quartic equation. Although the two layers TWI modal can reduce much computational burden than that of the three layers TWI modal, considering that the solving process must repeat 2M times for just a single image point in the imaging region, it is still very timeconsuming.
Secondly, to speed up the calculation of , a piecewise linear approximation method proposed in [10] as equation (4) shows is used here instead of solving the quartic equation. This approximation is accurate and fast for the calculation of refraction point for two layered ground penetrating radar (GPR).
Here the is the coordinate of the transmitter. Obviously, this can further reduce the computational burden than directly solving quartic equation. Since the two-layer TWI is quite similar to the two-layer GPR as showed above, we use (4) for TWI.  Thirdly, the idea of the redundant calculation elimination proposed in [10] is used to reduce the total times of the round-trip time delay calculation. Suppose that the imaging region is gridded into points ( 1 is an integer and the number of points in the range direction while N is the number of points in the depth direction). Since there are M antenna pair positions, the total of 2 inflection points and round-trip time delay must be calculated. As shown in Fig.2.(a), if we substitute and for and , respectively, in (4), the new solution of will be . That is to say, the relative position of to and is shift invariant with respect to . From this point of view, there are a lot of redundant round-trip time delay calculations in the classic BP algorithm. By using this shift invariant property, the total of 2 inflection points must be calculated to create the look-up table to record all the possible round-trip time delays. That is, computational burden of the round-trip time delay further reduces to 1/M.
As shown in Fig.2.(b), if the imaging points are symmetric with the antenna, the solution of the reflection points are also symmetric with the antenna. From this point of view, there are still a lot of redundant round-trip time delay calculations in CBP algorithm besides using the shift invariant property. By using this symmetric invariant property, the total of inflection points need to be calculated when calculate the possible round-trip time delay look-up table. That is, computational burden of the round-trip time delay further reduces 1/2.

Modified cross-correlated back-projection method
In [11], to suppress artifact, a cross-correlated back projection (CBP) algorithm was proposed instead of delay-and-sum beam forming in classic BP algorithm for the ground penetrating radar. That is: ∑ ∑ • (5) To save the computing, (5) is rewritten as: ∑ ∑ (6) For each (pixel), while in (5), 0.5 1 additional operations and0.5 1 multiply operations are needed, only 2 multiply operations and 2 2 multiply operations are used in (6). There are pixels, therefore much more computation will be saved by equation (6).

Simulation Results
To illustrate the proposed CBP imaging algorithm, the two dimensional imaging space of size 2m×2m containing two perfect electric conductor (PEC) targets behind a wall is considered. The imaging data are collected by GprMax2.0 software package [11] and the background subtraction method is used for imaging. In GprMax2.0, the parameters are set as follow: the imaging space is a rectangle zone with two coordinate points (0m,0m)and (2m,2m),one target is with the radius of 0.1m at P1(0.8m,0.8m) and the other is with the radius of 0.15m at P2(1.3m,1.3m).The wall is with the thickness of 0.3m and its relative permittivity is dielectric with 6.0 and conductivity 0.01 / .The length of the measurement line is 1.8 long and 0.3m before the wall. The distance between the transmitter and the receiver is 0.08m and the antenna pair moves along the measurement line with the step of 0.02m. The transmitter transmits ricker wave with the central frequency of 1GHz and the receiver receives the Ex component of the scattering field. The imaging space is discredited into 200×200 pixels. The simulations are done by matlab 2014a on win7 and 2.4GHz CPU.
To make a comparison, three imaging methods including the classic BP by equation (3), the CBP by equation (5) and the proposed CBP by equation (6) are used with the same data. The simulation result is shown in Fig.3. From Fig.3, it can be seen clearly that the two target can be imaged by all these methods. In Fig.3.(b), there are some artifacts, which may be regard as targets; in Fig.3.(c),the background of the image is very clear and the targets can be easily identified. The superiority in artifact suppression by the proposed CBP algorithm over the classic BP is obvious. Because the imaging method are the same,  Fig.4.
It can be seen clearly that for the dielectric targets, the CBP and the proposed CBP can suppress the artifacts much better than the classic BP. The simulation time for the classic BP algorithm is 214.8s, the CBP algorithm takes 6.10s and the proposed CBP takes 5.32s.
For TWI, the difference between CBP and the proposed CBP lies in two aspect: one is the time delay matrix computation. Except for the shift invariant property used in CBP, the symmetric invariant property is also used for the proposed CBP; the other is the imaging equation. Equation (6) is equal to equation (5), but the computation is much less. Therefore the image by CBP and the proposed CBP is the same, but the computation time is different.

5.Conclusion
This work deals with the imaging of the targets behind the wall by the back projection method based on cross correlation (CBP). By using the delay-multiply-sum idea and the square of the collected data, the artifacts, can be efficiently suppressed. Moreover, the much less time is used by the linear approximation when calculating the position of inflection point. The application of the proposed algorithm is practically valuable in the imaging of objects behind the wall with fast speed and high quality. Further work will exert on the application of CBP to the 3D TWI imaging with real data.