# Lowering latency and processing burden in computational imaging through dimensionality reduction of the sensing matrix

Thomas Fromenteze, O. Yurduseven, P. del Hougne, D.R. Smith

► To cite this version:

Thomas Fromenteze, O. Yurduseven, P. del Hougne, D.R. Smith. Lowering latency and processing burden in computational imaging through dimensionality reduction of the sensing matrix. Scientific Reports, 2021, 11 (1), pp.3545. 10.1038/s41598-021-83021-6. hal-03195840

HAL Id: hal-03195840

https://hal.science/hal-03195840

Submitted on 26 May 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
<figure><img src="image_2.png" /></figure>

Check for updates

OPEN

# Lowering latency and processing burden in computational imaging through dimensionality reduction of the sensing matrix

Thomas Fromentèze<fnref n="1" />,⁵✉, Okan Yurduseven²,⁵, Philipp del Hougne³ & David R. Smith⁴

Recent demonstrations have shown that frequency-diverse computational imaging systems can greatly simplify conventional architectures developed for imaging by transferring constraints into the digital layer. Here, in order to limit the latency and processing burden involved in image reconstruction, we propose to truncate insignificant principal components of the sensing matrix that links the measurements to the scene to be imaged. In contrast to recent work using principle component analysis to synthesize scene illuminations, our generic approach is fully unsupervised and is applied directly to the sensing matrix. We impose no restrictions on the type of imageable scene, no training data is required, and no actively reconfigurable radiating apertures are employed. This paper paves the way to the constitution of a new degree of freedom in image reconstructions, allowing one to place the performance emphasis either on image quality or latency and computational burden. The application of such relaxations will be essential for widespread deployment of computational microwave and millimeter wave imagers in scenarios such as security screening. We show in this specific context that it is possible to reduce both the processing time and memory consumption with a minor impact on the quality of the reconstructed images.

Imaging with microwaves brings about several advantages, such as all-weather operation, the use of non-ionizing radiation, the ability to see through most optically opaque materials and the availability of mature component technology. Due to these advantages, microwave imaging has proven to be useful in a wide range of applications ranging from security-screening<fnref n="1" />⁻⁴ to biomedical imaging⁵⁻⁸ and non-destructive testing⁹⁻¹². Despite these advantages, conventional microwave imaging systems suffer from a number of challenges. In general, microwave imaging requires the synthesis of composite apertures to raster-scan the scene, either mechanically, such as in synthetic aperture radar¹³⁻¹⁶, or electronically, such as in phased arrays¹⁷⁻²¹. Mechanical scanning is not desirable in that data acquisition time can be significant, posing a major challenge for real-time operation. All-electronic operation can be achieved using phased arrays. However, beam synthesis using phased arrays requires that each antenna within the composite aperture has a phase-control circuit, or a phase-shifter, and preferably a power amplifier (to compensate for the insertion losses of the phase shifters). As a result, using phased arrays, all-electronic operation can require a large number of phase-shifting circuits and power amplifiers, significantly increasing the complexity and power-consumption of such apertures. Applications requiring high depth resolution will also be particularly constrained by the fractional bandwidth limitations of the active components of each transmission and reception chain.

To address these challenges, unusual modalities have been proposed that leverage ideas related to computational imaging²²⁻²⁷. Numerous proofs of concept of frequency-diverse computational imaging have recently been proposed²⁸⁻³⁵, where the scene information is encoded onto a set of measurement modes that exhibit quasi-randomness across the operating frequency band. In other words, the scene information is sampled on quasi-random bases by stepping through a number of frequency points across the operating frequency band (frequency sweep)³⁶. Leveraging this quasi-random set of modes offered by engineered frequency-diverse systems

<footnote n="1">¹University of Limoges, CNRS, XLIM, UMR 7252, 87000 Limoges, France. ²Centre for Wireless Innovation (CWI), Institute of Electronics, Communications and Information Technology (ECIT), School of Electronics, Electrical Engineering and Computer Science (EEECS), Queen's University Belfast, Belfast BT3 9DT, UK. ³Univ Rennes, CNRS, IETR - UMR 6164, 35000 Rennes, France. ⁴Department of Electrical and Computer Engineering, Center for Metamaterials and Integrated Plasmonics, Duke University, Durham, NC 27708, USA. ⁵These authors contributed equally: Thomas Fromentèze and Okan Yurduseven. ✉email: thomas.fromenteze@unilim.fr</footnote>
eliminates the need for beam synthesis and requires a simple frequency sweep to achieve imaging. A three-dimensional (3D) microwave image of the scene can then be reconstructed by interacting the measurements of the scene with the transfer function of the frequency-diverse imaging system using computational reconstruction techniques, such as direct algorithms, including the matched-filter technique, or iterative algorithms, including the least-square technique³⁷.

Using such frequency-diverse schemes for computational imaging of large scenes at microwave frequencies, such as in security-screening, can be a computationally demanding task. For a 3D scene discretised at the resolution limit of a finite-size aperture, the number of voxels to be reconstructed can be significant. Several techniques can be adopted to reduce the computational complexity of the imaging problem. For example, additionally using data from optical sensors (sensor fusion), the field-of-view for a given (microwave) imaging problem can be constrained to a volume bounded by the surface of the imaged object³². It is reported in³⁸ that adopting such a prior constraint on the imaging problem can reduce the number of voxels to be reconstructed, and hence the complexity of the imaging problem, by more than 90%. Another approach includes the parallelization of the image reconstruction algorithm using general-purpose graphics processing units (GPGPUs)³⁹ and field-programmable-gate-array (FPGA) structures⁴. However, despite the significant reduction in the size of the computational problem enabled by applying the aforementioned techniques, the reconstruction time can be on the order of several seconds³².

In view of this, it is evident that additional techniques are needed to further reduce the computational complexity of the imaging problem for real-time operation. Inspired by well-known lossy compression schemes such as JPEG or MP3, we hypothesize that a moderate loss of information can be tolerated in order to substantially reduce the memory and processing burden. Principal component analysis (PCA) is a mathematical tool that can be used in this perspective to decompose a matrix on a set of orthogonal bases⁴⁰. In the considered imaging problem, the sensing matrix links the space to be imaged to the measured signals; by decomposing the sensing matrix into its principle components, we can identify independent structures of which the sensing matrix is composed and their respective contribution levels. We then select only the most significant principal components in order to compute a model of the sensing matrix with reduced dimensions.

