Article Polarimetric Emission of Rain Events: Simulation and Experimental Results at X-Band

Accurate models are used today for infrared and microwave satellite radiance simulations of the first two Stokes elements in the physical retrieval, data assimilation etc. of surface and atmospheric parameters. Although in the past a number of theoretical and experimental works have studied the polarimetric emission of some natural surfaces, specially the sea surface roughened by the wind (Windsat mission), very limited studies have been conducted on the polarimetric emission of rain cells or other natural surfaces. In this work, the polarimetric emission (four Stokes elements) of a rain cell is computed using the polarimetric radiative transfer equation assuming that raindrops are described by Pruppacher-Pitter shapes and that their size distribution follows the Laws-Parsons law. The Boundary Element Method (BEM) is used to compute the exact bistatic scattering coefficients for each raindrop shape and different canting angles. Numerical results are compared to the Rayleigh or Mie scattering coefficients, and to Oguchi’s ones, showing that above 1-2 mm raindrop size the exact formulation is required to model properly the scattering. Simulation results using BEM are then compared to the experimental data gathered with a X-band polarimetric radiometer. It is found that the depolarization of the radiation caused by the scattering of non-spherical raindrops induces a non-zero third Stokes parameter, and the differential phase of the scattering coefficients induces a non-zero fourth Stokes parameter.


