1. Introduction

Microwave scattering measurements in free space are a common tool for the analysis of the electromagnetic signature of macroscopic targets. In some cases however, the method is used to simulate experimentally light scattering experiments. The underlying theoretical principle on which such experiments are based is referred to as the Scale Invariance Rule (SIR). It states that for an arbitrary fixed aggregate embedded in an infinite, homogeneous, linear, isotropic, and non absorbing medium, a dimensionless amplitude scattering matrix can be defined [1]. In other words, it means that if one scales both the object size (effective radius r) and incident wavelength λ by the same factor, keeping fixed their size parameter x = 2πr/λ and complex refractive index, the scattering behaviour of the analog object is identical to the original one. This so-called Microwave Analogy principle has been applied since the fifties in order to analyze optical scattering properties of non-spherical particles and systems [2, 3]. Since in the optical region, systems of interest are generally very small (nanometers to micrometers scale), it is very difficult to control their refractive index, exact shape and orientation and hence validate approximate models. In the microwave region on the other hand, objects can be made macroscopic and consequently are easier to manufacture in a controlled fashion. The first experiments were realized by Greenberg and later by Giese [2, 3]. Both used a single wavelength (3.18 cm and 8 mm respectively) and measurements were performed in a single plane. Phase measurements were out of their reach however, except in the forward direction. This was a serious drawback since this information allows one to subtract in a straightforward manner the incident field and obtain the purely scattered signal. This technical difficulty was solved with the development of the University of Florida microwave scattering facility [4]. The equipment permitted broadband measurements (2.7–4 mm) and gave access to the phase. Nevertheless, this phase information was never exploited to its full potential and measurements were almost always restricted to co-polarized or total intensities [4, 5, 6, 7, 8]. Furthermore, the types of targets studied were predominantly either complex (but not controlled [4, 5, 7, 8]), or controlled (but regular [6, 9]). Hence it became desirable to extend those measurements to complex targets with fully known geometrical properties, and at the same time analyze the phase of the scattered fields for all (co- and cross-) polarizations. Thanks to the anechoic chamber at the Institut Fresnel, such experiments are now possible. This set-up, as will be described in more details in section 3, has an additional advantage over the previous ones, that it allows to move incident and scattered directions independently.

This paper is organized as follows. We start by describing the complex target and how it was built, then the technical details of the experimental set-up are explained. Before starting with the discussion of the results, we briefly describe the numerical procedures that were used to obtain the data presented in this work. The objective of this, is to make it as easy as possible for other people to use the data and compare with their computational methods. Once all those aspects are settled, we present the most relevant results and finally end with some conclusions and perspectives.

2. Model aggregate

2.1. Choice of the model

The aim of this study is to compare the microwave signal scattered by a complex aggregate with fully established electromagnetic and geometrical properties with the solutions generated by different computational methods. Therefore, we need a model of such a complex aggregate. Our choice was guided by the analogy with soot particles. These particles are usually modeled as an agglomeration of mono- or polydisperse spherical particles. The geometry of the agglomerate can be considered to a good degree of approximation as fractal. Hence, a population of such aggregates can be described by the well known power law:

where N is the number of monomers in the aggregate, k ₀ is the prefactor, R_g the gyration radius, a the monomer radius and D_f the fractal dimension.

To determine the optimal size of the constituent monomers of our model aggregate, we had to find the multiplication factor which realizes the scale translation from the optical to the microwave region of the electromagnetic spectrum as stated by the SIR. The choice of this factor is limited by different constraints set by the size of real soot aggregates, the frequency range of the microwave facility and practical considerations related to the construction of the aggregate. For soot particles, the diameters of the monomers ranges between 10 and 60 nm [10]. The detection range in the microwave region stretches from 2 to 20 GHz. A multiplication factor of 1.4 × 10⁵ gives a translated particle radius of a = 2.50 mm for a soot particle with radius = 18 nm. The same factor makes the correspondence between the optical [0.4 – 0.7] μm and the microwave [3.06 – 5.35] GHz (corresponding wavelength range is [56.04 – 97.97] mm) regions, the latter falling into the detection range of our facility. For the other numerical values of the fractal parameters we chose: N = 74, k ₀ = 2,D_f= 1.7 [11, 12].

