Microwave Near-Field Dynamical Tomography of Thorax at Pulmonary and Cardiovascular Activity

The developed near-field microwave diagnostics of dynamical lung tomography provide information about variations of air and blood content depth structure in the processes of breathing and heart beating that are unattainable for other available methods. The method of dynamical pulse 1D tomography (profiling) is based on solving the corresponding nonlinear ill-posed inverse problem in the extremely complicated case of the strongly absorbing frequency-dispersive layered medium with the dual regularization method—a new Lagrange approach in the theory of ill-posed problems. This method has been realized experimentally by data of bistatic measurements with two electrically small bow-tie antennas that provide a subwavelength resolution. The proposed methods of 3D lung tomography based on the multisensory pulse, multifrequency, or multi-base measurements are based on solving the corresponding integral equations in the Born approximation. The experimental 3D tomography of lung air content was obtained by the results of the multiple 1D pulse profiling by pulse measurements in several grid points over the planar square region of the thorax. Additionally, the possible applicability of multifrequency measurements of scattered harmonic signals in the monitoring of lungs was demonstrated by four-frequency measurements in the process of breathing. The results demonstrated the feasibility of the proposed control in the diagnosis of some lung diseases.


Introduction
The proposed electromagnetic tomography is based on the near-field microwave monitoring of lung parameters (depth profiles of air and the blood relative content) in the processes of pulmonary and cardiovascular activity by measured variations in the scattered pulse shape. To date, there are no alternative methods for the dynamical monitoring of these biomedical parameters; thus. the proposed tomography can be applied as a new biomedical diagnostic of certain lung diseases. Algorithms of this dynamical tomography are based on solving the inverse scattering problem, and unlike most of the available methods based on signal imaging, the proposed tomography allows for obtaining quantitative distributions of the physical parameters of the lungs.
By now, studies on microwave diagnostics have mainly been associated with methods of the multisensory tomography of women breast tumors and brain tumors or strokes by measuring signal transmission parameters (e.g., [1][2][3][4][5][6]). A large number of measurement schemes, methods of visualization of corresponding absorbing inhomogeneities, and various measurement systems have been worked out and studied on phantoms and real patients (see in the most competent reviews [7,8]). However, results demonstrate that the quality of the visualization gives little essentially new information compared to the available methods of diagnostics. The conclusion of these reviews stated that "there are very few commercially available and cost-effective microwave-based systems for breast cancer imaging or detection" and "there is still an enormous amount of research and development to be carried out in order to achieve the full capacity of this technology".
Previous studies on the remote microwave diagnostics of human breathing and heartbeat were initially associated with the task of the radar detection of a human behind the walls of a building or under snow-and landslides, and they were used in some security applications [9][10][11][12][13][14][15][16][17]. In these works, the diagnostics consisted of the visualization (imaging) of signal variations received in the far-field zone that were mainly formed due to refection from the surface of the human body and are primarily related to mechanical vibrations of the skin surface in the process of breathing or heart activity. At that, because of the low penetration at high frequencies, signals scattered from the inner body structure were too small to discern this structure-especially at the poor resolution in the far field. To improve the penetration depth and resolution, as far back as 1979, it was proposed to use near-field contact sensors and lower frequencies (740 MHz-1.5 GHz), and they applied measurements at 915 MHz for the diagnosis of pulmonary edema [18]. Later, with the improvement of sensors, experimental studies at frequencies of 0.5-10 GHz were carried out [19,20], where phase variations of the signal related to breathing and heart activity were observed at the frequency of 0.92 GHz. Recently, measurements using advanced techniques and matched sensors have demonstrated the difference in heart-related signals between healthy people and people after a heart attack, as well as the sensitivity of scattered signals to the pulmonary edema [21].
These results and our experiments stimulated us to develop a new and mathematically consistent diagnostic based on the solution of the inverse scattering problems based on data of bistatic near-field measurements of ultra-wideband (0.05-5 GHz) scattered pulse signals. To realize the high-enough depth resolution of the diagnostics, near-field measurements were used, which, as in near-field optics, provide the necessary subwavelength resolution beyond the Rayleigh limit for deep-penetrating probing signals and sufficiently low frequencies. Such a resolution at low frequencies could not be achieved for far-field measurements, while the use of higher frequencies could only be effective for skin diagnostics because of their strong absorption in living tissues.
In this paper, we applied pulse sounding; however, since in our case the pulse length is comparable to the probing depth, we applied not only the pulse spectrum analysis [31] but also a new method based on the analysis of the shape of the scattered pulses itself. Measurements were carried out with the transceiver system that includes two identical bowtie antennas. The sharp localization of currents on the emitting antenna produces a broad k-space spectrum of the probing field that provides a high-depth resolution of subsurface inhomogeneities in the near zone (see Figure A1 in Appendix A). The same antennas were used in our previous studies on the near-field subsurface diagnostics [28][29][30][31]33].
The tomography analysis based on the solution of the nonlinear ill-posed inverse scattering problem is especially complicated in the considered case of time-varying inhomogeneities in the strongly absorbing and frequency dispersive lung tissues under a multilayer thorax cover. In the developed dynamical 1D tomography (profiling) of the air and blood content of the lungs in the processes of breathing and cardiac activity, we developed algorithms based on the dual regularization method, a new powerful method in the Lagrange approach in the theory of nonlinear ill-posed problems [30] that is free from restriction of the perturbation theory. The developed method also includes the retrieval of the multilayer thorax structure that covers lung tissues, where it turned out to be possible to apply another, simpler algorithm. The paper presents the results of a numerical and experimental study of such dynamical tomography. For 3D tomography, where the dual regularization method appeared too resource-consuming, we proposed 3D integral equations obtained in the Born approximation. Since a multisensory system was not available, we only implemented 1D dynamical tomography (profiling) with the hope that this diagnostic would be further upgraded to 3D tomography. Nevertheless, the first experimental results of 3D lung tomography are also presented here. In this experiment, 3D tomography images were obtained from air content depth profiles in the phase of inhalation retrieved using the 1D algorithm at several points of the 2D thorax surface.
Additionally, the possibility of realizing multifrequency methods based on measuring the complex amplitudes of the scattered harmonic signals was proposed. It should be noted that such measurements could be more informative than measurements of scattered realvalued pulse signals, but it appeared impossible to implement sufficiently fast frequency commutation on the available equipment to obtain data at the required number of frequencies. However, the sensitivity of multifrequency data to variations of the lung structure during breathing is demonstrated here by measurements at four chosen frequencies. The proposed electromagnetic tomography is based on the near-field microwave monitoring of lung parameters (depth profiles of air-and blood-relative content) in the processes of pulmonary and cardiovascular activity by measured pulse shape variations of the scattered signals. For data, there are no alternative methods for the dynamic monitoring of these biomedical parameters; thus, the proposed tomography can be used as a new biomedical diagnostic of lungs. The algorithms of this dynamical tomography are based on solving the inverse scattering problem and, unlike most existing methods based on signal imaging, allowing for obtaining quantitative distributions of the physical parameters of lungs.