Introduction: Theoretical Formulation of the Polarimetric Emission by Rain
Accurate models are used today for infrared and microwave satellite radiance simulations of the first two Stokes elements in the physical retrieval, data assimilation etc. of surface and atmospheric parameters [1].Although in the past a number of theoretical and experimental works have studied the polarimetric emission of some natural surfaces, specially the sea surface roughened by the wind (Windsat mission) [2], very limited studies have been conducted on the polarimetric emission of rain cells or other natural surfaces [3][4][5][6].
Emission models at vertical and horizontal polarizations originated by rain cells have been intensively studied (see for example the intercomparison performed in [7]), but very few studies have been conducted to study the full polarimetric emission of rain cells [3][4][5].
The polarimetric emission of a rain cell can be computed by means of the Radiative Transfer Equation (RTE) [6]: where ( ) T z is the temperature distribution, the Stokes vector B T is related to the electric fields incident in the antenna at vertical and horizontal polarizations (E v and E h [V/m]) by [8]: and: • T v and T h are brightness temperatures at vertical and horizontal polarizations, respectively (Instead of T v and T h , the first and second Stokes parameters are sometimes defined as • T U and T V are the third and fourth Stokes elements, respectively, • λ is the electromagnetic wavelength, η is the wave impedance, k B is the Boltzmann's constant, and B is the radiometer's bandwidth.The extinction matrix e k in (1) is given by: 0 where n 0 is the number of particles per unit of volume, k is the wave number, f pq are the forward scattering amplitudes, and the operator < > stands for the average over the ensemble of scatterers in the volume (raindrops, each one with a canting angle).
The emission vector F is given by: and the phase matrix P is given by: ( ) In Equation (5) the f pq coefficients are the vector-scattering amplitude functions that provide the amplitude, phase and polarization information of the scattered field at q-polarization ( ) , ', , ' jkr q q p q E e f e r θ φ = = on top of the atmosphere.This condition means that external radiations (cosmic and galactic noise) are negligible in front of the radiation introduced by the rain cell itself.Section 2 describes the formulation of the computation of the f pq scattering coefficients using the Boundary Elements Method (BEM).Numerical results are presented and compared with previously reported values and Rayleigh and Mie scattering coefficients.Section 3 presents experimental data gathered with an X-band fully-polarimetric radiometer and compares it with the numerical predictions performed using the polarimetric RTE and the BEM.Finally, Section 4 summarizes the main conclusions.

Computation of the Scattering by Raindrops Using the Boundary Element Method
The scattering properties of spherical hydrometeors can be computed using [9]: • the Rayleigh approximation for small spheres in comparison with wavelength, • the Mie expression for spheres whose size is comparable to wavelength, or • optical approximation for spheres much larger than the wavelength.In practice, many hydrometeors are not spherical, but their shape can often be approximated by an ellipsoid or spheroid.Figure 1a shows the raindrop shapes for 13 different equivalent radii from 0.25 to 3.25 mm [10].As it can be appreciated, the larger the equivolumetric raindrop radius, due to air friction, the larger the shape distorsion from the spherical shape (This process is more complex since raindrops merge and break, and their shape oscillates around the equilibrium shape).Recent measurements of drop shapes by video disdrometer under calm air conditions seem to show that the Pruppacher-Pitter model slightly overestimates the real drop deformation [11][12][13] although the contribution of large drops to the total EM scattering properties is relatively small because of their small number density. Figure 1b shows the Laws-Parsons equivalent radii distribution as a function of the rain rate [14].If the ellipsoid size is small compared with a wavelength, a Rayleigh static solution can be used.For spheroidal obstacles the spheroidal function expansion method seems to be a good approximation [15].However, for large drops the shape is not that simple and it is necessary to apply a numerical method to determine the scattering coefficients f pq .A number of different methods has been developed: pointmatching, T-matrix, unimoment, and Fredholm Integral Equation [16,17].
In this paper the Boundary Element Method (BEM) [18] is applied to compute the electromagnetic fields inside and outside dielectric scatterers.The advantage in the Boundary Element Method (BEM) arises from the fact that only the boundary of the domain of the integral equation requires meshing.Thus, the dimension of the problem is effectively reduced by one, that is, the integral equation governing a three dimensional region is transformed into one over its surface.The discretization of only dielectric surfaces reduces user input and storage requirements.The method not only produces precise results with far less data as compared to the methods of finite differences and finite elements, but also works with open region problems without any truncation of the region.The fact that the elements have a parabolic shape makes that BEM models solve problem geometries accurately.

Formulation of the Boundary Element Method (BEM)
In the BEM the scattering problem is formulated in terms of a surface integral equation.The region V is source-free and it is filled with a homogeneous material bounded by a surface S (Figure 2a).The surface S is divided into curvilinear elements.Assuming a time dependence of the electromagnetic fields as ( ) exp j t ω , the integral equation for node 'i' located at r on S is given by [18]: where: ( ) ε is the dielectric permittivity, G i is the free space Green's function, and is the distance between the field point and the source point on the surface.The surface is discretized using the Finite Element Method with second order isoparametric elements and the electromagnetic fields on each element are approximated using the same weighting functions as for the surface's discretization.
Applying Equation (6) to each node 'i' successively, a system of equations is obtained, whose unknowns are the nodal values of the fields: where [D], [H], [C] and [E tan ] are called the complex sub-matrices of fields and coefficients [19] and the subscript 'tan' indicates the tangential components.

BEM Applied to Raindrop Scatterers
In order to apply Equation ( 8), the whole space is divided into two homogeneous regions (Figure 2a): region R 1 inside the raindrop is water with permittivity ε 1 and permeability μ 0 , and region R 0 outside the raindrop is free space.The surface S refers to the boundary between R 0 and R 1 .
The plane wave E e e − = incident on a drop located at the origin induces a field inside the drop and a scattered field S E outside the drops that in the direction 2 K is equal to: where k=2π/λ 0 is the free-space propagation constant, and ( ) , f K K is the vector scattering amplitude function that relates the amplitude, phase and polarization information of the scattered field relative to the incident one.This function is obtained from the solution of the boundary value problem at the surface of a raindrop.
The total fields are then the sum of the incident fields plus the scattered ones: , .
For regions R 0 and R 1 , Equation ( 8) can be written as: and applying the boundary conditions over the surface S, (11a) and (11b) can be cast as: Equation ( 12) is then solved using the conjugate-gradient method and the far fields are computed from the tangential field components over S (equivalent currents) and the vector scattering amplitudes.
Figure 2b presents the raindrop scattering geometry.The incident field is assumed to come from a direction forming an angle θ with respect to ẑ , and the drop axis is ˆp z ( ˆp z ≠ ẑ ).For θ = 0º, α is the drop canting angle (nutation) and β the angle between the azimuth of the direction of the incident electric field (direction of observation) and the wind speed (precession).angle between the azimuth of the direction of observation and the wind direction (precession).
As an example, Figure 3 presents the bistatic scattering patterns in the XZ plane (see Figure 2b) computed using the BEM for 1.4, 2.56, 10.68, 19, 22, 37 and 94 GHz, for Pruppacher-Pitter raindrop shapes with different equivolumetric radii (Figure 1a) [4] for the case of incidence in the Z-direction and X-polarization.Figure 4 shows an example of the bistatic power scattering patterns at 34.8 GHz, in excellent agreement with Oguchi's results [16,17].
Figures 4a and 4b correspond to the patterns for the copolar scattered field in the vertical and horizontal planes, respectively.It should be pointed out that forward scattering is dominant, and that the drop shape asymmetry produces a slight difference in the amplitude of the two copolar scattered fields in the XZ and XY planes.Figures 4c and 4d

Inter-Comparison between Scattering Methods
An excellent agreement has been found between simulation results obtained with the boundary element method and Oguchi's results [15] (Figure 5).At 13 GHz, forward scattering amplitudes computed using the Rayleigh approximation only agree with the BEM solution for drops up to 1 mm equivalent radius (Figure 5).Results computed at 13 GHz with Mie's formulation approach the BEM results (Figure 6a) better than the Rayleigh's ones (Figure 5), specially at large equivalent radius (> 2 mm).However, at higher frequencies (34.8 GHz), the scattering by the "kidney-bean" shaped raindrops deviates too much from Mie's solution and the real part of the scattering functions are underestimated (Figure 6b), while BEM and Oguchi's results keep showing an excellent agreement [17].Real and imaginary parts of the amplitude scattering functions vs. raindrop radii at 13 and 34.8 GHz.Direction of propagation of the incident wave: x , vertical polarization ( ẑ ).BEM solution (solid line) shows an excellent agreement with Oguchi 1977 [17] (not shown), while Mie solution (dashed line) agrees only up to 1 to 1.5 mm radii.

Polarimetric Emission of Rain Events
In the down-welling emission case, the boundary condition at the top of the atmosphere can be reasonably assumed to be [ ] ( , , ) 0 0 0 0 In order to solve Equation (1), the scattering functions have to be computed at arbitrary angles θ, α and β (Figure 2b), which depend on the observation geometry of the radiometer and the canting angles suffered by falling raindrops under a horizontal wind (In the presence of wind, raindrops migrate laterally within the air with negligible deformation and cant.For a constant shear flow in the horizontal direction, the axis of symmetry of the drop lies in the plane determined by the vertical direction and the wind speed vector, with a canting angle in the upwind direction.Depending on its size, the canting angle of the raindrop is estimated to be between 0º and 25º, but the average canting angle is just a few degrees).
The evaluation of equations (3)(4)(5) requires two steps: a) computation of the scattering functions at incident polarization q and scattered polarization p, for each raindrop size according to the observation geometry and canting angle, and b) averaging of the scattering functions over the ensemble of raindrop sizes according to its distribution: where a 0 is the equivalent radius, and N (a 0 ) is the Laws-Parsons drop size distribution [14].Experimental data acquired with an X-band polarimetric radiometer developed at UPC [5, Figure 5] at 30° elevation are presented.Rain rate data was acquired every 5 minutes from three rain gauge stations of the Barcelona city network lying in the ground projection of the antenna beam (Δθ E 3dB = 30°, Δθ H 3dB = 20°), and from the meteorological station of the Institut de Tècniques Energètiques located on Campus, local rain rate was acquired every 30 s, while wind speed and direction were acquired every 5 seconds.
A rain event corresponding to August 18, 1998 is presented.The corresponding meteorological data is plotted in Figure 7. Figure 8 presents the simulated (left) and measured (right) Stokes parameters, respectively.The numerical modeling included data from the on Campus meteorological station only, since it was the one for which wind speed and direction were available.However, since the radiometer collects radiation from all the raindrops within its antenna beam, some peaks in the measured Stokes vector elements do not appear in the simulation results (Figure 8).The limited availability of wind speed and direction data does not allow to simulate better beam filling effects.
As expected, due to the flattened raindrop shape, the attenuation is slightly higher at horizontal polarization than at vertical polarization and therefore T h is slightly higher than T v .Non-zero T U and T V parameters were also measured, although the values of T V are higher than expected, which can be attributed to a too simple modeling of the atmosphere or, probably, to the presence of ice particles that increase the differential phase shift between polarizations.In Figure 8, for high rain rates, changes in wind direction seem to be responsible for the amplitude modulation of T h and T v .

