A Forward-Backward Iterative Procedure for Improving the Resolution of Resonant Microwave Sensors

This paper sets out a method for improving the resolution of resonant microwave sensors. Usually, the frequency response of these devices is associated with a low quality factor, and consequently with a low resolution in terms of tracking capacity of the resonance frequency shift. Furthermore, since only a finite number of samples can be acquired during the measurement process, the “true” resonance frequency may not be included in the set of acquired data. In order to have an accurate estimate of the resonance frequency, high performance systems with very fine frequency sampling are thus required. To limit these drawbacks, an iterative algorithm is presented which aims to refine the response of resonant microwave sensors by means of a suitable post-processing. The algorithm evaluation is first carried out on synthetic data, and then applied on experimental data referring to a practical scenario, which is inherent to return loss measurements performed by a microwave patch antenna immersed in a water-glucose solution with different concentrations.


Introduction
Nowadays, microwave sensors are found in a variety of contexts, as their field of application is extremely broad. They have very interesting peculiarities, and among them microwave (resonant) sensors based on resonant devices are particularly attractive, as their frequency response is highly sensitive to the characteristics of the investigated medium in which they are inserted [1]. For example, referring to biomedical applications, in [2], the design of a resonant microwave sensor for applications related to blood glucose monitoring is proposed, showing that, in order to obtain a sensor suitable to provide an accurate estimation of the glucose concentration level, it is necessary to take into account the variations of both the real and the imaginary part of the complex permittivity of the medium under investigation as a function of frequency, yet discussed in [3]. In [4], still considering the problem of blood glucose monitoring, a portable prototype of a planar microwave sensor is shown, consisting of four distinct hexagonal-shaped complementary split ring resonators. Instead, in [5], microwave and millimeter wave dielectric spectroscopy performed at the cellular and molecular level is discussed, in the context of early cancer diagnostics.
In [6], a planar two-dimensional inductor-capacitor circuit is instead used as a resonant sensor for monitoring of temperature, humidity, and pressure. In [7], a chain of complementary split-ring resonators is used both for the measurement of dielectric constants as well as for faults detection. Resonant sensors are also used for angular displacement and velocity measurements in satellite applications [8], for detection of acetone vapour [9], and for the characterisation of aqueous solutions [10].
An important class of microwave resonant sensors are the frequency-variation sensors, for which the behavior of the frequency response depends on the variations of the measurand [1]. Typically, such sensors can have a low quality factor, which results in a relatively much wider than ideal frequency response pattern. This in turn results in a low resolution, and therefore in a low tracking ability of the resonance frequency shift [11]. This problem can be addressed from the circuit point of view, by means of some active microwave device which allows an increase in the quality factor of sensors [11,12]. However, another important aspect to be considered is related to the fact that the frequency response samples, acquired by the measurement instrumentation, represent a discrete set with finite extension, so there is a very high probability that the sample associated with the resonance frequency will not be captured correctly.
In order to limit the above drawbacks, an iterative algorithm is proposed in this work which leads to post-process the acquired samples of the frequency response, in order to significantly increase the quality factor and therefore achieve a more accurate estimate of the resonance frequency. This algorithm has high performance even starting from a relatively very low number of acquired samples of the sensor response within a frequency band with a high extension. In this way, even unsophisticated measuring equipment can be exploited.
An appropriate performance analysis of this algorithm is first carried out using synthetic data; then, the proposed method is tested on measured data collected by a microstrip patch antenna immersed in a water-glucose solution, which represents a simplified case in the context of a blood glucose monitoring problem.

