Full-Vectorial 3D Microwave Imaging of Sparse Scatterers through a Multi-Task Bayesian Compressive Sensing Approach

In this paper, the full-vectorial three-dimensional (3D) microwave imaging (MI) of sparse scatterers is dealt with. Towards this end, the inverse scattering (IS) problem is formulated within the contrast source inversion (CSI) framework and it is aimed at retrieving the sparsest and most probable distribution of the contrast source within the imaged volume. A customized multi-task Bayesian compressive sensing (MT-BCS) method is used to yield regularized solutions of the 3D-IS problem with a remarkable computational efficiency. Selected numerical results on representative benchmarks are presented and discussed to assess the effectiveness and the reliability of the proposed MT-BCS strategy in comparison with other competitive state-of-the-art approaches, as well.


Introduction
Microwave imaging (MI) techniques are aimed at inferring the complex permittivity distribution within an inaccessible investigation domain from the scattering interactions between the matter and probing electromagnetic (EM) waves [1].They have been successfully applied in several diagnostic scenarios including non-destructive testing and evaluation [2,3], through-wall imaging [4], subsurface prospecting [5][6][7][8][9][10], and structural health monitoring [11].Moreover, they represent a very appealing technology in many biomedical applications [12] such as, for instance, breast cancer detection [13][14][15][16][17][18][19][20][21][22] thanks to the use of non-ionizing radiations.To date, significant efforts have been mostly devoted to the development of two-dimensional (2D) MI algorithms, mainly based on transverse-magnetic (TM) [23][24][25][26][27][28] or transverse-electric (TE) polarized [29] configurations, rather than fully-vectorial three-dimensional (3D) ones [30].As a matter of fact, under the assumption that the EM properties of the unknown scattering scenario are invariant along a longitudinal direction, the arising inverse scattering (IS) problem can be recast to the solution of simplified scalar Helmholtz equations [1].On the other hand, 2D-MI approaches are prone to errors when finite-volume scatterers are under test [31] because of the over-simplified modelling of the scattering phenomena.The slower evolution of 3D-MI techniques has been mainly caused by the higher complexity of both data-collection/storage and image reconstruction processes with respect to the tomographic (2D) case.Moreover, a significantly larger number of unknowns has to be retrieved and it becomes very hard to manage when there is the need of high-resolution images (e.g., a realistic discretization of a human thorax needs millions of voxels for having a clinical significance [32]).Furthermore, the amount of non-redundant information on the investigation domain achievable from measurements is upper-bounded and the ratio between unknowns and scattering data turns out to be very high [33].Owing to such limitations, solving 3D-MI problems faces hard challenges and it requires the non-trivial implementation of effective countermeasures to both the non-linearity and the ill-posedness issues of the arising full-vectorial IS problem.
Dealing with 3D scenarios, synthetic aperture radar (SAR)-based methodologies such as confocal MI [34] and synthetic near-field focusing [5] have been proposed.They are based on the emission of wide-band pulses from multiple transmitting positions and the successive processing of the collected echoes.Only target detection and localization (i.e., a qualitative imaging of the investigation domain) is typically yielded, while quantitatively retrieving the distribution of the EM properties needs the numerical solution of the non-linear scattering equations.Towards this end, effective 3D-MI approaches have been presented in the scientific literature.They are based on the processing of time-domain [19,35] or frequency-domain [36,37] data and, due to the extremely-wide dimension of the solution space, which is usually proportional to the number of discretization domains, state-of-the-art methodologies are primarily based on deterministic methods (e.g., Gauss-Newton (GN) [36,38], Conjugate-Gradient (CG) [39], Level Set (LS) [16], and Inexact-Newton (IN) [40] methods, possibly formulated in different functional spaces [41]), even though they do not a priori guarantee to reach the actual solution (i.e., the global optimum of the cost function quantifying the mismatch between measured and estimated scattering data).To overcome such a drawback, either very efficient forward solvers [42,43] have been introduced or stochastic multi-agent inversion algorithms, suitably integrated with multi-resolution strategies [37] for a sustainable customization to the 3D case, have been successfully exploited.Alternatively, computationally-efficient approaches to the MI problem have been recently explored within the compressive sensing (CS) framework [23,29,[44][45][46][47][48][49].Despite the early stage [47,50] of their implementation in 3D cases, very interesting applicative examples are already available [20].
According to the CS theory, sparseness priors can be enforced to solve the IS problem and to yield a regularized solution provided that it admits a representation with few non-null coefficients in a suitably chosen basis [44,45].However, since available CS solvers generally deal with linear problems, many sparsity-promoting approaches have been formulated within Born-like approximations, their success being limited to weak scatterers or to specific applications where a qualitative imaging is enough [29].Alternatively, contrast source inversion (CSI)-based formulations of the IS problem can be successfully employed to yield accurate reconstructions also in the presence of scatterers with high EM contrast with respect to the surrounding medium [23].Following such a line of reasoning, this paper is aimed at presenting a novel computationally-efficient approach, based on a Bayesian CS (BCS) method, to solve the 3D-IS problem concerned with non-weak scatterers.Towards this end, the full-vectorial IS problem is formulated within a probabilistic CSI framework and then it is efficiently solved through a customized multi-task BCS (MT-BCS) strategy.
The outline of the paper is as follows.The formulation of the 3D-CSI MI problem is detailed in Section 2, while the proposed BCS-based solution strategy is described in Section 3. Selected numerical results, from representative test cases, are presented and compared with competitive state-of-the-art alternatives in Section 4. Finally, some concluding remarks are drawn (Section 5).

Mathematical Formulation
Let us consider a 3D isotropic non-magnetic [µ (r) = µ 0 ] scattering scenario characterized by a relative permittivity distribution, ε (r), and a conductivity profile, σ (r), r being the position vector defined as r = xu x + yu y + zu z , u p being the unit vector along the p-th (p = {x, y, z}) direction.The goal of the MI problem at hand is to estimate the contrast function [51] within the investigation domain D (i.e., ∀r ∈ D), ω = 2π f , ε 0 , and σ 0 being the angular frequency (A time dependency factor exp (jωt) is assumed and omitted hereinafter to simplify the notation, but without loss of generality in the mathematical formulation.), the background permittivity and conductivity (with σ 0 = 0 hereinafter), respectively, and standing for "defined as".Towards this end, a set of V monochromatic ( f being the working frequency) plane-waves impinging from known angular directions (θ v , ϕ v ), v = 1, ..., V, with known electric field is used to successively probe the unaccessible domain D ( Figure 1).
D Under such hypotheses and by adopting a CSI formulation [52] for the scattering problem, the electromagnetic interactions between the scatterers in D and the v-th (v = 1, ..., V) incident wave can be mathematically described through the following integral data equation where Ω is an observation domain external to D (Ω ∩ D = ∅-Figure 1), .., V) scattered field defined as the difference between the v-th (v = 1, ..., V) electric field with, E v (r) = ∑ p={x, y, z} E v,p (r) u p (i.e., the v-th total electric field), and without, E i v (r) (i.e., the v-th incident electric field), the scatterers in the background medium • stands for the scalar product, and is the dyadic Green's function for the homogeneous free-space background medium of dielectric and magnetic properties ε 0 and µ 0 , respectively, I being the unit tensor.In (3), J v [J v (r) = ∑ p={x, y, z} J v,p (r) u p ] is the v-th (v = 1, ..., V) unknown contrast current density induced in the investigation domain (r ∈ D) by the v-th probing field (i.e., the v-th illumination of D) that models the scattering profile of D and it is defined as the v-th (v = 1, ..., V) equivalent source radiating in the background medium an electromagnetic field equal to the v-th (v = 1, ..., V) scattered field E s v (r).Furthermore, the v-th (v = 1, ..., V) incident field E i v (r) complies with the so-called integral state equation within the investigation domain (r ∈ D).
To numerically deal with the inverse problem at hand, a set of N 3D rectangular pulse functions defined in a set of N cubic (also called voxels) sub-domains (D = ∪ N n=1 D (n) -Figure 1) is adopted for yielding the following piece-wise constant representation of the v-th (v = 1, ..., V) unknown contrast source where J (n) v,p = J v,p (r n ) (p = {x, y, z}) and r n denotes the barycentre of the n-th (n = 1, ..., N) voxel D (n) .By sampling the scattered field in M probing locations of Ω (r m ∈ Ω, m = 1, ..., M), it is possible to rewrite (3) in the following discrete form [1] where G pq (r m , r ) dr (p, q = {x, y, z}), or in compact matrix notation as follows where .T indicates the transpose operator, is the external 3D Green's operator (G ext ∈ C 3M×3N ) where the (m, n)-th (m = 1, ..., M; n = 1, ..., N) entry of the (p, q)-th (p, q = {x, y, z}) sub-matrix, G ext pq ∈ C M×N , is equal to G (mn) pq .Under the hypothesis that the contrast distribution in D is sparse with respect to the basis at hand (7), it also turns out that the v-th (v = 1, ..., V) contrast source (8)