Our proposal is distinct from recent ideas in the literature to optimize microwave imaging systems using PCA. One such approach, proposed in⁴¹ and⁴², is based on the radiation of illumination patterns specific to each type of scene to be imaged, necessarily requiring the use of reconfigurable systems for scene-dependent beam-synthesis. Provided a prior knowledge about the nature of the scene to be imaged is available, these authors have shown that it is possible to limit the number of sequential captures necessary for image reconstruction compared to the use of purely random patterns. Recent work in⁴³ took this idea even further by directly integrating a model of the physical layer into an artificial neural network in order to jointly learn optimal measurement and processing strategies based on a priori knowledge of scene, task and measurement constraints. Since this "learned sensing" strategy enables one to minimize the acquisition of task-irrelevant information, it is highly task-specific and requires a supervised learning technique. In contrast to these prior works, the approach we propose here applies PCA to the sensing matrix rather than the expected scene. Consequently, our approach is not scene dependent and does not require the use of sequential measurements relying on active reconfigurable antennas. We propose a system-dependent but scene-independent method relying on a frequency-sweep to generate a succession of random illumination patterns that interrogate the scene to be imaged; via PCA we can limit the dimensionality of the sensing matrix and thereby the computational complexity of the image reconstruction.

This paper demonstrates that our approach makes it possible to reduce computation times and memory consumption for imaging applications while having only a minor impact on the quality of the reconstructed images. Our approach allows each application to dispose of a new degree of freedom in image reconstruction by moving the constraints to either image quality or numerical performance. The outline of this paper is as follows: In "Operation principle" section, we explain the concept of frequency-diverse computational microwave imaging and the application of PCA to the sensing matrix. In "Results and discussion" section, we present the images of a human-sized object reconstructed with and without PCA-based dimensionality reduction; moreover, we provide a quantitative analysis of the reconstructed images in terms of imaging quality, mean-square-error (MSE), memory consumption, and reconstruction time. Finally, we provide concluding remarks in the last section.

## Operation principle

**Frequency-diverse computational microwave imaging.** In frequency-diverse computational imaging, a set of frequency-dependent spatially-varying field patterns are used to interrogate the scene to be imaged. Using the first Born approximation, the measured signal is linked to the imaged object by means of the sensing matrix (or the transfer function) of the imaging system as follows:

$$g_{i,j}(\omega) = \int_{\mathbf{r}} E_{T_x}^i(\mathbf{r}, \omega) E_{R_x}^j(\mathbf{r}, \omega) \rho(\mathbf{r}) d\mathbf{r} + n(\omega). \quad (1)$$

Equation (1) is the integral form of the forward model where *g* denotes the measurement signal, *ρ* is the imaged scene, *i* and *j* are the indices of the transmit and receive frequency-diverse antennas, and *E*<sub>Tx</sub> and *E*<sub>Rx</sub> are the radiated *E*-field patterns of the transmit and receive antennas propagated to the scene using the free-space Green's function⁴⁴. Finally, *n*(*ω*) represents additive Gaussian white noise. For the sake of clarity, the principle of dimensionality reduction proposed in this paper is implemented within the framework of a scalar wave approximation. An adaptation of the studied techniques can nevertheless be carried out to the reconstruction of complete susceptibility tensors by means of recent developments in polarimetric computational imaging³³,⁴⁵. Under the first Born approximation, the sensing matrix is proportional to the dot product of the transmit and receive frequency-diverse antenna fields: *H* ~ *E*<sub>Tx</sub>*E*<sub>Rx</sub> (Fig. 1).
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Illustration of a computational imaging system. The fields radiated in transmission and reception by a frequency diverse system interact with a target to be imaged. The scene $\rho$ is interrogated by a succession of quasi-random field distributions varying according to the frequency and measurement ports. By means of simple frequency sweeps, spatial information is quickly encoded into a set of measurements gathered in $\mathbf{g}$.</figcaption></figure>

Retrieving an estimate of the reflectivity of the scene $\hat{\rho}$ is an inverse problem that is solved numerically based on a discretized version of the forward model (1). For a frequency-diverse aperture consisting of $T$ transmit antennas and $R$ receive antennas, and sampling the frequency band at $F$ frequency points, the total number of measurement modes is $M = T \times R \times F$. The discretization of the scene is performed at the resolution limit of the aperture in the range and cross-range planes. The range resolution is governed by the imaging bandwidth $B = 9$ GHz in the K-band (17.5–26.5 GHz): $\delta_r = c/2B = 1.67$ cm; the cross-range resolution is governed by the aperture size $D$, target distance $d$, and wavelength $\lambda$: $\delta_{cr} = \lambda d/D$, resulting in $\delta_{cr} = 7$ mm cross-range resolution for imaging at $d = 1$ m distance. The $M \times N$ sensing matrix $\mathbf{H}$ links the $N$ voxels to the $M$ measurement modes. The discretized forward model is the following:

$$ \mathbf{g}_{M \times 1} = \mathbf{H}_{M \times N} \rho_{N \times 1} + \mathbf{n}_{M \times 1}. \qquad (2) $$

In (2), although the field definitions forming the sensing matrix are scalar, bold font is used to denote the vector and tensor quantities. The sensing matrix $\mathbf{H}$ does not have an exact inverse because it is not square ($M \neq N$) and full rank. Several techniques can be used to retrieve an estimate of the scene information. Among those, direct reconstruction techniques can retrieve the scene estimate without the need for additional iteration steps as follows:

$$ \hat{\rho} = \mathbf{H}^{+} \mathbf{g} = \mathbf{H}^{+} (\mathbf{H} \rho + \mathbf{n}) \approx \rho, \qquad (3) $$

where $\mathbf{H}^+$ denotes the reconstruction operator applied to the measurements in order to retrieve $\hat{\rho}$.

The simplest techniques, known as matched filtering, consists in using the conjugate transpose of the sensing matrix, denoted $\mathbf{H}^\dagger$, as reconstruction operator.

**Dimensionality reduction with principal component analysis.** As noted above, for imaging problems with large scenes, such as in security screening, it can be challenging to retrieve the estimate of the scene in a timely manner due to the large size of the sensing matrix. In fact, processing a big sensing matrix might not even be possible due to a lack of available computational resources. For instance, in security-screening applications at K-band frequencies, the size of the sensing matrix without the use of any a priori knowledge was shown to be as large as 90 Gb$^{38}$. Therefore, it is evident that additional techniques are needed to solve the inverse problem of Eq. (2) to retrieve $\hat{\rho}$, particularly for real-time applications.

The first step we take to lower the reconstruction problem's computational complexity is the technique at the heart of the present paper: we use PCA to remove insignificant dimensions of the sensing matrix $\mathbf{H}$. PCA
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Top: Schematic representation of a measurement performed by a computational imaging system. The scene to be imaged $\rho$ is multiplied by a pre-determined sensing matrix **H**. Linear combinations specific to the imaging system encode the spatial information into a set of measurements gathered in **g**. Using information from optical sensors (sensor fusion), the model can be reduced on the fly to select only a portion of the region of interest, represented by the red rectangle. Bottom: The SVD allows one to have a new perspective on this interaction. Each column of **V** represents a spatial mask whose degree of similarity to the scene $\rho$ is estimated by correlation. Each of these masks is multiplied by a singular value that makes it possible to determine their level of contribution to the measured vector **g**, through linear combinations of columns of **U**. Once again, sensor fusion makes it possible to reduce the dimensions of the **V** matrix by selecting on the fly only useful portions of the region of interest.</figcaption></figure>