Conclusions
This paper has presented the application of the Boundary Element Method (BEM) to the computation of the scattering functions of the fields scattered on Pruppacher Pitter raindrops with a Laws-Parsons distribution.Numerical simulations are in good agreement with the reported values up to 1-2 mm raindrop size.The accuracy of the method makes it suitable to evaluate the scattering functions of raindrops above 35 GHz, not reported in the literature.Simulation results show that: a) raindrop shape is responsible for T h > T v (at zero canting angle), b) under horizontal winds, T U is originated by the cross-polarization of the emitted/attenuated radiation by raindrops rotated a certain canting angle (size dependent), and c) T V is originated by the differential phase of the scattering amplitude functions between orthogonal polarizations.
The inter-comparison of simulation results with X-band experimental data is in good agreement except for the amplitude of the fourth Stokes vector element, higher than predicted, which can be probably originated by the phase shift difference originated by ice particles in the clouds.
as in Equation (3) the operator < > stands for the average over the ensemble of scatterers in the volume.The solution of Equation (1) entails the application of the following boundary conditions:

Figure 2a .
Figure 2a.Geometry of the scattering problem.Figure 2b.Scattering by a raindrop: Electric incident field propagates along ẑ , drop axis

Figure 2b .
Figure 2a.Geometry of the scattering problem.Figure 2b.Scattering by a raindrop: Electric incident field propagates along ẑ , drop axis along ˆp z , α: drop canting angle (nutation), and β: correspond to the cross-polar scattered field patterns, which are zero for spheroidal drops.