Materials and Methods
Let us consider a scattering region with the complex permittivity depth profile ε(z, T) = ε 0i (z) + ε 1 (z, T) that consists of the basic multilayer thorax structure of tissues that cover lungs ε 0i (z) and the time-dependent depth profiles of lung permittivity ε 1 (z, T) that are determined by variations of air-and blood-relative content profiles ρ air (z, T), ρ blood (z, T) in lungs in the process of breathing and in heart activity. The parameters of the lung are assumed to be unchangeable during the duration of the pulse ∆t. The shape of the pulse (the dependence of the received real-valued pulse signal on time t) s(t, T) = s[ε(z, T)](t) is determined by profiles ε(z, T). Calculations of the pulse shape by parameters of the lungs are based on the plane wave decomposition of the probing and scattering fields and are described in Appendix A. The shape of the pulse signal scattered from the probed medium is expressed as the inverse Fourier transform of its frequency spectrum s(ω, T): where ω is the cyclic frequency. Because it would be unsuitable to use such resourceconsuming numerical methods as FDTD in pulse calculations, we applied here the approach [28,31] based on the plane wave decomposition of fields radiated by currents generated by a pulse input signal on the surface of a transmitted antenna. It is worth noting here that some problems with such a representation may arise when one solves the inverse problems of restoring inhomogeneities with sharp boundaries (with a derived singularity) because of their smoothing or distortion by the Gibbs effect. The measured scattered pulse signal can be obtained from complex amplitudes of the scattered harmonic field components s ω (ω, T) generated by the corresponding spectral components of the currents. Using the formulas of this theory given in Appendix A, the received signal can be calculated from the spectrum s ω (ω, T) for any complex permittivity profile as: where Re(*) marks the real part of the integral expression. For brevity, above and below, we used the same symbols for the respective Fourier transforms; thus, they are defined by their arguments. The spectrum of the signal in (1) is expressed as The statement of the inverse problem of the reconstruction of the depth profile ε(z) can be formulated conditions for the equality of the calculated and measured pulse parameters (pulse spectrum or shape). These functional equations are strongly nonlinear and underdetermined; it is impossible to obtain both parts of the complex-valued function of the two variables ε(z, ω) to form parameters of the real-valued scattered pulse. However, it is possible to simplify the statement of the problem by using its specific features.
You can first solve the problem of determining the cover layered structure of the thorax ε 0i and then the problem of restoring the profile of the lung permittivity. To determine ε 0i , it is necessary to use measurements in the phase of complete exhalation and diastole (ρ air (z) = 0, ρ blood (z) = 0), when the signal variations associated with the heartbeat stop, and the deflated lungs can be considered a homogeneous half-space with the known dielectric parameters. This strongly nonlinear ill-posed problem can be solved using, as a first guess, the available a priory information about this structure, including the known permittivity of tissues of layers [34] and their approximate (or measured by other means) thickness. The depth positions of the interfaces for five layers (skin, fat, muscles, ribs, and lungs) were chosen as unknown parameters. Because the number of these parameters is rather small, this discrete problem is not so ill-posed-it can be solved by minimizing of the quadratic discrepancy between the calculated and measured values of the signal spectrum using the enumerative technique in the reasonable range of layer thicknesses around the first guess z i0 : D ≡ z 1i = ∆zi, i = 1, .., i max ; z 2j = ∆zj, j = i max , .., j max ; z 3k = ∆zk, k = j max , .., k max ; z 4l = ∆zl, l = k max , .., l max , where ∆z is the step of discretization. To improve the retrieval accuracy, the optimal band of the analysis ∆ω o is chosen outside the more high-frequency spectrum of random errors based on the results of numerical simulation. The frequency spectra of the permittivity of the main (unchangeable) layered structure of the thorax are known for all their tissues, including the permittivity of the cover tissues and lungs in phases of complete exhalation and inhalation for all frequencies of the spectrum of the received signal [34]. Then, using the obtained cover lung structure ε 0i , it is possible to solve the inverse problems of restoring the profiles of air and blood content in the lungs in the process of breathing and cardiac activity. At the same time, to determine complex-valued variations of the permittivity of the lung profile ε(z, ω, T) from variations of the real-valued parameters of the air ρ air (z) and blood ρ blood (z) content, we used Rayleigh formulas [35] for dielectric mixtures of lung tissues with air and blood: ε − ε lung ε + 2ε lung = ρ air,blood ε air,blood − ε lung ε air,blood + 2ε lung (5) which leads to the nonlinear expression for the variation of permittivity ε 1 (ρ air,blood ) = ε − ε lung in lungs: ε 1 (ρ air,blood ) = ε lung + 2ε lung ρ air,blood ε air,blood −ε lung ε air,blood +2ε lung 1 − ρ air,blood ε air,blood −ε lung ε air,blood +2ε lung − ε lung .
Variations in the scattered signal associated with breathing and cardiac activity are rather small and cannot be seen everywhere on the surface of the chest-especially variations associated with the heartbeat, which are most often recorded in the full inhalation phase, since the blood wave fades in lung capillaries in the exhalation phase. Periodic processes associated with breathing and cardiac activity have different periods; moreover, breathing can be controlled. Therefore, it is possible to isolate the contributions of air and blood to the variations in the permittivity ε 1 and solve the corresponding inverse problems separately using the cover structure ε 0i , previously obtained from (4).
Let us first consider the statement of these inverse problems based on the spectral approach. The equations of the corresponding inverse problem can be written (omitting the time variable T) as the equality condition: where the depth profile ρ(z) = ρ air (z, T) or ρ(z) = ρ blood (z, T), and s[ρ(z)](ω) is the spectrum of the scattered signal calculated from (2), (A18) by the permittivity profile ε(z, ω) obtained for the given profile ρ(z) from (4) and (6); s 0 (ω) is the spectrum of the measured scattered pulse. To suppress the effect of random errors, the problem is solved in the informative frequency band ∆ω o outside the more high-frequency region of the error spectrum, where the signal exceeds the error level. This band is determined from results of numerical simulation.
To solve this nonlinear ill-posed inverse problem, the dual regularization method [30] was adopted. It consists of minimizing the modified Lagrange function at the simultaneous maximization of the corresponding regularized functional of the dual problem. The modified Lagrange function for this method is written as where ρ 2 L 2 = 1 ∆z ∆z ρ(z) 2 dz, λ = (λ 1 , λ 2 ), ∆z is the region of analysis,µ > 0, · is a scalar product, λ = (λ 1 , λ 2 ), s = (Res, Ims) marks two-dimensional vectors, and µ > 0. When the free parameter µ is large enough, the minimum of the modified Lagrange function L[ρ](λ) over ρ(z) exists for certain at any λ. The regularized dual problem is the problem of maximizing the concave functional in the Hilbert space L 2 2 (ω 1 , ω 2 ) that is expressed as where the maximum is found over λ from the set The desired solution is obtained as the saddle point of this process of minimization of (8) with respect to ρ at a simultaneous maximization of (9) with respect to the dual variable λ. The numerical algorithm that implements this method is based on the usual gradient minimization of (8) beginning with the zero first guess when maximizing the functional (9) using a quite obvious explicit expression for its gradient [30].
This method has been studied in numerical simulations and successfully applied to experimental data. However, this algorithm turned out to be too time-consuming for multiple operations in the diagnostics of dynamic processes, in particular due to operations with 2D Lagrange coefficients. To avoid this difficulty, we developed and applied an algorithm based on a functional equation that solves the considered inverse problems in terms of the real-valued shape of the scattered pulse. In this case, solution algorithms with 1D Lagrange coefficients are much less time-consuming and are proven to be suitable for repeated use in diagnosing dynamic processes of breathing and cardiac activity. The corresponding inverse problem equation can be written (omitting the argument T) as the equality condition: where s[ρ(z)](t) and s 0 (t) are shapes of the modified pulses calculated by the corresponding spectra s[ρ(z)](ω), s 0 (ω) that, such as the spectra in (7), are obtained in the informative frequency band ω ∈ ∆ω 0 . When solving this nonlinear equation, the dual regularization method includes 1D Lagrange coefficients that make it much less time-consuming. The desired solution ρ(z) is obtained as a saddle point of the process of minimizing the following modified Lagrange functional: (11) over the functional parameter ρ(z) while simultaneously maximizing the corresponding regularized concave functional of the following dual problem: (12) with respect to the 1D Lagrange coefficient λ from the set Λ µ ≡ {λ ∈ L 2 : ||λ||≤ µ}, D = {ρ ∈ L 2 (z n , 0) : 0 ≤ ρ(z) ≤ ρ max }; ∆t is the pulse length. The algorithm of this method was studied in numerical simulation and applied in the experimental study of the proposed dynamical tomography.