originates from the general context of spectral theory: **H** can be diagonalized using the singular value decomposition (SVD) technique⁴⁶ that factorizes it into three sub-matrices,

$$ \mathbf{H} = \mathbf{U}\mathbf{\Sigma}\mathbf{V}^{\dagger}. \qquad (4) $$

In (4), **U** and **V** are orthonormal matrices and **Σ** is a diagonal matrix whose non-negative coefficients $\sigma_n$ are organized by convention in descending order. An intuitive illustration of this decomposition is given in Fig. 2. The distribution of singular values (also referred to as singular value spectrum) directly reflects the degree of correlation between the rows of the sensing matrix⁴⁷. An ideal imaging system will ensure that each measurement is a new (completely independent) source of information and that the sampling of the reconstructed image perfectly matches the resolution limit of this system. Under such ideal conditions, the SVD of the sensing matrix yields a decomposition with singular values of identical levels⁴⁸, reflecting that **H** carries an optimal amount of information in relation to its dimensions. In practice, space and frequency dimensions are often oversampled in order to ensure a certain immunity to measurement uncertainties and additive noise. SVD is often used to study and improve the conditioning in inverse problems, limiting the amplification of subspaces that are weakly contributing and particularly sensitive to the impact of additive noise³⁰,⁴⁹.

Since the number of reconstructed voxels *N* is lower than the number of measurement modes *M*, there are *N* non-zero singular values $\sigma_n$ of which $K < N$ are deemed significant. The interaction with the scene to be imaged can thus be approximated as follows⁵⁰:

$$ \mathbf{g} = \mathbf{U}\mathbf{\Sigma}\mathbf{V}^{\dagger}\rho + \mathbf{n} \qquad (5) $$

$$ = \sum_{n=1}^{N} \sigma_n \mathbf{u}_n \mathbf{v}_n^{\dagger} \rho + \mathbf{n} \qquad (6) $$

$$ \mathbf{g} \approx \sum_{n=1}^{K} \sigma_n \mathbf{u}_n \mathbf{v}_n^{\dagger} \rho, \qquad (7) $$

where the measurement is approximated in the last step by a succession of $K < N$ orthogonal masks applied to the scene, neglecting the contribution of the lowest singular values and their associated singular vectors $\mathbf{u}_n$ and $\mathbf{v}_n^\dagger$. Rather than interrogating *N* voxels of the scene, it is then possible to try to reconstruct only *K* coefficients corresponding to spatial correlations with the masks formed by the most significant singular vectors $\mathbf{v}_n^\dagger$.
We now define a new matrix formalism of the imaging problem by selecting only the first K singular vectors, forming $U_s$ and $V_s$, and by creating a new diagonal matrix $\Sigma_s$ with the K most significant singular values. The measurement then takes the following form, revealing the principal components $P_c = U_s\Sigma_s$ in the form of an operator interrogating $\alpha = V_s^\dagger\rho$, a vector of K coefficients corresponding to independent linear combinations of the information of the scene to be imaged:

$$g \approx U_s \Sigma_s V_s^\dagger \rho \qquad (8)$$

$$\approx P_c \alpha. \qquad (9)$$

The matrix $P_c$ has reduced dimensions $M \times K$ and inversion strategies similar to the ones already described for $H$ can be implemented to obtain an estimate of the vector $\hat{\alpha}$ as follows:

$$\hat{\alpha} = P_c^+ g. \qquad (10)$$

Finally, taking advantage of the orthonormal properties of the $V_s$ matrix ($\langle V_i | V_j \rangle = \delta_{ij}$), the scene estimation is straightforward:

$$\hat{\rho} = V_s \hat{\alpha} \qquad (11)$$

$$\approx V_s V_s^\dagger \rho \qquad (12)$$

$$\hat{\rho} \approx \rho. \qquad (13)$$

Rather than storing the matrix $H$ of dimensions $M \times N$ for the application of iterative techniques or its reconstruction operator $H^+$ for direct estimations, it is thus possible to store only two smaller matrices $P_c$ (dimensions: $M \times K$) and $V_s$ (dimensions: $N \times K$) or their respective reconstruction operators.

Having detailed the PCA-based dimensionality reduction, we now discuss how it can be combined with a priori knowledge from optical sensors to reduce the number of voxels to be reconstructed. Prior work$^{32,38,39,51}$ has shown that based on optical sensors the surface boundaries of the object under consideration can be detected to define a volume around this surface that constrains the domain for microwave imaging (Fig. 2). Thereby, the number of voxels constituting the scene is reduced via support detection to $N_o < N$ elements, determining a new formalism valid for a given position of the target:

$$g_{M \times 1} = H_{oM \times N_o} \rho_{oN_o \times 1} + n_{M \times 1}. \qquad (14)$$

It should be carefully noted that, as shown in Fig. 2, the pixel/voxel set $N_o$ selected based on information from optical sensors depends directly on the target position and must be continuously updated from the initial set of $N$ elements defining the entire imageable scene.

Only the $V_s$ matrix is impacted by this support detection, reducing the number of reconstructed pixels from $N$ to $N_o$. One can then estimate $\hat{\rho}_o$ in this reduced volume as follows:

$$\hat{\rho}_o = V_{so} \hat{\alpha}. \qquad (15)$$

Our approach not only greatly reduces the memory consumption of an imaging system, but it also accelerates image reconstruction by reducing the number of required computational operations. We note that our approach remains fully compatible with the use of a priori knowledge from optical sensors which allows one to reduce the reconstructed domain in real time. In the following section, we validate our technique by studying a body scanner in a context oriented toward the detection of concealed threats.

# Results and discussion

## Considered imaging setup.
Figure 3 shows the imaging scenario that we consider in simulation. This imaging system consists of 16 transmit and 16 receive metasurface antennas ($10 \times 10$ cm$^2$), arranged in a Mills-Cross system layout. The latter represents a thinned array architecture while ensuring a similar k-space support extent in comparison to populating the entire aperture$^{31,52}$. A human-sized object is placed at a distance of $d = 1$ m from the aperture. The frequency-diverse imager operates at K-band frequencies (17.5–26.5 GHz), sampling the 9 GHz bandwidth at 101 frequency points which corresponds to a frequency interval of 90 MHz. We choose to work with 101 frequency sampling points based on the upper bound on the number of useful measurement modes that a frequency-diverse antenna with a finite quality (Q) factor can provide$^{36}$, which is given by $QB/f_c$; for our system with antenna quality factor $Q = 330$, bandwidth $B = 9$ GHz and central operating frequency $f_c = 22$ GHz, the resulting upper bound is 135 modes. The synthesized frequency-diverse aperture thus produces $M = 25,856$ measurement modes. The initial region of interest contains $N = 12 \times 118 \times 251 = 355,416$ voxels; based on a priori knowledge from optical sensors, we narrow the number of voxels to be reconstructed down to $N_o = 16,292$.