Iterative Method to Improve the Resolution of Resonant Sensors
The algorithm presented in this section is based on an iterative procedure inspired by the Papoulis-Gerchberg algorithm [13,14]. However, while the objective of the Papoulis-Gerchberg algorithm is to estimate the Fourier transform of a signal starting from its truncated version, the aim of the proposed method is to progressively narrow the sensor response in such a way to reach, from the starting signal, a sinc function with a very narrow main lobe, whose point of maxima coincides with the resonance frequency.
Let us consider a signal as that depicted in Figure 1a, which is representative of the sensor return loss measured between the frequencies f min and f MAX , namely [15]: (1) The function Γ( f ) into Equation (1) represents the reflection coefficient at the microwave sensor input. An interesting aspect to be remarked is that here we are going to work with values expressed in dB for the magnitude of the reflection coefficients (regardless of the presence of the minus sign in the definition), and thence the computation of the Fourier transform of the return loss is as equivalent to a "sort" of complex cepstrum [16] of the inverse Fourier transformation of |Γ( f )|.
A brief consideration is relative to the presence of noise. In this paper, we assume the signal is smooth enough so that the noise can be considered almost full negligible. However, in practice this cannot always be ensured, and thence the size of the signal-to-noise ratio (SNR) must be taken into account in order to be able to adequately process the return loss. However, one solution is to carry out more acquisitions of the return loss, and then to perform a sample average of these acquisitions, in order to significantly improve the SNR, without distorting the "true" signal. Now, let us move on to introduce the various steps that make up the algorithm we intend to propose. Assume that RL( f min ) < RL( f MAX ) (similar arguments hold true when RL( f min ) > RL( f MAX )) and that the curve has a single maximum in [ f min , f MAX ]. The algorithm consists of the following steps: 1.
Shift RL( f ) down by a quantity equal to RL( f ) = RL( f MAX ), so that now the "new" return loss, RL( f ), is as shown in Figure 1d, i.e., it is bereft of end points.

4.
Set up a threshold 0 < α < 1 for RL( f )/ max{RL( f )} in order to identify the two values f min and f MAX such that ( Figure 1d): thus obtaining the narrower frequency observation interval [ f min , f MAX ].

5.
Compute the inverse Fourier transform of RL( f ) (for ν ∈ [ν min , ν MAX ]): in order to obtain the complex function of Figure 2a. As seen, rl(ν) exhibits an almost linear unwrapped phase φ(ν) = ∠rl(ν) (bottom panel of Figure 2a) in correspondence of the main lobe (top panel of Figure 2a), which actually encodes information about the peak of the resonance frequency. 6.