Numerical Simulation
As it is well-known, in the considered ill-posed inverse problems, it is impossible to establish a universal relationship between the error parameters (such as the signal-to-noise ratio) and the solution accuracy parameters, which strongly depend on the complexity of the restored function. Each algorithm for solving any ill-posed problem, developed for a specific application, must be studied, tested, and optimized in numerical simulations, taking into account the parameters of the medium under study, the available equipment, measurement errors, and typical and extremal functions expected under the given conditions to be restored. In Figure 1, the pulse shape and its spectrum are shown for the generator used in the experiments.
As it is well-known, in the considered ill-posed inverse problems, it is impossible to establish a universal relationship between the error parameters (such as the signal-tonoise ratio) and the solution accuracy parameters, which strongly depend on the complexity of the restored function. Each algorithm for solving any ill-posed problem, developed for a specific application, must be studied, tested, and optimized in numerical simulations, taking into account the parameters of the medium under study, the available equipment, measurement errors, and typical and extremal functions expected under the given conditions to be restored.
In Figure 1, the pulse shape and its spectrum are shown for the generator used in the experiments.  Figure 2 shows the optical thickness of the probing signal (integral of the absorption coefficient) in a five-layered living tissue (skin, fat, muscle, ribs, and lungs) for the inhalation and exhalation phases in the signal frequency band. It can be seen that the frequencydepth distribution of the attenuation of the probing field in the lungs is sensitive to the air content in the lungs, at least up to penetration depths of about 10 cm, which provides the informativity of the proposed method.   Figure 2 shows the optical thickness of the probing signal (integral of the absorption coefficient) in a five-layered living tissue (skin, fat, muscle, ribs, and lungs) for the inhalation and exhalation phases in the signal frequency band. It can be seen that the frequencydepth distribution of the attenuation of the probing field in the lungs is sensitive to the air content in the lungs, at least up to penetration depths of about 10 cm, which provides the informativity of the proposed method. establish a universal relationship between the error parameters (such as the signal-tonoise ratio) and the solution accuracy parameters, which strongly depend on the complexity of the restored function. Each algorithm for solving any ill-posed problem, developed for a specific application, must be studied, tested, and optimized in numerical simulations, taking into account the parameters of the medium under study, the available equipment, measurement errors, and typical and extremal functions expected under the given conditions to be restored.
In Figure 1, the pulse shape and its spectrum are shown for the generator used in the experiments.  Figure 2 shows the optical thickness of the probing signal (integral of the absorption coefficient) in a five-layered living tissue (skin, fat, muscle, ribs, and lungs) for the inhalation and exhalation phases in the signal frequency band. It can be seen that the frequencydepth distribution of the attenuation of the probing field in the lungs is sensitive to the air content in the lungs, at least up to penetration depths of about 10 cm, which provides the informativity of the proposed method.  The numerical simulation of the cover structure retrieval included: for the simulated five-layer structure in the exhalation phase; (b) Calculation of the corresponding scattered pulse s(t); (c) Adding simulated uncorrelated normally distributed random errors with rms that corresponded to that in the measured data; (d) Recalculating the pulse spectrum with errors; (e) Allocation of the informative analysis band ∆ω o outside the band of the more highfrequency spectrum of random errors; (f) Solving the inverse scattering problem (4) and comparison of the results with the preset structure z 1 , z 2 , z 3 , z 4 .
The numerical simulation of the ρ air (z) profile retrieval using algorithms (10)-(12) included the calculation of the signal spectrum s[ρ(z)](ω) and pulse shape s[ρ(z)](t) with the cover layer structure retrieved from (4), adding random errors, the recalculation of the spectrum with errors in the analysis band ∆ω o , the calculation of the modified pulse s[ρ(z)](t), the solution of the inverse problem using algorithms (10)- (12), and a comparison of the results with the preset profiles. Figure 3 shows the scattered pulse signals and their frequency spectra calculated from the known parameters of the complex permittivity for the inflated and deflated lung tissues. In these calculations, lungs were considered a half-space with constant permittivity under the four-layered cover tissues (skin, fat, muscle, ribs).
(e) Allocation of the informative analysis band o   outside the band of the more highfrequency spectrum of random errors; (f) Solving the inverse scattering problem (4) and comparison of the results with the preset structure The numerical simulation of the () air z  profile retrieval using algorithms (10)- (12) included the calculation of the signal spectrum , the solution of the inverse problem using algorithms (10)- (12), and a comparison of the results with the preset profiles. Figure 3 shows the scattered pulse signals and their frequency spectra calculated from the known parameters of the complex permittivity for the inflated and deflated lung tissues. In these calculations, lungs were considered a half-space with constant permittivity under the four-layered cover tissues (skin, fat, muscle, ribs). As seen in Figure 3, the differences between the pulse shapes calculated for inflated and deflated lungs are rather small, while the corresponding differences in the spectra are more pronounced. Our simulation showed the approximate frequencies at which the spectral contribution of the modeled inhomogeneities disappeared. It turns out that the spectrum of random errors is predominantly located at higher frequencies, which makes it possible to limit the spectral bandwidth of the analysis, eliminating the influence of the corresponding components of these errors. As seen in Figure 3, the informative part of their spectra lies at f < 1.5 GHz. The possible effect of the experimental random is investigated by calculating the spectrum of the scattered signal with added random errors at an rms value of 1% (relative to the pulse amplitude) that corresponds to the estimated error level of the measurement noise. It was obtained that the spectral errors increased at frequencies f > 1 GHz; thus, the informative band appeared even narrower. The bandwidth of 0.05-0.9 GHz was found and used to be optimal in the further processing. To suppress As seen in Figure 3, the differences between the pulse shapes calculated for inflated and deflated lungs are rather small, while the corresponding differences in the spectra are more pronounced. Our simulation showed the approximate frequencies at which the spectral contribution of the modeled inhomogeneities disappeared. It turns out that the spectrum of random errors is predominantly located at higher frequencies, which makes it possible to limit the spectral bandwidth of the analysis, eliminating the influence of the corresponding components of these errors. As seen in Figure 3, the informative part of their spectra lies at f < 1.5 GHz. The possible effect of the experimental random is investigated by calculating the spectrum of the scattered signal with added random errors at an rms value of 1% (relative to the pulse amplitude) that corresponds to the estimated error level of the measurement noise. It was obtained that the spectral errors increased at frequencies f > 1 GHz; thus, the informative band appeared even narrower. The bandwidth of 0.05-0.9 GHz was found and used to be optimal in the further processing. To suppress the noise effect during further processing by the pulse signal, it was proposed to recalculate it using the inverse Fourier transform by the informative part of its frequency spectrum and to use the resulting modified pulse s(t) in solving the corresponding inverse problems (4), (10)- (12).
In Figure 4, one can see the results of the numerical simulation of air content profiling (10)- (12) in combination with the reconstruction (4) of the five-layer structure of the permittivity of cover tissues at the rms of random errors δ s = 1%.
The numerical simulation of the 1D tomography presented in Figure 4 was carried out for two preset profiles of air content ρ air (z) in the lungs. The monotonous exponential profile corresponds to the most likely simple profiles during breathing; the Gaussian profile with a sharp maximum simulated a layer of air accumulated under the pleura during pneumothorax. Algorithms (10)-(12) used modified pulses (obtained by recalculating the informative part of their spectrum outside the more high-frequency spectrum of random errors) as the input "measured data". late it using the inverse Fourier transform by the informative part of its frequency spectrum and to use the resulting modified pulse () st in solving the corresponding inverse problems (4), (10)- (12).
In Figure 4, one can see the results of the numerical simulation of air content profiling (10)- (12) in combination with the reconstruction (4) of the five-layer structure of the permittivity of cover tissues at the rms of random errors s  = 1%.   Figure 4d). The middle row shows simulated (lines 1) and retrieved (lines 2) lung air content profiles ρ air (z): (c)-Gaussian preset, and retrieved from the data in Figure 4a (line 1) profiles; (d)-exponential preset and retrieved from the data in Figure 4a (line 2) profiles. The lower row shows the preset (lines 1) and retrieved (lines 2) profiles of the real part of the thorax permittivity ε (z) at a frequency f = 250 MHz, which combine the structure of the chest cover reconstructed from (4) and the permittivity profiles in the lungs: (e)-profiles correspond to the Gaussian profiles in Figure 4c; (f)-profiles correspond to the exponential profiles in Figure 4d.
In this simulation, a preset structure of the cover tissues ε 0i (z) was used with the inflated lungs (with full inhalation). The same structure was used to calculate the signal parameters in Figure 3. The retrieved structure of the cover tissues was retrieved from the solution of the inverse problem using algorithm (4) from the informative part of the spectrum of the simulated modified signal. In Figure 4e,f, the retrieved and preset real parts of the thorax permittivity profiles ε(z) = ε 0i (z) + ε 1 [ρ air (z)] are shown for the frequency f = 250 MHz. In these profiles, the cover structure ε 0i (z) obtained from (4) was combined with lung permittivity ε 1 [ρ air (z)] for the Gaussian and exponential profiles ρ air (z) (lines 2 in Figure 4c,d) retrieved from the solution of the inverse problems (10)- (12).
From this simulation, it was found that the rms relative value of respiration-related variations in the scattered signal (informative part) was about 5%; thus, the simulated measurements with errors of 1% are quite informative, especially considering that the applied spectral cleaning significantly suppresses noise effect and, consequently, increases the information content of the input data. The errors in retrieved positions of the boundaries z 1 , z 2 , z 3 , z 4 determined by (4) ranged from 0.5 mm for z 1 to δz = 1.5 mm for z 4 . Such restoration accuracy provides a sufficiently high quality of restoration both for a monotonic profile and for a profile with a sharp maximum, as shown in Figure 4.
Based on the presented results, it can be concluded that the proposed algorithms can be effective for lung diagnostics.

