Abstract

Multitaper methods for a scan-free spectrum estimation that uses a rotational shear interferometer are investigated. Before source spectra can be estimated the sources must be detected. A source detection algorithm based upon the multitaper F-test is proposed. The algorithm is simulated, with additive, white Gaussian detector noise. A source with a signal-to-noise ratio (SNR) of 0.71 is detected 2.9° from a source with a SNR of 70.1, with a significance level of 104, ∼4 orders of magnitude more significant than the source detection obtained with a standard detection algorithm. Interpolation and the use of prewhitening filters are investigated in the context of rotational shear interferometer (RSI) source spectra estimation. Finally, a multitaper spectrum estimator is proposed, simulated, and compared with untapered estimates. The multitaper estimate is found via simulation to distinguish a spectral feature with a SNR of 1.6 near a large spectral feature. The SNR of 1.6 spectral feature is not distinguished by the untapered spectrum estimate. The findings are consistent with the strong capability of the multitaper estimate to reduce out-of-band spectral leakage.

© 2006 Optical Society of America

Full Article  |  PDF Article

References

  • View by:
  • |
  • |
  • |

  1. S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
    [CrossRef]
  2. S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
    [CrossRef]
  3. D. J. Thomson, "Spectrum estimation and harmonic analysis," Proc. IEEE 70, 1055-1096 (1982).
    [CrossRef]
  4. L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge U. Press, 1995).
  5. D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
    [CrossRef] [PubMed]
  6. D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty—I," Bell Syst. Tech. J. 40, 43-63 (1961).
  7. A. Hanssen, "Higher-dimensional spectral estimation by means of multitapering," in Proceedings of 1997 International Conference on Information, Communications and Signal Processing (IEEE, 1997), Vol. 3, pp. 1603-1607.
  8. L. L. Scharf.Statistical Signal Processing: Detection, Estimation, and Time Series Analysis (Addison-Wesley, 1991).
  9. E. W. Weisstein, "Fast fourier transform," from mathworld—a wolfram web resource, 1998, http://mathworld.wolfram.com/FastFourierTransform.html.
  10. D. Percival and A. Walden, Spectral Analysis for Physical Applications (Cambridge U. Press, 1993).
    [CrossRef]

2003 (1)

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

2002 (1)

S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
[CrossRef]

1999 (1)

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

1982 (1)

D. J. Thomson, "Spectrum estimation and harmonic analysis," Proc. IEEE 70, 1055-1096 (1982).
[CrossRef]

1961 (1)

D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty—I," Bell Syst. Tech. J. 40, 43-63 (1961).

Basty, S.

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

Brady, D.

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

Brady, D. J.

S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
[CrossRef]

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Brady, R. B.

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Gallicchio, J. R.

S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
[CrossRef]

Hanssen, A.

A. Hanssen, "Higher-dimensional spectral estimation by means of multitapering," in Proceedings of 1997 International Conference on Information, Communications and Signal Processing (IEEE, 1997), Vol. 3, pp. 1603-1607.

Kraut, S.

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
[CrossRef]

Mandel, L.

L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge U. Press, 1995).

Marks, D. L.

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Munson, D. C.

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Neifeld, M. A.

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

Percival, D.

D. Percival and A. Walden, Spectral Analysis for Physical Applications (Cambridge U. Press, 1993).
[CrossRef]

Pollak, H. O.

D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty—I," Bell Syst. Tech. J. 40, 43-63 (1961).

Scharf., L. L.

L. L. Scharf.Statistical Signal Processing: Detection, Estimation, and Time Series Analysis (Addison-Wesley, 1991).

Slepian, D.

D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty—I," Bell Syst. Tech. J. 40, 43-63 (1961).

Stack, R. A.

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Thomson, D. J.

D. J. Thomson, "Spectrum estimation and harmonic analysis," Proc. IEEE 70, 1055-1096 (1982).
[CrossRef]