7.
Determine the two valuesν min andν MAX for |rl(ν)|/ max{|rl(ν)|} in such a way ( Figure 2a): in order to pick out the interval [ν min ,ν MAX ] inside of which there is a greater reliability of the phase information. 8.
Determine the magnitude of the Fourier transform of rl(ν) (for f ∈ [ f min , f MAX ]): which returns a new return loss (blue curve of Figure 3) that is narrower than the original one (of Figure 1a). By calculating the new return loss only in the interval [ f min , f MAX ], we are doing a sort of magnification in the region in which the resonance frequency is assumed to be. 10. Compute the inverse Fourier transform of the squared current return loss (i.e., of the square of the blue curve of Figure 3): so that the magnitude and the unwrapped phase of rl(ν) are now like the blue curves sketched in the top and bottom panels of Figure 4a, respectively. It should be noted that here we are implicitly assuming that RL( f ) is nonzero only in [ f min , f MAX ]. By comparing the previous and the current version of rl(ν), now the magnitude main lobe is wider (this is due to squaring RL( f ) in (7)) and the region where the instantaneous phase is linear wider. In fact, in correspondence with the threshold β, two new values of ν are identified, i.e., |rl that ν min <ν min and ν MAX >ν MAX . By way of example, if from this step one returned to step 7 and the same step 7 were performed together with steps 8 and 9, the new return loss would be given by the red curve of Figure 5. 11. Repeat points 7 to 10 until the function rl(ν) is such that |rl(ν)|/ max{|rl(ν)|} ≥ β ∀ ν ∈ [ν min , ν MAX ] and, in the end, compute the final return loss by Equation (6), whose point of maxima gives the estimated resonance frequency. For greater clarity, it is worth remembering that ν min and ν MAX are the extremes of the whole observation interval of the function rl(ν).
For a better understanding of the algorithm process, a flowchart description is reported in Figure 6.
It is worthwhile to remark that the different Fourier transformations involved in the previous steps are actually achieved by the discrete time Fourier transform (DTFT) procedure, which allows to control in a very flexible way both the number and the positioning (e.g., periodic or aperiodic) of the observation points of the various transformations in specific intervals. This entails that even starting with N f number of points, which correspond to the frequencies actually collected by the sensor, the number of points in the ν-domain, N ν , and subsequently in the f -domain, N f , can be chosen on the basis of convenience. This implies a very advantageous practical aspect, since it may happen that the frequency step allowed by the acquisition system is not sufficiently fine for the actual resonance peak to be captured well. The proposed method, on the other hand, allows one to reduce the difference between the frequencies of adjacent samples, and therefore a better estimate of the resonance frequency is expected.    Now, before moving on to the assessment of the algorithm under consideration, it is worth investigating further some practical aspects in relation to how α is determined and the implementation aspects of the Fourier transforms related to the functions rl(ν) and RL( f ). In particular, Equation (3) would require that, strictly speaking, α be determined with respect to the maximum value of the return loss, which however is actually unknown in real scenarios. However, from a practical point of view, in reality α is calculated with respect to the maximum value of the acquired return loss, that is, max{RL( f )} is the highest value among those acquired during the measurement process. As for computation of rl(ν) and RL( f ), Equation (4) is implemented by sampling in N ν points of the interval [ν min , ν MAX ] the inverse DTFT of the acquired N f points in [ f min , f MAX ] of the initial return loss; Equation (6) is implemented by sampling in N f points of the interval [ f min , f MAX ] the DTFT of the rl(ν) samples (of the total N ν samples) that fall in [ ν min , ν MAX ]. Obviously, it is possible to change parameters N ν , N f , the positioning of the samples of rl(ν) and RL( f ) and also the thresholds α and β even during the execution of the algorithm, in order to perform a kind of "dynamic" performance control. Furthermore, we choose to keep the extension of the temporal observation interval constant, thus fixing ab initio in the way that is explained below the values of ν min and ν MAX , even if, in general, the value of these parameters could also be changed while running the algorithm.
Another important consideration concerns the setting of the acquisition band, i.e., the choice of f min and f MAX . Due to the algorithm structure, in order to obtain an accurate estimate of the resonance frequency, it is required that f min and f MAX must be such that the "acquired" return loss has an increasing trend up to the point of maxima (or to the points of maxima in case of multiple acquired maximum values), inside [ f min , f MAX ], and then after the same point is decreasing. Thence, during the acquisition stage, a simple procedure for analysing the trend of the initial return loss could also be set up in order to verify whether the above condition is satisfied.
Finally, the frequency and time steps should be also considered. In fact, if the various functions are sampled uniformly (which we assume below), then there is an initial frequency step given by ∆ f = ( f MAX − f min )/(N f − 1) (assuming that any sampling errors on the part of the measurement system are negligible), a subsequent frequency step given by

and a time step given by
In the following, the proposed algorithm is applied both to synthetic and experimental signals by resorting to the environment for numerical computation MATLAB [17].