Methods of 3D Tomography
The development of 3D subsurface tomography for diagnosing the 3D dynamics of air and blood content in the lungs is a much more difficult task. Its statement can be based on the three-dimensional nonlinear integral relation (A19) between the scattered electric field, the probing field, and the scattered inhomogeneity of the complex permittivity, which was obtained from Maxwell's equations. In general, even its solution as a direct problem of calculating the scattered field is too resource-intensive, which makes it impossible to use (A19) when solving a 3D inverse problem of retrieving the permittivity with the double regularization method described above. In Appendix A, we propose methods for 3D permittivity tomography based on the Born approximation (A20) of Equation (A19), which, of course, leads to serious limitations of possible applications. This approach can be applied to tomography 3D variations in blood content, since the corresponding variations in the dielectric constant of the lungs are quite small; but in the case of 3D diagnostics of the air composition, the Born approximation could be valid only for a sufficiently small inhalation of the lungs.
Unfortunately, even in the Born approximation (A20), we encountered difficulties in the numerical solution of this 3D integral equation. To overcome these difficulties, we proposed and applied the method based on 2D scanning data with a bistatic transceiver system at a fixed source-receiver offset (base) δr or at δr = 0. Under this condition, Equation (A20) can be represented by the convolution equation with respect to transversal coordinates (A22) that, using Fourier transform with respect to these coordinates, is reduced to the onedimensional Fredholm integral Equation (A23) for the scattered field in the k-space. Using this equation, the corresponding integral equations were obtained for the multifrequency tomography (A24), for the pulse tomography (A27), and for the single-frequency multi-base tomography based on measurements of the scattered signal with the variable parameters δx, δy, δz of the source-receiver base δr (A29). These equations must be multifold solved for each pair of k-spectrum components k x , k y ; then, the 3D complex permittivity is obtained in Cartesian coordinates by the inverse 2D Fourier transform.
It is clear that scanning methods in the 3D tomography of the time-dependent media are inapplicable; thus, it can be realized only with a multisensory system of small enough sensors, where data should be properly measured to reproduce data of scanning at the fixed source-receiver base δr.
The condition ε 1 << ε can be met for small lung inhalation ρ air << 1, which gives us the following linear expression: Note that, in this case, it is necessary to use the value of ε lung for the lungs in the full exhalation phase [34].
In the diagnostics of variations related to the blood content, where ε blood ≈ ε lung , one has ε + 2ε lung ≈ 3ε lung , ε blood + 2ε lung ≈ 3ε lung , and we obtained the required linear expression: where it is necessary to use the value of ε lung for lungs in the full exhalation phase [34]. Then, using the Born approximations (A24)-(A29) of the nonlinear integral Equation (A19), we obtained the tomography equations for the methods mentioned above, omitting the time variable T. The k-space integral equation of the multifrequency tomography was obtained from (A24), (14), (15) as: where K air m f = − 3 2 ε lung K ω , K blood m f = (ε blood − ε lung )K ω , and K ω is expressed by (A23). The k-space integral equation of the pulse spectrum tomography was obtained from (A27), (13), (14) as: where K air pulse = − 3 2 ε lung K, K blood pulse = (ε blood − ε lung )K, and K is expressed by (A28). For the air content tomography, it may be more efficient to use a more general Equation (13) that gives us the following equations for multifrequency and pulse tomography: One can also propose a method of multi-based single-frequency tomography, in which frequency independence is the main advantage. In this case, the integral equation in the k-space is defined simply by Formula (A29): with the kernel (A30). Functions s ω and K ω are the same as in the multifrequency tomography (15) but at ω = const. The solving algorithm can be based on the signal dependence on any of the base parameters δx, δy, δz. The most serious problem of this method is the difficulty of placing the required number of sensors on the surface of the chest.
Equations (16)- (20) are ill-posed Fredholm integral equations of the first kind that must be multiplied for each pair k x , k y , for example, by using the algorithm of Tikhonov's method of generalized discrepancy in the Hilbert space W 1 2 developed in [27]. The desired 3D structure of the air-or blood-relative content in lungs in Cartesian coordinates should be obtained from the k-space solution of corresponding equations using the 2D inverse Fourier transform. For methods (16), (17), it was obtained from: ρ air,blood (x, y, z) = ρ air,blood (k x , k y , z) exp(ik x x + ik y y)dk x dk y .