The gains achievable in numerical performance by applying the dimensionality reduction of the sensing matrix that we propose inherently depend on the correlations between the utilized measurement modes. The considered imaging hardware has been optimized to generate a quasi-random set of measurement modes, meaning that an effort was made to keep the correlations low. This is evidenced by the singular value spectrum of the sensing matrix shown in Fig. 4 which is relatively flat. Therefore, gains in numerical performance achievable via dimensionality reduction of the sensing matrix are somewhat limited; nonetheless, even under these (for our
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Synthesized composite frequency-diverse Mills-Cross aperture for imaging a human-sized object in the scene. The 16 metasurface antennas used in transmission (reception) are shown in orange (blue). The rectangular-shaped red volumetric frame represents the complete domain imaged by the **H**-matrix interrogating a set of *N* = 355,416 voxels. Based on a priori knowledge, *N*₀ = 16,292 are retained (depicted by red cubes forming the human-sized object) in a defined volume encompassing the front surface of the target.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Singular value spectrum of the full **H**-matrix. Imaging results are presented in this section by selecting four subsets of principal components, represented by cases A, B, C, and D, corresponding to the 2500, 5000, 7500, and 10,000 largest singular values, respectively.</figcaption></figure>

proposal unfavourable) conditions, it is possible to illustrate the benefits of the proposed method. We anticipate that systems with lower performance and consequently a higher level of correlation between each measurement may benefit more strongly from a greater model reduction with a minor impact on the quality of the reconstructed images. The reconstructions in the following are performed on a computer with 128 Gb of RAM and an Intel Xeon processor with 10 cores running at 2.4 GHz (E5-2640 v4). To provide characteristic computation times that do not directly depend on the specifications of the utilized machine, we calculate ratios with respect to reference computation times evaluated with the full sensing matrix **H**.

In the following, two data processing scenarios are studied to illustrate the impact of the proposed technique on the quality of the reconstructed images and on the associated computation times. First, we consider a matched filter as example of a direct reconstruction technique; second, we consider the generalized minimal residual method (GMRES) as example of an iterative reconstruction technique.
<figure><img src="image_7.png" /><figcaption>Figure 5. Reference scene reconstructions with matched filtering without the proposed PCA-based dimensionality reduction of the sensing matrix. Left: Scene estimate $\hat{\rho}$ computed for the entire region of interest (N = 355,416 voxels). Right: Scene estimate $\hat{\rho}_o$ computed only for the voxels selected from optical sensors (N$_o$ = 16,292 voxels). The colormap is given in decibels with normalized data.</figcaption></figure>

**Case 1: direct reconstruction by matched filtering.** The first approach using a matched filter is based on a simple phase compensation of the sensing matrix. In the absence of any dimensionality reduction, image reconstruction by matched filtering takes the following form:

$$ \hat{\rho}_{(\text{MF})} = \mathbf{H}^{\dagger}\mathbf{g}. \qquad (16) $$

For reference, we show reconstruction results with matched filtering without the proposed PCA-based dimensionality reduction in Fig. 5, illustrating the difference between reconstructing all N voxels or only reconstructing the selected subset of N$_o$ voxels.

Matched filtering not only has the advantage of reducing the dimensions of the matrices manipulated for image reconstruction (in comparison to other techniques), but also of limiting the impact of clutter related to imperfect estimations outside of the reduced region of interest. Using this first simple reconstruction method as a benchmark, the advantages and disadvantages associated with PCA-based dimensionality reduction are now investigated.

A selection of four reconstruction examples with matched filtering is presented in Fig. 6, computed following Eqs. (10) and (11). Reconstructions are carried out by calculating $\hat{\rho}_{o(\text{PCA,MF})} = \mathbf{V}_{so} \mathbf{P}_c^{\dagger} \mathbf{g}$. We vary the number of principal components (PC) in steps of 500. In each case, the $\mathbf{P}_c$ and $\mathbf{V}_s$ matrices can be pre-calculated and stored in memory to be used on different measurements g; this reduces the computational burden during operation since for a given scenario it is not necessary to change the number of selected principle components. The reconstructed images are then compared to the reconstructions achieved with $\mathbf{H}_o^{\dagger}$. It is therefore intuitive that beyond a certain number of principal components, the intermediate factorization step ultimately requires more calculation than for a unique multiplication by the conjugate-transpose of the sensing matrix. This phenomenon is seen in Fig. 6 when more than 5700 PCs (22%) are selected, leading to an average computation time overcom ing the 196 ms threshold required to evaluate Eq. (16).

Example B in Fig. 6 requires a reconstruction time slightly below that of the reference case; upon visual inspection, the quality of the reconstructed image appears to be comparable to cases C and D which only use 29% and 38.7% of the principal components. We evaluate the MSE by taking the reconstruction with the $\mathbf{H}_o$ matrix as reference in order to assess the loss in image quality. Indeed, the gradient of the MSE curve seen in Fig. 6 is very low at the beginning. Thus, by moderately sacrificing the accuracy of the reconstructions it is possible to improve the reconstruction speed (and hence refresh rate) with our proposed approach.

An analysis of the achieved reduction of the computational burden for the considered matched-filtering technique is shown in Fig. 7 as a function of the number of retained principal components. The $S_{MF}$ curve corresponds to the ratio between the raw image-computation time using PCA-based dimensionality reduction of the sensing matrix relative to that using the full $\mathbf{H}_o$-matrix. We also plot a second curve, $S_{MFwmi}$, that in addition to the raw computation times accounts for the time it takes to extract the $\mathbf{H}_o$ or $\mathbf{V}_{so}$ sub-matrices from $\mathbf{H}$ or $\mathbf{V}_s$, respectively. As seen in the figure, for the four considered truncations A to D, the gains in computation time are comparable with and without accounting for the time that the matrix-indexing takes. Matrix-indexing only begins to adversely affect the latency gain once a very large number of singular values is selected. Yet, as seen in Fig. 6, the absolute extraction times of the sub-matrices in our study for each new position of the target in the
<figure><img src="image_8.png" /></figure>

<figure><img src="image_9.png" /></figure>

<figure><img src="image_10.png" /></figure>

<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /></figure>

<figure><img src="image_13.png" /></figure>