Walden, A.

D. Percival and A. Walden, Spectral Analysis for Physical Applications (Cambridge U. Press, 1993).
[CrossRef]

Weisstein, E. W.

E. W. Weisstein, "Fast fourier transform," from mathworld—a wolfram web resource, 1998, http://mathworld.wolfram.com/FastFourierTransform.html.

Wolf, E.

L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge U. Press, 1995).

Bell Syst. Tech. J. (1)

D. Slepian and H. O. Pollak, "Prolate spheroidal wave functions, Fourier analysis and uncertainty—I," Bell Syst. Tech. J. 40, 43-63 (1961).

Opt. Commun. (1)

S. Basty, M. A. Neifeld, D. Brady, and S. Kraut, "Nonlinear estimation for interferometric imaging," Opt. Commun. 228, 249-261 (2003).
[CrossRef]

Proc. IEEE (1)

D. J. Thomson, "Spectrum estimation and harmonic analysis," Proc. IEEE 70, 1055-1096 (1982).
[CrossRef]

Proc. SPIE (1)

S. Kraut, J. R. Gallicchio, and D. J. Brady, "High-resolution direction finding and scan-free spectrum estimation with rotational-shear interferometric sensor arrays," in Algorithms and Systems for Optical Information ProcessingVI, B. Javidi and D. Psaltis, eds., Proc. SPIE 4789, 267-275 (2002).
[CrossRef]

Science (1)

D. L. Marks, R. A. Stack, D. J. Brady, D. C. Munson, and R. B. Brady, "Visible cone-beam tomography with a lensless interferometric camera," Science 284, 2164-2166 (1999).
[CrossRef] [PubMed]

Other (5)

L. Mandel and E. Wolf, Optical Coherence and Quantum Optics (Cambridge U. Press, 1995).

A. Hanssen, "Higher-dimensional spectral estimation by means of multitapering," in Proceedings of 1997 International Conference on Information, Communications and Signal Processing (IEEE, 1997), Vol. 3, pp. 1603-1607.

L. L. Scharf.Statistical Signal Processing: Detection, Estimation, and Time Series Analysis (Addison-Wesley, 1991).

E. W. Weisstein, "Fast fourier transform," from mathworld—a wolfram web resource, 1998, http://mathworld.wolfram.com/FastFourierTransform.html.

D. Percival and A. Walden, Spectral Analysis for Physical Applications (Cambridge U. Press, 1993).
[CrossRef]

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 (11)

Fig. 1
Fig. 1

(Left) RSI pattern produced by two equal-intensity point sources separated by Δφ = 8 deg and Δγ = 8.7 deg. (Right) Magnitude of the Fourier transform. Note the double-peaked Lorentzian spectra along the radial lines containing power.

Fig. 2
Fig. 2

Estimate of RSI fringe-pattern Fourier transform (magnitude). Interferogram produced by two equal-intensity point sources with identical, double-peaked Lorentzian spectra. The sources are separated in roll and elevation, thus rotating their corresponding radial lines of maximum power in the Fourier transform with respect to each other as well as radially shifting their spectral features. The interferogram is zero padded to four times its original size before DFT computation. This estimate is the 2D analog of Eq. (9). Note the leakage, suppressed 10 dB from the spectra peaks, spanning as much as half of the spatial frequencies along the vertical and horizontal lines.

Fig. 3
Fig. 3

Power spectral estimates of a sinusoid with a frequency of 0.1 cycles∕sample. The x is the true spectrum, the periodogram (thin curve, narrow mainlobe) is the square of the magnitude of Eq. (9), and the Slepian-computed spectral estimate (thick curve, wide mainlobe) is the square of the magnitude of Eq. (11), where w[n] is the zeroth-order Slepian with NW = 4 and N = 101. The Slepian estimate is biased by approximately 30 dB less outside the frequency interval (−0.1 − 4∕101, −0.1 + 4∕101) ∪ (0.1 − 4∕101, 0.1 + 4∕101). Note that the resolution in the Slepian-computed estimate is traded, relative to the frequency resolution of the periodogram, for a reduction in sidelobe level. Finally, note the approximately 2 dB of inband bias present in the Slepain-computed estimate.

