Abstract

We present a sequential algorithm for estimating both concentration dependence on range and time and backscatter coefficient spectral dependence of optically thin localized atmospheric aerosols using data from rapidly tuned lidar. The range dependence of the aerosol is modeled as an expansion of the concentration in an orthonormal basis set whose coefficients carry the time dependence. Two estimators are run in parallel: a Kalman filter for the concentration range and time dependence and a maximum-likelihood estimator for the aerosol backscatter wavelength and time dependence. These two estimators exchange information continuously over the data-processing stream. The state model parameters of the Kalman filter are also estimated sequentially together with the concentration and backscatter. Lidar data collected prior to the aerosol release are used to estimate the ambient lidar return. The approach is illustrated on atmospheric backscatter long-wave infrared (CO2) lidar data.

© 2008 Optical Society of America

Full Article  |  PDF Article

References

  • View by:
  • |
  • |
  • |

  1. The text Elastic Lidar, by V. A. Kovalev and W. E. Eichinger (Wiley, 2004) has an extensive treatment and bibliography of lidar analysis methods for aerosols distributed continuously over extended paths.
    [CrossRef]
  2. We are assuming a single aerosol material here. Generalization to more than one material present at a given time and range is currently being investigated.
  3. A note on units. From Eqs. we see that the physical units of αC must be those of β, 1/(m sr), for example. Since concentration is typically expressed in units of mg/m3, this requires the units of α to be m2/(mg sr). Setting g=1 of course means that these units no longer apply and, consequently, only relative backscatter and concentration can be estimated by this choice. This is described in more detail in Section .
  4. B. R. Frieden, Probability, Statistical Optics, and Data Testing (Springer-Verlag, 1983).
  5. A Gaussian model for Sg is not the only choice, of course. Using a step function with width rg would give Sg(f)=σg2rgsinc⁡(frg). This model has many more ripples than the Gaussian and provides less noise damping for a given rg.
  6. Even with the choice of maximum likelihood, there is the decision of whether the estimates should be computed independently for each time step or be processed by a recursive algorithm, such as recursive least squares. We have adopted the former choice here since the estimates are to be used to train and test classifiers for which independent training samples are required. Actual deployment may be better served by a recursive algorithm.
  7. A related concern is the huge amount of background data needed to construct suitable sample estimates of the range-covariance matrices. This is not compatible with the limited amount of background data typically available.
  8. See S. Mallat, Wavelet Tour of Signal Processing, Section 9.1.3 and Appendix A6, for an excellent discussion of the KL expansion (Academic, 1998).
  9. If our covariance function, given by Eq. , had circular periodicity, Mallat shows that the Fourier basis would be the KL basis. We cannot claim this, but for target aerosols not close to either end of the sampling range, there is little error in making this assumption.
  10. See, for example, F. L. Lewis, Optimal Estimation (Wiley, 1986).
  11. A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).
  12. R. H. Shumway and D. S. Stoffer, Time Series Analysis and Its Applications (Springer-Verlag, 2000), Chap. 4.
  13. S. M. Kay, Fundamentals of Statistical Signal Processing: Estimation Theory (Prentice Hall, 1993) has very readable treatments of ML and related topics.

1977 (1)

A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).

Dempster, A. P.

A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).

Frieden, B. R.

B. R. Frieden, Probability, Statistical Optics, and Data Testing (Springer-Verlag, 1983).

Kay, S. M.

S. M. Kay, Fundamentals of Statistical Signal Processing: Estimation Theory (Prentice Hall, 1993) has very readable treatments of ML and related topics.

Laird, N. M.

A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).

Lewis, F. L.

See, for example, F. L. Lewis, Optimal Estimation (Wiley, 1986).

Mallat, S.

See S. Mallat, Wavelet Tour of Signal Processing, Section 9.1.3 and Appendix A6, for an excellent discussion of the KL expansion (Academic, 1998).

Rubin, D. B.

A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).

Shumway, R. H.

R. H. Shumway and D. S. Stoffer, Time Series Analysis and Its Applications (Springer-Verlag, 2000), Chap. 4.

Stoffer, D. S.

R. H. Shumway and D. S. Stoffer, Time Series Analysis and Its Applications (Springer-Verlag, 2000), Chap. 4.

J. R. Stat. Soc. Ser. B. Methodol. (1)