In order for the SIR to be completely satisfied, the refractive index should also remain constant under the scale translation operation. In other words, the refractive index should be identical in both the optical and in the microwave region. Soot has a high absorbtion in the visible spectrum. Finding a material possessing the same order of absorption in microwaves is not an easy task. Our material choice was polyacetal, which although exhibiting almost no absorption, does have a real part of the refractive index quite similar to that of soot. For the experiments, we concentrated on the interval [15 – 20] GHz instead of the complete [2 – 20] GHz one. The reason for doing this is that because the scattering pattern under 10 GHz does not exhibit many detailed features (the size parameter x of the single spheres varies from 0.16 to 0.52 in the frequency interval [3 – 10] GHz), thus making it harder to carry out meaningful comparisons with the numerical methods. Over the interval [15 – 20] GHz, the refractive index can be considered constant and we took n = √ε = 1.688. This value is slightly different from the one we used in [12]. For the present study we reduced the uncertainty on the refractive index to approximately 1% (instead of 5% previously) by refining the EpsiMu facility result with diffraction measurements of a single sphere made of the same material.

2.2. Target construction method

Once the properties of the model were chosen, we could begin the design phase. We first created a computer generated model (see Fig. 1(a)). The algorithm we adopted generates the aggregate by starting with a single sphere to which a second one is stuck randomly to the first one. For the next steps, the algorithm adds a new sphere randomly and checks whether the new aggregate still satisfies Eq. (1) until the predetermined number of spheres is attained.

Once the list of the particles coordinates was generated, we could start the actual fabrication. Thereto we used a 3D micro milling machine (Isel, Flatcom 20) adapted to manipulate spheres of 5.0 mm in diameter. They were brought into position and then stuck together in a two-step procedure. The first step consisted in slightly heating the polyacetal spheres surface at the contact point. This was necessary to give extra grip to the glue and hence additional strength to the structure. The second step is the application of the glue and its hardening. In Fig. 1, a computer generated representation and two viewpoints of the finished aggregate are shown.

 

Fig. 1. Computer generated representation of the aggregate using the output of the algorithm (a). Front (b) and side (c) views of the constructed aggregate. In pictures (b) and (c), one can appreciate the polystyrene supports, transparent to microwaves, on which the aggregate rests. The plexiglass spheres are removed after the aggregate has been aligned.

Download Full Size | PDF

3. Experimental set-up: microwave scattering facility

The experimental set-up developed by members of the Institut Fresnel in the anechoic chamber of the ‘Centre Commun de Ressources Microondes’ (already described in [13] and also used for inverse problem applications [14]) had to be optimized once again to allow the measurement of the scattering phenomena obtained with the proposed scale translation needed to move from optical waves to microwaves. In fact, since the considered aggregate was the smallest target ever measured with this equipment, the scattering measurements presented here required shortening the cable lengths, optimization of the signal levels and measurement parameters to exploit the full dynamic range of this setup. Another advantage of the microwave domain is that both the amplitude and phase of the electric fields can be measured using a network analyzer (HP8510 here).

This equipment allows the 3D characterization of antennas with various polarizations and to measure the scattering from targets in almost all directions from a 4 meter diameter sphere, with working frequencies between 2 and 20 GHz. As can be seen from Fig. 2 (a) and (b), the emitter can be moved along the vertical arch. This introduces an extra degree of freedom, which allows making measurements outside the azimuthal plane (xy plane). For highly non-spherical and static targets, this arch permits the testing of computational methods to a higher degree of precision since the object doesn’t need not be moved which, depending on the set-up, is less precise than displacing the emitter and receiver.