Dynamical Pulse 1D Tomography (Profiling)
In our experiments, we used the equipment that includes processing, a transceiver system with two identical bow-tie antennas with the length of arms 3.8 cm and the width of 5.4 cm that were used in [28][29][30]33], an oscilloscope GZ10E, broadband (0.05-5 GHz) pulse generators GZ1120ME-01,03 with the pulse time 50 ns at 1024 sweep points (a resolution of about 0.05 ns) and a repetition time of 10 µs, and the control computer for recording signals. The pickup time added up to 10 s (500 pulses). The safety of these measurements in the band of cellular communications was confirmed by the conclusion of the Ethics Committee that was based on estimations based on data of measurements on biological tissues, defining that the upper limit of the mobile phone radiation determined in Europe by the SAR parameter (Specific Adsorption Rate) as 2 W/kg (calculated for 10 g of tissues) more than an order of magnitude exceeds the SAR value estimated for our system radiation.
Variations in the profile of the air content associated with breathing and cardiac activity were rather small and could be observed with our equipment only in some places on the chest; variations in the blood content profile were even smaller and were observed only in the phase of the full inhalation. They stopped in the exhalation phase, when the pulse wave stopped in the capillaries of the lungs. It can be noted that the blood wave propagated only in the capillaries of the lungs and not in the capillaries of different tissues. Of course, the sensitivity of measurements may also depend on the individual depths of the cover tissues. Figures 5 and 6 show the measurement scheme and experimental setup.
The results of the experiments are presented in Figure 7. These results showed that variations of the received signals in the processes of respiration and cardiac activity were clearly visible; moreover, in more detailed images (Figure 7c,d), these variations significantly exceed random errors, which are also clearly visible. The calculated level of these errors was about 1%; the same level was assumed in the numerical simulation shown in Figure 4. It can also be seen that the changes in signal amplitude during inhalation and exhalation in Figure 7a are similar to those observed in the simulated signals in Figure 3a. The effect of heart activity shown in Figure 7b,d is smaller, but as can be seen in Figure 7d  The results of the experiments are presented in Figure 7.  The results of the experiments are presented in Figure 7. The images in Figure 7e,f represent variations in air and blood content and were obtained from the measured signal variations in Figure 7a,b using algorithms (4), (10)- (12). They visually demonstrate the processes of the inflation and blood filling of the lungs during breathing and cardiovascular activity (they can be seen even more clearly in the animation given in the animation.pptx file included in the Supplementary Materials). Once again, we note that this analysis provides a quantitative diagnosis of the considered physical parameters of the lungs. Thus, the obtained results demonstrate the feasibility of the proposed diagnostics.

Tomography Based on Multi-Position Profiling
In the absence of a multi-sensor system and in the approximation of a smoothly inhomogeneous ρ air (x, y, z) in x and y directions, we tried to implement the 3D tomography of the lung air content ρ air (x, y, z), applying one-dimensional pulse profiling at 9 grid points over a square x-y area with sizes 6 × 6 cm on the chest with the step ∆x = ∆y = 3 cm. Measurements were taken at each point, with breath holding in the inhalation phase. To increase the depth resolution in the transceiver system, smaller bow-tie antennas with an arm length of 2.3 cm and a width of 3.0 cm were used (see in Figure 8).
Due to the small number of point measurements, we used a less resource-consuming 1D algorithm (4, 7-9) for solving the inverse problem from the spectrum of the received pulse signal to obtain the 3D distribution of air content ρ air (x, y, z) in the lungs. Images of such pulse tomography are shown in Figure 9. These results showed that variations of the received signals in the processes of respiration and cardiac activity were clearly visible; moreover, in more detailed images ( Figure  7c,d), these variations significantly exceed random errors, which are also clearly visible. The calculated level of these errors was about 1%; the same level was assumed in the numerical simulation shown in Figure 4. It can also be seen that the changes in signal amplitude during inhalation and exhalation in Figure 7a are similar to those observed in the simulated signals in Figure 3a. The effect of heart activity shown in Figure 7b,d is smaller, but as can be seen in Figure 7d, it significantly exceeded the level of measurement errors.
The images in Figure 7e,f represent variations in air and blood content and were obtained from the measured signal variations in Figure 7a,b using algorithms (4), (10)- (12). They visually demonstrate the processes of the inflation and blood filling of the lungs during breathing and cardiovascular activity (they can be seen even more clearly in the animation given in the animation.pptx file included in the Supplementary Materials). Once again, we note that this analysis provides a quantitative diagnosis of the considered physical parameters of the lungs. Thus, the obtained results demonstrate the feasibility of the proposed diagnostics.