Inversion Method
The matrix relation in (10) is representative of a severely ill-posed problem being (i) ill-defined since its solution is not-unique due to the presence of non-radiating sources in D that afford a null/not-measurable field in Ω and (ii) ill-conditioned since the condition number of G ext , η, is large (η 1).To counteract such a negative feature for yielding stable reconstructions of the EM properties of the investigation domain from the solution of (10) when processing noisy data, as well, enforcing the a-priori knowledge of the sparseness of the unknown equivalent source with respect to a suitable basis turns out to be an effective regularization strategy.Towards this end, the CS formulation based on the Bayesian theory is adopted to re-formulate the 3D-CSI problem at hand in a probabilistic sense.Such a choice is done to avoid the need of numerically checking the fulfillment of the restricted isometry property (RIP) of the 3D Green's operator (11) as required by deterministic CS solvers.Indeed, whether such a compliance test is already very hard from a computational viewpoint for small-scale problems, it becomes computationally unfeasible in 3D scattering scenarios since exponentially heavier than for 2D cases.
In order to apply computationally-efficient BCS solvers, ( 10) is first rewritten as a larger real-valued linear problem where R ( . ) and I ( . ) being the real and the imaginary part function, respectively.According to such a description, the original inverse scattering problem can be then formulated as follows 3D-CSI BCS-Based Problem Formulation-Starting from the measurement of E v,p (v = 1, ..., V, p = {x, y, z}), and the knowledge of G p (p = {x, y, z}), determine the sparsest provided that the support of In order to enforce the physical correlation existing among the V solutions of (14) as well as between the unknown contrast sources and the vectorial components of the known scattered field, a customized version of the MT-BCS solver [53] is adopted by setting the number of parallel "tasks" to H = (3 × V) (owing to the 3D nature of the problem).More in detail, ( 14) is firstly reformulated in the Bayesian MT framework as [29,53] starting from the definition of a shared prior P J v that statistically links the 3 × V parallel problem unknowns as [29,53] where ψ (ψ ψ (t) ; t = 1, ..., 6 × N ) is the set of shared-among the V-views and the three Cartesian scattered field components-hyper-parameters.By assuming a hierarchical Gaussian model for P J v ψ and substituting ( 16) in (15) [29,53], the MT-BCS v-th (v = 1, ..., V) source term is then given after simple manipulations [29] by the following expression where A = diag ψ and ψ ψ(t) ; t = 1, ..., 6 × N .These latter are determined with a computationally-efficient relevant vector machine (RVM) [54] by solving the following maximum L ψ |α, β being the logarithmic marginal likelihood function for the fully-vectorial problem given by where I is the identity matrix, while α and β are user-defined control parameters.Once J v (v = 1, ..., V) has been estimated through (17), the contrast function, τ (r) = ∑ N n=1 τ (n) Ψ (n) (r), is retrieved by computing the corresponding expansion coefficient vector τ = τ (n) ∈ C; n = 1, ..., N , whose n-th (n = 1, ..., N) entry is given by where while the coefficient 17) according to the following mapping