Our technique is based on the measurements of the total field (the field in presence of the target) and secondly the measurement of the incident field (field without a target) which is subtracted from the total field to obtain the scattered field. Small scattering fields, like those studied here, require high sensitivity in the measurements which consequently increases the impact of noise on the complex subtraction of the total and incident fields (a special drift compensation procedure was developed [15] to compensate for drifts of 0.1% in amplitude and 0.01 radians in phase). A trade off between increasing the scattered field level and saturating the receiver in the forward diffraction zone had to be found. To increase the scattered field level, a high gain horn antenna (ARA MWH-1826B, gain 23.4 dBi) was chosen as the source in order to put a maximum of energy in the object region. This however increases the signal level when the receiving antenna is in front of the source and can saturate the receiver which results in an imperfect subtraction in the forward scattering direction [4]. Thus, for this study, the source signal level has been optimized to be just below the saturation limit when the measured signal is the highest, in order to exploit the total dynamic range of our system using a 1024 averaging factor (dynamic range of about 89 dB with the external mixers 85320A [16]) and allowing precise field measurements in the forward cone. The different polarization cases are obtained through the rotation of the linearly polarized antennas and some of the differences visible on the cross-polarization curves in [12] were eliminated thanks to a better alignment of the source-target couple.

 

Fig. 2. On the left, a picture of the anechoic chamber and on the right the name conventions for the polarization components of the electric field in the experimental set-up. Here the red dashed line corresponds with the vertical arch and the green dashed one with the receiver positions in our configuration.

Download Full Size | PDF

The polarization name convention is depicted in Fig. 2(b). It is imposed by the set-up architecture. The polarizations of the fields are defined as the projections:

with k = i, s which stand for incident and scattered respectively. E ^k is the incident/scattered field, and (e _θ,e _ϕ) are the basis vectors of a spherical coordinate system, where the somewhat unusual convention for the polar and azimuthal angle notations (θ, ϕ) is the one defined in [17].

4. Numerical analysis

The numerical analysis for this paper is developed on different levels. For one, we will use two basic types of electromagnetic solvers for the modeling of the scattering phenomena. Since these are well known methods, we will outline only briefly their respective mathematical frameworks in order to make the differences clear. However, we do try to specify the technical details necessary for reproducing our calculations. The main reason to do this is because the data presented here will be made public. Hence, it is of high importance to know the datails for others to be able to reproduce our simulations and data analysis.

4.1. Electromagnetic scattering codes: a short description

We will start by describing briefly the electromagnetic scattering solvers. We employed four different codes. Two of them are based on the Volume Integral Formulation. The remaining two are of the multiple Mie sphere solution type, which we shall refer to as T-matrix method.

For the volume discretization codes, we used ddscat7.0 developed by Draine [18, 19] and a second one by Eyraud [20]. The first one (ddscat7.0) is an implementation of the Discrete Dipole Approximation (DDA), while the second one can be classified under the Method of Moments (MoM). The T-matrix codes are the ones developed by Mackowski and Mishchenko [21, 22] and Stout and Auger [23].

4.1.1. Volume Integral methods

The scattered field is obtained with a 3D volume integral formulation which allows the computation of the scattering by a 3D inhomogeneous target. This object is enclosed within a zone of volume Ω considered to be non homogeneous. Indeed, for each point of this zone, the permittivity can be written as:

where ε ₀ is the vacuum permittivity, ε_r(r′) is the real part of the permittivity and ε_i(r′) is the imaginary part at point r′. The scattered field E^s on the receiver positions r ∈ Γ is obtained with the observation equation [24]:

where G_0r is the free-space dyadic Green function between r′ ∈ Ω and r ∈ Γ, χ is the contrast (χ(r′) = k ²(r′) - k ² ₀) with k the wavenumber in the zone Ω, k ₀ the wavenumber in vacuum and E is the field inside the zone Ω.

The field inside the target verifies the equation (r′ ∈ Ω):

where E ⁰ is the incident field and G₀₀ is the free-space dyadic Green function.

In both methods (ddscat7.0 and MoM), Eq. (6) is numerically computed using a FFT-method coupled with a BiConjugate Gradient Stabilized solver [25]. For a detailed discussion about the difference between DDA and MoM, we refer to [26]. As the computation complexity is then O(NlogN) and the memory requirement is O(N), these methods allows the calculation of the scattered field for large objects. In our case, the aggregate is a relatively large target at the frequency of 20 GHz since the enclosing sphere has a volume of ≈ 150λ ³. Cubic cells with lattice spacing parameter d are used in both methods for the computations presented here.

4.1.2. T-matrix method