**Figure 6.** Study of the impact of dimensional reductions by PCA for direct reconstructions by matched filtering. The upper left graph corresponds to the raw image computation time and the upper right one includes, in addition to the computation time, the time required to extract the rows corresponding to the reduced region of interest from the $V_s$-matrix according to the target position. The MSE is calculated with each reconstructed spatial reflectivity by considering the reconstruction with $H_o$ as benchmark. The averaged computation times are represented for each reconstruction and the general trend is highlighted by means of a regression. The dashed red line represents the reconstruction time using the $H_o^\dagger$ matrix in both cases. An animation presenting the successive reconstructions carried out by gradually increasing the number of principal components is available in the Supplemental Materials: https://bit.ly/3a7GFFy.

<figure><img src="image_14.png" /></figure>

**Figure 7.** Reduction of the computational burden through PCA-based dimensionality reduction of the sensing matrix in direct reconstructions (matched filtering). We plot the reduction in computation time ($S_{MF}$, red) and memory usage ($G_{MF}$, blue). The red dashed curve denoted by $S_{MFwmi}$ additionally accounts for the impact of $H$ and $V_s$ matrix indexing times (wmi stands for "with matrix indexing").
region of interest are not compatible with real-time applications to a moving scene. The reason for this is that we performed our study on a CPU with MATLAB which suffers from notoriously poor matrix-indexing efficiency. However, space-partitioning strategies associated with the use of GPUs and low-level programming can achieve a latency compatible with real-time applications³².

For direct reconstruction techniques like matched filtering, the PCA-based dimensionality reduction also yields a much more direct benefit related to memory usage. Real-time applications require the ability to store all matrices needed for computations in a fast memory; body-scanning applications are then particularly constrained by the volume of data used to calculate each image. The memory saving factor is defined by the following equation:

$$G_{\text{MF}} = \frac{\text{size}(\mathbf{H})}{\text{size}(\mathbf{P}_c) + \text{size}(\mathbf{V}_s)} = \frac{M \times N}{(M + N) \times K} \quad (17)$$

As a reminder, *M* corresponds to the number of utilized measurement modes, *N* to the total number of voxels and *K* to the number of selected principal components. We emphasize that the total number of voxels *N* is taken into account in Eq. (17) so as to not restrict the study to a target specific sub-domain. Based on Eq. (17), we can determine a limiting value of *K* beyond which the memory consumption will exceed the one in case of using the entire sensing matrix, that is when $G_{\text{MF}} < 1$. The latter happens for $K > M \times N/(M + N) = 24,102$ in the considered imaging scenario, equivalent 93% of the total number of principle components in the studied case. We present the evolution of the memory saving as a function of the number of retained principle components in Fig. 7.

The application of our PCA-based dimensionality reduction enabled memory-consumption saving factors of 9.6, 4.8, 3.2 and 2.4 in the considered cases A, B, C, and D, respectively. We used complex values coded in simple precision on 8 bytes such that *H* occupies more than 68 Gb, which justifies the usefulness of the proposed approach. In conclusion of this first study on a direct reconstruction technique, we have evidenced the possibility to slightly accelerate the image reconstruction (on average about 148 ms vs 196 ms for the reference) while simultaneously significantly reducing the RAM consumption by a factor of 4.8—at the price of a minor image degradation that does not hamper potential security-screening applications.

**Case 2: iterative reconstruction with GMRES.** Reconstructions by matched filtering have the advantage of being fast; however, they do not allow for the compensation of magnitude terms. More advanced iterative techniques that improve the reconstruction quality can be particularly slow, especially when the dimensions of the problems to be treated are large. The PCA-based dimensionality reduction is now applied to this new context of iterative techniques to study the possible savings in computation time and memory consumption. Specifically, the generalized minimal residual method (GMRES) is applied, having yielded good performances in previous work³⁰,³². Since this technique requires the use of a square sensing matrix, the $\mathbf{P}_c^\dagger$ matrix is multiplied to each member of (8):

$$\mathbf{P}_c^\dagger \mathbf{g} \approx \mathbf{P}_c^\dagger \mathbf{P}_c \alpha. \quad (18)$$

By minimizing $||\mathbf{P}_c^\dagger \mathbf{g} - \mathbf{P}_c^\dagger \mathbf{P}_c \alpha||_2$, this iterative approach enables the reconstruction of the vector $\hat{\alpha}$ starting from a matched-filtering estimate $\mathbf{P}_c^\dagger \mathbf{g}$ and with use of the correlation matrix of the principal components $\mathbf{P}_c^\dagger \mathbf{P}_c$. In a similar way to the first study, the performances are studied for a set of reconstructions carried out using different amounts of leading singular values. Four examples are more carefully investigated by displaying the reconstructed images (Fig. 8).

Compared to reconstructions performed only by matched filtering, magnitude-term compensation provides more accurate estimates of the imposed constant reflectivity function over the entire target. However, the specular nature of the interaction of waves with the target tends to create shadows that cannot be removed unless the target is illuminated from other incident angles. This physical limitation, also affecting conventional imaging systems, can nevertheless be mitigated by merging a set of reconstructions from different perspectives³²,⁵³.

As in the previous section, we plot the convergence of the reconstruction quality in terms of the MSE for the use of different amounts of leading principle components in Fig. 8. The result obtained with only 6.8% of the principal components (case A) makes it possible to distinguish the first main structures of the image but still remains corrupted by the lack of information. The addition of information to reach 22.7% of the principal components allows the reconstruction of a much more exploitable image revealing significantly more detail. Once again, it is difficult to discern upon visual inspection the difference between the reconstructions carried out with 43% and 100% of the principal components.

The effect of the PCA-based dimensionality reduction on the calculation time depends again on the number of selected principal components. Although the results obtained with 100% of the principal components are identical to those obtained by applying the same technique with the *H* matrix, the PCA-based approach still requires in this case an additional calculation step justifying a longer time beyond a certain number of principal components.

As before, we also compare the memory consumption with and without dimensionality reduction. For this purpose, the ratio of matrix sizes required to reconstruct images from frequency measurements alone is calculated for the iterative GMRES technique. If the full sensing matrix is used by GMRES, the quantity $||\mathbf{H}^\dagger \mathbf{g} - \mathbf{H}^\dagger \mathbf{H} \rho||_2$ is minimized, which requires storage of both *H* and additionally the correlation matrix $\mathbf{H}^\dagger \mathbf{H}$. The reduction in memory usage achieved by PCA-based dimensionality reduction is thus given by

$$G_{\text{GMRES}} = \frac{\text{size}(\mathbf{H}^\dagger \mathbf{H}) + \text{size}(\mathbf{H})}{\text{size}(\mathbf{P}_c^\dagger \mathbf{P}_c) + \text{size}(\mathbf{P}_c) + \text{size}(\mathbf{V}_s)} \quad (19)$$
<figure><img src="image_15.png" /><figcaption>**Figure 8.** Study of the impact of dimensional reductions by PCA for iterative reconstructions by the general minimum residual method (GMRES). The upper left graph corresponds to the raw image computation time and the upper right one includes in addition the time required to extract the rows corresponding to the reduced region of interest from the $V_s$-matrix according to the target position. The MSE is calculated with each reconstructed spatial reflectivity by considering the GMRES reconstruction with the matrix $H_o$ as benchmark. The averaged computation times are represented for each reconstruction and the general trend is highlighted by means of a regression. The dashed red line represents the reconstruction time with the $H_o$ matrix in both cases. An animation presenting the successive reconstructions carried out by gradually increasing the number of principal components is available in the Supplemental Materials: https://bit.ly/2PuiWYl.</figcaption></figure>