Iterative Algorithm Application to Synthetic Data
In this section, we turn to check the presented algorithm by employing a simple signal that mimics the response of a resonant sensor. The main aim is to better clarify the key features of the method. To this end, we initially consider the following values for parameters: From the values of the these parameters, it follows that the initial frequency step is equal i.e., the combination of two "half" Gaussian functions, both with maximum values equal to 15 and mean (point of maxima) equal to −200.5 MHz. However, the first one has a standard deviation of 70 MHz while the second one of 100 MHz. A curve of this type can be close to the actual frequency response of a resonant sensor. It is worth pointing out that, for how the values of the above parameters are set, the frequency −200.5 MHz is not included in the set of initial return loss points. Therefore, here we propose both to tighten the initial return loss and to more accurately identify the frequency shift. Figure 7 shows the behavior of the percentage error, PE = |( f shi f t − f shi f t )/ f shi f t | × 100%, as a function of the two thresholds α and β, with f shi f t and f shi f t being the estimated and actual frequency shift (i.e., f shi f t = −200.50 MHz), respectively. As can be seen, these thresholds have a relatively strong impact on PE values, being able to notice that they are non-linearly related to PE. In particular, it can be observed that for the considered cases, the minimum percentage error value is for (α, β) = (0.85, 0.05). Figure 8 shows the results of applying the algorithm specifically in the case in which (α, β) = (0.85, 0.05). Figure 8a shows both the (normalised) initial return loss (blue curve), given by Equation (8), and the normalised one (red curve) obtained by applying steps 1 to 3 of the algorithm, which has no end points. The magnitudes of the (normalised) inverse DTFTs of these signals are reported in Figure 8b, which shows that they differ significantly from zero over a much smaller interval than the entire observation one [ν min , ν MAX ] = [−0.5, 0.5] µs. These are actually the parts that must be retained. Indeed, it is clear from Figure 8c that in this smaller interval the phase of rl(ν) shows a different slope as compared to the regions where the magnitude of rl(ν)/ max{|rl(ν)|} is negligible. Therefore, this justifies why we focus on the main lobe of rl(ν)/ max{|rl(ν)|}. Furthermore, it should be noted that if the end points are removed, the function rl(ν) has a different behaviour than that relating to the case in which this removal is not performed, as can be clearly recognised by observing the blue and red curves of Figure 8b,c. In particular, by looking at the inset of Figure 8b, it can be seen that removing the end points allows for a larger interval [ ν min , ν MAX ] (which is positive for narrowing the RL). In the end, Figure 8d shows the comparison between the initial return loss and the final return losses obtained at the end of the full application of the algorithm, although it must be taken into account that the blue dot sinc-like curve is relative to the case in which steps 1 to 3 are not applied to the initial return loss (blue curve of Figure 8a). In fact, although the final curves turns out to be much tighter than the initial ones; to the blue dot curve corresponds PE ≈ 9.23% while to the red one corresponds PE ≈ 0.09%. Therefore, this shows that the removal of the end points turns out to be beneficial from a performance point of view. Now, let us focus on the resolution capability of the algorithm under consideration. Here, we refer to the Rayleigh criterion for the diffraction limit to resolution, therefore going to verify if it is possible to distinguish two frequency shifts if their distance is at least equal to half the width of the main lobe of the "correct" final return loss. So, referring to the red curve of Figure 8d, we have that the width (beamwidth) of its main lobe is approximately equal to BW = 2.01 MHz. Consequently, for what is said above, let us now consider the case in which the initial return loss has an actual frequency shift equal to f shi f t + BW/2 ≈ −199.50 MHz (remember that f shi f t = −200.50 MHz), that is, it is shifted to the right with respect to the blue curve of Figure 8a of BW/2 (refer to Figure 9a). As can be seen from Figure 9b, in which now the blue curve coincides with the red curve of Figure 8d whilst the red-dashed curve with the final return loss of the present case, the performance in terms of resolution is satisfactory, as the point of maxima of the above red dashed curve has a high correspondence with −199.50 MHz. In fact, even in this case the percentage error is PE ≈ 0.09%. Obviously, also in this case, it has been placed (α, β) = (0.85, 0.05), Finally, it is worth asking about the performance behaviour of the algorithm in question when changing the number of initial points, N f . To do this, let us consider