In the T-matrix approach, fields are developed on basis sets of vector spherical wave functions (VSWFs) of the type encountered in Mie theory [23, 27, 28]). The T-matrix of an individual particle in a N-particle system, (labeled by j = 1,..N), completely characterizes its electromagnetic response by relating a column ‘matrix’, f ^(j), composed of the coefficients of the field scattered by a particle to a column ‘matrix’ , e ^(j), containing the field coefficients of the field incident on the particle (called the excitation field). In this convenient matrix notation [27], the T-matrix for the particle has the property, f ^(j) = T ^(j) e ^(j). For a sphere, this matrix is diagonal and its matrix elements can be expressed analytically in terms of the ‘Mie’ coefficients [29, 30].

In a multiple scattering context [31], the ‘excitation’ field associated with a particle is the superposition of the incident field with the field scattered by all the other particles in the system (excluding the scattered field of the particle itself). Upon invoking the analytic translation-addition theorem formulas [28], the multiple equations take the form of a set of coupled matrix equations. A frequently cited advantage of the above matrix formulation is that it is formally exact (for spheres) up to a necessary truncation of the basis set of the VSWFs.

Numerical implementation of the multiple scattering T-matrix algorithms encountered considerable difficulties up until the late 90s, but these problems have now been largely resolved by invoking a variety of techniques like iteration employing conjugate gradient methods [32]. For our simulations, we adopted methods that directly generate multiple scattering solutions, either through stabilized recursive algorithms or direct matrix inversions using analytic matrix balancing [23].

An advantage of the multiple scattering methods for future simulations is that they allow one to readily extract a variety of near field information including local absorption and optical forces [33] as well as offering the possibility of carrying out analytically procedures such as orientation averaging that are rather difficult to perform numerically.

As we said previously, we used two implementations of the same formulation. Differences arise only from slight variations in numeric implementations or choice of multipole cutoff. The Mackowski-Mishchenko codes use iteration, while the method of Stout-Auger relies on analytic matrix inversions in a recursive scheme. It has been previously demonstrated that the available T-matrix codes yield the same predictions when exploited to the same level of accuracy [34].

4.2. Calibration of experimental and simulated data

To compare the different simulated scattered fields between each other as well as with experimental data, we have to choose the same reference for all. The incident field is then taken with an amplitude equal to unity and a phase null at the centre of the target. To calibrate the experimental fields, one complex-coefficient C is determined using a reference target in the ϕ co-polarization case [14, 35]:

where E^s_c,i (resp. E^s_m,i) is the scattered field calculated (resp. measured), N_r the number of receiver positions and x̄ is the complex conjugated of x. A metallic sphere was chosen as a reference target as it can be accurately approximated by a perfect conductor in the microwave domain and as its scattered field can also be readily computed using the Mie series. All the experimental fields are then calibrated using this multiplicative coefficient C.

With the method described above, all experimental scattered fields are normalized with respect to the criterion: amplitude of the incident field equal to unity and phase null at the centre of the target. All the measurements presented in the following are calibrated and thus can directly be compared quantitatively to any simulation considering an incident field of amplitude one and phase null at the centre of the referential.

4.3. Convention for the incident and scattered polarization vectors

As explained in section 3, the emitter can be moved along the vertical arch shown in Fig. 2(a), while the receiver always remains in the so-called azimuthal plane (which corresponds to the xy-plane in the scheme of Fig. 2(b)). Since the emitter moves along the circular vertical arch, we shall define its position using the angle ϕ_E. We say that ϕ_E = 90° when it is located in the azimuthal plane and ϕ_E = 0° when it is at the zenith point of the arch. When the emitter is placed at ϕ_E = 90°, the azimuthal plane corresponds to the scattering plane (defined by incident and scattered wavevectors). We will call this configuration the ‘In Plane’ scattering configuration (IP). When the emitter is moved along the vertical arch out of the azimuthal plane (-11° ≤ ϕ_E ≤ 179°), we have a scattering plane which is different from the azimuthal plane. We will call this situation ‘Out of Plane’ scattering configuration (OP). In this case, it should be noted that the receiver position no longer corresponds to the scattering angle. By scattering angle we mean the angle 6S defined by cos θ_S = k̂_i · k̂ _s where k̂ _i and k̂ _s are the unit vectors corresponding with incident and scattered directions respectively. We adopt the notations θ_S for scattering angle and θ_R for the receiver position angle. We nevertheless represent all the experimental data with respect to the receiver positions (θ_R). In Fig. 3, we represented the scattering angle as a function of receiver position. For the IP configuration both are identical, but for the OP situation, we see that the scattering angle only spans the interval [30°, 123.83°] instead of [0°, 130°].

 