Figure 3 .Figure 4 .
Figure 3. Scattering power patterns in the XZ plane for incidence from Z-direction (indicated by the arrow) at X-polarization, and different equivolumetric drop radii from 0.5 (inner curves) to 3.0 mm (outer curves) and frequencies: a) 1.4 GHz, b) 2.56 GHz, c) 10.68 GHz, d) 19.0 GHz, e) 22.0 GHz, f) 37.0 GHz, and g) 94.0 GHz.Plots normalized and in decibels, maximum value shown in each figure.

Figure 5 .
Figure 5. Real and imaginary parts of the amplitude scattering functions vs. raindrop radii at 13 GHz.Direction of propagation of the incident wave: x , horizontal ( ŷ ) and vertical polarizations ( ẑ ).BEM solution (solid line) agrees well with Oguchi [15] (dotted line) for raindrops up to 3.25 mm radius.Rayleigh solution (dashed line) agrees with BEM solution up to 1 mm.

Figure 6 .
Figure 6.Real and imaginary parts of the amplitude scattering functions vs. raindrop radii at 13 and 34.8 GHz.Direction of propagation of the incident wave: x , vertical polarization ( ẑ ).BEM solution (solid line) shows an excellent agreement with Oguchi 1977[17] (not shown), while Mie solution (dashed line) agrees only up to 1 to 1.5 mm radii.