Fmcw Sparse Array Imaging and Restoration for Microwave Gauging

The application of imaging radar to microwave level gauging represents a prospect of increasing the reliability of target detection. The aperture size of the used sensor determines the underlying azimuthal resolution. In consequence , when FMCW-based multistatic radar (FMCW: frequency modulated continuous wave) is used, the number of antennas dictates this essential property of an imaging system. The application of a sparse array leads to an improvement of the azimuthal resolution by keeping the number of array elements constant with the cost of increased side lobe level. Therefore, ambiguities occur within the imaging process. This problem can be modelled by a point spread function (PSF) which is common in image processing. Hence, an inverse system to the imaging system is needed to restore unique information of existing targets within the observed radar scenario. In general, the process of imaging is of ill-conditioned nature and therefore appropriate algorithms have to be applied. The present paper first develops the degradation model, namely PSF, of an imaging system based on a uniform linear array in time domain. As a result, range and azimuth dimensions are interdependent and the process of imaging has to be reformulated in one dimension. Matrix-based approaches can be adopted in this way. The second part applies two computational methods to the given inverse problem, namely quadratic and non-quadratic regular-ization. Notably, the second one exhibits an ability to suppress ambiguities. This can be demonstrated with the results of both, simulations and measurements, and enables sparse array imaging to localize point targets more unambiguously.


Introduction
Industrial level gauging by means of microwaves offers a reliable and accurate method to determine the fill level of a reservoir.In most cases the microwave level gauge is a monostatic radar system.This becomes a problem when metallic fixtures like pipes or the enclosure impede detecting backscattering from volume surface.In order to improve the robustness of such systems introducing imaging by beamsteering is very promising.In addition, using more than one measurement dimension offers the possibility for estimating levels of non-planar surfaces.
Due to practical aspects, a switched approach using digital beam-forming (DBF) in conjunction with an antenna array will be preferred versus one mechanical panned antenna.In the following this method will be referred to as multistatic radar in which transmitting and receiving elements are colocated.Aiming to continue using the established concept of frequency modulated continuous wave (FMCW) with advantages, as for instance a straightforward RF hardware and a good range resolution (Stove, 1992), transmission of a single waveform at each antenna is supposed.
The ability to distinguish between two closely spaced targets is an essential characteristic of an imaging radar and depends on the aperture size used.In order to achieve a high angular resolution the array aperture must be large.Mostly, the number of antenna elements is a small number, meaning spacing of the antenna array becomes comparatively wide.As a result, such a sparse antenna array has a high side lobe response generating artefacts in the angular domain.
In general, such an imaging system with DBF can be represented by a convolution of the beam-pattern h ant and the underlying scene of targets x over azimuth ϕ (Gething, 1991).Therefore, the output of such an imaging system y is described by the well-known convolution expression in Eq. ( 1) inside the observed area .
In consequence, impulse response h ant is the origin of ambiguities and loss of resolution within a linear imaging process.A comparable issue is known in theory of image processing in which degradation of a recorded image is described by a point spread function (PSF) (e.g.Bertero et al., 1998).In order to improve the quality of an image, the concept of deconvolution as an inverse process of imaging is used.This was already applied to enhance angular resolution within an imaging system with mechanical beam-steering (Kolb et al., 2010).The present paper utilizes existing methods of image processing to reduce degenerating and ambiguous effects of a sparse array imaging system.The general method is first to express the forward problem which describes how to get the measurable data of an observed radar scene (Bertero et al., 1998).This is shown in Sect. 2 with respect to an antenna array with DBF.The second step deals with the corresponding inverse problem and its solution by regularization methods to restore unambiguously information of a radar scene.Two essential methods are applied to simulation data in Sect. 3 and measurement data in Sect. 4.