$$ = \frac{(N + M) \times N}{(K + M + N) \times K} \qquad (20) $$

Our dimensionality-reduction technique hence enables us to limit memory usage as long as $K < (\sqrt{M^2 + 5N^2 + 6MN} - (M + N))/2$, or 223,914 principal components in this study (impossible case here since we only dispose 25,856 singular values). The memory gain $G_{GMRES}$, which is thus always greater than 1, is plotted in Fig. 9.

Having considered a more computationally intensive image reconstruction based on an iterative technique in this section, the utility of the proposed approach was even more evident. For example C (see Figs. 8 and 9), the calculation time is divided by a factor of 3.5 and the RAM consumption by a factor of 46.5—yielding an image with a MSE of -20.3 dB. The lowest achievable MSE of -24.4 dB is only achieved if all principle components are used.

Finally, we use the GMRES results to illustrate the impact of the non-orthogonality of $V_{so}$. In order to avoid a dependence on the position of the target within the region of interest, the singular vectors $V_{so}$ are extracted from the singular vectors $V_s$ corresponding to the full sensing matrix $H$ for the entire region of interest. By limiting the support of the reconstruction, the orthogonality of the $V_{so}$ sub-matrix is not guaranteed, especially since the number of reconstructed pixels may be greater than the number of selected principal components. This limitation is illustrated in Fig. 10 by comparing the result obtained with GMRES using the matrix $H_o$ and the result obtained with GMRES using all principal components in the PCA-based approach. To mitigate the negative effect of this spatial filtering, one could compute the pseudo-inverse of $V_{so}$ for each new target position. However, this approach would require an additional calculation step for each reconstructed image, limiting the achievable refresh rate of the system under consideration.

**Effects of measurement noise.** This study concludes with an analysis of the impact of thermal noise usually associated with the use of active receivers. Previous work has highlighted the usefulness of PCA in the identification⁵⁴ and filtering⁵⁵ of noise components of signals in various imaging applications such as X-ray
<figure><img src="image_16.png" /><figcaption>**Figure 9.** Reduction of the computational burden through PCA-based dimensionality reduction of the sensing matrix in iterative reconstructions (GMRES). We plot the reduction in computation time ($S_{GMRES}$, red) and memory usage ($G_{GMRES}$, blue). The red dashed curve denoted by $S_{GMRESwmi}$ additionally accounts for the impact of **H** and $V_s$ matrix indexing times (wmi stands for "with matrix indexing").</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>**Figure 10.** Comparison of GMRES-based reconstructions performed using the PCA-based technique with all principal components (left) and using the $H_0$ matrix (right).</figcaption></figure>

tomography⁵⁶ and magnetic resonance imaging⁵⁷. The decomposition of linear operators by PCA facilitates the study of noise corruption, revealing its impact on the weakest principal components in a degree directly dependent on the signal-to-noise ratio (SNR). The elimination of these weakest components thus leads to the suppression of noise subspaces, implicitly realizing a regularization by truncated singular value decomposition⁵⁸.

The previous sections have highlighted the beneficial effects of dimensionality reduction of sensing matrices on image computation times and memory consumption. As a continuation of these efforts, it is finally proposed to study the effects of additive measurement noise on the proposed PCA-based technique. Following the model given by Eq. (2), thermal noise is added to the measurements, following a Gaussian distribution of zero mean and standard deviation $\sigma_n$ such that $n \sim N(0, \sigma_n^2)$. The noise standard deviation is determined to reach an objective SNR from the measured signal power. It will be considered in the following that the sensing matrix can be perfectly determined and will thus be noise-free. This assumption seems reasonable insofar as this matrix is measured in a preliminary stage without any specific time constraint, opening the possibility of averaging measurements and of increasing integration time of each acquired sample.

For a series of SNRs of noise-corrupted measured signals, a matched-filtering reconstruction is performed. The figure of merit is again defined as the MSE for each estimate, considering the noise-free reconstruction with $H_0^1$ as a reference. For the different considered SNRs, 10 measurements are simulated to compute the mean and standard deviation of the MSE in each case. Figure 11 displays the average and standard deviation of the MSE for four considered SNRs as a function of the number of retained principal components. For reference, the MSE achieved with the full sensing matrix $H_0$ is also indicated.
<figure><img src="image_18.png" /><figcaption>**Figure 11.** Evolution of the MSE of reconstructions computed by matched filtering (MF) in the presence of Gaussian additive noise. The MSE is calculated for each estimate with the noise-free reconstruction carried out with the $H_0^\dagger$ matrix chosen as a reference. The continuous lines represent the reconstructions performed using the proposed PCA-based technique while the dotted lines correspond to the reconstructions with the $H_0$ matrix. These results correspond to the means computed in each case with 10 simulations carried out with different noise realizations satisfying objective SNR, while the associated colored envelopes correspond to a variation of plus or minus one standard deviation around these means.</figcaption></figure>

As expected, the quality of the PCA-based reconstructions eventually converges to that of reconstructions based on the full sensing matrix $H_0$ once a sufficient number of principal components is retained. The latter depends on the SNR of the measurements. Indeed, as the noise level increases, a growing number of subspaces are corrupted, reducing the overall number of uncorrupted subspaces in $H_0$ that can contribute new information to the reconstruction. These results lead to the conclusion that the proposed PCA-based technique is particularly suitable for the reconstruction of images from measurements with unfavorable SNR. Under such conditions, PCA-based reconstructions can indeed be faster and less memory-intensive with little to no impact on the quality of the resulting images compared to results obtained from a full sensing matrix.

## Conclusion