Fig. 3. Scattering angle as a function of receiver position. The curve for the OP configuration was obtained for ϕ_E = 60°.

Download Full Size | PDF

The specificities of the OP configuration also have an impact on the definitions of polarization components. The Bohren and Huffman (BH) definition of the amplitude scattering matrix [30] is the most widely used in literature. It transforms the incident electric fields into the scattered ones in the far field as:

where ∥ and ⊥ stand for parallel and perpendicular to the scattering plane respectively. Here we will call the diagonal elements S ₁ and S ₂, the co-polarization elements and S ₃ and S ₄, the cross-polarizations. The BH convention differs considerably from the one that is used in our experimental set-up. In the above conventions, incident and scattered polarization vectors are determined by the scattering plane as clearly expressed by the previous definition. For the In Plane measurements, the scattering plane remains the same for different scattering angles. For the Out of Plane measurements however, at each scattering angle, we have a different scattering plane. In order to measure the matrix elements, incident and scattered polarization should be changed accordingly. This is very unpractical with our set-up and would be a new source of measurement errors. It is for this reason that in the experimental set-up, fixed polarizations of emitter and receiver are chosen. Hence, the measured matrix elements are no longer identical to those given by Eq. (8). In order to distinguish between these two conventions, we use another notation for the experimentally obtained matrix elements:

where the indices θ and ϕ indicate the polarization state as defined in section 3.

The two representations are related by a transform of the form:

where M _TP is the amplitude scattering matrix in the ‘Theta-Phi’ convention as used for the raw measurements and M _BH is the one in ‘Bohren/Huffman’ convention. The matrices T _i and T _s are the rotation matrices which transform the incident and scattered field components from the TP- to the BH-convention. The TP-convention has the disadvantage that the off-diagonal terms include the contributions of the S ₁ and S ₂ matrix elements, which means that S_ϕθ and S_θϕ can’t be interpreted as the pure cross-polarization signal. The S ₃ and S ₄ elements on the the other hand, do offer such an interpretation. For this reason, and because the light scattering community is familiar with it, we will present the results in the BH-convention. For the In Plane measurements, both transformation matrices in Eq. (10) reduce to the identity matrix.

5. Results

As explained in the previous section, the results we present are the amplitude scattering matrix elements as defined by Eq. (8) and calibrated as described in section 4.2.

Presentation of the results is divided in three parts. First we analyze two (S ₁ and S ₃) of the four matrix elements (both amplitude and phase) as a function of frequency for a single position (position (c) in Fig. 6). We concentrate efforts on those two elements since the discussion for the S ₂ and S ₄ elements is completely analogous. We then consider a single frequency and analyze the scattering matrix as a function of the aggregate orientation. Finally, we compare the different computational methods for one frequency and one orientation.

5.1. Measurements as a function of frequency for a single position

The amplitude of the S ₁ matrix element is depicted in Fig. 4(a). For each frequency, we show the measured amplitudes for both In and Out of azimuthal Plane configurations. Furthermore, in this series of results, we compared the experimental amplitudes with those obtained by Mackowski’s T-matrix code. Here we chose only to compare with one numerical method in order to make interpretation easier. In Fig. 4(b) we did the same for the phase of S ₁. We will start by discussing the amplitude. As can be seen, the In Plane configuration (IP) simulations compare very well to the experiment. For the Out of Plane configuration (OP), differences can be observed especially when the receiver position is close to 60° for 16, 17 and 18 GHz. Those differences are most probably due to tiny geometrical inconsistencies between the input coordinates of the spheres used by the code and the real positions of the spheres in the realized aggregate.

The biggest difference between IP and OP configurations, is the forward diffraction lobe. For the OP configuration, it is almost absent. This can be explained easily by using a scaling approach as applied by Sorensen [36]. It considers the aggregate as an assembly of regions with sizes q ^-1 where q is the amplitude of the scattering wavevector defined by:

 