Direct problem: modelling FMCW sparse array imaging
Due to the architecture of a particular RF frontend, more or less DBF setups are feasible.For reasons of specification, a measurement matrix S can be defined that relates transmitting channels a and receiving channels b = S • a.  -only one element is able to transmit respectively receive (SIMO/MISO), all others do the opposite.
Without loss of generality, the present paper concentrates on full measurements (Reflection and Transmission) but all considerations could also be used for other modes of measurements.

Point spread function
Mostly, far-field conditions in gauging situations are fulfilled.In consequence, a plane wavefront (WF) approach, according to Fig. 1, is used to develop an expression of the PSF.Starting from an arbitrary transmitted spectrum X(f ) at each single element, a wavefront orthogonal to the beamdirection is generated.With ϕ t being the direction of a wavefront to and from a potential target, the contribution X m (f ) of the m-th element is given by: where p y,m denotes the y-coordinate of the m-th element and c 0 equals the speed of light, respectively.Looking at an incoming wavefront from transmitting antenna m, the resulting spectrum Y nm (f ) at receiving antenna n becomes By defining a transfer function as the ratio of received spectrum Y nm (f ) over transmitted spectrum X(f ) and combining Eqs.
(2), (3), the following expression is found: It is desired to generate a wavefront in direction ϕ t .Therefore, usually complex weights w nm (f,ϕ) are chosen complex conjugate of H nm (f ), to realize a delay-and-sum beamformer.
Assuming a linear array geometry of N antennas with equidistant elements (spacing d) and positioning the phase centre in the middle of the array (y = 0), the output of the beam-former becomes: respectively where the substitutions u = sinϕ and u t = sinϕ t are used to achieve a shift-invariant behaviour.Therefore, Eq. ( 6) formulates the two-dimensional array factor (AF) of an intended point target in direction u t = 0, which represents the relating transfer function.Consequently, a subsequent inverse Fourier transform af (r,u) = IFT{AF(f,u)} complies to the impulse response of a target in direction u t = 0 and range r t = 0.In case of Reflection and Transmission-variant and a point target at direction u t = 0 the expression results in: It is obvious that range-relating domain f and angularrelating domain u are not separable, meaning both dimensions are interdependent.This needs to be considered in the process of describing the image process.
In order to take account of all relevant frequency dependent effects, remaining components corresponding to Fig. 2 have to be considered.In detail, there are two-way free-space loss α ∝ 1 f 2 and a f -domain window H win (f ), possibly used to suppress ringing in range-domain.According to this description, the PSF in f -domain results in: The PSF in range-domain is finally obtained by an inverse Fourier transform psf(r,u) = IFT{PSF(f,u)}.In contrast to Eq. ( 1), the imaging process depends on range r in addition to the angular domain and has to be completed with the replaced two-dimensional impulse response psf(r,u) of the imaging system to the following expression: To get an idea of the PSF's effects, Fig. 3 illustrates psf(r,u) of an X-band (8.2 to 12.4 GHz) FMCW sparse array imaging system with N = 4 single elements and element spacing of d = 40 mm; seriously violating the half-wavelength rule ( λ center 2 ≈ 14.6 mm).In consequence of this spatial undersampling, ambiguities against azimuth occur.The applied frequency window is of a Hanning characteristic.

Matrix-vector formulation
For reasons of computational aspects, a finite discrete representation of Eq. ( 9) is needed.In this context there is a need to take respect of the interdependent character of both dimensions within the PSF.Hence, in image processing a one-dimensional representation (lexicographic ordering) of two-dimensional data is utilized (e.g.Hunt, 1973), because successively processing of occurring two-dimensional convolution is illegitimate.
The vector representation of targets versus range and azimuth x(r,u) ∈ C N r ×N u is reordered according to where N r is the length along range r and N u is the length along azimuth-relating u, respectively.Relating to the channel model in Fig. 2 an additive noise source n∈ C N r to Eq. ( 10) gives the output of the system.Finally, a suitable one-dimensional reformulation of the convolution kernel of psf(r,u) to K∈ C N r •N u ×N r •N u enables to write the twodimensional convolution of Eq. ( 9) by a simple matrix-vector multiplication: It is important to note that convolution kernel K is not of block-circulant character, meaning this convolution cannot be expressed in a corresponding frequency domain.Formally stated, the underlying convolution of Eq. ( 9) is not of a periodic character and there is also no other way to define any other boundaries like zero-or reflexive-boundary conditions (see e.g.Nagy et al., 2004).This can also be recognized in the example PSF in Fig. 3 because psf(r,u) does not become 0 for u → ±∞.