A. P. Dempster, N. M. Laird, and D. B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm (with discussion),” J. R. Stat. Soc. Ser. B. Methodol. 39, 1-38 (1977).

Other (12)

R. H. Shumway and D. S. Stoffer, Time Series Analysis and Its Applications (Springer-Verlag, 2000), Chap. 4.

S. M. Kay, Fundamentals of Statistical Signal Processing: Estimation Theory (Prentice Hall, 1993) has very readable treatments of ML and related topics.

The text Elastic Lidar, by V. A. Kovalev and W. E. Eichinger (Wiley, 2004) has an extensive treatment and bibliography of lidar analysis methods for aerosols distributed continuously over extended paths.
[CrossRef]

We are assuming a single aerosol material here. Generalization to more than one material present at a given time and range is currently being investigated.

A note on units. From Eqs. we see that the physical units of αC must be those of β, 1/(m sr), for example. Since concentration is typically expressed in units of mg/m3, this requires the units of α to be m2/(mg sr). Setting g=1 of course means that these units no longer apply and, consequently, only relative backscatter and concentration can be estimated by this choice. This is described in more detail in Section .

B. R. Frieden, Probability, Statistical Optics, and Data Testing (Springer-Verlag, 1983).

A Gaussian model for Sg is not the only choice, of course. Using a step function with width rg would give Sg(f)=σg2rgsinc⁡(frg). This model has many more ripples than the Gaussian and provides less noise damping for a given rg.

Even with the choice of maximum likelihood, there is the decision of whether the estimates should be computed independently for each time step or be processed by a recursive algorithm, such as recursive least squares. We have adopted the former choice here since the estimates are to be used to train and test classifiers for which independent training samples are required. Actual deployment may be better served by a recursive algorithm.

A related concern is the huge amount of background data needed to construct suitable sample estimates of the range-covariance matrices. This is not compatible with the limited amount of background data typically available.

See S. Mallat, Wavelet Tour of Signal Processing, Section 9.1.3 and Appendix A6, for an excellent discussion of the KL expansion (Academic, 1998).

If our covariance function, given by Eq. , had circular periodicity, Mallat shows that the Fourier basis would be the KL basis. We cannot claim this, but for target aerosols not close to either end of the sampling range, there is little error in making this assumption.

See, for example, F. L. Lewis, Optimal Estimation (Wiley, 1986).

Cited By

OSA participates in CrossRef's Cited-By Linking service. Citing articles from OSA journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (10)

Fig. 1
Fig. 1

Transmitted (dashed) and received (solid) waveforms from a single time step, line 10R38.

Fig. 2
Fig. 2

Received waveform from a single time step, line 10R38, after background subtraction.

Fig. 3
Fig. 3

Processing flow of the pulse-deconvolution filter.

Fig. 4
Fig. 4

Pulse-deconvolved waveform from a single time step, line 10R38.

Fig. 5
Fig. 5

Processing flow of the dual aerosol filters. The numbers refer to equations in the text.

Fig. 6
Fig. 6

Estimated relative concentration versus range and time step.

Fig. 7
Fig. 7

Estimated relative concentration versus range at time step 120.

Fig. 8
Fig. 8

Estimated non-instrument-corrected backscatter versus wavelength index and time step.

Fig. 9
Fig. 9

Estimated non-instrument-corrected backscatter mean and standard error estimates.

Fig. 10
Fig. 10

Estimates of (a) transition functions Φ n and (b) process variance Q n .

Tables (1)

Tables Icon

Table 1 Variables in Range- and KL-Transform-Space

Equations (48)

Equations on this page are rendered with MathJax. Learn more.