Tomography Based on Multi-Position Profiling
In the absence of a multi-sensor system and in the approximation of a smoothly inhomogeneous ( , , ) air x y z  in x and y directions, we tried to implement the 3D tomography of the lung air content ( , , ) air x y z  , applying one-dimensional pulse profiling at 9 grid points over a square x-y area with sizes 6  6 cm on the chest with the step 3 cm xy  =  = . Measurements were taken at each point, with breath holding in the inhalation phase. To increase the depth resolution in the transceiver system, smaller bowtie antennas with an arm length of 2.3 cm and a width of 3.0 cm were used (see in Figure  8). Due to the small number of point measurements, we used a less resource-consuming 1D algorithm (4, 7-9) for solving the inverse problem from the spectrum of the received pulse signal to obtain the 3D distribution of air content ( , , ) air x y z  in the lungs. Images of such pulse tomography are shown in Figure 9.  Due to the small number of point measurements, we used a less resource-consuming 1D algorithm (4,(7)(8)(9) for solving the inverse problem from the spectrum of the received pulse signal to obtain the 3D distribution of air content ( , , ) air x y z  in the lungs. Images of such pulse tomography are shown in Figure 9. One can see in Figure 9 that the observed spatial variations in the structure of the air content ( , , ) air x y z  look realistic. It can be noted that these variations in the x-and y- One can see in Figure 9 that the observed spatial variations in the structure of the air content ρ air (x, y, z) look realistic. It can be noted that these variations in the x-and y-directions (along and across the direction of the ribs) are comparable; we do not see noticeable effects of anisotropy. This confirms the reliability of the results of dynamical tomography presented in Figure 4. It is important to note that this method, if the conditions of its applicability are met, could be more efficient than the methods described above based on solving three-dimensional equations, since the one-dimensional algorithm based on the powerful dual regularization method is free from the limitations associated with the Born approximation and the nonlinearity of the problem.

Multifrequency Observations of Breathing
As mentioned above, multifrequency measurements can be more informative than pulse sounding. With multi-frequency respiration monitoring, it is necessary to ensure a fast commutation between frequencies that lung parameters could be considered unchangeable during the commutation. Such data can be used in the future when solving the corresponding inverse problem for the implementation of the dynamic tomography of air and blood lung contents.
Unfortunately, it turned out to be impossible with the available equipment to provide sufficiently fast frequency commutation to increase the number of frequencies to the level necessary for solving the corresponding inverse problems of tomography. However, to demonstrate the possibility of such diagnostics, we implemented such measurements of the complex amplitudes of scattered harmonic signals during breathing at four frequencies (2, 3, 4, 5 GHz), which were carried out with the same transceiver system. Sufficiently fast four-frequency commutation was implemented using a variable-frequency oscillator. Data were registered for 10 s with a time step of 10 −5 s −1 and a single-frequency commutation time of 0.01 s so that the total commutation time was 0.04 s. Figure 10 shows the results of measurements during breathing.
The results show that the measurements are sensitive to lung transformation during respiration; moreover, this transformation is very uneven at different frequencies. This sensitivity indicates the possibility of using these measurements for subsurface multifrequency lung diagnostics in 1D profiling based on (4) and (8) or in 3D tomography (16) and (18). demonstrate the possibility of such diagnostics, we implemented such measurements of the complex amplitudes of scattered harmonic signals during breathing at four frequencies (2,3,4,5 GHz), which were carried out with the same transceiver system. Sufficiently fast four-frequency commutation was implemented using a variable-frequency oscillator. Data were registered for 10 s with a time step of 10 -5 s -1 and a single-frequency commutation time of 0.01 s so that the total commutation time was 0.04 s. Figure 10 shows the results of measurements during breathing. The results show that the measurements are sensitive to lung transformation during respiration; moreover, this transformation is very uneven at different frequencies. This sensitivity indicates the possibility of using these measurements for subsurface multifrequency lung diagnostics in 1D profiling based on (4) and (8) or in 3D tomography (16) and (18).

Discussion
Let us compare the developed dynamic tomography with the existing methods of subsurface diagnostics and consider the areas of possible biomedical application of this microwave tomography. Available methods for diagnosing pathologies in various

Discussion
Let us compare the developed dynamic tomography with the existing methods of subsurface diagnostics and consider the areas of possible biomedical application of this microwave tomography. Available methods for diagnosing pathologies in various diseases include magnetic resonance imaging (MRI), radiography, fluoroscopy, X-ray computed tomography (CT), and various types of ultrasound computed tomography (USCT) [36]. These methods are based on the use of various waves (electromagnetic, acoustic) capable of penetrating the human body and being transformed by the inhomogeneities of its tissues. The measured parameters of transmitted, reflected, or scattered waves provide valuable information about these inhomogeneities. However, all available methods have fundamental limitations.
X-ray images visualize the transverse 2D distributions of ray absorption in the human body (integral of the absorption coefficient along the ray). In these images, the shadows from the near-surface layers of human tissues are superimposed on the shadows from deeper layers, making them difficult to interpret. Three-dimensional CT images showed the extinction coefficient distribution (depending mainly on tissue density) in each cross section of the body, which visualizes the structure of living tissues. This method is much more informative than the method of simple X-ray diagnostics; however, it can damage DNA, is technically very complex, and is not free from artifacts. Unlike X-ray methods, our microwave diagnostics are not capable of damaging DNA; a probing signal with a spectrum in the cellular frequency band has a much lower thermal effect than mobile phones. Additionally, this ultra-wideband pulsed signal is not capable of causing any resonant effects in the DNA of living tissues. The safety of the measurements was confirmed in the supplied conclusion of the Ethics Committee.
Ultrasonic computed tomography (USCT) methods use high-frequency sound reflections from tissue interfaces and scattering from density inhomogeneities in the human body. The ability to select the depth in the reflection signal allows you to obtain the tomographic visualization of the tissue structure with a resolution determined by the choice of frequency. In addition, Doppler measurements make it possible to observe and measure currents in blood vessels. The information in conventional USCT data is mostly qualitative and is more difficult to express in physical terms than in other tomography methods; thus, its interpretation requires qualified service personnel.
Summing up, we can say that all of the above methods are based mainly on signal visualization and provide qualitative information about the tissues of the human body, which requires experienced highly qualified medical personnel for interpretation, as well as additional data (symptoms, blood tests, etc.). Our diagnostics provided quantitative information about the physical parameters of the lungs, which is unattainable for these methods.
In addition, few of the available methods are suitable for monitoring lung parameters, such as the recently developed real-time MRI (magnetic resonance imaging) [37]. The MRI method, based on the measurement of radio emission during the spin-orbit transition of hydrogen atoms in a magnetic field, visualizes the distribution of hydrogen atoms (mainly in water molecules) in living tissues, especially in aqueous fluids such as blood and lymph. In Figure 11, two screenshots of thorax from the computer playback of the web site video [38] with the animation of sectional images in the real-time MRI (Figure 11a,b) were given in comparison with the screenshots from the animation of our lung tomography (Figure 11c,d) (see also in the animation.pptx file in Supplementary Materials). As is easily seen in the MRI images (Figure 11a,b), the only visible objects in lungs are small white spots that correspond to the cross-sections of the largest blood vessels. They appeared in the systole (Figure 11a) and disappeared in the diastole (Figure 11b), where it is also possible to see the shift of the diaphragm position at breathing. In Figure 11c,d, where screenshots of our dynamical tomography are demonstrated, one can see that our tomography gives lung parameters that are unattainable by the real-time MRI. Here, it should be emphasized once again that, unlike MRI and other methods based on tomography (signal imaging), our tomography gives quantitative values of the physical parameters of the lungs.
Let us consider the pathologies in which the developed tomography can be applied. One of them is pneumothorax (the concentration of air and other gases in the pleural cavity), which can occur spontaneously (primary) or as a consequence of lung disease or chest trauma (secondary). The possibility of the proposed microwave diagnostics of such air layers in the lungs is demonstrated by numerical simulation in Figure 4c,d. Additionally, our method can be used to control operations of artificial pneumothorax, the introduction of air into the pleural cavity to obtain the collapse (atelectasis) of the affected lung.
Atelectasis (partial or complete pulmonary collapse of the lung) may be associated with a variety of lung pathologies. Primary atelectasis is neonatal atelectasis, in which the lungs remain unexpanded; secondary atelectasis may be associated with various lung diseases such as pneumonia, tumors, pulmonary shocks, pleural emphysema, lung injury (pneumothorax, hemothorax), hydrothorax, foreign body aspiration, or other pathologies. Corresponding strong distortions of air exchange profiles and their dynamics can be detected by the developed microwave method.
Some pulmonary pathologies are associated with the concentration of fluid in the lungs during pulmonary edema (blood plasma transudate from capillaries) associated with certain cardiovascular diseases of the lungs or with traumatic injuries [39]. Very similar symptoms can be associated with hydrothorax-the accumulation of fluid (transudate) in the pleural cavity. This may appear with cardiac decompensation. Pleural effusion can also be associated with hemothorax (blood concentration) or after cardiac surgery [39]. The rapid development of such pathologies requires constant monitoring, which can be provided by the proposed microwave diagnostics of the dynamics of air and blood content in the lungs.
Of great interest may be research on microwave diagnosis and the monitoring of COVID-19 pulmonary pathologies and their consequences. During the pandemic, CT imaging proved to be the most effective in diagnosing relevant lung pathologies that appear on CT images as "frosted glass" spots. However, long-term use of this method could be associated with exposure to X-rays. MRI diagnostics turned out to be significantly less effective and are used only as a background study. In this regard, it seems appropriate to study the possibility of using the developed method in the diagnosis of pathologies associated with COVID-19 [40,41].

Conclusions
In this work, the method of dynamical lung tomography was developed and studied both in numerical simulations and experimentally. The proposed near-field microwave pulse diagnostics made it possible to obtain information about the depth profiles of the relative content of air and blood in the lungs during breathing and the heart beating. This information is unattainable by other available methods based mainly on signal visualization (imaging), while the proposed tomographic analysis gives the physical parameters of the lungs. The physical novelty of this study lies in the recovery of the desired parameters from the shape of the scattered pulse from the solution of a nonlinear ill-posed inverse problem in the extremely complicated case of a highly absorbing and frequency-dispersive layered medium. The solution algorithms are based on the double regularization method in the Lagrangian approach of the theory of ill-posed problems, a new method that is free from the limitations of perturbation theory. They involve the shape of the scattered pulse or its frequency spectrum as input into the analysis. The proposed near-field tomography was experimentally implemented based on bistatic measurements of scattered UWB (0.05-5 GHz) pulsed signals using a transceiver system based on electrically small bow-tie antennas providing subwavelength resolution. The retrieved dynamics of the air and blood content profiles during breathing and cardiac activity were demonstrated.
In the Born approximation, integral equations were obtained that generalized the theory of three-dimensional microwave tomography to the cases of pulse and multifrequency measurements. In the absence of a multi-sensor system, we obtained and demonstrated in this paper 3D lung tomography derived from air content profiles obtained by measurements at multiple points on a square chest area in the inhalation phase.
Additionally, the possibility of implementing such tomography based on the data of multifrequency measurements of the complex amplitudes of scattered harmonic signals during four-frequency measurements of breathing was experimentally studied. The results show that measurements are sensitive to the lung transformation; therefore, such measurements could be applied to the subsurface multifrequency diagnostics of lungs. To perform this, it is necessary to increase the number of frequency channels to the level required to solve the corresponding inverse problems of tomography. Unfortunately, we failed to provide sufficiently fast frequency commutation on the existing equipment.
The possibility of using the proposed tomography for the medical diagnosis of certain lung diseases was discussed. Of course, the results shown in this article are very preliminary and are still very far from real applications in medical diagnostics. Individual variability of the structure and permittivity of tissues should be studied and taken into account (see, for example, [42]). We also hope that the measurement schemes, instrumentation, and sensitivity will be substantially optimized. However, these first results are encouraging and should stimulate further research.
Summarizing, we conclude that the presented results of modeling and experimenting demonstrate the possibility of using the developed dynamical tomography in the biomedical diagnostics of lung parameters that are unattainable in diagnostics by other known methods. These parameters are very informative in the diagnosis of pulmonary edema, pneumo-and hydrothorax, and atelectasis pneumonia and its consequences. The possible applicability of this new tomography can also be considered for the diagnosis of cardiopulmonary pathologies of COVID-19.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

Appendix A Appendix A.1. One-Dimensional Diagnostics
The shape (dependence on time) of the received real-valued impulse signal can be expressed by complex amplitudes of harmonic signals calculated at individual frequencies of the input signal spectrum (3) [23], which, omitting the time variable T, is expressed as: Let us consider the three-dimensional electromagnetic problem for a transceiver system that was used in the experiments. This system includes two identical planar bow-tie antennas placed on the surface of the body. The scattered harmonic signal s f (x r , y r , z r , ω) received by the antenna attached to the point x = x r , y = y r , z = z r of Cartesian coordi-nates is determined as a convolution of the transversal distribution of the scattered vector electrical field E 1 (x, y, z = z r , ω) and the receiver antenna function F: s ω (x r , y r , z r , ω) = E 1 (x, y, z = z r , ω)F(x r − x, y r − y, ω)dxdy. (A2) The two-dimensional Fourier transform (A2) in transverse coordinates gives the relationship in k-space between the signal and the plane wave decomposition in the form of a scalar product: where k x and k y are wavevector components, and F(k x , k y , ω) represents the k-space pattern of the receiving antenna.
To obtain the k-space pattern of the receiving antenna F, we obtained the k-space pattern of the transmitting antenna F tr (k x , k y , z, ω) through its field E 0 (x, y, z, ω) radiated by the antenna excited by the input pulse signal V(t) with the frequency spectrum V(ω). Spectral components of V(ω) produce the corresponding current on the antenna aperture with the spectral distribution j(r , ω) that was calculated in our measuring system by the built-in algorithm. The spectral components of the probing field E 0 (r, ω) in a surrounding media with the permittivity ε are determined by the current distributions j(x , y , z ) = j(x , y )δ(z − z tr ) of the antenna attached to the point x = 0, y = 0, z = z tr by the convolution with the Green function (tensor) The corresponding transversal field spectrum (plane wave decomposition) is: that is represented by the exact expression: where j(k x , k y , ω) = 1 j(x, y, ω)e −ik x x−ik y y dxdy. For field spectral components, the expression (A6) can be written in a compact form as: where the tensor g ji is the Green function ↔ g (k x , k y , ω) = ↔ G(k x , k y , z tr , z tr , ω) (subscribes i and j and denotes x, y, z components; the convenient summation over repeated indices is assumed); the double sign in (7) is plus at z − z tr > 0 and it is minus at z − z tr < 0; Cartesian coordinates, and ε = 1 (for antenna in the free space). As is seen in (A7), the field spectrum can be determined by its distribution on the antenna aperture as: On the other hand, by definition, this spectrum is expressed as: From (A8) and (A9), the expression for the transmitting antenna k-space pattern can be obtained: Then, taking into account that the antenna system under consideration is linear and necessary conditions of reciprocity are fulfilled, the desired k-space pattern of the receiving antenna F(k x , k y , ω) is obtained from the k-space pattern of the antenna in the transmitting mode (A10) calculated on the plane of the antenna z = z tr : where the constant is obtained by the calibration (may be frequency-dependent, taking into account the transfer function of the antenna-receiver communication line).
In the near-field zone, the antenna k-space pattern depends on the vertical distance z, whereas in the far-field zone, where near-field (evanescent) components with are fading, the field is mainly determined by far-field (propagating) components with k 2 ⊥ < Re(k 2 ) and the antenna radiation pattern does not depend on the distance z. The far-field components also determine the angular pattern of the antenna in the far zone that does not depend on the distance. In the probing field of electrically small (with subwavelength sizes) antennas, evanescent components dominate in the near-field zone, providing the sensitivity of the near-field diagnostics to subwavelength-size inhomogeneities in this zone. The resolution decreases with distance from the transceiver to the Rayleigh far-field resolution limit.
In Figure A1, the distribution of the current y-component j y (x, y, ω) on the electrically small bow-tie antenna and its k-space spectrum j y (k x , k y , ω) are demonstrated at frequency 1.5 GHz.
where the constant is obtained by the calibration (may be frequency-dependent, taking into account the transfer function of the antenna-receiver communication line).
In the near-field zone, the antenna k-space pattern depends on the vertical distance z, whereas in the far-field zone, where near-field (evanescent) components with Re( ) kk ⊥  and the antenna radiation pattern does not depend on the distance z. The far-field components also determine the angular pattern of the antenna in the far zone that does not depend on the distance. In the probing field of electrically small (with subwavelength sizes) antennas, evanescent components dominate in the near-field zone, providing the sensitivity of the near-field diagnostics to subwavelength-size inhomogeneities in this zone. The resolution decreases with distance from the transceiver to the Rayleigh far-field resolution limit.
In Figure A1, the distribution of the current y-component ( , , ) y j x y  on the electrically small bow-tie antenna and its k-space spectrum ( , , )  In Figure A1a, one can see a sharp localization of the current distribution with the very broad k-spectrum (approximately up to 25 cm −1 ) that is demonstrated in Figure A1b. The corresponding small-scale current components create a small-scale structure of the damp field components in the near zone and provide the required subwavelength resolution. To obtain the k-space spectrum of the received signal (A1), it is necessary to obtain the k-space of the scattered field 1 ( , , , ) In Figure A1a, one can see a sharp localization of the current distribution with the very broad k-spectrum (approximately up to 25 cm −1 ) that is demonstrated in Figure A1b. The corresponding small-scale current components create a small-scale structure of the damp field components in the near zone and provide the required subwavelength resolution. To obtain the k-space spectrum of the received signal (A1), it is necessary to obtain the k-space of the scattered field E 1 (k x , k y , z r , ω) emitted by the antenna attached to the point x tr = x r + δ x , y tr , y r + δ y , z tr = z r that produces the incident field: E 0 (r, ω) = j(x − x tr , y − y tr , ω) ↔ G(x − x , y − y , z, z tr , ω)dx dy (A12) with the k-space spectrum E 0i (k x , k y , z = 0, ω) = 1 4π 2 e −ik x x−ik y y dx 1 dy 1 j(x − x tr , y − y tr , ω) ↔ G(x − x , y − y , z tr , z = 0, ω)dxdy = − 2π ωε e −ik x x tr −ik y y tr +i k 2 −k 2 x −k 2 y z tr −ik x (x r +δr)−ik y (y r +δr)+i k 2 −k 2 x −k 2 y z tr j j (k x , k y , ω)[g ji (k x , k y , ω) + g ⊥ ji (k x , k y , ω)]. (A13) In the considered case of one-dimensional inhomogeneity under the layered structure of the lung cover tissues, this field can be determined by the reflection coefficients that can be calculated for an arbitrary 1D structure of permittivity: ωε j j (k x , k y , ω)[R (k x , k y , ω)g ji (k x , k y , ω) + R ⊥ (k x , k y , ω)g ⊥ ji (k x , k y , ω), ω)] ×e −ik x x tr −ik y y tr exp[i k 2 − k 2 x − k 2 y (z tr + z r )], where j x , j y (k x , k y , ω) are k-space components of the transverse spatial spectrum of the distributions of the current components on the antenna calculated by the built-in software, R ,⊥ (k x , k y , ω) are k-space components of the reflection coefficients from a probed medium at two polarizations TH( ) and TE(⊥) that can be calculated for an arbitrary onedimensional inhomogeneity ε(z), and δx, δy, δz are components of the offset δr between antennas of the transceiver system. The harmonic components of the scattered field (A14) are determined by the components of the plane wave decomposition of the incident field E 0i (k x , k y , z, ω) in (A13) at z = 0: E 0i (k x , k y , 0, ω) = 4π 2 j j (k x , k y , ω)G ji (k x , k y , z, z, ω)e −ik x (x r +δx)−ik y (y r +δy) = 4π 2 j j (k x , k y , ω)g ji (k x , k y , ω) exp[ik z (z r + δz)] = 4π 2 j j (k x , k y , ω)[g ji (k x , k y , ω) + g ⊥ ji (k x , k y , ω)]e −ik x (x r +δx)−ik y (y r +δy)+i k 2 −k 2 x −k 2 y (z r +δz) . (A15) Components of the k-space pattern in (A3) can be written in a compact form as: Finally, we express the frequency spectrum of the received scattered signal (A3) as the inverse Fourier transform of (A3) with the dot product of (A14) and (A16) as: s ω (x r , y r , ω) = const V(ω) exp −ik x δx − ik y δy + i k 2 − κ 2 x − κ 2 y (z tr + z r ) {R [ε(z, ω)](k x , k y ) j j (k x , k y , ω)g ji 2 + R ⊥ [ε(z, ω)](κ x , κ y ) j j (k x , k y , ω)g ⊥ ji 2 }dk x dk y .
where the signal s(x r , y r , ω) = s(ω) appears, of course, independent of transversal coordinates x r , y r . The value of const in (A15) (that can be frequency-dependent) is determined by fitting the calculated signal to the measurement data at the full exhalation phase, when the tissue structure is assumed to be known. It is also possible to represent the functional (A17) as s ω [ε(z, ω)](ω) = const(ω) V(ω) e −ik x δx−ik y δy+ k 2 −k 2 x −k 2 y (z tr +z r ) [j 2 x f x (k x , k y , ω) + j 2 y f y (k x , k y , ω)]dk x dk y , f x (k x , k y , ω) = R (k x , k y , ω)[( k 2 that gives us the desired spectrum of the measured pulse shape in (1) from the expression (A1) as s(ω) = 1 2 [s ω (ω) + s ω (−ω) * ]. The necessary relationships between the permittivity and the determined parameters of the air and blood content in the lungs can be based on (6), obtained from the Rayleigh formula (5) in the mixing theory [28].