Inverse problem: image restoration
The declaration of Eq. ( 11) gives a description of the imaging behaviour of a FMCW sparse array radar system.With the aim to compute an estimation x of the original scene x from the observation z an inverse system is needed.First of all, this can be formulated by a matrix inversion: which can also be expressed by a least square (LS) expression (Bertero et al., 1998): The LS expression solves the minimization problem but will fail due to ignoring the noise source n.In addition, it is common knowledge in image processing that inverse mapping in terms of Eq. ( 11) is of ill-conditioned nature.This means that there is no unique solution and more appropriate algorithms have to be used to get physically acceptable solutions for the approximation of targets x.
Generally, constraining the solution provides an opportunity to improve the approximation.It is common practice to apply constraints by a functional (x) of the solution to the minimization term Eq. ( 13) and to use a regularization parameter µ > 0 to balance between the residual norm and the solution constraint: The solution constraint is also referred to as solution norm, because a norm is frequently used as a constraining criterion, as described in the following methods.

Constrained Least Squares restoration
One way to arrive a constrained solution to Eq. ( 13) is the Tikhonov regularization (Bertero et al., 1998) which is often referred to as constrained least square (CLS) solution, because residual and solution norm are of a squared character (Euclidean norm).The specific formulation of Eq. ( 14) with a quadratic functional is as follows: where L is usually a differential operator, e.g.Laplacian to preserve edges in image processing.Due to the fact, that here point targets are of main interest, the identity matrix L=I is used.By applying a Lagrange multiplier an explicit solution is found and corresponds to a Wiener filter approach (Murli et al., 1999), whose regularization parameter µ complies with inversely proportional value of the occurring signal-tonoise ratio (SNR).Performing a parameter sweep of µ is the subsequent procedure to get a reasonable approximation of targets x.An appropriate value of µ can be selected by the L-curve criterion for instance.Therein a parametric log-log plot of residual norm versus solution norm is utilized to detect a pounce within the graph (Hansen, 2010).
With the intention to verify performance of the CLS method, a noisy channel (SNR ≈ 50 dB) with two pointtargets was simulated by a simple time-delay model of Reflection and Transmission-variant (see Table 1) in Matlab.A linear uniform array with equidistant spacing d = 20 mm and N = 7 channels (in plots referred to as full array) is used to form up an antenna array that nearly fits the half-wavelength rule when a frequency range of X-band (see above) is applied.The sparse array raw data (d = 40 mm and N = 4) is obtained by selecting every second antenna; corresponding PSF was shown in Fig. 3.The resulting Gaussian noise added output of the system z(r,ϕ) is illustrated in Fig. 4a, whereas factual position of point-targets are shown by green circles.Ambiguities of real existing targets are located quite plainly.
In conjunction with a parameter study of µ the CLS method is applied and a suitable approximation xCLS (r,ϕ) for µ=10 −2.5 with a logarithmic scaling is shown in Fig. 4b.Obviously, a slight ambiguity suppression can be recognized.An azimuthal sectional view of a selected range cell (highlighted by gray lines) is given in Fig. 5 and enables stating a suppression of approximately 5 dB in contrast to sparse array data.