Fig. 4
Fig. 4

Estimates of RSI interference-pattern Fourier transform (magnitudes). Left, estimate using the Dirichlet kernel; right, estimate using the zeroth-order 2D Slepian taper. There is an approximate 30 dB reduction in sidelobes between the left and the right estimates of the Fourier-transform magnitude. Note that the weak source in the left-hand image is barely visible but is clearly present in the right-hand image. Again, the resolution of the Slepian-tapered estimate is traded for a reduction in sidelobe level.

Fig. 5
Fig. 5

Magnitude of the discrete Fourier transform of an image with a row of ones. The image has 501 rows and 501 columns. Note the magnitude of the 1D Dirichlet kernel centered along the rows in the Fourier transform. Before transformation, the image is zero padded to four times its size in both dimensions. The DFT is not normalized. By normalizing the DFT by [(4 × 501)2]1∕2 to obtain a unitary version of the DFT, and by summing the peak of the Dirichlet kernel along the rows, one approximately obtains 501. This is the sum of the row of ones in the image, and it agrees with Eq. (18).

Fig. 6
Fig. 6

Two-dimensional Slepian tapers for detecting a source with tan−1 β j ∕α j = 135°. These tapers are formed from outer products with 1D tapers with NW = 500 and NW = 3. The tapers are computed aligned with the vertical and horizontal axes and are then rotated. The rotation is performed by rotating a square grid of points and then interpolating the tapers onto the rotated grid. In principle, the tapers can be computed off-line and stored in memory.

Fig. 7
Fig. 7

Source detection simulation. (Upper left) Simulated RSI interferogram produced by two sources with identical, double-peaked Lorentzian spectra, with an intensity ratio of 1:100. The sources are separated by 2.9° of roll about the RSI optical axis and by 1.2° of elevation from the optical axis. The SNR of the weak source is 1 / 2 W∕m2. The dim source at 135° is detected by the F-ratio of Eq. (19) with a significance of 10−4, and the F-ratio ratio formed by radial summing detects the same source at a significance level of 0.25.

Fig. 8
Fig. 8

(All plots) Thick curve, actual source spectrum. (Top plot) Thin curve, spectrum estimated by interpolating the 2D DFT onto the radial–azimuthal grid by using linear interpolation. The line of medium thickness is DFT interpolated by using nearest-neighbor interpolations. SNR = 1 × 106 (Middle plot) DFTs interpolated by using cubic and spline interpolations. Note that these 2D interpolation algorithms are those implemented in matlab Version 7.0.0.19901 (Rev. 14). (Bottom plot) DFT directly estimated along the detected angle at points evenly spaced in radial frequency. SNR = 1 × 104.

Fig. 9
Fig. 9

Same as in Fig. 8 with the exception that the SNR is 2∕2.

Fig. 10
Fig. 10

Effect of prewhitening. All plots, interference pattern pixel noise standard deviation of 10−4 W∕m2; top plot, nonprewhitened bright source spectrum estimates (angle, 137.9°; SNR, 1 × 106); second plot; nonprewhitened dim source spectrum estimates (angle, 135°; SNR, 1 × 104); third plot, prewhitened bright source spectrum estimates (angle, 137.9°; SNR, 1 × 106); bottom plot, prewhitened dim source spectrum estimates (angle, 135°; SNR, 1 × 104). By reducing the dynamic range, prewhitening reduces the leakage from the bright source into spectral estimates for the dim source. This is evidenced by the approach of the untapered, prewhitened spectral estimate to the multitaper estimate, which has superior leakage protection.

Fig. 11
Fig. 11