Fig. 4. Amplitude (a) and phase (b) of S ₁ as a function of frequency for a fixed aggregate position. (() IP simulated, () OP simulated, () IP experimental, () OP experimental).

Download Full Size | PDF

where θ_S is the scattering angle and λ is the wavelength of the incident radiation. All scatterers inside a region with radius equal to q ^-1 can be considered to scatter in phase.

When averaged over all aggregate orientations it can be shown that for qR_g < 1:

where I(q) is the co-polarized intensity. This equation is often called the Guinier equation [36] and describes the region of the forward scattering lobe. It depends only on the size of the object and is valid for all types of aggregates (not only fractal ones).

For the OP configuration, q is never zero since k _i and k _s are never colinear. For this aggregate, the minimum value of qR_g is attained for the scattering direction which lies in the plane defined by the vertical arch. In the frequency range [15–20] GHz, we have 6.82 < qR_g < 9.09 which clearly does not satisfy the Guinier condition, and hence, no forward scattering lobe is observed.

For the analysis of the phase, the scaling arguments are also quite useful. In Fig. 4(b) we represented the phase of S ₁ as a function of frequency. In the region corresponding with the forward scattering lobe, we see a constant value of the phase. As we have written before, we have here qR_g < 1 such that the whole aggregate lies in one q-region and hence all constituent spheres scatter in phase. When the scattering angle is increased, higher q-values are scanned, which corresponds to smaller spatial regions of the object and we are able to resolve smaller details of the object. For our aggregate, this implies that we gain insight of its internal anisotropy. We are sure this internal structure is not due to the excitation of higher multipole moments of the single particles, since in the considered frequency range, the scattering pattern is very close to the dipolar one except for a slight asymmetry between forward- and backscattering hemispheres.

Since for the OP measurements, qR_g is always larger than one, no constant phase region is observed around θR ≈ 0°.

We now turn to the off-diagonal elements of the scattering amplitude matrix. As for the co-polarized case, we only consider one of the two cross-polarizations, since the discussion for each is very similar. In Fig. 5(a), the amplitude of the S ₃ matrix element is represented, again as a function of frequency.

The interpretation of the cross-polarized elements is completely different from the co-polarized ones for both amplitude and phase. Since the aggregate is entirely composed of spherical particles, the cross-polarized signal is due solely to multiple scattering effects. While in the real aggregate there is some geometrical anisotropy due to the presence of glue between the spheres, we consider that this effect can be largely neglected, considering the very good accordance with the T-matrix results in which no such anisotropy is taken into account.

Since S ₃ (and S ₄) are different from zero due to multiple scattering, no forward scattering lobe is observed and the scattered cross-polarized amplitude is almost isotropic [36]. The presence of a forward lobe is a single scattering effect and gives direct size information. Since for the cross-polarizations, no such lobe is present, a model should be built in order to extract the relevant information. While for the co-polarized amplitude, it is quite easy to develop such a model, it is a much harder task for the cross-polarized intensities. This is an unsatisfactory situation since those polarizations contain information about the system shape and anisotropy.

In Fig. 5(a) we see that the results from the T-matrix computations compare extremely well with the experimental measurements. This is quite an achievement considering that the cross-polarized intensities are much more sensitive to orientation and misalignment than the co-polarized ones. To our knowledge, it is the first time that cross-polarized amplitudes are measured to such a high degree of precision with all the minima being well resolved. We think those results can give us important clues to understand better the multiple scattering phenomena in such aggregates. Finally we can say that the agreement with experiment for the IP and OP results are of the same quality.

 

Fig. 5. Amplitude (a) and phase (b) of S ₃ as a function of frequency for a fixed aggregate position. (() IP simulated, () OP simulated, () IP experimental, () OP experimental).

Download Full Size | PDF

 

Fig. 6. Axial rotations of the aggregate seen from above. Figure (c) corresponds with the position of the object in the experiments of section 5.1.

Download Full Size | PDF