Iterative Algorithm Application to Experimental Data
In this section, the algorithm previously discussed is applied to return loss data captured by the resonant microwave sensor described in [2], in the context of blood glucose monitoring. This sensor consists of a standard inset-fed microstrip patch antenna working in the Industrial, Scientific, Medical (ISM) band around the frequency f 0 = 2.4 GHz. The substrate has a relative dielectric constant equal to ε r = 10.
The accurate design of a typical microwave sensor able to capture dielectric variations of human body parameters should be performed by starting from a proper definition of a multilayer model for the human body portion assumed as the radiation medium. Specifically, an accurate body stratification should be assumed, by considering the proper behaviour of complex permittivity versus frequency (i.e., the dispersive behaviour) for each layer. The effectiveness of the dielectric models assumed for the various layers and the performance of the methods used to fit those models play a crucial role in the accuracy of the design process of the microwave sensor [2,3,18,19].
In this paper, without affecting the conceptual validity of the algorithm under investigation, we consider a simplified scenario that consists in estimating the resonance frequency of the sensor immersed in a water-glucose solution as the glucose concentration, GC, inside it varies. Consequently, consistent with the foregoing, for the design of the above sensor the properties of a water-glucose solution were first investigated, thus characterising the behaviour of its complex permittivity (the dielectric constant as well as the loss tangent) as a function of the concentration of glucose levels. Subsequently, data obtained from the above analysis were used to proceed to the actual design of the sensor through an optimisation performed using the Ansys software, in such a way as to have optimum matching conditions at the prescribed frequency f 0 = 2.4 GHz. In Figure 11, a photograph of the fabricated microwave sensor is reported, which is immersed into water-glucose solutions to perform reflection measurements at the Microwave Laboratory of University of Calabria. Figure 12 shows the experimental (normalised) return loss measured for the waterglucose solution with different values of the glucose concentration where it can be clearly appreciated how the resonance frequency shifts with GC.
To run the sharpening algorithm, the data measured within the frequency band [2,3] GHz are first shifted at the band [−0.5, 0.5] GHz. Accordingly, as mentioned above, this implies that the resonance frequency estimated by the algorithm returns the frequency shift, f shi f t , with respect to the middle frequency of the measurement band, i.e., with respect to 2.5 GHz. For each measured curve, the algorithm is run for different values of threshold N f and always for (α, β) = (0.85, 0.05) and N ν = N f = 4096. The results are summarised in Table 1, whilst Figure 13 shows the comparisons between the curves of Figure 12 and the curves representing the final return loss obtained after applying the algorithm presented above, for N f = 1001, being able to observe well that the curves relating to the final return loss are much narrower than the initial ones and the main lobes are clearly distinguish-able from each other. As previously ascertained, since the frequency shift estimation is slightly dependent on N f , assuming the above is also valid in this case, it is appropriate to assume as the most reliable frequency shift, for each value of the glucose concentration (that is, for each case), the one corresponding to N f = 1001 and (α, β) = (0.85, 0.05). However, referring to Figure 14, which shows the behaviour of the sample mean and sample standard deviation curves as functions of glucose concentration, derived from the values in Table 1 (e.g., the three values −64.62 MHz, −64.24 MHz, and −65.02 MHz in Table 1, for GC = 100 mg/dL, have µ shi f t ≈ −64.63 MHz and σ shi f t ≈ 0.39 MHz, as can be appreciated qualitatively in Figure 14), in this case it is possible to note that the values of the sample standard deviation are all "sufficiently" much smaller than the respective values of the sample mean, thence indicating a relatively small dispersion and therefore a high precision, with respect to the N f variation, for each glucose concentration.

Conclusions
In this work, an iterative algorithm was presented by which to significantly increase the resolution of resonant microwave sensors, in such a way as to be able to track the shift of the resonance frequency with great accuracy. A test procedure was first carried out on synthetic signals, generated ad hoc, in order to highlight its peculiar characteristics. Subsequently, the algorithm was successfully applied to a simplified scenario related to the problem of blood glucose monitoring. The above algorithm revealed a high flexibility from the implementation point of view, and the achieved results proved to be satisfactory. More complex test scenarios will be considered in future works.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.