Numerical Assessment
In this Section, representative results from a wide numerical analysis are presented and discussed to assess the reconstruction capabilities and the robustness of the proposed 3D-MI approach.Moreover, some practical guidelines for the optimal setting of its key calibration parameters will be provided for helping the interested users.Towards this end, the following reference scenario has been considered.
A cubic volume of side L = 1.25 [λ], λ being the free-space wavelength, has been chosen as investigation domain D and it has been probed by V = 48 plane-waves impinging from the V angular directions (θ where M ( • being the floor operator) and ρ = 3.2 [λ].To avoid the so-called inverse crime [1], two different voxel-based discretization of D have been considered for the inverse (N = 10 × 10 × 10) and the forward (N f wd = 20 × 20 × 20) problem, respectively.As for a quantitative evaluation of the reconstructions, the following error metric has been computed for the whole investigation domain (R = tot ⇒ N tot = N), the scatterers support (R = int ⇒ N int = O), and its complementary region (It is worthwhile to remark that the scatterer support is not an a-priori information exploited during the inversion, but it is rather employed in the post-processing phase only to compute the error metrics (24).
and τ (n) being the actual and the estimated contrast of the n-th (n = 1, ..., N) voxel in D, respectively.The first test case is concerned with the retrieval of a single (K = 1, K being the number of disconnected scattering regions lying in D) off-centered scatterer composed by a single-voxel (O = 1) of side x = y = z = 0.125 [λ] with contrast τ = 1.0 (Figure 2a).First, a sensitivity analysis for setting the optimal trade-off values of the control parameters of the MT-BCS solver [i.e., α and β in (19)] has been carried out.More specifically, the reconstruction error (24) has been computed by varying the values of the calibration coefficients within the ranges α ∈ 1, 10 2 and β ∈ 10 −5 , 10 −2 for different signal-to-noise ratios (SNRs).Figure 2b shows that the reconstruction "quality" is almost insensitive to the choice of α when SNR > 5 [dB] and it mainly depends on the noise level.
As it can be observed, the object support is always faithfully detected regardless of the amount of noise.There are only slight deviations from the actual contrast value (e.g., τ SNR=50 [dB] = 0.81-Figure 3a; τ SNR=5 [dB] = 0.77-Figure 3d) as quantitatively indicated by the corresponding values of the error index ξ tot being equal to ξ tot SNR=50 [dB] = 9.58 × 10 −5 (Figure 3a) and 3d) at the highest and lowest SNRs, respectively.
In order to prove that such results are not due to a customization of the MT-BCS solver to the specific scenario at hand, a set of W = 100 inversions has been performed by randomly changing the position of the same target within the investigation domain.The results of such a statistical analysis are summarized in Table 1 where the minimum (ξ min tot min w=1, ..., W {ξ tot,w }), the maximum (ξ max tot max w=1, ..., W {ξ tot,w }), the average (ξ Whatever the SNR, the scatterer retrieval is very accurate since on average ξ avg tot < 1.0 × 10 −4 (Table 1) and the variance of the reconstruction error, ξ var tot , is very small (i.e., ξ var tot < 5.0 × 10 −11 -Table 1).
Such a positive outcome holds true when dealing with stronger scatterers, as well.Thanks to the CSI formulation, which avoids any physical assumption/approximation in the application of the CS to the 3D scattering equations (Section 2), faithful data inversions are yielded also when increasing the actual contrast well-beyond the value of the first example.As a proof, the plots of ξ tot (Figure 4a), ξ int (Figure 4b), and ξ ext (Figure 4c) versus τ and the noise level indicate that D can be carefully imaged when the scatterer is at least up to τ = 4.0.
As expected, the reconstruction progressively degrades when higher and higher contrasts are at hand especially when processing highly blurred data.Figure 5 5c).Nevertheless, it is worth pointing out that these inaccuracies correspond to voxels with a very low contrast (i.e., τ < 3 × 10 −4 -Figure 5c) that can be easily filtered out by a simple (i.e., the result is not so sensitive to the choice of the threshold value τ th ) thresholding (Figure 5d-τ th = 10 −3 ) and the scatterer is always correctly localized (Figure 5c as well as Figure 5d vs. Figure 5a).The previous examples have dealt with real-valued target contrasts.However, the imaginary part of τ is not enforced to be zero during the inversion (see Section 3).In order to assess the reliability of the MT-BCS also when complex contrast values are at hand, the retrieval of a K = 1, O = 1 target with τ = 1.0 − 0.6j (Figure 6a,b) has been considered next.The analysis of the retrieved profiles (Figure 6) shows that the location, size, and contrast of the scatterer is correctly retrieved by the proposed imaging strategy both in low (e.g., SNR = 50 [dB]-Figure 6c,d) and in high noise conditions (e.g., SNR = 5 [dB]-Figure 6e,f), as shown by the corresponding error indexes (e.g., ξ tot SNR=50 [dB] = 8.63 × 10 −5 -Figure 6c,d 6e,f).
Being assessed the effectiveness of the proposed approach in reconstructing the sparsest (O = 1) actual profile, let us now analyze its performance for scenarios exhibiting a lower sparsity order with respect to the selected voxel basis (7).Towards this end, a second single-voxel scatterer has been added to the configuration of the first test case (Figure 2a), but in a different position in D (K = O = 2, τ = 1.0-Figure 7a).As it can be inferred from the plots in Figure 7b-e, the MT-BCS is always able to correctly retrieve both scatterers even from very low-noise data (e.g., ξ tot SNR=5 [dB] = 9.80 × 10 −4 -Figure 7e).However, the same statistical analysis of the single-scatterer/single-voxel case (i.e., randomly changing the locations of the scatterers) here results in wider variations of the integral error (i.e., ξ var tot ≥ 4.5 × 10 −7 (Table 2) vs. ξ var tot ≥ 1.11 × 10 −11 (Table 1)).Therefore, to better understand how the inversion accuracy is affected by the positions of the scatterers within the investigation domain, a further set of experiments has been performed by deterministically changing the relative locations of the K = 2 scatterers, d being the distance of each one of them from the origin.Some representative test cases are reported in Figure 8: 8c,d), and d 8e,f).As expected, the reconstruction worsens when the two objects get closer (Figure 9) until some artifacts appear 2.41-Figure 8b vs.
)-Behavior of the total error, ξ tot , versus the distance of the scatterers from the origin, d.
Further lowering the sparseness of the actual profile causes a decrement of the CS performance as proven by the results of the third test case when O = 4 and the K = 4 single-voxel scatterers are randomly placed in D (Figure 10a).Unavoidably, the MT-BCS gets worse and undesired artifacts occur also when processing low-noisy data (e.g., ξ tot SNR=50 [dB] = 2.37 × 10 −3 -Figure 10b), but still without preventing the possibility to correctly identify K = 4 separate objects (Figure 10b-e).In addition, in this case, changing the distance d of the objects with respect to the origin influences the quality of the retrieval as pictorially shown in Figure 11  Indeed, it turns out that more accurate images of D are yielded when the scatterers are far (i.e., d ↑ ⇒ ξ tot ↓) as confirmed by the plot of ξ tot vs. d (Figure 12a).However, it cannot be neglected that a simple filtering (Figure 12b-τ th = 10 −3 ) allows one to clearly resolve the scatterer support even in the most critical case (Figure 11a).The next numerical test is devoted to validate the MT-BCS in a more complex 3D imaging scenario concerned with non-uniformly shaped scatterers.More specifically, K = 3 objects sized x , 2 y , 2 z = [0.125, 0.25, 0.125] [λ], and 3 x , 3 y , 3 z = [0.25, 0.25, 0.125] [λ] (O = 7) with contrast τ = 1.0 have been imaged (Figure 13a).
Despite the increased complexity and the reduced intrinsic-sparsity order of the scattering configuration, faithful reconstructions of the 3D contrast distribution are yielded (Figure 13b-e) even though there is an over-estimation of the scatterers supports when highly-blurred data are at hand (e.g., Figure 13b vs. Figure 13e) and, consequently, an increase of the reconstruction error with the noise level (e.g., ξ tot SNR=50 [dB] ≈ 2.04).Moreover, the previous considerations regarding the fidelity of the proposed strategy with respect to the target contrast are confirmed also when dealing with non-uniformly shaped scatterers (Figure 14).
Finally, it is interesting to underline the advantage of using the MT extension of the BCS-based approach for solving the 3D-CSI inversion problem instead of its "naive" single-task implementation (ST-BCS [23]) that does not impose any physical correlation in solving (14).Towards this end, let us consider the benchmark O = 6 voxel arrangement in Figure 15a with two (K = 2) "L-shaped" homogeneous (τ = 2.0) objects.For illustrative purposes, the reconstructions with the MT-BCS (Figure 15b,c) and the ST-BCS (Figure 15d,e) when processing two different sets of noisy data [SNR = 20 [dB]-Figure 15b,d; SNR = 10 [dB]-Figure 15c,e] are shown.As it can be visually inferred, the MT strategy turns out to be more effective than the ST one in both locating and shaping the non-connected scattering regions as well as in estimating the actual contrast value.Such an inference is quantitatively confirmed by the error indexes since 15b vs. Figure 15d) and 15c vs. Figure 15e) (Table 3).Similar outcomes can also be drawn when changing the contrast of the scatterers as indicated by the behaviour of the reconstruction error ξ tot versus τ in Figure 16. ) )-Behavior of the total error, ξ tot , as a function of the object contrast, τ, when processing the scattering data with the MT-BCS, the ST-BCS, and the CG methods.
To conclude the numerical assessment of the reconstruction capabilities of the MT-BCS, it has been compared with a competitive non-CS state-of-the-art approach.Towards this end, a deterministic CG-based inversion tool-still based on a CSI formulation of the scattering problem-has been applied to the same scenario in Figure 15a.The retrieved images of the dielectric profile of the investigation domain are shown in Figure 17a  Without imposing sparseness priors, only the presence of K = 2 scatterers lying in D can be deduced, but their contrasts are strongly under-estimated and their supports/shapes are unreliably predicted.Comparatively, the MT-BCS enables a reduction of the reconstruction error of about    3) times.Such conclusions are not limited to the contrast τ = 1.0, but they also arise for stronger scatterers as detailed by the plot of ξ tot versus τ in Figure 16.
As for the computational efficiency, the MT-constrained exploitation of a sparseness promoting inversion technique allows a non-negligible reduction of the computational time, ∆t, when processing 3D scattering data (For the sake of fairness, all the computational times refer to non-optimized Matlab codes executed on a single-core laptop running at 2.20 GHz).Indeed, the MT-BCS is not only faster than the ST-BCS, thanks to the joint processing of the data (i.e., ∆t MT−BCS ∆t ST−BCS ≈ 0.52-Table 3), but it also overcomes the CG speed (i.e., ∆t MT−BCS ∆t CG ≈ 0.032-Table 3).