We turn now to the phase. As expected, the correspondence is not as good as in the co-polarized case. For the cross-terms, the signal is much lower (less than -50 dB globally and less than -70 dB in the minima) than for the co-terms. For such low signal levels, the phase is almost undetermined and is therefore very sensitive to all kinds of error sources: systematic, alignment/geometric or random noise. Nevertheless, the overlap is still quite good considering the technical difficulties associated with the measurement of S ₃ and S ₄. Before continuing our analysis, we recall that the phase for these elements does not have the same interpretation as for the diagonal elements. The scaling approach is no longer useful here. If we consider the previous situation of qR_g < 1, then the scattered waves from each scatterer are by definition not in phase since they are the result of multiple scattering. In this sense every scattering direction is equivalent. Hence, the interpretation of q as a measure of the resolution is no longer valid. This explains why both In and Out of Plane results of amplitude and phase as a function of the receiver position do not exhibit qualitative differences between them.

5.2. Measurements for different orientations at a single frequency

Here we will analyze how the scattering diagrams evolve when the orientation of the target is changed. As explained in section 2.2, the aggregate is fixed to the support which allows us only to rotate the object around its vertical axis. The configurations that were chosen can be seen in Fig. 6. While measurements have been made in the range [14.5 – 20] GHz in steps of 0.5 GHz, for brevity we restrict ourselves to the 20 GHz results.

We start by describing general aspects of the scattering profile. We see immediately that the width of the forward diffraction lobe shows strong variations with orientation. The diffraction lobe, as we already discussed, is inversely proportional to target size. In the results of Fig. 4 and Fig. 5, the aggregate was facing the emitter with its “flat” part (position (c) in Fig. 6). When it is placed in the positions corresponding with the other angles as represented in Fig. 6, the surface facing the emitter gets smaller and the diffraction lobe broadens.

It should be noted that since the object is not axisymmetric, positions 90° and -90° (or 45° and -45°) are not equivalent.

In the analysis of the results as a function of frequency we encountered some discrepan- cies between simulations and experimental results for the OP configuration. We can now see whether those differences remain. However, as mentioned above, we anticipated an effect related with position. This is confirmed, since now (see Fig. 7) no such differences are observed for the S ₁ element. One could object this sudden accordance is because we are looking at only one frequency. However we confirm that the same observation can be made for the other measured frequencies. This is a quite convincing argument in favor of the hypothesis of possible geometrical inconsistencies between the real model and the virtual one.

5.3. Comparison of numerical methods at a single frequency

5.3.1. Extremal curves

Before starting with the comparison of different models with the measurements, we shall discuss the experimental errors. In Fig. 9(a), we plotted the S ₁ amplitude as a function of the receiver position at 20 GHz. We also plotted the so-called extremal curves. These are obtained by taking the mean field value (taken from 1024 measurements) and adding numerically a noise with experimentally determined properties [37]. For each receiver position, 500 noisy samples are randomly generated and the maximum and minimum values are kept. These extremal curves give an estimation of the measurement error exclusively related with the instrumental and environmental noise. We excluded possible errors related with the positioning of the target. As can be seen from Fig. 9(a), the noise does not exceed 1 dB. This observation will be useful when comparing the different methods.

5.3.2. Comparison of the methods

In Fig. 10 and Fig. 11 we represented the experimental data obtained at 20 GHz, together with the simulated data generated by the four mentioned codes: T-matrix (D. Mackowski), ddscat7.0 (B.T. Draine), Method of Moments (C. Eyraud), T-matrix (B. Stout)

Obviously, all four codes are performing very well in order to reproduce the experimental results. The amplitude of the S ₁ element is reproduced almost within experimental errors. For the phase, differences cannot be explained by the experimental noise only. In order to estimate correctly the degree of confidence with which the computational methods reproduce the experimental results and the importance of each error source in the final result, a much more detailed study should be realized. With the presented results, we can roughly estimate the positioning error and indicate the presence of the geometrical inconsistencies in the model. It is nevertheless possible and useful here to make a more quantitative analysis. Therefore we shall use an error function defined as:

where we considered the array of scattered complex fields for each detector position as a single vector ((E^s _m,1,…,E^s_m,Nr) = E ^s_m and analogously for the simulated fields) and applied the definition of the quadratic vector norm. The values of f ^err for the four matrix elements and for each method are given in table 1. We see immediately that the values of the error function for the elements S ₁ and S ₂ are an order of magnitude lower than the ones for S ₃ and S ₄. This is simply due to the fact that the signal level is 10 times lower for S ₃ and S ₄ than for S ₁ and S ₂. Hence for a constant noise level (∥E ^exp - E ^sim∥² is approximately constant for the co- and cross-terms), we have an error function which is ten times larger. On the whole we see that the f ^err -values are of the same order of magnitude for each of the methods and matrix elements.