Since computational imaging systems typically cannot ensure perfect spatial and spectral orthogonality of the utilized measurement modes, there is usually a non-negligible level of correlation between the obtained measurements. These correlations imply that upon a change of basis via singular value decomposition, many components of the sensing matrix contribute only weakly; neglecting such insignificant components thus entails only a minor deterioration of the imaging quality but offers significant reductions in memory usage and enhancements in computation speed. Although the suppression of the least contributing principal components generally results in a loss of quality of the reconstructed images, this study also highlighted that many of these suppressed components cannot contribute to the reconstruction under realistic operating conditions due to the presence of noise anyway. Our study illustrated the proposed dimensionality reduction of the sensing matrix with a case study of a specific computational imaging system based on a Mills-Cross aperture composed of metasurface antennas and a complex target (a human), considering both direct and iterative image reconstruction techniques. The finite spatial extent and intrinsic losses of the considered metasurface antennas result in notable correlations between the utilized measurement modes. Our results offer the possibility to trade-off memory needs and computation time with the image reconstruction quality, offering a new degree of freedom in computational imaging to place the emphasis on either reconstruction quality or numerical performance. Our method can readily be applied to other imaging systems and contexts that rely on sensing matrices with inevitable correlations³⁶,⁵⁹. Finally, we note that the reconstruction techniques used in this paper are agnostic to the scene to be imaged, thus imposing no restrictions on the nature or on the sparsity of the region of interest. Future efforts may investigate whether the presented dimensionality reduction technique can be combined with structure-aware imaging algorithms⁴¹⁻⁴³,⁶⁰.

Received: 21 April 2020; Accepted: 1 January 2021

Published online: 11 February 2021

## References

1. Sheen, D. M., McMakin, D. L. & Hall, T. E. Three-dimensional millimeter-wave imaging for concealed weapon detection. *IEEE Trans. Microwave Theory Tech.* **49**, 1581–1592 (2001).
2. Salmon, N. et al. Digital beam-forming for passive millimetre wave security imaging (2007).

3. Nanzer, J. A. *Microwave and Millimeter-Wave Remote Sensing for Security Applications* (Artech House, Norwood, 2012).

4. Ahmed, S. S. *Electronic Microwave Imaging with Planar Multistatic Arrays* (Logos Verlag Berlin GmbH, Berlin, 2014).

5. Klemm, M. et al. Microwave radar-based breast cancer detection: Imaging in inhomogeneous breast phantoms. *IEEE Antennas WIREL. PROPAG. Lett.* **8**, 1349–1352 (2009).

6. Nikolova, N. K. Microwave imaging for breast cancer. *IEEE Microwave Mag.* **12**, 78–94 (2011).

7. Grzegorczyk, T. M. et al. Fast 3-d tomographic microwave imaging for breast cancer detection. *IEEE Trans. Med. Imaging* **31**, 1584–1592 (2012).

8. Fear, E. C. et al. Microwave breast imaging with a monostatic radar-based system: A study of application to patients. *IEEE Trans. Microwave Theory Tech.* **61**, 2119–2128 (2013).

9. Schajer, G. S. & Orhan, F. B. Microwave non-destructive testing of wood and similar orthotropic materials. *Subsurface Sens. Technol. Appl.* **6**, 293 (2005).

10. Kharkovsky, S. & Zoughi, R. Microwave and millimeter wave nondestructive testing and evaluation—overview and recent advances. *IEEE Instrum. Meas. Mag.* **10**, 26–38 (2007).

11. Hasar, U. Non-destructive testing of hardened cement specimens at microwave frequencies using a simple free-space method. *NDT & E Int.* **42**, 550–557 (2009).

12. Marks, D. L., Yurduseven, O. & Smith, D. R. Sparse blind deconvolution for imaging through layered media. *Optica* **4**, 1514–1521 (2017).

13. Amineh, R. K., McCombe, J. & Nikolova, N. K. Microwave holographic imaging using the antenna phaseless radiation pattern. *IEEE Antennas WIREL. PROPAG. Lett.* **11**, 1529–1532 (2012).

14. Martinez-Lorenzo, J. A., Quivira, F. & Rappaport, C. M. Sar imaging of suicide bombers wearing concealed explosive threats. *Prog. Electromagn. Res.* **125**, 255–272 (2012).

15. Moreira, A. et al. A tutorial on synthetic aperture radar. *IEEE Geosci. Remote Sens. Mag.* **1**, 6–43 (2013).

16. Laviada, J., Arboleya-Arboleya, A., Alvarez-Lopez, Y., Garcia-Gonzalez, C. & Las-Heras, F. Phaseless synthetic aperture radar with efficient sampling for broadband near-field imaging: Theory and validation. *IEEE Trans. Antennas Propag.* **63**, 573–584 (2014).

17. Ku, B.-H. et al. A 77–81-ghz 16-element phased-array receiver with ±50° beam scanning for advanced automotive radars. *IEEE Trans. Microwave Theory Tech.* **62**, 2823–2832 (2014).

18. Ralston, T. S., Charvat, G. L. & Peabody, J. E. Real-time through-wall imaging using an ultrawideband multiple-input multiple-output (mimo) phased array radar system. In *2010 IEEE International Symposium on Phased Array Systems and Technology*, 551–558 (IEEE, New York, 2010).

19. Withington, S., Saklatvala, G. & Hobson, M. P. Partially coherent analysis of imaging and interferometric phased arrays: Noise, correlations, and fluctuations. *JOSA A* **23**, 1340–1348 (2006).

20. Charvat, G. L., Kempel, L. C., Rothwell, E. J., Coleman, C. M. & Mokole, E. L. An ultrawideband (ubw) switched-antenna-array radar imaging system. In *2010 IEEE International Symposium on Phased Array Systems and Technology*, 543–550 (IEEE, New York, 2010).

21. Fenn, A. J., Temme, D. H., Delaney, W. P. & Courtney, W. E. The development of phased-array radar technology. *Lincoln Lab. J.* **12**, 321–340 (2000).

22. Brady, D. J. *Optical Imaging and Spectroscopy* (Wiley, New York, 2009).

23. Cossairt, O. S., Miau, D. & Nayar, S. K. Scaling law for computational imaging using spherical optics. *JOSA A* **28**, 2540–2553 (2011).

24. Shrekenhamer, D., Watts, C. M. & Padilla, W. J. Terahertz single pixel imaging with an optically controlled dynamic spatial light modulator. *Opt. Express* **21**, 12507–12518 (2013).

25. Welsh, S. S. et al. Fast full-color computational imaging with single-pixel detectors. *Opt. Express* **21**, 23068–23074 (2013).

26. Shin, D., Kirmani, A., Goyal, V. K. & Shapiro, J. H. Photon-efficient computational 3-d and reflectivity imaging with single-photon detectors. *IEEE Trans. Comput. Imaging* **1**, 112–125 (2015).

27. Fromenteze, T., Decroze, C., Abid, S. & Yurduseven, O. Sparsity-driven reconstruction technique for microwave/millimeter-wave computational imaging. *Sensors* **18**, 1536 (2018).

28. Hunt, J. et al. Metamaterial apertures for computational imaging. *Science* **339**, 310–313 (2013).

29. Fromenteze, T., Decroze, C. & Carsenat, D. Waveform coding for passive multiplexing: Application to microwave imaging. *IEEE Trans. Antennas Propag.* **63**, 593–600 (2014).

30. Fromenteze, T. et al. Computational imaging using a mode-mixing cavity at microwave frequencies. *Appl. Phys. Lett.* **106**, 194104 (2015).

