772
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
Multilag Correlation Estimators for Polarimetric Radar Measurements in the Presence of Noise LEI LEI,*,# GUIFU ZHANG,*,1,# RICHARD J. DOVIAK,@,*,1 ROBERT PALMER,*,1,# BOON LENG CHEONG,# MING XUE,1 QING CAO,# AND YINGUANG LI*,# * School of Electrical and Computer Engineering, University of Oklahoma, Norman, Oklahoma 1 School of Meteorology, University of Oklahoma, Norman, Oklahoma # Atmospheric Radar Research Center, University of Oklahoma, Norman, Oklahoma @ NOAA/National Severe Storms Laboratory, Norman, Oklahoma (Manuscript received 11 January 2011, in final form 18 January 2012) ABSTRACT The quality of polarimetric radar data degrades as the signal-to-noise ratio (SNR) decreases. This substantially limits the usage of collected polarimetric radar data to high SNR regions. To improve data quality at low SNRs, multilag correlation estimators are introduced. The performance of the multilag estimators for spectral moments and polarimetric parameters is examined through a theoretical analysis and by the use of simulated data. The biases and standard deviations of the estimates are calculated and compared with those estimates obtained using the conventional method.
1. Introduction In addition to reflectivity factor (Z), radial velocity (y r), and spectrum width (sy), a polarimetric radar also measures the polarimetric parameters: differential reflectivity (ZDR), copolar cross-correlation coefficient (rhy where h and y stand for horizontal and vertical, respectively), and differential phase (fDP) and its derivative-specific differential phase (KDP; Doviak and Zrnic´ 2006; Bringi and Chandrasekar 2001). The polarimetric radar measurements provide extra information about targeted media, such as phase composition and shape of hydrometeors, and allow better hydrometeor classification and quantitative estimation of the physical states of precipitation (Zrnic´ and Ryzhkov 1999; Zhang et al. 2001; Brandes et al. 2002). The quality of polarimetric radar data (PRD; i.e., the polarimetric parameters), however, degrades when signal-to-noise (SNR) decreases, which limits the usage of PRD to the high SNR region. For example, by using the conventional estimator (Park et al. 2009) ZDR and rhy are used only in the region of SNR that is larger than 5 dB in most hydrometeor classification
Corresponding author address: Guifu Zhang, Suite 5900, 120 David L. Boren Blvd., School of Meteorology, University of Oklahoma, Norman, OK 73072. E-mail:
[email protected] DOI: 10.1175/JTECH-D-11-00010.1 Ó 2012 American Meteorological Society
algorithms and for rain estimation. However, by using the multilag rhy estimator, it is shown herein that the PRD can be used to the region of SNR as low as 0 dB. Data quality depends on the signal processing method used for parameter estimation as well as on radar characteristics. In conventional autocovariance processing, Z and ZDR are estimated by subtracting noise power from lag 0 of the measured autocorrelation functions (ACF); rhy is estimated from the cross-correlation function (CCF) of horizontal (H) and vertical (V) copolar signals at lag 0 and normalized by the noise-corrected ACFs, also at lag 0; yr is estimated from lag 1; and sy is estimated from lag 0 and lag 1 of the ACF (Doviak and Zrnic´ 2006). Because noise power is measured after a volume scan when the transmitter is turned off and the beam is at a high elevation angle and directed at clear skies or light precipitation, the true noise power at the angle of weather data differs from the measured noise power at the higher angle (Fang et al. 2004). This causes biases in spectral moment and PRD estimates when the conventional estimator is used. Because noise power measurements use many more samples than M, the number used to estimate moments and PRD, the noise power is measured with negligible variance, but it can have significant bias. Thus, polarimetric parameter estimates can have asymptotic bias (i.e., bias independent of the number M of samples processed), especially when
JUNE 2012
LEI ET AL.
the SNR is low (Melnikov and Zrnic´, 2007). For example, if proper noise corrections are not made, then the copolar correlation coefficient rhy will have a negative asymptotic bias because the normalizing factors should be signal power, not signal plus noise power. Another source of bias of rhy is bias originating from the limited number of samples used in the estimates; this is always positive as will be shown (Fig. 11a). There have been some attempts to improve the quality of PRD by reducing the effects of noise and freeing it from asymptotic bias. For example, an LDR estimator introduced by Hubbert et al. (2003) is effective in low SNR regions. It uses the cross-to-cross covariances in contrast to using just the autocovariance of the crosspolar time series to calculate cross-polar power. In addition, one-lag estimators were introduced to estimate ZDR and rhy from lag-1 correlations (Melnikov and Zrnic´ 2007). The one-lag estimators avoid lag-0 data to estimate the polarimetric parameters ZDR and rhy, and hence the one-lag estimates of these two parameters are unbiased by noise. Furthermore, the standard deviations of one-lag estimators of ZDR and rhy are practically the same as the conventional estimates when spectral width is less than 6 m s21 and SNRs are larger than 5 dB; at wider spectral widths, the standard deviations of one-lag estimates are larger than the ones from the conventional algorithm (Doviak and Zrnic´ 2006, section 8.5.2). Because one-lag estimators only use lag 1, some useful data at other lags within the correlation time (i.e., data from 2t c to t c) are discarded. That is, because data from contiguous lags are correlated the added information is likely to be small for high SNR signals, but significant at low SNR (Zhang et al. 2003, 2004). This added information and its impact on improving estimates at low SNR are quantified in this paper. Correlation estimators that do not use zero lag produce moment and PRD estimates that are free from asymptotic bias, which is independent of the number M of samples and results when the measured noise power deviates from the true noise power. All useful ACF lag correlations except for the ones at lag 0 are used to fit a Gaussian function. This research is an extension of the multilag correlation estimator developed for cross-beam wind measurement using spaced antenna interferometry (Zhang et al. 2004). It is true that some observed Doppler spectra do not follow the Gaussian function shape well (Moisseev et al. 2008). These include those caused by second-trip echoes, nonuniform wind shear, and nonuniform reflectivity in regions of uniform shear and turbulence. The secondtrip echoes can be detected and excluded before applying the multilag processing. The main purpose of this work is to improve spectral moment and PRD estimation
773
for weak weather echoes, such as those from stratiform cloud–precipitation, and for power spectra that are mostly Gaussian (Doviak and Zrnic´ 2006, section 5.2). The noise effects on the nonzero lags of ACF and CCF at low SNR can cause the non-Gaussian shape of the observed correlation functions, but these noise effects can be reduced by fitting the Gaussian model at multiple lags. The window effect can be neglected if the number of samples is large enough (Doviak and Zrnic´ 2006, section 5.1.4); therefore, we still use a Gaussian model to fit ACF or CCF data in this research. The multilag estimators for polarimetric parameters are described, and the corresponding expressions are presented in section 2. Biases and standard deviations of the multilag estimates are derived and plotted in section 3. The performance of the multilag estimators are examined and compared with conventional pulse-pair processing (PPP; Doviak and Zrnic´ 2006). Theoretical results of the estimation errors are provided and verified by numerical simulations. Conclusions and discussions are provided in last section.
2. Multilag estimators The idea of a multilag estimator is to use many available and informative lag correlation estimates to fit a Gaussian function, and hence to obtain more accurate estimates of spectral moments and PRD at low SNR than those obtained using conventional estimators. Weather signals from regions of uniform reflectivity, shear, and turbulence should have ACFs closely resembling a Gaussian function (Doviak and Zrnic´ 2006). Depending on the correlation time (inversely proportional to spectrum width) of the weather signal, the number of sequential correlation lags, excluding zero lag, can be two, three, four, etc. They are correspondingly called two-lag, threelag, four-lag, and x-lag estimators (Lei 2009, Lei et al. 2009a). To illustrate the idea, a numerical simulation is shown in Fig. 1 for a weather signal that is generated using the spectrum method of Zrnic´ (1975). In this case, the frequency is 3 GHz, the number of pulses is 128, and pulse repetition time (PRT) is 0.001 s; these three radar parameters are used throughout this work. Figure 1 shows simulated estimates and fitted Gaussian functions if SNR is 3 dB and spectrum width is 3 m s21. The raw ACF estimates are shown as solid dots connected with a dotted line. The fitted Gaussian functions for two- and four-lag estimators are shown as dashed and solid lines, respectively. Because the ACF at lag 0 is excluded from the fittings, the resulting fitted Gaussian functions provide better estimates of the spectral moments than those obtained using the conventional PPP in which an
774
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
radial velocity, Nh,y is the expected value of white noise power in the H or V channels, dm 5 1 for m 5 0, and zero otherwise; and rhy is the copolar correlation coefficient at lag 0 and fdp is differential phase. If the beams are matched, then the Doppler mean velocities y and the correlation times tc are the same for weather signals from the H and V channels (e.g., y h 5 y y 5 y; Melnikov and Zrnic´ 2007; Sachidananda and Zrnic´ 1985, 1986; Bringi and Chandrasekar 2001). Taking the natural logarithm of the magnitude of both sides of Eq. (1), the expected ACF is rewritten as ym 5 ln[jRh,y (mTs )j] 5 am2 Ts2 1 b
FIG. 1. An example of multilag estimated ACF. Data from two or four lags are used to determine Gaussian parameters to estimate the expected ACF magnitude.
for m 5 1, 2, 3 . . . , where 21/2t 2c 5 a, ln(Sh,y ) 5 b. Then, the estimated ym is given by ln[jR^h,y (mTs )j] 5 y^m , and the caret ( ) denotes the estimated value. Using both the expected ^ is ACF and estimated ACF the merit function F(^ a, b)
^
^ 5 F(^ a, b) independent measurement of noise power is used to estimate the signal power from the ACF datum at zero lag. The fitted Gaussian function is then used to calculate moments (i.e., power, velocity, spectrum width, etc.) and PRD. Detailed fitting procedure and radar moment estimators are in the following subsections.
a. General expressions for fitted ACF and CCF The expected ACF of weather radar signals is mostly represented by the Gaussian form (Janssen and Van Der Spek 1985) and is given by Doviak and Zrnic´ (2006, p. 125), which is extended to include the autocovariance function for the horizontally H and vertically V polarized signals as well as their cross-covariance functions. The expected ACF is jpmy 1 Nh,y dm , Rh,y (mTs ) 5 Sh,y r(mTs ) exp 2 yN
(1)
and the expected CCF is given by qffiffiffiffiffiffiffiffiffiffi jpmy 1 jfdp , Chy (mTs ) 5 Sh Sy rhy r(mTs ) exp 2 yN
N
å
m51
(m2 a^Ts2 1 b^ 2 y^m )2
5 a^2 Ts4
N
N
N
^ 2 å m2 2 2^ aTs2 å m2 y^m å m4 1 2^abT s
m51
1 N b^2 2 2b^
m51
N
å
m51
m51
N
y^m 1
å
m51
y^2m ,
(4)
where N is the number of lags used and N $ 2 (e.g., N 5 ^ 3 indicates that data at lags 1, 2, and 3 are used); a^ (or b) is the estimate of a (or b) obtained from the fitted Gaussian correlation function using N lags, and y^m is the estimated value of ym using M 2 1 2 m lags, where M is the number of weather signal samples. According to the definition of the least squares fit, the merit function Eq. (4) reaches its minimum value when 8 ^ ›F(^ a, b) > > 50 < ›^ a . ^ > a, b) > :›F(^ 50 ›b^ These two equations are solved simultaneously for a^ ^ and after manipulation the solution can be put in and b, the reduced form
(2) where the subscripts h, y, and hy mean the parameters are calculated by using signals from the H or V channels, and from both the H and V channels [i.e., for Chy (mTs ) ]. The m 5 0, 1, 2 . . . , N is the lag number, Ts is the pulse repetition time, and S is the signal power; r(mTs ) 5 exp[2(mTs )2 /2t 2c ] is the correlation coefficient of the weather signals and tc 5 l/4psy , where sy is the spectrum width and y N is the Nyquist velocity; y is the mean
(3)
N
30 a^ 5
å
[6m2 2 (N 1 1)(2N 1 1)]^ ym
m51 , Ts2 N(N 2 1)(N 1 1)(2N 1 1)(8N 1 11)
(5a)
N
6 b^ 5
å
m51
(3N 2 1 3N 2 1 2 5m2 )^ ym N(N 2 1)(8N 1 11)
where N $ 2.
,
(5b)
JUNE 2012
^ From the above equations and the definition of a^, b, and y^m , the general expression for multilag power and 8 > < (N ) (N ) S^h,y 5 jR^h,y (0)j 5 exp> :
(N )
s ^ h,y
N
6
å
f(3N 2
m51
(6c)
where qm are integers to unwrap the phase for lag m.
*
(N ) ^ Chy (0) 5 exp
N
3
å
m52N
(6a)
(6b)
To Gaussian fit the CCF from a set of measurements, a similar procedure as that used to estimate the ACF is applied to Eq. (2). The only difference is that negative lags are used together with zero and positive lags, that is, m 5 2N,2(N21), . . . , 21, 0, 1, . . . , N. As a result, the four-lag estimator, for example, actually uses nine lags of data. The reason why more lags can be used in fitting the CCF than the ACF is that the CCF is not symmetric and lag 0 is not biased by noise. After applying the least squares fit and similar calculations used to estimate the ACF, the result of CCF fitting is
+
f(3N 2 1 3N 2 1 2 5m2 ) ln[jC^hy (mTs )j]g (2N 2 1)(2N 1 1)(2N 1 3)
For polarimetric radar parameters, the expressions for the differential reflectivity and correlation coefficient are " (N ) # S^ (N ) ^ (8) ZDR 5 10 log10 h(N ) , S^y (N )
5m2 )
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u N u u 30 å f[6m2 2 (N 1 1)(2N 1 1)] ln[jR^ (mT )j]g pffiffiffiffiffiffiffiffiffi pffiffiffi u h,y s l 22^ a l 2t 5 , 5 2 m51 2 4p 4p Ts N(N 2 1)(N 1 1)(2N 1 1)(8N 1 11)
N ny o N ^ (mT ) 1 2pq ] , ^y (N ) 2 1 å [arg R h,y h,y s m N m51 mp
r^hy
spectrum width for H or V channel signals can be written as
9 > = 1 3N 2 1 2 ln[jR^h,y (mTs )j]g > ;, N(N 2 1)(8N 1 11)
where the superscript (N) means an N-lag-fitted estimate. The multilag method can also be applied to obtain an estimate of y [i.e., the unwrapped phase angle of ACFs at multiple lags can be fitted by a linear line to obtain y (Lee 1978; May et al. 1989)]. If measurements are made at several lags (i.e., N lags) and the ambiguities are resolved, a multilag estimator for the mean Doppler velocity is
(N ) ^ Chy (0) 5 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi , (N ) (N ) S^h S^y
(9)
(N ) (N ) (N ) where S^h , S^y , and jC^hy (0)j are calculated using Eqs. (6a) and (7). Differential phase fdp can also be obtained from multilag data of the CCF angle. CCF data at positive and negative lags are multiplied to cancel the Doppler velocity; that is,
^ (N ) 5 f dp
775
LEI ET AL.
N 1 arg[C^hy (mTs )C^hy (2mTs )]. å 2(N 1 1) m50
(10)
.
(7)
Expressions of the statistical characteristics for the Doppler velocity and differential phase estimates are lengthy and will not be presented in this paper. Equations (6)–(9) are the general expressions of power, spectrum width, differential reflectivity, and copolar cross-correlation coefficient magnitude estimates obtained using multilag estimators. For example, if N 5 4, it means that lags 1, 2, 3, and 4 of the ACF data and lag 0, lag 61, lag 62, lag 63, and lag 64 of CCF data are used to estimate polarimetric parameters. The value of the chosen N depends on the SNR, number of pulses, spectrum width, noise type, etc. The following sections provide detailed description for N 5 2, 3, and 4.
b. Specific estimators TWO-LAG ESTIMATOR The two-lag estimator uses lags 1 and 2 of ACF and lags 0, 61, and 62 of CCF to estimate polarimetric parameters. Substituting N 5 2 into Eqs. (6)–(9), the
776
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
formulas to estimate spectrum width, differential reflectivity, and correlation coefficient are given by (2) S^h,y 5
jR^h,y (Ts )j4/3 , jR^ (2T )j1/3 h,y
l (2) s ^ h,y 5 pffiffiffiffiffi 24pTs
(11)
s
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi lnjR^h,y (Ts )j 2 lnjR^h,y (2Ts )j, (12)
VOLUME 29
such as power, differential reflectivity, and correlation coefficient, the statistical performance of the two-lag estimator is poor based on the analysis shown in section 3. This is so because the exponent 4/ 3 in the numerator of Eq. (11) is larger than 1; this increases the noise effect at low SNRs. To decrease the exponent and improve the statistical performance, more lags are needed (e.g., three or four lags), as shown in appendix A.
! jR^h (Ts )j4/3 jR^y (2Ts )j1/3 (2) ^ , ZDR 5 10 log10 jR^ (2T )j1/3 jR^ (T )j4/3
(13)
3. Performance of the estimators
(2) [jR^ (2T )jjR^ (2T )j]1/6 (2) s y s r^hy 5 C^hy (0) h , [jR^h (Ts )jjR^y (Ts )j]2/3
(14)
In this section, the performance of the multilag estimators is examined through perturbation analysis and compared with the performance of conventional estimators. Theoretical statistical biases and standard deviations of power, spectrum width, and differential reflectivity estimates are calculated and verified by simulations.
h
s
y
s
(2 )
where jC^hy (0)j is obtained by substituting N 5 2 into Eq. (7). For spectrum width estimates, the two-lag estimator (i.e., using lags 1 and 2) has been recommended to be used in place of the conventional estimator (i.e., one using lag 0 and lag 1) if SNR is low (Doviak and Zrnic´ 2006). The two-lag estimator has better performance than the conventional estimator at narrower spectrum widths, but poorer performance at wider spectrum widths (Srivastava et al. 1979). However, for other parameters
a. General expression of statistical analysis for multilag estimators To calculate the bias and standard deviation of multilag estimators, perturbation analysis is used (Zhang et al. 2004; Melnikov and Zrnic´ 2007) and terms to second order of the Taylor expansion are retained. For example, the following is the Taylor expansion of signal power in ^ ^ ^ )j, jR(2T several variables from jR(T s s )j, to jR(NTs )j:
8 ^ 9 nN ^ [jR(Ts )j 2 jR(Ts )j]n1 . . . [jR(NT > > s )j 2 jR(NTs )j] > > > > = ‘ ‘ < n ! . . . n ! 1 N (N ) ^ ^ " # . S [jR(Ts )j, . . . , jR(NT )j] 5 . . . å å s > ›n111nN S(N ) n1 50 nN 50 > > > > > : nN [jR(Ts )j, . . . , jR(NTs )j] ; ^ )jn1 . . . ›jR(NT ^ ›jR(T s s )j The R(Ts) is the ACF and S(N) is the signal power estimated by using N lags of ACF; n1, . . . , nN are the number of derivatives for each variable in the Taylor expansion. Similar expansion can also apply to other radar
(15)
parameters. If an N lag estimator is used, Sh,y and sy both have N variables, but ZDR has 2N variables and rhy has 4N 1 1 variables. The difference between the estimated and true S(N) is expressed by
(N ) ^ )j, . . . , jR(NT ^ dS(N ) 5 S(N ) [jR(T [jR(Ts )j, . . . , jR(NTs )j] s s )j] 2 S 8 ^ 9 n1 nN ^ [jR(Ts )j 2 jR(Ts )j] . . . [jR(NT > > s )jj 2 jR(NTs )j] > > > > > < = ‘ ‘ > n1 ! . . . nN ! # " 2 S(N ) [jR(Ts )j, . . . , jR(NTs )j] 5 å ... å n 11n > > ( N ) 1 N n1 50 nN 50> > › S > > > [jR(Ts )j, . . . , jR(NTs )j] > : ; ^ )jn1 . . . ›jR(NT ^ ›jR(T )jnN
s
s
9 8 ^ nN ^ [jR(Ts )j 2 jR(Ts )j]n1 . . . [jR(NT > > s )jj 2 jR(NTs )j] > > > = < ‘ ‘ > n ! . . . n ! 1 N # " ; 5 å ... å n 11n ( N ) NS > › 1 n1 50 nN 50> > > > > [jR(Ts )j, . . . , jR(NTs )j] ; : ^ )jn1 . . . ›jR(NT ^ ›jR(T )jnN s
s
(16)
JUNE 2012
777
LEI ET AL.
if n1, n2, . . . , nN, the number of derivatives of the Taylor expansion are not simultaneously equal to zero. According to the definition of bias and variance, statistical bias and variance of S are defined as Bias[S(N ) ] 5 hdS(N ) i,
* Bias[S(N ) ] 5
(17)
Var[S(N ) ] 5 h[dS(N ) ]2 i.
Putting Eq. (16) into Eqs. (17) and (18), we arrive at an expression for the bias
+
8 ^ 9 nN ^ [jR(Ts )j 2 jR(Ts )j]n1 . . . [jR(NT > > s )j 2 jR(NTs )j] > > > < = ‘ ‘ > n1 ! . . . nN ! # " . . . , å å > > ›n1 11nN S(N ) n1 50 nN 50 > > > > : ; nN [jR(Ts )j, . . . , jR(NTs )j] ^ )jn1 . . . ›jR(NT ^ ›jR(T s s )j
and the following expression for the variance:
*0
Var[S(N ) ] 5
(18)
(19)
+
8 ^ nN 912 ^ [jR(Ts )j 2 jR(Ts )j]n1 . . . [jR(NT s )j 2 jR(NTs )j] > > > > > B ‘ < =C ‘ > n1 ! . . . nN ! B C B å ... å C . # " B C > ›n1 11nN S(N ) nN 50 > > > @n1 50 A > > : ; [jR(T )j, . . . ,jR(NT )j] s s ^ )jn1 . . . ›jR(NT ^ ›jR(T )jnN s
In Eqs. (19) and (20), n1, n2, . . . , nN are not simultaneously equal to zero. Following the same procedure used to obtain the bias and variance of the signal power (N ) (N ) estimates, the bias and variance of s ^ (yN ) , Z^DR , and r^hy can be obtained. Detailed calculations for the two- and four-lag estimators are in appendixes B and C. To simplify the calculations, the Taylor expansion is kept to second derivative term (i.e., n1 1 1 nN # 2).
b. Statistical analysis for power, spectrum width, differential reflectivity, and correlation coefficient 1) SIGNAL POWER For the conventional estimator, asymptotic bias comes from errors in estimating noise power. Noise power for the Weather Surveillance Radar-1988 Doppler (WSR88D) is measured at the highest elevation angle of each volume scan and while the transmitter is shut off. Measured noise power deviates from true noise power because the measured noise power is measured at different times and directions. Noise power originates within the receiver chain from lossy [e.g., waveguides, radio frequency (RF) filters, etc.], and active components (e.g., low noise amplifiers etc.) as well as from external sources (e.g., the earth, cloud/precipitation, sun, etc.). Noise power level depends on beam direction, temperatures of the radar’s lossy components, and external sources, etc. (Fang et al. 2004). Noise power variations of 1 dB are observed frequently with WSR-88D, although most noise bias is within 0.5 dB of the mean value (Melnikov and
(20)
s
Zrnic´ 2007). Sometimes, exceptionally strong electrical emissions cause noise to increase more than 8 dB above the receiver noise level (Fang et al. 2004). Because of the measured noise power deviating from true noise power, the conventional estimator does not perform well at low SNR. The one-lag and multilag estimators are efficient ways to mitigate this problem. The following five power estimators: conventional and one, two, three, and four lag are expressed as S^h,y 5 R^h,y (0) 2 N h,y , (conventional estimator) (21) (1 ) S^h,y 5 jR^h,y (Ts )j, (one-lag estimator)
(22)
and the two-, three-, and four-lag signal power estimators are given by Eqs. (11), (A1), and (A5). To compare these estimators, the bias and standard derivation are calculated for each.
(i) Analytical expressions for conventional power estimate bias and variance In (21), N h,y is an average of noise power over many more than M samples and the variance of N h,y should be minimal (Fang et al. 2004). Therefore, we assume N h,y is measured with negligible variance, but it can have bias, and thus S^h,y bias depends on how much N h,y deviates from its true value. If noise power has negative bias, then the signal power will have positive bias, and vice versa. The variance of the conventional estimate S^h,y can be found in Melnikov and Zrnic [2004, Eq. (A6)] as
778
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
" # 2(SNR)h,y 1 1 1 2 2 ^ ^ , Var(Sh,y ) 5 h(dPh,y ) i 5 Sh,y 1 MI M(SNR)2h,y (23) " MI 5 M 1 1 2
#21
M21
å
m51
2
(1 2 m/M) j r(mTs )j
, (24)
in which M is the number of samples used to estimate S h,y, and MI is the equivalent number of independent samples of the signal power having a correlation coefficient r(mTs) (section 6.3.1.2, Doviak and Zrnic´ 2006). For large M and sy /2yN 1, MI ’ 4MTs p1/2 sy /l.
(ii) Analytical expressions for one-lag power estimate bias and variance The one-lag signal power estimator [i.e., Eq. (22)] is a biased estimator because the ACF at lag 1 underestimates signal power unless the correlation time is very large (or the spectrum width is very small; see Fig. 2). The negative bias is given by
FIG. 2. Fractional bias of power estimates M 5 128, Ts 5 0.001 s, l 5 0.1 m (these radar parameters apply to all figures). For the one-lag estimator, bias is independent of SNR. (1)
Bias[S^h,y ] 5 hjR^h,y (Ts )j 2 Sh,y i 5 Sh,y r(Ts ) 2 Sh,y , (25) and the variance of the one-lag signal power estimates is
(1 )
Var[S^h,y ] 5 h[jR^h,y (Ts )j 2 jRh,y (Ts )j]2 i 5 hjR^h,y (Ts )j2 i 1 S2h,y r2 (Ts ) 2 2Sh,y r(Ts )hjR^h,y (Ts )ji 5 hjR^h,y (Ts )j2 i 2 S2h,y r2 (Ts ).
(26a)
Using Eq. (C5) in Melnikov and Zrnic´ (2004) to approximate the hjR^h,y (Ts )j2 i in Eq. (26a) and then combine the (1) result with Eq. (1), the Var[S^h,y ] is " # 2(SNR)h,y 1 1 1 (1 ) 2 2 S2h,y r2 (Ts ) 1 2 Var[S^h,y )] ’ jRh,y (Ts )j 1 1 r (Ts )MI (M 2 1)r2 (Ts ) (SNR)2h,y " # 2(SNR)h,y 1 1 1 2 1 2 5 jRh,y (Ts )j r (Ts )MI (M 2 1)r2 (Ts ) (SNR)2h,y # " 2(SNR)h,y 1 1 1 2 2 1 2 5 Sh,y r (Ts ) r (Ts )MI (M 2 1)r2 (Ts ) (SNR)2h,y " # 2(SNR)h,y 1 1 1 2 . (26b) 1 ’ Sh,y MI M(SNR)2h,y
Note that this equation is the same as Eq. (23). This is because we have assumed Ts is much smaller than the correlation time, resulting in signal estimates at lag 0 and 1 that are nearly equivalent.
(iii) The performance comparison of power estimators For the two-, three-, and four-lag estimators of signal power [i.e., computed from Eqs. (11), (A1), and (A5)],
JUNE 2012
779
LEI ET AL.
FIG. 3. Normalized standard deviation of signal power estimates.
the bias and variance can be calculated using Eqs. (19) and (20). Detailed calculations are shown in appendix B. Fractional biases of power estimates using different estimators (i.e., one, two, and four lag) are shown in Fig. 2. Results presented in this paper are focused on weather types (e.g., stratiform) that produce weather signals with spectrum widths less than 4 m s21; the median value of spectrum widths for all types of weather other than squall lines is less than 2.5 m s21 (Fang et al. 2004). For the conventional estimator, the bias is asymptotic and depends on how much the measured noise value deviates from the truth. Bias is the function of SNR, M (number of samples), and spectrum width. Increasing M and SNR decreases the bias. In Fig. 2, the one-lag estimator is a strongly biased estimator and has low bias only when the spectrum width is very small. The four-lag estimator has less bias than the two-lag estimator for a spectrum width less than 3.3 m s21 when the SNR is less than 10 dB, although the biases for all the estimators, except for the one-lag estimator, are small. As the SNR increases, the two- and four-lag estimator bias approaches zero. The theoretically derived standard deviation of power is shown in Fig. 3 and is verified with simulations in Fig. 4. In the simulation, H and V channel signals are generated by the spectrum method (Zrnic´ 1975). Then, the copolar correlation coefficient and power differences are generated for the polarimetric radar parameters (Galati and Pavan 1995). The conventional estimator has the same standard deviation as the one-lag estimator because it is calculated here assuming that the noise power is accurately measured. If not, the standard deviation of the conventional estimator will be larger. The two-lag estimator does not have significantly better performance than the conventional or one-lag estimator, and the reason is given in last paragraph in section 2.
FIG. 4. Comparison of theory and simulation for the standard deviation of power estimates at SNR 5 5 dB.
The four-lag estimator is slightly better than other estimators when the SNR is under 10 dB and the spectrum width is less than 3.5 m s21 with the given radar parameters. In conclusion, the multilag estimators can decrease the standard deviation of the signal power estimates only by small amounts (Fig. 3). In Fig. 4, the theoretical results and simulations differ at larger spectrum widths for the four-lag estimator. This is caused by the undersampling rate at larger spectrum widths. In other words, when the spectrum width is large, the ACF is narrow. Thus, PRT needs to be short enough to obtain accurate ACF estimates, especially when using more lags. By decreasing the PRT, the difference between theory and simulation at larger spectrum width diminished significantly. Another reason for the difference between theory and simulation comes from the Taylor expansion. In the calculation of bias and standard deviation, the Taylor expansion is used. However, this expansion is only valid when the estimated value is close to the true value. Thus, the theory results will be valid only when the bias and standard deviation of the spectral moments and PRD estimates are not large (within 10% of the true value). If the number of samples or SNR increases, then the differences between the theoretical and simulated results are smaller.
2) SPECTRUM WIDTH Four spectrum width estimators: conventional (lag 0 and 1) and two, three, and four lags are expressed, respectively, by l s ^ h,y 5 pffiffiffi 2 2pTs
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ln(jS^h,y j) 2 ln[jR^h,y (Ts )j],
(conventional estimator)
(27)
780
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
FIG. 6. Standard deviation of spectrum width estimates.
FIG. 5. Bias of spectrum width estimates resulting from (a) sy estimation error and (b) equal measured H and V noise power biases at SNR 5 5 dB.
and by Eq. (12), (A2), and (A6) for the two-, three-, and four-lag estimators. The standard deviation of the conventionally derived estimates can be found in Zrnic´ (1977, 1979). The bias and variance for the two- and fourlag estimators can be calculated from Eqs. (12), (A2), and (A6) using the perturbation method through Eqs. (15)– (20). The detail calculations can be found in appendix B. The biases of spectrum width estimators are compared in Fig. 5a. Because the estimation error of s ^ y is not much smaller than its expected value s ^ y (i.e., if the spectrum width is less than 1 m s21; Fig. 6), the condition for Taylor expansion up to the second order is not satisfied. Hence, the results at very narrow spectra are not reliable, but those at sy .1 m s21 are accurate. At low SNR, the four-lag estimator provides notable decreases of spectrum width bias. On the other hand, the bias of sy estimates, obtained using the conventional estimator, is strongly dependent on the measured noise bias; if noise bias is negative (e.g., 20.5- and 21-dB
differences between true noise and measured noise in Fig. 5b), s ^ y bias is positive and large (Fig. 5b). From Fig. 6, it is seen that the standard derivation of the spectrum width estimates strongly depends on the spectrum width. In general, spectrum width estimated using more lags has smaller standard deviation at small spectrum widths. However, when the spectrum width becomes large, the multilag fitting method has poorer performance than the other estimators. This is so because data at large lags introduce more uncertainty than information. In addition, the standard deviation of the conventional estimator (i.e., the one using lags 0 and 1) is shown in Fig. 6.6 of Doviak and Zrnic´ (2006), and the standard deviation of the two-lag estimator, which uses lags 1 and 2, is shown in their Fig. 6.7. The good agreement of our results with those in the cited reference supports our discussion of Fig. 6. In Fig. 7, the theoretical analysis of standard deviation is verified with simulations. They match reasonably well over the spectrum width interval from about 0.3 to 3.6 m s21 for the four-lag estimator, and for spectrum widths larger than about 1 m s21 for the two-lag estimator. The poorer performance at small spectrum widths is caused by not having enough independent samples to yield small relative errors required for the Taylor expansion to be satisfactory. For the four-lag estimator, the difference between theory and simulation at large spectrum widths is due to the fact that the PRT is not very small compared with the correlation time.
3) DIFFERENTIAL REFLECTIVITY The differential reflectivity is a ratio of the reflected horizontal and vertical signal powers. Five estimators [conventional (lag 0) and one, two, three, and four lag) are examined. They are explicitly written as
JUNE 2012
LEI ET AL.
781
FIG. 7. Comparison of theory and simulation for the standard deviation of spectrum width estimates at SNR 5 5 dB.
! S^h ^ , (conventional estimator) ZDR 5 10 log10 S^y
(28)
^ jRh (Ts )j (1) ^ , (one-lag estimator) ZDR 5 10 log10 jR^y (Ts )j
(29)
whereas Eqs. (13), (A3), and (A7) are the two-, three-, and four-lag estimators. The one-lag estimator is introduced in Melnikov and Zrnic´ (2004, 45–54, 2007). (1 ) The bias and variance of Z^DR and Z^DR can also be found in the report. The conventional estimator calculated herein assumes that the noise level is accurately measured. If not, the bias of the conventional estimator (Fig. 8b) will be larger than the results shown in Fig. 8a. For Eqs. (13), (A3), and (A7), the bias and variance can be calculated using the perturbation method through general Eqs. (15)–(20). The detailed calculations can be found in appendix C. It is seen (Fig. 8a) that as more lags are used, ZDR bias decreases although the bias is very small for all estimators. However, if there is bias in the measured noise power, then the improvement of multilag estimators over the conventional estimator is significant (Fig. 8b), where the same noise power bias is assumed to be equal in both H and V channels. On the other hand, the H and V noise powers can differ and the difference can be as much as 1 dB (Melnikov and Zrnic´ 2007). In this case the ZDR bias of conventional estimators would be larger than that shown in Fig. 8b. The multilag estimator can mitigate this kind of bias as the one-lag estimator does. Figure 9 shows the standard derivations of different estimates as a function of the spectrum width. If the spectrum width is small, then a higher-order multilag estimator can be used to provide better performance.
FIG. 8. Bias of differential reflectivity estimates resulting from (a) estimation error and (b) equal measured H and V noise power biases of 20.5 and 21.0 dB at SNR 5 5 dB, rhy 5 0.97, and ZDR 5 1 dB.
For example, although the bias of differential reflectivity can be improved by only 0.01 dB when SNR 5 0 dB, the standard deviation can be improved by about 0.1 dB at spectrum widths less than 3 m s21 if a four-lag estimator is used. However if the spectrum width is large, then the four-lag fitting method has poorer performance than the other estimators. In Fig. 10, the theoretical analysis of the standard deviation SD(Z^DR ) is compared with simulations; the comparison shows they match reasonably well except for the four-lag estimator at large spectrum widths. The improvement of multilag estimator forZDR is very limited.
4) CORRELATION COEFFICIENT The correlation coefficient is the correlation between the horizontally and vertically copolarized weather echoes. Both ACF and CCF are used to estimate the correlation coefficient. For the ACF, lags 1, 2, 3. . . are used; for the CCF, lags. . .23, 22, 21, 0, 1, 2, 3. . . are used.
782
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
FIG. 9. Standard deviation of differential reflectivity estimates of rhy 5 0.97 and ZDR 5 1 dB.
The four types of estimators are conventional (lag 0), and one, two, three, and four lag, respectively,
r^hy (0) 5
jC^hy (0)j , (S^ S^ )1/2
(conventional estimator)
(30)
h y
(1 )
r^hy (0) 5
jC^hy (2Ts )j 1 jC^hy (Ts )j , 2 jR^h (Ts )R^y (Ts )j]1/2
(one-lag estimator) (31)
whereas Eqs. (14), (A4), and (A8) are the two-, three-, and four-lag estimators. In Eqs. (14), (A4), and (A8), we use both positive and nonpositive lags to estimate (N ) C^hy (0). Figure 11a shows that the four-lag estimator produces significantly lower bias estimates than the other estimators do at low SNR (e.g., , 5 dB) and for spectrum widths less than about 3.5 m s21. However, if there is bias in the measured noise power, then the improvement of multilag estimators over the conventional estimator is even larger (Fig. 11b). Furthermore, the SD(^ rhy ) for the four-lag estimator is significantly lower at low SNR (Fig. 12). However, these results are a bit suspect because the SD(^ rhy ) obtained with the four-lag estimator does not cross that obtained with the two-lag estimator as we expected; this crossover can be seen more clearly in the simulated data of Fig. 13. The reason for the crossover is that the four-lag estimator has better performance than the two-lag estimator at small spectrum width (i.e., when the spectrum width is small more lags are within the correlation time), while the two-lag estimator has better performance than the four-lag estimator at large spectrum width (i.e., the lags outside the
VOLUME 29
FIG. 10. Comparison of theory and simulation for the standard deviation of differential reflectivity estimates. SNR 5 5 dB, rhy 5 0.97 and ZDR 5 1 dB.
correlation time can cause estimation errors). The lack of a crossover point in the theoretical curves might be due to truncating terms in the Taylor expansion to second order. Higher-order expansion is needed to obtain more accurate results. In Fig. 13, simulation and theory for the standard derivations of different estimates are compared, and the simulation does show that the SD(^ rhy ) obtained with the four-lag estimator does exhibit a crossover with results obtained with the two-lag estimator. The SD(^ rhy ) result obtained from simulation with the two-lag estimator is worse than the SD(^ rhy ) obtained with the one-lag estimator (Fig. 13). This small discrepancy, which increases at lower SNR (Fig. 12), could also be related to limitations of the perturbation method. For example, the estimated value cannot be far away from true value, and it is likely that SNR cannot be too small (e.g., less than 0 dB) in order to satisfy this condition. In Fig. 14, a dual-polarization radar simulation verifies that the multilag estimator outperforms other estimators for rhy. The input of this radar simulation is from the Advanced Regional Prediction System (ARPS) model (Lei 2009; Lei et al. 2009b; Xue et al. 2000, 2001). The prognostic state variables include the three wind components, potential temperature, pressure, turbulent kinetic energy, mixing ratios for water vapor, rainwater, cloud water, cloud ice, snow, and hail. The ground truth of rhy, SNR, and spectrum width are calculated from ARPS and shown in the right column (Figs. 14d–f). The model-generated ground truth correlation coefficient is larger than 0.96 in most places, even though SNR is as small as 0 dB. Estimates made with the four-lag estimator are closer to the ground truth than those made with the one-lag and conventional estimators. For this
JUNE 2012
LEI ET AL.
783
FIG. 12. Standard deviation of the copolar cross-correlation coefficient estimates; rhy 5 0.97 and ZDR 5 1 dB.
FIG. 11. Bias of copolar cross-correlation coefficient estimates resulting from (a) estimation error and (b) equal measured H and V noise biases at SNR 5 5 dB; rhy50.97 and ZDR 5 1 dB.
than 4 m s21); the median value of spectrum widths for all types of weather other than squall lines is less than 2.5 m s21 (Fang et al. 2004). The multilag estimator produces meteorological parameter estimates with smaller bias and standard deviation than conventional estimators when the spectrum width is small and echoes have low SNR. The multilag estimators are also immune from the bias error in noise estimation. Performances for each of the four variables (i.e., signal power, spectrum width, differential reflectivity, and copolar correlation coefficient) are listed separately in terms of 1) bias reduction and 2) standard deviation reduction.
simulation, noise power is not calculated because there is not a standard procedure for noise correction to all radar parameters, and thus the noise power is not subtracted to estimate the signal power needed in Eq. (9) for the conventional estimate of rhy (Fig. 14a). The omission of noise power from the calculation of signal power causes significant bias when compared with the ground truth (Fig. 14d); the bias is largest where SNR (Fig. 14e) is weakest. Figure 14c shows that ‘estimation improves significantly when the four-lag estimator is used.
4. Summary A multilag estimator has been developed to improve the 10-cm wavelength radar estimation of polarimetric parameters using weak signals principally from clouds/ light precipitation where turbulence and shear are weak (i.e., SNR is less than 5 dB, and spectrum widths are less
FIG. 13. Comparison of theory and simulations for the copolar correlation coefficient estimates. SNR 5 5 dB, rhy 5 0.97, and ZDR 5 1 dB.
784
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
FIG. 14. Comparison of copolar correlation estimates made using (a) the conventional, (b) one-lag, and (c) four-lag estimators. Also plotted are (d) ground truth of rhy, (e) SNR, and (f) sy.
d
d
d
Signal power: Signal power bias decreases when using the multilag estimators (Fig. 2). The standard deviation of signal power shows no significant decrease when using multilag estimators (Fig. 3). Spectrum width: Spectrum width bias decreases significantly using the four-lag estimator (Fig. 5). If H and V noise power estimates have equal biases (e.g., 20.5 and 21 dB), then the improvement in the spectrum width bias by using a multilag estimator versus using the conventional one can vary from 0.5 to 3 m s21 (Fig. 5b); the improvement is even larger if the H and V noise power biases are not equal. The standard deviation of spectrum width decreases significantly when using multilag estimators (Fig. 6). Differential reflectivity: Differential reflectivity bias decreases when using multilag estimators (Fig. 8). If H and V noise power estimates have equal biases (e.g., 20.5 and 21 dB), then the improvement of differential reflectivity bias by using a multilag estimator versus using the conventional one is around 0.035 and 0.06 dB, respectively (Fig. 8b); the improvement is even larger if the H and V noise power biases are not equal. The standard deviation of differential reflectivity shows no significant decrease when using multilag estimators (Fig. 9).
d
Copolar correlation coefficient: Copolar correlation coefficient bias decreases using the multilag estimators (Fig. 11). If the H and V noise power estimates have biases (e.g., 20.5 and 21 dB), then the improvement of copolar correlation coefficient bias by using a multilag estimator versus using the conventional one is around 0.03 and 0.06, respectively (Fig. 11b); the improvement is even larger if the H and V noise power biases are not equal. The standard deviation of copolar correlation coefficient decreases significantly when using multilag estimators (Fig. 12).
In this study, equal weights have been used in Gaussian fitting to derive the multilag estimator. In practice, weights should be dependent on the merit and reliability of the estimates. It is expected that adjustable weights and an adjustable number of lags would produce even more accurate estimates. Hence, an adaptive multilag estimator should be developed, which can automatically choose how many lags to use and how to weight them in Eq. (4) according to initial estimates of the spectrum width and SNR. With adaptive processing, the N-lag estimator that gives the best performance can be chosen. The adaptive multilag estimator is under development
JUNE 2012
785
LEI ET AL.
APPENDIX A
and being tested with the University of Oklahoma Polarimetric Radar for Innovations in Meteorology and Engineering (OU-PRIME) data (Cao et al. 2010). In addition, to improve the data quality non-Gaussian spectra will be studied in future.
a. Three-lag estimator
Acknowledgments. The work was supported by NOAA Grant NA08OAR4320904 and NSF Grant AGS-1046171. We also thank the anonymous reviewers for their careful reviews and suggestive comments.
The three-lag estimator uses lags 1, 2, and 3 of ACF and lags 0, 61, 62, and 63 of CCF to estimate radar parameters. This is a case for N 5 3 in the general formulas of Eqs. (6)–(9). Substituting N 5 3 into the general expressions, the formulas to estimate spectrum width, differential reflectivity, and correlation coefficient are given by
(3) S^h,y 5
jR^h,y (Ts )j6/7 jR^h,y (2Ts )j3/7 , jR^ (3T )j2/7 h,y
(3 )
s ^ h,y 5
Three-Lag and Four-Lag Estimators
l 28pTs
(A1)
s
qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 11 3 lnjR^h,y (Ts )j 1 2 3 lnjR^h,y (2Ts )j 2 13 3 lnjR^h,y (3Ts )j,
jR^ (T )j6/7 jR^ (2Ts )j3/7 jR^y (3Ts )j2/7 (3) Z^DR 5 10 log10 h s 2/7 h , jR^h (3Ts )j jR^y (Ts )j6/7 jR^y (2Ts )j3/7 [jR^h (3Ts )jjR^y (3Ts )j]1/7 (3) (3) r^hy 5 jC^hy (0)j . [jR^h (Ts )jjR^y (Ts )j]3/7 [jR^h (2Ts )jjR^y (2Ts )j]3/14 In Eq. (A3), the exponential factors are 6/ 7, 3/ 7, and 2/ 7, which are much smaller than the exponential factor of the two-lag estimator, which is as large as 4/ 3. We expect that the statistical performance of the threelag estimator is better than the two-lag estimator and
(4 ) S^h,y 5
(4 )
s ^ h,y
(A3)
(A4)
conventional estimator at low SNR and narrow spectrum width.
b. Four-lag estimator The four-lag estimators for the spectral moments and PRD are
jR^h,y (Ts )j54/86 jR^h,y (2Ts )j39/86 jjR^h,y (3Ts )j14/86 jR^ (4T )j21/86 h,y
(A2)
(A5)
s
vffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi u ^ ^ u l t 13 3 lnjRh,y (Ts )j 1 7 3 lnjRh,y (2Ts )j , 5 pffiffiffiffiffiffiffiffi ^ 4 129pTs 23 3 lnjRh,y (3Ts )j 2 17 3 lnjR^h,y (4Ts )j
" # jR^h (Ts )j54/86 jR^h (2Ts )j39/86 jR^h (3Ts )j14/86 jR^y (4Ts )j21/86 (4 ) ^ , ZDR 5 10 3 log10 jR^h (4Ts )j21/86 jR^y (Ts )j54/86 jR^y (2Ts )j39/86 jR^y (3Ts )j14/86 [jR^h (4Ts )jjR^y (4Ts )j]21/172 (4 ) (4 ) r^hy 5 jC^hy (0)j . [jR^h (Ts )jjR^y (Ts )j]27/86 [jR^h (2Ts )jjR^y (2Ts )j]39/172 [jR^h (3Ts )jjR^y (3Ts )j]7/86
(A6)
(A7)
(A8)
786
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
In Eqs. (A5), (A7), and (A8), the exponents are smaller than the corresponding exponents for the twoand three-lag estimators. Therefore, we expect that the performance of the four-lag estimator is better at low SNR and narrow spectrum width (i.e., at long correlation time). The higher-order multilag (e.g., four lag) estimators give better estimates at narrow spectrum widths. However, the maximum number of lags to be used is determined by the correlation time (inverse proportional to spectrum width). If the lag is larger than the correlation time, using this lag will bring error because the noise contaminates ACFs and CCFs at larger lags. There is a limit as to how many lags to use in multilag estimators; this limit depends on the number of pulses, the SNR, and the correlation time. For example,
(2)
dS^h ’
the absolute value of ACF is used to weight the multilag polypulse-pair estimate for Doppler velocity (Lee 1978; May et al. 1989). The number of lags to be used can be determined according to measured radar parameters.
APPENDIX B Bias and Variance of the Spectral Moment Estimates (Sh and sy) a. Power (two-lag estimator) Using Eqs. (19) and (20) to calculate bias and variance for the two-lag power estimates, and keeping terms to second order, we obtain
(2 ) (2 ) (2) ›S^h ›S^h 1 ›2 S^h [djR^h (Ts )j]2 djR^h (Ts )j 1 djR^h (2Ts )j 1 2 ›jR^ (T )j2 ›jR^h (Ts )j ›jR^h (2Ts )j h s (2)
1
’
VOLUME 29
(2)
1 ›2 S^h ›2 S^h ^ (2T )j]2 1 [dj R [djR^h (Ts )j][djR^h (2Ts )j] h s ^ (T )j›jR^ (2T )j 2 ›jR^ (2T )j2 ›j R h s h s h s
Sh r(Ts )4/3 4 1 2 2 2 2 4 A 2 B 1 A B AB , 1 2 3 9 9 9 r(2T )1/3 3
(B1)
s
(2) with a similar expression for dS^y , and where
A5
B5
djR^h (Ts )j jR^h (Ts )j 2 jRh (Ts )j 5 , jRh (Ts )j jRh (Ts )j djR^h (2Ts )j jR^ (2Ts )j 2 jRh (2Ts )j 5 h . jRh (2Ts )j jRh (2Ts )j
We shall use the following equation: ^ 5 jC 1 DCj 5 (C 1 DC)1/2 (C* 1 DC*)1/2 jCj ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi s ffisffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi DC DC* 5 jCj 11 11 , C C* given by Zhang et al. [2003, Eq. (A4)], in which C is Rh(Ts) [or Chy(Ts) in appendix C when we calculate rhy ]. To calculate A and B the approximation pffiffiffiffiffiffiffiffiffiffiffiffiffi 1 1 1 1 x ’ 1 1 x 2 x2 2 8 is used to obtain
2 2 2 ^ ’ jCj 1 1 DC 1 DC* 1 jDCj 2 DC 2 DC* jCj 2 2 2 2C 2C* 8C 8C* 4jCj
0
!
^ 2 jCj DC jCj DC* jDCj2 DC2 DC*2 ’ 1 1 2 2 , jCj 2C 2C* 8C2 8C*2 4jCj2 (B2)
where, DC 5 C^ 2 C. Therefore, DR^h (Ts ) DR^h*(Ts ) jDR^h (Ts )j2 1 1 A’ 2Rh (Ts ) 2Rh*(Ts ) 4jRh (Ts )j2 2
B’
[DR^h (Ts )]2 8[Rh (Ts )]2
2
[DR^h* (Ts )]2 8[Rh* (Ts )]2
DR^h (2Ts ) DR^h*(2Ts ) jDR^h (2Ts )j2 1 1 2Rh (2Ts ) 2Rh*(2Ts ) 4jRh (2Ts )j2 2
[DR^h (2Ts )]2 8[Rh (2Ts )]2
2
[DR^h*(2Ts )]2 8[Rh* (2Ts )]2
Bias is the expectation of dSh, hence,
.
(B3)
JUNE 2012
787
LEI ET AL. (2 )
(2 )
r(Ts )4/3 Bias[S^h ] hdS^h i 5 ’ Sh Sh r(2Ts )1/3
4 1 2 2 4 A 2 B 1 A2 1 B2 2 AB 3 3 9 9 9
;
(B4)
variance of the Sh estimates is given by 2
(2 ) (2 ) r(Ts )8/3 r(Ts )8/3 Var[S^h ] h[dS^h ]2 i 4 1 2 2 2 2 4 16 2 1 2 8 A 2 B 1 A B AB A B AB . ’ 5 ’ 1 2 1 2 3 9 9 9 9 9 9 S2h S2h r(2Ts )2/3 3 r(2Ts )2/3 (B5)
To arrive at the above expressions, we have retained terms to second order in A and B, whereas analytical formulas for the expected value of other terms in the expressions for A and B are in list B at the end of this appendix.
(4 )
dS^h ’
b. Power (four-lag estimator) A similar approach is used to calculate the bias and variance of the four-lag estimates. Because the estimate bias and variance for the three-lag estimator is between that for the two- and four-lag estimators, we omit results obtained for the three-lag estimator. Following is the statistical analysis result for the four-lag estimator:
Sh r(Ts )54/86 r(2Ts )39/86 r(3Ts )14/86 54 39 14 21 1 54 32 2 1 39 47 2 1 14 72 2 A1 B1 C2 D2 A 2 B 2 C 86 86 86 86 2 86 86 2 86 86 2 86 86 r(4Ts )21/86 1 21 107 2 54 39 54 14 54 21 39 14 39 21 14 21 D 1 AB 1 AC 2 AD 1 BC 2 BD 2 CD , 1 2 86 86 86 86 86 86 86 86 86 86 86 86 86 86
(B6)
where A and B can be found in Eq. (B3). Here, C and D are given by djR^h (3Ts )j jR^h (3Ts )j 2 jRh (3Ts )j 5 jRh (3Ts )j jRh (3Ts )j DR^h (3Ts ) DR^h*(3Ts ) jDR^h (3Ts )j2 [DR^h (3Ts )]2 [DR^h* (3Ts )]2 1 1 2 2 ’ 2Rh (3Ts ) 2Rh*(3Ts ) 8[Rh (3Ts )]2 8[Rh* (3Ts )]2 4jRh (3Ts )j2 djR^h (4Ts )j jR^h (4Ts )j 2 jRh (4Ts )j 5 D5 jRh (3Ts )j jRh (4Ts )j ^ ^ DRh (4Ts ) DRh*(4Ts ) jDR^h (4Ts )j2 [DR^h (4Ts )]2 [DR^h* (4Ts )]2 1 1 2 2 . ’ 2 2 2Rh (4Ts ) 2Rh*(4Ts ) 8[Rh (4Ts )] 8[Rh* (4Ts )]2 4jRh (4Ts )j C5
(B7)
Bias is the expectation of dS^h ; hence, (4 ) (4 ) Bias[S^h ] h[dS^h ]i 5 Sh Sh
’
r(Ts )54/86 r(2Ts )39/86 r(3Ts )14/86
54 39 14 21 1 54 32 2 1 39 47 2 1 14 72 2 A1 B1 C2 D2 A 2 B 2 C 86 86 86 86 2 86 86 2 86 86 2 86 86 r(4Ts )
1 21 107 2 54 39 54 14 54 21 14 39 21 39 14 21 D 1 AB 1 AC 2 AD 1 BC 2 BD 2 CD . (B8) 1 2 86 86 86 86 86 86 86 86 86 86 86 86 86 86 21/86
Variance is " #2
(4 ) (4) r(Ts )54/86 r(2Ts )39/86 r(3Ts )14/86 Var[S^h ] h[dS^h ]2 i 54 39 14 21 2 5 ’ . A1 B1 C2 D 86 86 86 86 S2h S2h r(4Ts )21/86
(B9)
788
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
The expressions of A, B, C, and D are shown in Eqs. (B3) and (B7), in which the subterms can be looked up in list B.
c. Spectrum width (two-lag estimator) For the two-lag estimator bias and variance are calculated by using Eqs. (15)–(20) and procedures are similar to signal power to obtain jR (T )j 21/2 l 1 1 1 (2 ) d^ sh ’ pffiffiffiffiffi A 2 B 2 A 2 1 B2 ln h s 2 2 24pTs 2 jRh (2Ts )j 23/2 jR (T )j l 1 1 1 ln h s 2 A2 2 B2 1 AB , 1 pffiffiffiffiffi 2 2 24pTs 4 jRh (2Ts )j
(B10)
where A and B are given by Eq. (B3). Hence, we find that the two-lag spectrum width bias is (2)
(2 )
sh ]i Bias[^ sh ] 5 h[d^
jRh (Ts )j 21/2 l 1 1 2 1 2 ln A2B2 A 1 B ’ pffiffiffiffiffi 2 2 24pTs 2 jRh (2Ts )j 23/2
jR (T )j l 1 1 2 1 2 ln h s 1 pffiffiffiffiffi 2 A 2 B 1 AB , 2 2 24pTs 4 jRh (2Ts )j
(B11)
and the two-lag variance of the spectrum width estimates is l Var[^ sh ] 5 h[d^ sh ] i ’ pffiffiffiffiffi 24pTs (2 )
(2 ) 2
!2 jRh (Ts )j 21 1 ln h(A2 1 B2 2 2AB)i. 4 jRh (2Ts )j
(B12)
d. Spectrum width (four-lag estimator) Likewise, following similar procedures used to obtain the perturbations for the two-lag estimator, we can calculate the perturbations for the four-lag estimator
(4)
d^ sh
" # 13 7 21/2 l 1 jRh (Ts )j jRh (2Ts )j 13 2 7 2 3 2 17 2 ln ’ pffiffiffiffiffiffiffiffi 13A 1 7B 2 3C 2 17D 2 A 2 B 1 C 1 D 2 2 2 2 4 129pTs 2 Rh (3Ts )3 Rh (4Ts )17 " #23/2 13 7 l 1 jRh (Ts )j jRh (2Ts )j 132 2 72 32 172 2 ln 2 1 pffiffiffiffiffiffiffiffi A 2 B2 2 C 2 2 D 2 13 3 7 3 AB 3 17 2 2 2 2 4 129pTs 4 Rh (3Ts ) Rh (4Ts ) ! 1 13 3 3 3 AC 1 13 3 17 3 AD 1 7 3 3 3 BC 1 7 3 17 3 BD 2 3 3 17 3 CD , (4 )
(B13)
(4 )
Bias[^ sh ] 5 h[d^ sh ]i
" #21/2
jRh (Ts )j13 jRh (2Ts )j7 l 1 13 2 7 2 3 2 17 2 ’ pffiffiffiffiffiffiffiffi 13A 1 7B 2 3C 2 17D 2 2 1 1 ln A B C D 2 2 2 2 4 129pTs 2 jRh (3Ts )j3 jRh (4Ts )j17 " #23/2* jR (T )j13 jR (2Ts )j7 l 1 132 2 72 2 32 2 172 2 ln h s 3 h 2 A B C D 2 13 3 7 3 AB 1 pffiffiffiffiffiffiffiffi 2 2 2 2 2 2 2 4 129pTs 4 jRh (3Ts )j jRh (4Ts )j17 !+ 1 13 3 3 3 AC 1 13 3 17 3 AD 1 7 3 3 3 BC 1 7 3 17 3 BD 2 3 3 17 3 CD
,
(B14)
JUNE 2012 (4)
789
LEI ET AL. (4 )
Var[^ sh ] 5 h[d^ sh ]2 i l ’ pffiffiffiffiffiffiffiffi 8 129pTs
!" 2
jRh (Ts )j13 jjRh (2Ts )j7
#21
ln jRh (3Ts )j3 jRh (4Ts )j17
(169A2 1 49B2 1 9C2 1 289D2 1 182AB 2 78AC
2 442AD 2 42BC 2 238BD 1 102CD),
where A and B are given by Eq. (B3), and C and D are given by Eq. (B7).
e. List B [Eqs. (B16)–(B19)] Note that in calculating the expected values of A, B, etc., we can use the property that DC and DC* have hD2 R^h (mTs )i [Rh (mTs )]2 hjDR^h (mTs )j2 i 2
jRh (mTs )j
5
’
hD2 R^h*(mTs )i [Rh* (mTs )]2
’
(B15)
zero expected value. The analytical formulas for other terms in Eqs. (B3) and (B7) are given by the equations in this list. For the terms having the ACF computed from H polarized signals (the forms are exactly the same for the ACF computed from V polarized signals):
1 2 1 1 r2 (mTs ) , MI M (SNR)h
(B16)
2(SNR)h 1 1 1 1 , r2 (mTs )MI Mr2 (mTs )(SNR)2h
(B17)
# " hDR^h (m1 Ts )DR^h*(m2 Ts )i hDR^h*(m1 Ts )DR^h (m2 Ts )i (m1 1 m2 )2 Ts2 r[(m2 2 m1 )Ts ] 1 2 5 ’ , 1 exp MI M(SNR)h r(m1 Ts )r(m2 Ts ) Rh (m1 Ts )Rh*(m2 Ts ) Rh*(m1 Ts )Rh (m2 Ts ) 4t 2c (B18) # " hDR^h (m1 Ts )DR^h (m2 Ts )i hDR^h*(m1 Ts )DR^h*(m2 Ts )i (m2 2 m1 )2 Ts2 r[(m2 1 m1 )Ts ] 1 2 1 5 ’ . exp MI M(SNR)h r(m1 Ts )r(m2 Ts ) Rh (m1 Ts )Rh (m2 Ts ) Rh*(m1 Ts )Rh*(m2 Ts ) 4t 2c (B19)
The solutions for the V channel signals have the same form. To demonstrate the derivation procedure used to obtain the four approximate analytical formulas, we provide, in what follows, an example derivation for one [i.e., (B18)] of these formulas. Equation (B18) is
hDR^h (m1 Ts )DR^h*(m2 Ts )i , Rh (m1 Ts )Rh*(m2 Ts ) where m1 and m2 are positive integers. The expected value of the numerator can be expressed as
hDR^h (m1 Ts )DR^h* (m2 Ts )i 5 h[R^h (m1 Ts ) 2 Rh (m1 Ts )][R^h* (m2 Ts ) 2 Rh* (m2 Ts )]i 5 hR^h (m1 Ts )R^h* (m2 Ts )i 2 Rh (m1 Ts )Rh* (m2 Ts ) 8 * +9 M2m1 M2m2 [ES ( p) 1 EN ( p)][ES* ( p 1 m ) 1 EN * ( p 1 m )] = < 1 h h h h 1 1 1 2 Rh (m1 Ts )Rh* (m2 Ts ) 5 å å :M 2 m1 M 2 m2 p51 p951 [ES* ( p9) 1 EN *( p9)][ES ( p9 1 m ) 1 EN ( p9 1 m )] ; h h h h 2 2 5
1 1 M 2 m1 M 2 m2
M2m1 M2m2
å
p51
å
p951
hhESh ( p)ESh* ( p 1 m1 )ESh* ( p9)ESh ( p9 1 m2 )i
S* N* S S N* S* N 1hEN h ( p)Eh ( p 1 m1)Eh ( p9)Eh ( p91 m2 )i 1 hEh ( p)Eh ( p 1 m1 )Eh ( p9)Eh ( p91 m2 )ii
2 Rh(m1Ts)Rh* (m2Ts )
790
5
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
1 1 M 2 m1 M 2 m2
M2m1 M2m2
å
p51
å
p951
VOLUME 29
[hESh ( p)ESh* ( p9)ihESh* ( p 1 m1 )ESh ( p9 1 m2 )i
N* S S* N* N 1 hESh ( p9 1 m2 )ESh* ( p 1 m1 )ihEN h (p)Eh ( p9)i 1 hEh ( p)Eh ( p9)ihEh ( p 1 m1 )Eh ( p9 1 m2 )i]
5
1 1 M 2 m1 M 2 m2
M2m1 M2m2
å
p51
å
p951
fRh [(p9 2 p)Ts ]Rh [(p 1 m1 2 p9 2 m2 )Ts ]
1 Rh [(p 1 m1 2 p9 2 m2 )Ts ]Nh d( p92p) 1 Rh [(p9 2 p)Ts ]Nh d( p1m 2p92m ) g 1
5
1 1 M 2 m1 M 2 m2
å
p51
å
p951
fRh [(p9 2 p)Ts ]Rh [[(p 1 m1 2 p9 2 m2 )Ts ]g 1
5 Rh (m1 Ts )Rh* (m2 Ts )
1 Td
# (m2 2 m1 )Ts t m1 m2 Ts2 t2 dt, 3 exp 2 2 2 1 tc t 2c t 2c "
For the two-lag estimator bias and variance are calculated by using Eqs. (15)–(20), and procedures are similar to signal power to obtain
(B21) if Ts is much shorter than the correlation time and a dwell time is much larger than the correlation time (Zhang et al. 2003; Appendix A in Zhang et al. 2004). Using the equation ! rffiffiffiffi ð p b2 2 4ac 2 , a.0 exp(2ax 2 bx 2 c) dx 5 exp a 4a Ð we obtain: exp[2(t2 /t 2c ) 2 [f(m2 2 m1 )Ts tg/t 2c ] 1 pffiffiffiffi 2 2 [(m1 m2 Ts )/t c ] dt 5 pt c exp[(m1 1 m2 )2 Ts2 /4t2c ] and pffiffiffiffi MI 5 MTs / pt c . Therefore, (B21) can be written as # " (m1 1 m2 )2 Ts2 1 , and R (m T )R * (m T ) exp 5 MI h 1 s h 2 s 4t 2c # " hDR^h (m1 Ts )DR^h*(m2 Ts )i (m1 1 m2 )2 Ts2 1 ’ exp MI Rh (m1 Ts )Rh*(m2 Ts ) 4t 2c 1
r[(m2 2 m1 )Ts ] 2 . M(SNR)h r(m1 Ts)r(m2 Ts)
2Nh Rh [(m1 2 m2 )Ts ] . (B20) M 2 m1
10 4 1 (2 ) lnjR^h (Ts )j 1 lnjR^y (2Ts )j Z^DR 5 ln10 3 3 1 4 2 lnjR^h (2Ts )j 2 lnjR^y (Ts )j . (C1) 3 3
It can be shown that the first term in the above equation can be approximated by ð 1 Rh (t)Rh* [t 1 (m2 2 m1 )Ts ] dt ’ Td
ð
2
M2m1 M2m2
(2)
dZ^DR ’
10 4 1 1 4 A1 B2 C2 D ln10 3 3 3 3
2 1 1 2 2 A2 2 B2 1 C2 1 D2 , 3 6 6 3
where
A5
djR^h (Ts )j , jRh (Ts )j
C5
djR^h (2Ts )j , jRh (2Ts )j
A’
DR^h (Ts ) DR^h*(Ts ) jDR^h (Ts )j2 1 1 2Rh (Ts ) 2Rh*(Ts ) 4jRh (Ts )j2 2
B5
8[Rh (Ts )]
djR^y (2Ts )j , jRy (2Ts )j
and D 5
[DR^h (Ts )]2 2
2
(C2)
djR^y (Ts )j jRy (Ts )j
[DR^h* (Ts )]2 8[Rh* (Ts )]2
.
Similar approximations apply to B, C, and D. Hence, the bias for the two-lag estimator is
APPENDIX C (2 )
Bias and Variance of the Polarimetric Parameter Estimates (ZDR and r hy) a. Differential reflectivity (two-lag estimator) The two-lag estimator of differential reflectivity is given as
(2 )
Bias[Z^DR ] 5 h[dZ^DR ]i 10 4 1 1 4 ’ A1 B2 C2 D ln10 3 3 3 3
2 2 1 2 1 2 2 2 , 2 A 2 B 1 C 1 D 3 6 6 3
(C3)
JUNE 2012
791
LEI ET AL.
(2 )
(2 )
Var[Z^DR ] 5 h[dZ^DR ]2 i + 10 2 4 1 1 4 2 2 1 2 1 2 2 2 2 A1 B2 C2 D2 A 2 B 1 C 1 D ’ ln10 3 3 3 3 3 6 6 3 2
10 16 2 1 2 1 2 16 2 8 8 32 2 8 8 A 1 B 1 C 1 D 1 AB 2 AC 2 AD 2 BC 2 BD 1 CD . (C4) ’ ln10 9 9 9 9 9 9 9 9 9 9
To arrive at the above expressions, we have retained terms to second order in A, B, C, and D, whereas analytical formulas for the expected value of other terms in the expressions for A–D are in lists B and C.
bias and variance for the three-lag estimator is between that for the two-lag estimator and four-lag estimators, we omit results obtained for the three-lag estimator. Following is the statistical analysis result for the four-lag estimator:
b. Differential reflectivity (four-lag estimator) A similar approach is used to calculate the bias and variance of the four-lag estimates. Because the estimate
(4 )
dZ^DR ’
10 54 39 14 21 54 39 14 21 A1 B1 C2 D2 E2 F 2 G1 H ln10 86 86 86 86 86 86 86 86 27 2 39 2 7 2 21 2 27 2 39 2 7 2 21 2 B 2 C 1 D 1 E 1 F 1 G 2 H , 2 A 2 86 172 86 172 86 172 86 172
(C5)
where djR^h (Ts )j , jRh (Ts )j djR^y (Ts )j E5 , jRy (Ts )j
A5
djR^h (2Ts )j , jRh (2Ts )j djR^y (2Ts )j F 5 , jRy (2Ts )j
djR^h (3Ts )j , jRh (3Ts )j djR^y (3Ts )j G5 , jjRy (3Ts )j
B5
C5
djR^h (4Ts )j jRh (4Ts )j djR^y (4Ts )j H5 . jRy (4Ts )j D5
Approximation is shown as A’
DR^h (Ts ) DR^h*(Ts ) jDR^h (Ts )j2 [DR^h (Ts )]2 [DR^h* (Ts )]2 1 1 2 2 ; 2Rh (Ts ) 2Rh*(Ts ) 8[(Rh (Ts )]2 8[Rh* (Ts )]2 4jRh (Ts )j2
similar approximations apply to B, C, D, E, F, G, and H as (4 )
(4 )
Bias[Z^DR ] 5 h[dZ^DR ]i 10 54 39 14 21 54 39 14 21 A 1 B 1 C 2 D 2 E 2 F2 G 1 H ’ ln10 86 86 86 86 86 86 86 86
27 2 39 2 7 2 21 2 27 2 39 2 7 2 21 2 2 A 2 B 2 C 1 D 1 E 1 F 1 G 2 H , 86 172 86 172 86 172 86 172 2 2 + 10 54 39 14 21 54 39 14 21 (4 ) (4 ) 2 A1 B1 C2 D2 E2 F2 G1 H . Var[Z^DR ] 5 h[dZ^DR ] i ’ ln10 86 86 86 86 86 86 86 86
(C6)
(C7)
c. Correlation coefficient (two-lag estimator) The two-lag estimator of correlation coefficient is given by (2 ) r^hy 5 jC^hy (22Ts )j23/35 jC^hy (2Ts )j12/35 jC^hy (0)j17/35 jC^hy (Ts )j12/35 jC^hy (2Ts )j23/35 3 jR^ (2T )j1/6 jR^ (2T )j1/6 jR^ (T )j22/3 jR^ (T )j22/3 .
h
s
y
s
h
s
y
s
(C8)
792
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
VOLUME 29
For the two-lag estimator bias and variance are calculated by using Eqs. (15)–(20) and procedures are similar to the signal power to obtain 3 12 17 12 3 1 1 2 2 (2) d^ rhy ’ rhy r(2Ts )[(1/3)2(6/35)] r(Ts )[(24/35)2(4/3)] 2 A 1 B 1 C 1 D 2 E 1 F 1 G 2 H 2 I 35 35 35 35 35 6 6 3 3 1 3 3 1 12 12 1 17 17 1 12 12 2 2 1 A2 1 2 1 B2 1 2 1 C2 1 2 1 D2 1 2 2 35 35 2 35 35 2 35 35 2 35 35 1 3 3 1 1 1 1 1 1 1 2 2 2 2 2 2 21 E 1 21 F 1 21 G 1 2 2 2 1 H2 1 2 2 35 35 2 6 6 2 6 6 2 3 3 1 2 2 3 12 3 17 3 12 2 2 2 1 I2 1 2 AB 1 2 AC 1 2 AD 1 2 3 3 35 35 35 35 35 35 3 3 3 1 3 1 3 2 3 2 1 2 2 AE 1 2 AF 1 2 AG 1 2 2 AH 1 2 2 AI 35 35 35 6 35 6 35 3 35 3 12 17 12 12 12 3 12 1 12 1 BC 1 BD 1 2 BE 1 BF 1 BG 1 35 35 35 35 35 35 35 6 35 6 12 2 12 2 17 12 17 3 17 1 2 BH 1 2 BI 1 CD 1 2 CE 1 CF 1 35 3 35 3 35 35 35 35 35 6 17 1 17 2 17 2 12 3 12 1 CG 1 2 CH 1 2 CI 1 2 DE 1 DF 1 35 6 35 3 35 3 35 35 35 6 12 1 12 2 12 2 3 1 3 1 DG 1 2 DH 1 2 DI 1 2 EF 1 2 EG 1 35 6 35 3 35 3 35 6 35 6 3 2 3 2 1 1 1 2 1 2 1 2 2 EH 1 2 2 EI 1 FG 1 2 FH 1 2 FI 35 3 35 3 6 6 6 3 6 3 1 2 1 2 2 2 2 GH 1 2 GI 1 2 2 HI , 1 6 3 6 3 3 3
(C9)
where A5
djC^hy (22Ts )j , jChy (22Ts )j
F 5
djR^h (2Ts )j , jRh (2Ts )j
A’
* (22T )]2 * (22Ts ) DC^hy (22Ts ) DC^hy jDC^hy (22Ts )j2 [DC^hy (22Ts )]2 [DC^hy s 1 1 2 2 . 2Chy (22Ts ) 2Chy * (22Ts ) * (22Ts )]2 8[Chy (22Ts )]2 8[Chy 4jChy (22Ts )j2
B5
G5
djC^hy (2Ts )j , jChy (2Ts )j
djR^y (2Ts )j , jRy (2Ts )j
C5
H5
djR^h (Ts )j , jRh (Ts )j
Similar approximations apply to B, C, D, . . . , to I. Hence, we find that the two-lag correlation coefficient bias is (2 )
(2 )
rhy ]i, Bias[^ rhy ] 5 h[d^
djC^hy (0)j , jChy (0)j
D5
I5
djC^hy (Ts )j , jChy (Ts )j
E5
djC^hy (2Ts )j jChy (2Ts )j
djR^y (Ts )j jRy (Ts )j
Because of the lengthy derivations, calculations for the performance of the four-lag estimator are not shown here. However, they follow procedures used for the twolag estimator.
(C10)
d. List C [Eqs. (C12)–(C21)]
and the variance is (2 )
(2)
rhy ]2 i. Var[^ rhy ] 5 h[d^
(C11)
The analytical formulas for the terms in A, B, . . . , I are given by the equations in this list, as well as equations in list B:
JUNE 2012
LEI ET AL.
793
# " hDR^h (m1 Ts )DR^y (m2 Ts )i hDR^h*(m1 Ts )DR^y*(m2 Ts )i r2hy (m2 2 m1 )2 Ts2 , 5 ’ exp Rh (m1 Ts )Ry (m2 Ts ) Rh*(m1 Ts )Ry*(m2 Ts ) MI 4t 2c
(C12)
# " hDR^h*(m1 Ts )DR^y (m2 Ts )i hDR^h (m1 Ts )DR^y*(m2 Ts )i r2hy (m1 1 m2 )2 Ts2 , 5 ’ exp Rh*(m1 Ts )Ry (m2 Ts ) Rh (m1 Ts )Ry*(m2 Ts ) MI 4t 2c
(C13)
* (mTs )i hD2 C^hy (mTs )i hD2 C^hy 1 ’ , 5 2 2 M Chy (mTs ) * (mTs )] [Chy I hjDC^hy (mTs )j2 i jChy (mTs )j2
’
(C14)
(SNR)h 1 (SNR)y 1 1 1 1 , MI r2hy r2 (mTs ) M(SNR)h (SNR)y r2hy r2 (mTs )
* (m1 Ts )DC^hy (m2 Ts )i * (m2 Ts )i hDC^hy hDC^hy (m1 Ts )DC^hy 5 Chy Chy (m1 Ts )Chy * (m1 Ts )Chy (m2 Ts ) * (m2 Ts ) " # (m1 1 m2 )2 Ts2 1 1 r[(m2 2 m1 )Ts ] 1 1 ’ exp 1 , 1 (SNR)y 4t 2c MI r2hy Mr2hy r(m1 T)r(m2 Ts ) (SNR)h # " * (m1 Ts )DC^hy * (m2 Ts )i hDC^hy (m1 Ts )DC^hy (m2 Ts )i hDC^hy (m1 2 m2 )2 Ts2 1 , 5 ’ exp Chy (m1 Ts )Chy (m2 Ts ) Chy MI 4t2c * (m1 Ts )Chy * (m2 Ts )
(C15)
(C16)
(C17)
# " * (m1 Ts )DR^h*(m2 Ts )i hDC^hy (m1 Ts )DR^h (m2 Ts )i hDC^hy (m1 2 m2 )2 Ts2 1 1 r[(m2 1 m1 )Ts ] 1 1 5 ’ exp , 2 M M (SNR) Chy (m1 Ts )Rh (m2 Ts ) Chy r(m T )R T ) T )r(m T ) (m 4t c * 1 s h*(m2 s I 1 s 2 s h (C18) # " * (m1 Ts )DR^h (m2 Ts )i hDC^hy (m1 Ts )DR^h*(m2 Ts )i hDC^hy (m1 1 m2 )2 Ts2 1 1 r[(m2 2 m1 )Ts ] 1 1 5 ’ exp , 2 MI M r(m1 Ts )r(m2 Ts ) (SNR)h Chy (m1 Ts )Rh*(m2 Ts ) Chy 4t c * (m1 Ts )Rh (m2 Ts ) (C19) # " * (m1 Ts )DR^y*(m2 Ts )i hDC^hy (m1 Ts )DR^y (m2 Ts )i hDC^hy (m1 2 m2 )2 Ts2 1 1 r[(m2 1 m1 )Ts ] 1 1 5 ’ exp , 2 MI M r(m1 Ts )r(m2 Ts ) (SNR)y Chy (m1 Ts )Ry (m2 Ts ) Chy 4t c * (m1 Ts )Ry*(m2 Ts ) (C20) # " * (m1 Ts )DR^y (m2 Ts )i hDC^hy (m1 Ts )DR^y*(m2 Ts )i hDC^hy (m1 1 m2 )2 Ts2 1 1 r[(m2 2 m1 )Ts ] 1 1 5 ’ exp . MI M r(m1 Ts )r(m2 Ts ) (SNR)y Chy (m1 Ts )Ry*(m2 Ts ) Chv 4t2c * (m1 Ts )Ry (m2 Ts ) (C21)
To demonstrate the derivation procedure used to obtain the 10 approximate analytical formulas, we provide, in
what follows, an example derivation for one [i.e., (C12)] of these formulas. Equation (C12) is
794
JOURNAL OF ATMOSPHERIC AND OCEANIC TECHNOLOGY
hDR^h (m1 Ts )DR^y (m2 Ts )i , Rh (m1 Ts )Ry (m2 Ts )
VOLUME 29
where m1 and m2 are positive integers. The expected value of the numerator can be expressed as
hDR^h (m1 Ts )DR^y (m2 Ts )i 5 h[R^h (m1 Ts ) 2 Rh (m1 Ts )][R^y (m2 Ts ) 2 Ry (m2 Ts )]i 5 hR^h (m1 Ts )R^y (m2 Ts )i 2 Rh (m1 Ts )Ry (m2 Ts ) 8 * +9 M2m1 M2m2 S* (q 1 m ) 1 EN * (q 1 m )] < 1 = [ESh (q) 1 EN (q)][E 1 h h h 1 1 5 2 Rh (m1 Ts )Ry (m2 Ts ) å å S* N* :M 2 m1 M 2 m2 q51 q951 [ESy (q9) 1 EN ; y (q9)][Ey (q9 1 m2 ) 1 Ey (q9 1 m2 )] 5
5
1 1 M 2 m1 M 2 m2
1 1 M 2 m1 M 2 m2
M2m1 M2m2
å
q51
å
q951
M2m1 M2m2
å
q51
å
q951
hESh (q)ESh* (q 1 m1 )ESy (q9)ESy * (q9 1 m2 )i 2 Rh (m1 Ts )Ry (m2 Ts )
hESh (q)ESy * (q9 1 m2 )ihESh* (q 1 m1 )ESy (q9)i
M2m M2m
1 2 1 1 hChy [(q9 1 m2 2 q)Ts ]Chy [(q 1 m1 2 q9)Ts ]i å å M 2 m1 M 2 m2 q51 q951 ð 1 Chy (t 1 m2 Ts )Chy (m1 Ts 2 t) dt ’ Td 2 ð (m1 2 m2 )Ts t 1 t dt. exp 2 2 1 5 Chy (m2 Ts )Chy (m1 Ts ) Td tc t 2c
5
Using the equation ð
Nh,y
! rffiffiffiffi p b2 2 4ac , a.0 exp(2ax 2 bx 2 c) dx 5 exp a 4a 2
we obtain # " (m1 2 m2 )2 Ts2 1 2 5 . R (m T )R (m T )r exp MI h 1 s y 2 s y 4t 2c Therefore, # " hDR^h (m1 Ts )DR^y (m2 Ts )i (m2 2 m1 )2 Ts2 r2hy . ’ exp Rh (m1 Ts )Ry (m2 Ts ) MI 4t 2c
APPENDIX D List of Symbols Chy (mTs ): KDP : M MI N
Cross-correlation function. Specific differential phase Number of signal samples Number of independent samples Number of lags used for multilag estimator
R(mTs ) S SNR PRD Ts Td yr yN Z ZDR rhy r(mTs ) sy tc
White noise power in the horizontal (h) and vertical channels (y) Autocorrelation function Signal power Signal-to-noise ratio Polarimetric radar data (e.g., ZDR, etc.) Pulse repetition time Dwell time Doppler velocity Nyquist velocity Reflectivity factor Differential reflectivity Copolar correlation coefficient Temporal correlation coefficient Spectrum width Correlation time
REFERENCES Brandes, E., G. Zhang, and J. Vivekanandan, 2002: Experiments in rainfall estimation with a polarimetric radar in a subtropical environment. J. Appl. Meteor., 41, 674–685. Bringi, V. N., and V. Chandrasekar, 2001: Polarimetric Doppler Weather Radar: Principles and Applications. Cambridge University Press, 636 pp. Cao, Q., G. Zhang, and R. D. Palmer, 2010: Real-time implementation of enhanced moment estimation on the C-band
JUNE 2012
LEI ET AL.
polarimetric research radar OU-PRIME. ERAD 2010: Sixth European Conf. on Radar Meteorology and Hydrology, Sibiu, Romania, ERAD. [Available online at http://www.erad2010. org/pdf/POSTER/02_Advances/06_ERAD2010_0053_modificat. pdf.] Doviak, R. J., and D. S. Zrnic´, 2006: Doppler Radar and Weather Observations. 2nd ed. Dover Publication, 562 pp. Fang, M., R. J. Doviak, and V. Melnikov, 2004: Spectrum width measured by WSR-88D radar: Error sources and statistics of various weather phenomena. J. Atmos. Oceanic Technol., 21, 888–904. Galati, G., and G. Pavan, 1995: Computer simulation of weather radar signals. Simul. Pract. Theory, 3, 17–44. Hubbert, J. C., V. N. Bringi, and D. Brunkow, 2003: Studies of the polarimetric covariance matrix. Part I: Calibration methodology. J. Atmos. Oceanic Technol., 20, 696–706. Janssen, L. H., and A. Van Der Spek, 1985: The shape of Doppler spectra from precipitation. IEEE Trans. Aerosp. Electron. Syst., 21, 208–219. Lee, R. W., 1978: Performance of the poly-pulse-pair Doppler estimator. Lassen Research Company Tech. Memo. 78-03, 28 pp. Lei, L., 2009: Simulation and processing of polarimetric radar signals based on numerical weather prediction model output. M.S. thesis, Dept. of Electrical and Computer Engineering, University of Oklahoma, 75 pp. ——, G. Zhang, R. D. Palmer, B. L. Cheong, M. Xue, and Q. Cao, 2009a: A multi-lag correlation estimator for polarimetric radar parameters in the presence of noise. Preprints, 34th Conf. on Radar Meteorology, Williamsburg, VA, Amer. Meteor. Soc., 13.25. [Available online at http://arrc.ou.edu/;boonleng/pdf/ lei1109_multilag.pdf.] ——, ——, B. L. Cheong, R. D. Palmer, and M. Xue, 2009b: Simulations of polarimetric radar signals based on numerical weather prediction output. Preprints, 25th Conf. on IIPS, Phoenix, AZ, Amer. Meteor. Soc., 2.8. [Available online at http://ams.confex.com/ams/pdfpapers/145884.pdf.] May, P. T., T. Sato, M. Yamamoto, S. Kato, T. Tsuda, and S. Fukao, 1989: Errors in the determination of wind speed by Doppler radar. J. Atmos. Oceanic Technol., 6, 235–242. Melnikov, V. M., and D. S. Zrnic´, 2004: Simultaneous transmission mode for the polarimetric WSR-88D: Statistical biases and standard deviations of polarimetric variables. NOAA/NSSL Rep., 84 pp. [Available online at http://cimms.ou.edu/;schuur/ jpole/SHV_statistics.pdf.]
795
——, and ——, 2007: Autocorrelation and cross-correlation estimators of polarimetric variables. J. Atmos. Oceanic Technol., 24, 1337–1350. Moisseev, D., C. Nguyen, and V. Chandrasekar, 2008: Clutter suppression for staggered PRT waveforms. J. Atmos. Oceanic Technol., 25, 2209–2218. Park, H. S., A. V. Ryzhkov, D. S. Zrnic´, and K. E. Kim, 2009: The hydrometeor classification algorithm for the polarimetric WSR-88D: Description and application to an MCS. Wea. Forecasting, 24, 730–748. Sachidananda, M., and D. S. Zrnic´, 1985: Z measurement considerations for a fast scan capability radar. Radio Sci., 20, 907–922. ——, and ——, 1986: Differential propagation phase shift and rainfall estimation. Radio Sci., 21, 235–247. Srivastava, R. C., A. R. Jameson, and P. H. Hildebrand, 1979: Time-domain computation of mean and variance of Doppler spectra. J. Appl. Meteor., 18, 189–194. Xue, M., K. K. Droegemeier, and V. Wong, 2000: The Advanced Regional Prediction System (ARPS)—A multi-scale nonhydrostatic atmospheric simulation and prediction tool. Part I: Model dynamics and verification. Meteor. Atmos. Phys., 75, 161–193. ——, and Coauthors, 2001: The Advanced Regional Prediction System (ARPS)—A multi-scale nonhydrostatic atmospheric simulation and prediction tool. Part II: Model physics and applications. Meteor. Atmos. Phys., 76, 143–165. Zhang, G., J. Vivekanandan, and E. Brandes, 2001: A method for estimating rain rate and drop size distribution form polarimetric radar measurements. IEEE Trans. Geosci. Remote Sens., 39, 831–841. ——, R. J. Doviak, J. Vivekanandan, W. O. J. Brown, and S. A. Cohn, 2003: Cross-correlation ratio method to estimate crossbeam wind and comparison with a full correlation analysis. Radio Sci., 38, 8052, doi:10.1029/2002RS002682. ——, ——, ——, ——, and ——, 2004: Performance of correlation estimators for spaced-antenna wind measurement in the presence of noise. Radio Sci., 39, RS3017, doi:10.1029/ 2003RS003022. Zrnic´, D. S., 1975: Simulation of weatherlike Doppler spectra and signals. J. Appl. Meteor., 14, 619–620. ——, 1977: Spectral moment estimates from correlated pulse pairs. IEEE Trans. Aerosp. Electron. Syst., AES-13, 344–354. ——, 1979: Estimation of spectral moments for weather echoes. IEEE Trans. Geosci. Electron., GE-17, 113–128. ——, and A. V. Ryzhkov, 1999: Polarimetry for weather surveillance radars. Bull. Amer. Meteor. Soc., 80, 389–406.