The accuracy obtained for each method depends obviously on the parameters chosen to perform the computations. For ddscat7.0 and MoM, it depends essentially on the lattice parameter d. Its value is imposed by the size of the scatterer and the refractive index. The most commonly used discretization parameter that takes into account both physical variables is y = |n|kd where n is the refractive index and k the magnitude of the wavevector in vacuum [38, 39]. It is considered that accurate results (accuracy is usually tested by comparison with Mie theory) of the differential scattering cross-sections (within one percent) are obtained for |n|kd ≲ 0.5 [40]. At 20 GHz this sets an upper limit for d around 0.71 mm. In order not to be just on the upper limit, we chose d = 0.6 mm for both ddscat7.0 and MoM calculations.

 

Fig. 7. Amplitude (a) and phase (b) of S ₁ for different positions, at a single frequency = 20 GHz. (() IP simulated, () OP simulated, () IP experimental, () OP experimental).

Download Full Size | PDF

 

Fig. 8. Amplitude (a) and phase (b) of S ₃ for different positions, at a single frequency = 20 GHz. (() IP simulated, () OP simulated, () IP experimental, () OP experimental).

Download Full Size | PDF

 

Fig. 9. Amplitude (a) and phase (b) of S ₁ at 20 GHz in IP configuration. () Mean value, () Extremal values.

Download Full Size | PDF

 

Fig. 10. Amplitude and phase of S ₁ and S ₂ matrix elements at 20 GHz in IP configuration. (() Experiment, () T-Matrix Mackowski, () ddscat7.0, () MoM, () T-Matrix Stout).

Download Full Size | PDF

 

Fig. 11. Amplitude and phase of S ₃ and S ₄ matrix elements at 20 GHz in IP configuration. (() Experiment, () T-Matrix Mackowski, () ddscat7.0, () MoM, () T-Matrix Stout).

Download Full Size | PDF

Table 1. The values of the error function f ^err for each matrix element and computational method.

View Table

Furthermore, for DDA-type methods, a polarizability prescription has to be assigned. Depending on this choice, better accuracies can be obtained [39]. For the computations presented here, we used the Gutkowicz-Krusin & Draine formulation [41] which presents some subtle corrections to the Lattice Dispersion Relation (LDR) formulation [40].

For the T-matrix computations, the accuracy depends on the maximum order of the cluster-centered expansion. That is why the codes try to estimate the sufficient number of orders for which the results converge and are accurate within preset limits. As a first estimate, one can apply Wiscombe’s Lorentz/Mie truncation criterion [42]:

where x_e is the size parameter of the circumscribing sphere of the cluster. However, both T-matrix codes have more sophisticated tests in order to determine N_O,c. Since all those estimations rely on the internal convergence tests, it would be highly desirable to have another independent criterion based on comparison with experiment. Such a criterion would have a more objective value but is clearly harder to establish. Indeed, it is much easier to cover a wide parameter space by means of numerical simulations than by experiments. Moreover, different kinds of errors have to be corrected before the experiments can be used properly. However, in our opinion, experiments remain the only real test of a theory. Hence, the attainable precision of experimental studies gives a quantitative estimation of the precision limit to which simulations have to be compared.

6. Conclusion

In this study we presented experimental measurements of the full Jones matrix. These results were obtained for a fully characterized complex aggregate (both geometrical and electromagnetic properties). Results for amplitude and phase were presented for different frequencies of the incident radiation. Different orientations of the target were also considered. Finally, a comparison was made with four different computational methods. It was shown, that for the chosen configuration, all methods deliver results with comparable precision. A list with the coordinates of the aggregate’s spheres, the experimental results and the simulations obtained with Mackowskis T-matrix code can be downloaded at (http://www.fresnel.fr/3Ddirect/database.php). Furthermore, such studies are of special relevance to the inversion community, since the original complex target is known and an efficient evaluation of the inversion method can be realized.