31. Yurduseven, O., Gollub, J. N., Marks, D. L. & Smith, D. R. Frequency-diverse microwave imaging using planar mills-cross cavity apertures. *Opt. Express* **24**, 8907–8925 (2016).

32. Gollub, J. et al. Large metasurface aperture for millimeter wave computational imaging at the human-scale. *Sci. Rep.* **7**, 42650 (2017).

33. Fromenteze, T. et al. Computational polarimetric microwave imaging. *Opt. Express* **25**, 27488–27505 (2017).

34. Tondo Yoya, A. C., Fuchs, B. & Davy, M. Computational passive imaging of thermal sources with a leaky chaotic cavity. *Appl. Phys. Lett.* **111**, 193501 (2017).

35. Jouade, A., Meric, S., Lafond, O., Himdi, M. & Ferro-Famil, L. A passive compressive device associated with a Luneburg lens for multibeam radar at millimeter wave. *IEEE Antennas WIREL. PROPAG. Lett.* **17**, 938–941 (2018).

36. Marks, D. L., Gollub, J. & Smith, D. R. Spatially resolving antenna arrays using frequency diversity. *JOSA A* **33**, 899–912 (2016).

37. Yurduseven, O. et al. Resolution of the frequency diverse metamaterial aperture imager. *Prog. Electromagn. Res.* **150**, 97–107 (2015).

38. Yurduseven, O., Gollub, J. N., Rose, A., Marks, D. L. & Smith, D. R. Design and simulation of a frequency-diverse aperture for imaging of human-scale targets. *IEEE Access* **4**, 5436–5451 (2016).

39. Marks, D. L., Yurduseven, O. & Smith, D. R. Fourier accelerated multistatic imaging: A fast reconstruction algorithm for multiple-input-multiple-output radar imaging. *IEEE Access* **5**, 1796–1809 (2017).

40. Shlens, J. A tutorial on principal component analysis. arXiv preprint arXiv:1404.1100 (2014).

41. Liang, M., Li, Y., Meng, H., Neifeld, M. A. & Xin, H. Reconfigurable array design to realize principal component analysis (pca)-based microwave compressive sensing imaging system. *IEEE Antennas WIREL. PROPAG. Lett.* **14**, 1039–1042 (2015).

42. Li, L. et al. Machine-learning reprogrammable metasurface imager. *Nat. Commun.* **10**, 1082 (2019).

43. del Hougne, P., Imani, M. F., Diebold, A. V., Horstmeyer, R. & Smith, D. R. Learned integrated sensing pipeline: Reconfigurable metasurface transceivers as trainable physical layer in an artificial neural network. *Adv. Sci.* **7**, 1901913 (2020).

44. Lipworth, G. et al. Metamaterial apertures for coherent computational imaging on the physical layer. *JOSA A* **30**, 1603–1612 (2013).

45. Peng, R., Yurduseven, O., Fromenteze, T. & Smith, D. R. Advanced processing of 3d computational microwave polarimetry using a near-field frequency-diverse antenna. *IEEE Access* **8**, 166261–166272 (2020).

46. Fang, Q., Meaney, P. M. & Paulsen, K. D. Singular value analysis of the Jacobian matrix in microwave image reconstruction. *IEEE Trans. Antennas Propag.* **54**, 2371–2380 (2006).

47. Roy, O. & Vetterli, M. The effective rank: A measure of effective dimensionality. In *2007 15th European Signal Processing Conference*, 606–610 (IEEE, New York, 2007).

48. del Hougne, P., Davy, M. & Kuhl, U. Optimal multiplexing of spatially encoded information across custom-tailored configurations of a metasurface-tunable chaotic cavity. *Phys. Rev. Appl.* **13**, 041004 (2020).
49. Fromenteze, T. et al. Unification of compressed imaging techniques in the microwave range and deconvolution strategy. In *2015 European Radar Conference (EuRAD)*, 161–164 (IEEE, New York, 2015).

50. Lipworth, G. et al. Comprehensive simulation platform for a metamaterial imaging system. *Appl. Opt.* **54**, 9343–9353 (2015).

51. Marks, D. L., Yurduseven, O. & Smith, D. R. Near-field multistatic radar reconstruction with stretched-phase Fourier accelerated multistatic imaging. *IET Radar Sonar Navig.* **11**, 1718–1729 (2017).

52. Fromenteze, T., Kpré, E. L., Carsenat, D., Decroze, C. & Sakamoto, T. Single-shot compressive multiple-inputs multiple-outputs radar imaging using a two-port passive device. *IEEE Access* **4**, 1050–1060 (2016).

53. Gonzalez-Valdes, B. et al. Millimeter wave imaging architecture for on-the-move whole body imaging. *IEEE Trans. Antennas Propag.* **64**, 2328–2338 (2016).

54. Pyatykh, S., Hesser, J. & Zheng, L. Image noise level estimation by principal component analysis. *IEEE Trans. Image Process.* **22**, 687–699 (2012).

55. Chen, G. & Qian, S.-E. Denoising of hyperspectral imagery using principal component analysis and wavelet shrinkage. *IEEE Trans. Geosci. Remote Sens.* **49**, 973–980 (2010).

56. Burdet, P., Croxall, S. A. & Midgley, P. Enhanced quantification for 3D SEM-EDS: Using the full set of available X-ray lines. *Ultra-microscopy* **148**, 158–167 (2015).

57. Manjón, J. V., Coupé, P. & Buades, A. MRI noise estimation and denoising using non-local PCA. *Med. Image Anal.* **22**, 35–47 (2015).

58. Jensen, S. H., Hansen, P. C., Hansen, S. D. & Sorensen, J. A. Reduction of broad-band noise in speech by truncated GSVD. *IEEE Trans. Speech Audio Process.* **3**, 439–448 (1995).

59. Erdős, L. Universality of Wigner random matrices: A survey of recent results. *Russ. Math. Surv.* **66**, 507 (2011).

60. Theodoridis, S., Kopsinis, Y. & Slavakis, K. Sparsity-aware learning and compressed sensing: An overview. *arXiv preprint arXiv:1211.5231* (2012).

## Acknowledgements

T.F. and O.Y. would like to thank the Hubert Curien Partnership (PHC Alliance) for their support. The work of O.Y. was funded by the Leverhulme Trust under Research Leadership Award RL-2019-019.

## Author contributions

T.F. and O.Y. developed the theoretical principle, carried out the numerical studies and wrote a first draft. T.F., O.Y., P.d.O, and D.S. participated in its edition and in the associated scientific discussions.

## Competing interests

The authors declare no competing interests.

## Additional information

**Supplementary Information** The online version contains supplementary material available at https://doi.org/10.1038/s41598-021-83021-6.

**Correspondence** and requests for materials should be addressed to T.F.

**Reprints and permissions information** is available at www.nature.com/reprints.

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_19.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2021