P r ( i , j , k ) = k = 1 k G ( i , j , k ) P t ( i , j , k k + 1 ) + ε r ( i , j , k ) ,
G ( i , j , k ) = G 0 ( i , j , k ) + G a ( i , j , k ) .
G 0 ( i , j , k ) = K j O ( k ) r k 2 T ( j , k ) β ( j , k ) ,
G a ( i , j , k ) = K j O ( k ) r k 2 T ( j , k ) α ( i , j ) C ( i , k ) ,
G a ( i , j , k ) K j O ( k 0 ) r k 0 2 T ( j , k 0 ) α ( i , j ) C ( i , k ) g ( j , k 0 ) α ( i , j ) C ( i , k ) ,
P r ( i , j , k ) = P r 0 ( i , j , k ) + k = 1 k G a ( i , j , k ) P t ( i , j , k k + 1 ) + ε r ( i , j , k ) ,
P r 0 ( i , j , k ) = k = 1 k G 0 ( i , j , k ) P t ( i , j , k k + 1 ) , G a ( i , j , k ) = α ( j , k ) C ( i , k ) .
P ¯ r 0 ( j , k ) = 1 T b i = 1 T b P r 0 ( i , j , k ) ,
G ^ a ( i , j , k ) = IFFT ( Y ( i , j , f ) · p r ( i , j , f ) ) , Y ( i , j , f ) S g ( f ) p t ( i , j , f ) * S g ( f ) | p t ( i , j , f ) | 2 + S n ( f ) = SNR ( f ) p t ( i , j , f ) * SNR ( f ) | p t ( i , j , f ) | 2 + 1 , SNR ( f ) S g ( f ) S n ( f ) ,
p r ( i , j , f ) = FFT ( P r ( i , j , k ) P ¯ r 0 ( j , k ) ) ( f ) , p t ( i , j , f ) = FFT ( P t ( i , j , k ) ) ( f ) ,
S g ( f ) = σ g 2 d z exp [ z 2 / r g 2 ] exp [ i 2 π f z ] = σ g 2 r g π exp [ ( π r g f ) 2 ] ,
j = 1 M [ α ( i , j ) ] 2 = 1 .
G ^ a ( i , j , k ) = α ( i , j ) C ( i , k ) + ε g ( i , j , k ) ,
E ( ε g ( i , j , k ) ε g ( i , j , k ) ) = Λ j ( k k ) δ i i δ j j = δ i i δ j j n = 1 L μ n ( j ) ϕ n ( k ) ϕ n ( k ) ,
ϕ n ( k ) 2 N sin π n ( k 1 ) N .
f 0 ( G ^ a ) = ( 2 π ) T b M N / 2 j = 1 M | Λ j | T b / 2 exp [ 1 2 i = 1 T b j = 1 M k , k = 1 N G ^ a ( i , j , k ) Λ j 1 ( k k ) G ^ a ( i , j , k ) ] ,
Λ j 1 ( k k ) = n = 1 L 1 μ n ( j ) ϕ n ( k ) ϕ n ( k ) .
l 0 = 1 2 i = 1 T b j = 1 M k , k = 1 N G ^ a ( i , j , k ) Λ j 1 ( k k ) G ^ a ( i , j , k ) T b 2 j ln | Λ j | .
l 0 = 1 2 i = 1 T b j = 1 M n = 1 L [ s n ( i , j ) ] 2 μ n ( j ) T b 2 n , j ln μ n ( j ) , s n ( i , j ) k = 1 N G ^ a ( i , j , k ) ϕ n ( k ) .
μ ^ n ( j ) = 1 T b i = 1 T b [ s n ( i , j ) ] 2 .
C ( i , k ) = n = 1 L a n ( i ) ϕ n ( k ) ,
l = 1 2 i = 1 T j = 1 M n = 1 L [ s n ( i , j ) α ( i , j ) a n ( i ) ] 2 μ ^ n ( j ) T 2 n , j ln μ ^ n ( j ) .
l α ( i , j ) = n = 1 L 1 μ ^ n ( j ) [ s n ( i , j ) α ( i , j ) a n ( i ) ] a n ( i ) .
α ^ ( i , j ) = n 1 μ ^ n ( j ) s n ( i , j ) a ^ n ( i ) n [ a ^ n ( i ) ] 2 μ ^ n ( j ) .
α ^ ( i , j ) | α ^ ( i , j ) | α ^ ( i , · ) , α ^ ( i , · ) [ j = 1 M α ^ ( i , j ) 2 ] 0.5 .
V ( α ^ ( i , j ) α ^ ( i , j ) ) E ( α ^ ( i , j ) α ( i , j ) ) ( α ^ ( i , j ) α ( i , j ) ) = δ i i δ j j n [ a ^ n ( i ) ] 2 / μ ^ n ( j ) .
s n ( i , j ) = α ^ ( i 1 , j ) a n ( i ) + ν n ( i , j ) ,
ν n ( i , j ) = k = 1 N ε g ( i , j , k ) ϕ n ( k ) .
E ( ν n ( i , j ) ν n ( i , j ) ) = δ i i δ j j δ n n μ n ( j ) .
a n ( i ) = Φ n a n ( i 1 ) + q n ( i ) ,
Γ n ( i ) = Φ n 2 Γ n ( i 1 ) + Q n .
a ^ ( i ) = Φ a ^ ( i 1 ) , Γ ^ ( i ) = Φ 2 Γ ^ ( i 1 ) + Q , Γ ^ ( i ) = [ ( Γ ^ ( i ) ) 1 + α ^ ( i 1 ) T μ ^ 1 α ^ ( i 1 ) ] 1 , K ( i ) = Γ ^ ( i ) α ^ ( i 1 ) T μ ^ 1 , a ^ ( i ) = a ^ ( i ) + K ( i ) [ s ( i ) α ^ ( i 1 ) a ^ ( i ) ] , C ^ ( i ) = a ^ ( i ) ϕ T .
C ( i , k ) = n = 1 L a n ( i ) ϕ n ( k ) ,
a n ( i ) = Φ n a n ( i 1 ) + q n ( i ) ,
s n ( i , j ) = k G ^ a ( i , j , k ) ϕ n ( k ) .
H ( Φ , Q | Φ ^ ( m 1 ) , Q ^ ( m 1 ) ) E ( 2 ln L a , s ( Φ , Q ) | s , Φ ^ ( m 1 ) , Q ^ ( m 1 ) ) ,
2 ln L a , s ( Φ , Q ) = T n = 1 L ln ( Q n ) + i = 1 T n = 1 L ( a n ( i ) Φ n a n ( i 1 ) ) 2 / Q n .
E ( a n ( i ) ) 2 = a ^ n 2 ( i ) + Γ ^ n ( i ) , E ( a n ( i ) a n ( i 1 ) ) = a ^ n ( i ) a ^ n ( i 1 ) + Φ ^ n ( m 1 ) Γ ^ n ( i 1 ) , E ( a n ( i 1 ) ) 2 = a ^ n 2 ( i 1 ) + Γ ^ n ( i 1 ) ,
H ( Φ , Q | Φ ^ ( m 1 ) , Q ^ ( m 1 ) ) = T n = 1 L ln ( Q n ) + n = 1 L Q n 1 [ S 11 , n 2 S 10 , n Φ n + S 00 , n Φ n 2 ] , S 11 , n = i = 1 T ( a ^ n 2 ( i ) + Γ ^ n ( i ) ) , S 10 , n = i = 1 T ( a ^ n ( i ) a ^ n ( i 1 ) + Φ ^ n ( m 1 ) Γ ^ n ( i 1 ) ) , S 00 , n = i = 1 T ( a ^ n 2 ( i 1 ) + Γ ^ n ( i 1 ) ) .
Φ ^ n ( m ) = S 10 , n / S 00 , n , Q ^ n ( m ) = 1 T [ S 11 , n S 10 , n 2 / S 00 , n ] .
S 11 , i n = γ S 11 , i 1 n + a ^ n 2 ( i ) + Γ ^ n ( i ) , S 10 , i n = γ S 10 , i 1 n + a ^ n ( i ) a ^ n ( i 1 ) + Φ ^ n ( i 1 ) Γ ^ n ( i 1 ) , S 00 , i n = γ S 00 , i 1 n + a ^ n 2 ( i 1 ) + Γ ^ n ( i 1 ) , Φ ^ n ( i ) = S 10 , i n / S 00 , i n , g i = γ g i 1 + 1 , Q ^ n ( i ) = [ S 11 , i n S 10 , i n 2 / S 00 , i n ] / g i .
f ( x | θ ) = i = 1 n f ( x i | θ ) .
θ ^ ML = arg max θ L ( θ | x ) .
θ ^ ML = arg max θ i = 1 n ln f ( x i | θ ) .
S j ( θ | x ) l ( θ | x ) θ j = i = 1 n 1 f ( x i | θ ) f ( x i | θ ) θ j ,
lim n E ( θ ^ ML ) = θ ¯ ,
lim n Pr [ θ ^ ML θ ¯ > ε ] = 0
[ Λ θ ] i j E ( ( θ i θ ¯ i ) ( θ j θ ¯ j ) ) = lim n { [ E ( S i ( θ ^ ML | x ) S j ( θ ^ ML | x ) ) ] 1 } .

Metrics