Sparse restoration
The method of CLS regularization yields an insufficient side lobe suppression.Therefore it seems appropriate to put more information in the regularization process to improve suppression of ambiguities.In relation to the preface of this section, this can be performed by the functional (x) in Eq. ( 14).Relating to the gauging situation, this can be achieved by uprating of non-sparse solutions of x, which fits to the fact that mainly a few isolated point targets occur.From a mathematical point of view the p norm solves this particular problem: In image processing the utilization of the p norm is referred }. (17) The rationale behind the p norm is illustrated in Fig. 6 by applying CLS-, 1-and 0.5-norm to a two component vector As an example, the marker belongs to The 2-norm yields the absolute value of a vector, no matter if the vector is sparse or not.The p norm yields similar values if the vector is sparse but much greater values if the vector is dense.In connection with the minimization term Eq. ( 17) and translated into the inverse problem of imaging this means that a solution with a few number of point targets is preferred to any smooth solution.
In order to test the procedure of NQ regularization to raw data of a sparse array imaging system, same raw data as for CLS regularization is processed for p = 1 and p = 0.5, respectively.The utilized algorithm that already handles complex-valued problems is described in (Cetin, 2001).Again, a parameter sweep of µ is performed and a manual selection of regularization parameter (µ p=1 ≈ 0.6, µ p=0.5 ≈ 0.3) is done.The results are shown in Fig. 4c and d for both values of p, respectively; sectional view is found in Fig. 5 again.As a result, the side lobe suppression is considerably improved for both choices of p in comparison with CLS method.Especially, in case of p = 0.5 the all over side lobe level reaches a value below the level of the full array as a consequence of the sparse regularization approach.A side lobe suppression of approximately 45 dB is achieved and therefore this approach is well suited to enable a FMCW sparse array imaging system to locate point targets more unambiguously.

Measurements
Measurements were performed to verify practical suitability of the presented regularization methods to suppress ambiguities of sparse array imaging.The RF frontend was realized by a vector network analyser and a coaxial switching matrix   to realize a switched multiport measurement system.With regard to the simulation example, 7 antipodal Vivaldi antennas (see Fig. 7) were used as single elements of the X-band linear uniform array with d = 20 mm to form the full array.Two corner reflectors served as targets and were positioned in front of the array.The distance between sensor and targets is approximately 3.25 m with a slight offset among each other.
The raw data (Fig. 8a) is processed in the same way as simulation data above with CLS and sparse restoration.Figure 9 again illustrates an azimuthal sectional view of a selected range cell view of all restoration methods.Similar to the simulation results applying sparsity to the restoration problem obviously results in good side lobe suppression (p = 0.5, Fig. 8b).An attenuation of approximately 45 dB is achieved and is therefore lower than the overall side lobe level of the full array.

Conclusions
The reliability of a microwave gauging system benefits from the concept of imaging radar.However, an adequate azimuthal resolution of a FMCW imaging system with a limited number of antennas can only be achieved by applying a sparse antenna array.Consequently, by violating the half-wavelength rule, a high side lobe arises and ambiguities within the imaging process occur.By modelling the process of FMCW sparse array imaging using a common wavefront approach, the corresponding PSF was achieved, which has to be inverted to get unambiguously information of real targets.Due to this way of proceeding, the present paper was able to apply two relevant regularization methods from image processing for the task of restoration.Furthermore, significant improvements in side lobe suppression have been achieved by utilizing sparse restoration for both simulated and measured data of a sparse antenna array.

Fig. 3 .
Fig. 3. PSF of an X-band FMCW sparse array imaging system with N=4 single elements and element spacing of d=40 mm.

Fig. 5 .
Fig. 5. Sectional view of one range cell of simulated data.

Fig. 9 .
Fig. 9. Sectional view of one range cell of measured data.

Table 1 .
Some significant variants for multistatic radar.

339, 2012 S. Kolb and R. Stolle: FMCW sparse array imaging and restoration for microwave gauging
•N u ×1 is introduced to give a more realistic model of the process of DBF FMCW array imaging.The noisy observation of the channel z∈ C N r •N u ×1 with analogically mapping www.adv-radio-sci.net/10/333/2012/Adv.Radio Sci., 10, 333-