Conclusions
An innovative approach to efficiently solve the full-vectorial 3D-IS problem has been presented.The retrieval of the volumetric contrast distribution of sparse non-weak scatterers has been tackled as a probabilistic CSI-based problem, which has been efficiently solved through a customized MT-BCS approach.The numerical analysis has pointed out the following key features of the proposed technique:

•
Reliable 3D reconstructions of the EM properties of the imaged domain are yielded processing scattering data also blurred with a non-negligible amount of additive noise;

•
The inversion accuracy of the proposed CS-based approach depends on the degree of sparseness of the actual scenario with respect to the expansion basis at hand.However, it can be fruitfully and profitably applied when other/different (non-voxel) representations of the contrast source/contrast function are chosen [46]; The MT implementation of the BCS-based inversion remarkably overcomes its single-task (ST-BCS) counterpart thanks to the profitable exploitation of the existing correlations between the V views and the scattered field components;

•
The MT-BCS positively compares with other state-of-the-art approaches, also deterministic and non-CS, in terms of both reconstruction accuracy and computational efficiency.
Moreover, the methodological advancements of this work with respect to the state-of-the-art on the topic [47,48] include (i) the generalization of the MT-BCS strategy to handle 3D-IS problems, differently from previous customizations of such an inversion paradigm which deal only with two-dimensional formulations [48], (ii) the derivation of a BCS-based imaging approach able to retrieve 3D target contrast information, unlike state-of-the-art Bayesian CS contributions only dealing with the reconstruction of equivalent sources [47], and (iii) the analysis and validation of suitable operative guidelines for the optimal setting of the key calibration parameters of the introduced methodology.Future works will be aimed at extending the capabilities of the proposed approach to effectively deal with non voxel-sparse targets as well as with other applicative scenarios of great interest (e.g., subsurface imaging) including the processing of multi-frequency data [55].

Figure 1 .
Figure 1.Geometry of the 3D-MI microwave imaging problem.

avg tot 1 W 2 ]
∑ W w=1 ξ tot,w ), and the variance [ξ var tot 1 W ∑ W w=1 ξ tot,w − ξ avg tot values of the total error among the W experiments are reported.

Figure 15c-Table
Figure 15c-Table3) times.Such conclusions are not limited to the contrast τ = 1.0, but they also arise for stronger scatterers as detailed by the plot of ξ tot versus τ in Figure16.As for the computational efficiency, the MT-constrained exploitation of a sparseness promoting inversion technique allows a non-negligible reduction of the computational time, ∆t, when processing 3D scattering data (For the sake of fairness, all the computational times refer to non-optimized Matlab codes executed on a single-core laptop running at 2.20 GHz).Indeed, the MT-BCS is not only faster
SN R =