(Top) Unprewhitened dim source spectrum estimates (angle, 135°; SNR, 28.28). Note the substantial sidelobe levels present adjacent to the spectral peak in the untapered estimate. The multitaper estimate successfully distinguishes the small peak at 1.5 m−1. This feature contains 0.1% of the intensity of the 40 W∕m2 source and has a SNR of 1.6. (Bottom), AR-2 prewhitened spectral estimates (angle, 135°; SNR, 28.28). Sidelobes are visible next to strong spectral features despite prewhitening. The approximate impulse response of the NW = 4, 8 taper multitaper estimate is centered about the strong spectral features. The spacing of frequency samples is 6.19 × 103 m−1. The prewhitened multitaper estimate distinguishes the peak at 1.5 m−1 from the noise floor. The prewhitened untapered estimate does not.

Equations (42)

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

I D ( x , y ) = 2 I o + 2 Re { Γ ( x , y , Δ z ) } ,
Γ ( x , y , Δ z ) = I ( α , β ) dαdβ 0 exp ( i 2 π ν Δ z c ) × S ( α , β , ν ) exp [ i 2 πν c 2 sin × θ ( x β y α ) ] d ν
I ( α , β ) = j = 1 N I j δ ( α α j , β β j ) ,
Γ ( x , y , 0 ) = j = 1 N I j ( α j , β j ) 0 S ( α j , β j , ν ) × exp [ i 2 πν c 2 sin θ ( x β j y α j ) ] d ν .
x m [ n ] = x a ( n T s ) = x a ( n ) .
x ^ e ( f ) = n = 0 N 1 x m [ n ] exp ( i 2 π f n ) .
x m [ n ] = 1 / 2 1 / 2 x ^ a ( f ) exp ( i 2 π f n ) d f
x ^ e ( f ) = n = 0 N 1 1 / 2 1 / 2 x ^ a ( f ) exp ( i 2 π f n ) d f exp ( i 2 π f n ).
x ^ e ( f ) = 1 / 2 1 / 2 x ^ a ( f ) exp [ i π ( f f ) ( N 1 ) ] × sin [ N π ( f f ) ] sin [ π ( f f ) ]  d f
= x ^ a ( f ) D ( f ) ,
x ^ e ( f ) = n = 0 N 1 w [ n ] x m [ n ] exp ( i 2 π f n )
= x ^ a ( f ) w ^ ( f ) .
w ^ ( f ) = n = 0 N 1 w [ n ] exp ( i 2 π f n ) .
= r = 0 N 1 c = 0 N 1 I D [ r , c ] δ K ( r , m c + b ) ,
= r = 0 N 1 c = 0 N 1 { 1 / 2 1 / 2 1 / 2 1 / 2 I ^ D actual ( f x , f y ) × exp ( i 2 π f x c ) exp ( i 2 π f y r )  d f x d f y } × δ K ( r , m c + b )
= 1 / 2 1 / 2 1 / 2 1 / 2 I ^ D actual ( f x , f y ) { c = 0 N 1 exp [ i 2 π c ( f x + m f y ) ] exp ( i 2 π b f y ) } d f x d f y
= 1 / 2 1 / 2 1 / 2 1 / 2 I ^ D actual ( f x , f y ) D ( f x + m f y ) × exp ( i 2 π b f y ) d f x d f y
Δ D ( 0 ) 1 / 2 1 / 2 1 / 2 1 / 2 I ^ D actual ( f x , f y ) × exp ( i 2 π b f y ) δ ( f y , f x / m ) d f x d f y ,
F = ( K 1 ) y ̃ ( P 1 ̃ ) N × N y ̃ y ̃ ( P 1 ̃ ) N × N y ̃ ,
x m [ j , k ] = r = 1 P c = 1 P a [ r , c ] x m [ j r , k c ] + n [ j , k ] ,
z [ j , k ] = x d [ j , k ] r = 1 P c = 1 P a [ r , c ] x d [ j r , k c ]
z ^ ( f j , f k ) = x ^ d ( f j , f k ) { 1 r = 1 P c = 1 P a [ r , c ] exp ( i 2 π f j r ) × exp ( i 2 π f k c ) } .
x ^ e ( f j , f k ) =
z ^ e ( f j , f k ) 1 r = 1 P c = 1 P a [ r , c ] exp ( i 2 π f j r ) exp ( i 2 π f k c ) .
S m ( ν ) = 1 2 W k = 1 K D k | y ^ k [ ν cos ( α m ) , ν sin ( α m ) ] | k = 1 K D k     2 ,
y ^ k ( f y , f x ) = r = 0 N 1 c = 0 N 1 t k [ r , c ] y [ r , c ] exp ( i 2 π f y r ) × exp ( i 2 π f x c )
Var   S m ( ν ) = 1 4 W 2 1 ( k = 1 K D k     2 ) 2 × k = 1 K D k     2   Var | y ^ k [ ν   cos ( α m ) , ν sin ( α m ) ] | ,
Var | y ^ k | = 4 π 2 σ 2 ,
Var   S m ( ν ) = 4 π 8 W 2 σ 2 k = 1 K D k     2 .
Var   S m ( ν ) = 3.85 × 10 15 W 2 / m 4 / m 2 ,
I ( x , y , α j , β j ) = Re { 0 S ( α j , β j ) exp [ i 2 π g ( ν ) ( α j y β j x ) ] } d ν ,
Q ( α j , β j , ν ) = S ( α j , β j , | ν | ) 2 ,
I ( x , y , α j , β j ) = Q ( α j , β j , ν ) exp [ i 2 π g ( ν ) ( α j y β j x ) ] d ν .
I ^ e ( f x , f y , α j , β j ) = 2 t ^ ( 0 , 0 ) r , c = 0 N 1 t [ r , c ] × I ( c Δ x , r Δ y , α j , β j ) × exp ( i 2 π f x c Δ x ) × exp ( i 2 π f y r Δ y ) ,
I ^ e ( f x , f y , α j , β j ) = 2 t ^ ( 0 , 0 ) r , c = 0 N 1 t [ r , c ] Q ( α j , β j , ν ) × exp [ i 2 π g ( ν ) ( α j r Δ y β j c Δ x ) ] exp ( i 2 π f x c Δ x ) × exp ( i 2 π f y r Δ y ) .
I ^ c ( f x , f y , α j , β j ) = 2 t ^ ( 0 , 0 ) Q ( α j , β j , ν ) r , c = 0 N 1 t [ r , c ] × exp { i 2 π c Δ x [ f x + g ( ν ) β j ] } × exp { i 2 π r Δ y [ f y g ( ν ) α j ] }
= 2 t ^ ( 0 , 0 ) Q ( α j , β j , ν ) t ^ ( f x + g ( ν ) β j , f y g ( ν ) α j ) .
t ^ ( u , v ) { t ^ ( 0 , 0 ) , u [ W , W ] , v [ W , W ] 0 , otherwise ,
A ( f x , α j , β j ) = { ν : ( f x + g ( ν ) β j [ W , W ] )       &  ( α j β j f x g ( ν ) α j [ W , W ] ) } .
I ^ c ( f x , f y , α j , β j ) { 2 t ( 0 , 0 ) A Q ( α j , β j , v ) t ^ ( f x + g ( ν ) β j , f y g ( ν ) α j ) d ν , f y = α i / β j 0 , otherwise .
I ^ e ( f x , f y , α j , β j ) { 2 t ( 0 , 0 ) [ t ^ ( 0 , 0 ) 2 W Q ( α j , β j , g 1 ( f x / β j ) ) ] , f y = α j / β j 0 , otherwise .
S j ( ν , α j , β j ) I ^ e ( β j g ( ν ) , α j g ( ν ) ) 2 W .

Metrics