Fizeau Fourier transform imaging spectroscopy: missing data reconstruction

doi:10.1364/OE.16.006631

« Show journal navigation

Fizeau Fourier transform imaging spectroscopy: missing data reconstruction

Samuel T. Thurman and James R. Fienup »View Author Affiliations

Optics Express, Vol. 16, Issue 9, pp. 6631-6645 (2008)
http://dx.doi.org/10.1364/OE.16.006631

View Full Text Article

Acrobat PDF (660 KB)

Journal Search
Article Lookup

Article Tools

Citations

Alert me when this article is cited
Export Citation/Save

Full Text
Article Info
References (22)
Cited By
Figures (11)
Multimedia (1)
Metrics

Abstract

Fizeau Fourier transform imaging spectroscopy yields both spatial and spectral information about an object. Spectral information, however, is not obtained for a finite area of low spatial frequencies. A nonlinear reconstruction algorithm based on a gray-world approximation is presented. Reconstruction results from simulated data agree well with ideal Michelson interferometer-based spectral imagery. This result implies that segmented-aperture telescopes and multiple telescope arrays designed for conventional imaging can be used to gather useful spectral data through Fizeau FTIS without the need for additional hardware.

1. Introduction

Imaging spectroscopy is the process of collecting of both spatial and spectral information for a scene. Several grating- or prism-based imaging spectrometers have been deployed for remote sensing applications. An alternative to dispersive spectroscopy is the Fourier transform technique [1

J. Kauppinen and J. Partanen, Fourier Transforms in Spectroscopy , (Wiley-VCH, Berlin, 2001). [CrossRef]

]. The conventional design for such a system is based on a Michelson imaging interferometer [2-4

N. J. E. Johnson, “Spectral imaging with the Michelson interferometer,” in Infrared Imaging Systems Technology, W. L Wolfe and J. Zimmerman, eds., Proc. SPIE 226, 2–9 (1980).

]. Panchromatic intensity measurements are made with various optical path differences (OPDs) between the arms of the interferometer, and spectral information is recovered through post-processing. Alternatively, Fourier transform imaging spectroscopy (FTIS) can be performed with a Fizeau imaging interferometer by introducing OPDs between subapertures of a multiple telescope array (MTA) or a segmented-aperture system [5-7

M. Frayman and J. A. Jamieson, “Scene imaging and spectroscopy using a spatial spectral interferometer,” in Amplitude and Intensity Spatial Interferometry, J. B. Breckingridge, ed., Proc. SPIE 1237, 585–603 (1990). [CrossRef]

]. Kendrick et al. [6

R. L. Kendrick, E. H. Smith, and A. L. Duncan, “Imaging Fourier transform spectrometry with a Fizeau interferometer,” in Interferometry in Space, M. Shao, ed., Proc. SPIE 4852, 657–662 (2003). [CrossRef]

] demonstrated the Fizeau FTIS concept with a two-telescope array system. An advantage of Fizeau FTIS is that it can be implemented with a MTA or segmented-aperture system, using existing optical delay lines or subaperture piston actuators that are normally used to phase up such systems for conventional imaging. Thus, such systems do not require any additional hardware or modifications to implement Fizeau FTIS, whereas the implementation of Michelson FTIS would typically require the addition of a Michelson interferometer science instrument to an existing imaging system.

In a previous paper [7

S. T. Thurman and J. R. Fienup, “Multi-aperture Fourier transform imaging spectroscopy: theory and imaging properties,” Opt. Express 13, 2160–2175 (2005). [CrossRef] [PubMed]

], we analyzed the imaging properties of Fizeau FTIS. A unique aspect of Fizeau FTIS is that spectral information is not obtained for a finite bandwidth of low spatial frequencies. Section 2 reviews the relevant equations and discusses this problem in more detail. Nonlinear reconstruction algorithms are required for reconstructing the missing information to yield spectral imagery that is more easily interpretable. In Sec. 3, we describe an algorithm that uses a gray-world model to reconstruct the missing data. The need for nonlinear reconstruction algorithms suggests that the Fizeau configuration is suboptimum for performing FTIS. Simulation results presented in Sec. 4, however, indicate that spectral imagery of moderate quality can be obtained with Fizeau FTIS. This result demonstrates that a MTA or segmented-aperture system designed for conventional panchromatic imaging has the added value of being able to fulfill some spectral imaging tasks without additional instrumentation. Section 5 is a discussion and summary.

2. Imaging model

The imaging properties of Fizeau FTIS are discussed in detail in Ref. [7

S. T. Thurman and J. R. Fienup, “Multi-aperture Fourier transform imaging spectroscopy: theory and imaging properties,” Opt. Express 13, 2160–2175 (2005). [CrossRef] [PubMed]

]. Here we only present equations necessary for formulating the image reconstruction algorithm. Furthermore, we only consider the case of introducing OPDs between two groups of subapertures for performing Fizeau FTIS, i.e., one group of subapertures has a fixed optical path length while the optical path length through the other group is varied. The generalized pupil function of the system can be written as

T (ξ, η, v, τ) = T_{1} (ξ, η, v) + T_{2} (ξ, η, v) \exp (i 2 π v τ),

(1)

where T ₁(ξ,η,ν) and T ₂(ξ,η,ν) are the respective generalized pupil functions for the two subaperture groups, (ξ,η) are pupil plane coordinates, ν is the optical frequency, τ=OPD/c is the time-delay introduced between the subaperture groups, and c is the speed of light. The coherent impulse response [8

J. Goodman, Introduction to Fourier Optics 2nd ed. , (McGraw-Hill, New York, 1996).

] corresponding to each T_q (ξ,η,ν) is

t_{q} (x, y, v) = \frac{1}{λ^{2} f^{2}} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} T_{q} (ξ, η, v) \exp [- i 2 π (\frac{x}{λ f} ξ + \frac{y}{λ f} η)] d ξ d η,

(2)

where q=1 and 2 is an index for the subaperture groups, λ=c/ν is the optical wavelength, f is the system focal length, and (x,y) are focal plane coordinates. It is convenient to define the spectral point spread functions (SPSFs), h _p,q(x,y,ν), as

h_{p, q} (x, y, v) = t_{p}^{*} (x, y, v) t_{q} (x, y, v),

(3)

for p and q=1 and 2, such that the incoherent point spread function (PSF) for the entire system can be written as a function of the time delay, τ, between the subaperture groups, as

h (x, y, v, τ) = h_{1,1} (x, y, v) + h_{1,2} (x, y, v) \exp (- i 2 π v τ)

+ h_{2,2} (x, y, v) + h_{2,1} (x, y, v) \exp (i 2 π v τ) .

(4)

The image intensity as a function of τ is given by

I (x, y, τ) = \int_{0}^{\infty} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} κ S_{o} (x', y', v) h (x - x', y - y', v, τ) d x' d y' d v,

(5)

where κ=λ ²/π, S _o(x,y,ν) is the spectral density of a spatially incoherent object, and the system is assumed to have unit magnification for notational convenience.

In performing Fizeau FTIS, a number of images are recorded by a focal plane array with different time delays between the subapertures. Following the normal Fourier transform spectroscopy processing steps, the spectral image, S _i(x,y,ν), is then computed as the Fourier transform of the measured intensity modulation at each image point, i.e.,

S_{i} (x, y, v) = \int_{- \infty}^{\infty} [I (x, y, τ) - I_{avg} (x, y)] \exp (- i 2 π τ v) d τ,

(6)

where

I_{avg} (x, y) = \underset{T \to \infty}{Lim} \frac{1}{2 T} \int_{- T}^{T} I (x, y, τ) d τ

(7)

is the average image intensity at each image point. Reference [7

S. T. Thurman and J. R. Fienup, “Multi-aperture Fourier transform imaging spectroscopy: theory and imaging properties,” Opt. Express 13, 2160–2175 (2005). [CrossRef] [PubMed]

] discusses the linear relationship between S _i(x,y,ν) and S _o(x,y,ν) involving the SPSFs, which enables the real part of the complex-valued spectral image S _i(x,y,ν) to be written as

S_{i}^{(Re)} (x, y, v) = \frac{κ}{2} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} S_{o} (x', y', v) [h_{1, 2} (x - x', y - y', v)

+ h_{2, 1} (x - x', y - y', v) d x' d y',

(8)

for ν>0. The problem of the missing low spatial frequencies can be understood by taking the spatial Fourier transform of S ^(Re) _i (x, y,ν) to yield [9

Note that the additional factor of A/(λ ² f ² _i ) at the beginning of Eq. (9) was omitted from Eqs. (9), (16), and (17) of Ref. [6] in error.

]

G_{i}^{(Re)} (f_{x}, f_{y}, v) = \frac{κ A_{pup}}{2 λ^{2} f^{2}} {[H_{1, 2} (f_{x}, f_{y}, v) + H_{2, 1} (f_{x}, f_{y}, v)] G}_{o} (f_{x}, f_{y}, v),

(9)

where A _pup is the combined area of all of the subapertures, G _o(f_x ,f_y ,ν) is the spatial Fourier transform (FFT) of S _o(x,y,ν), and the terms H_p,q (f_x ,f_y ,ν) are spectral optical transfer functions (SOTFs) defined as the normalized Fourier transform of the corresponding SPSFs, i.e.,

H_{p, q} (f_{x}, f_{y}, v) = \frac{1}{A_{pup}} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} T_{p} (ξ, η, v) T_{q}^{*} (ξ + λ {ff}_{x}, η + λ {ff}_{y}, v) d ξ d η .

(10)

Note that the SOTFs are given by cross correlations between subaperture group pupil functions. Thus, the terms H _1,2(f_x ,f_y ,ν) and H _2,1(f_x ,f_y ,ν) vanish necessarily in some finite region around the dc spatial frequency (f_x ,f_y )=(0,0) for non-overlapping subaperture groups, as is the case for a MTA or segmented aperture system. Low spatial-frequency data in this region are missing from Fizeau FTIS spectral imagery.

The measured image intensity I(x,y,τ) does, however, contain panchromatic information about the missing low spatial frequencies that can be used to facilitate the reconstruction of this missing data. The spatial Fourier transform of I _avg(x,y) is given by

G_{avg} (f_{x}, f_{y}) = \int_{0}^{\infty} \frac{κ A_{pup}}{λ^{2} f^{2}} [H_{1, 1} (f_{x}, f_{y}, v) + H_{2, 2} (f_{x}, f_{y}, v)] G_{o} (f_{x}, f_{y}, v) d v .

(11)

Note that the SOTF terms H _1,1(f_x ,f_y ,ν) and H _2,2(f_x ,f_y ,ν) are given by the autocorrelation of the pupil functions T ₁(ξ,η,ν) and T ₂(ξ,η,ν), respectively, and are non-zero in the region around the dc spatial frequency. Thus, G _avg(f_x ,f_y ) contains spectrally-integrated information about the missing data.

3. Reconstruction algorithm

We wish to reconstruct an object spectral density estimate Ŝ _o(x,y,ν) from the data and any available a priori knowledge about the object. The algorithm we describe for doing this uses a Wiener filter [10

C. W. Helstrom, “Image restoration by the method of least squares,” J. Opt. Soc. Am. 57, 297–303 (1967). [CrossRef]

] to reconstruct spatial frequencies inside the support of the SOTFs H _1,2(f_x ,f_y ,ν) and H _2,1(f_x ,f_y ,ν) and reconstructs the remaining spatial-frequency data inside the support of the conventional optical transfer function (OTF) from the panchromatic intensity measurements using a gray-world approximation. The gray-world approximation assumes that to first order S _o(x,y,ν) can be approximated as

S_{o} (x, y, v) \approx f (x, y) ψ (v),

(12)

where f(x,y) is the total flux at each image point (x,y) and ψ(ν) is a gray-world spectrum. In practice, the gray-world approximation can be fairly good for natural scenery within certain spectral intervals. For example, the gray-world approximation yields a root mean square (RMSE) error of 0.099 [W·m⁻²·sr⁻¹·µm⁻¹] and a normalized RMSE of 6.5% for the non-gray spectral imagery used for the simulation in Section 4.

In practice, the FTIS data is discretely sampled, whereas the equations in Section 2 are in terms of continuous variables. For convenience, we use the same notation for both the continuous and discretely sampled quantities. For example, while I(x,y,τ) is a function of the continuous variables x, y, and τ in Eq. (5), we also use I(x,y,τ) to represent the image intensity measured at discrete image points and time delays, such that x ∊{−Δ_x N_x /2, −Δ_x(N_x −2)/2, …, Δ_x(N_x −2)/2}, y ∊{−Δ_y N_y /2, −Δ_y(N_y −2)/2, …, Δ_y(N_y −2)/2}, and τ ∊{−Δ_τ N_τ /2, −Δ_τ (N_τ −2)/2, …, Δ_τ(N_τ −2)/2}, where Δ_x and Δ_y are the sample spacings in the focal plane, i.e., the detector pixel pitch, Δ_τ is the time delay sample spacing, and N_x , N_y , and N_τ are the number of samples for each variable. Additionally, the unprocessed spectral image S _i(x,y,ν) is calculated using a fast Fourier transform (FFT), and is discretely sampled along the spectral dimension with a sample spacing of Δν=1/(N_τ Δ_τ). Likewise, the sample spacings for spatial-frequency domain quantities calculated using FFTs are Δ_fx=1/(N_x Δ_x) and Δ_fy=1/(N_y Δ_y).

A gray-world spectrum ψ(ν) for the reconstruction can be estimated from the unprocessed spectral imagery by minimizing the following cost function

E_{ψ} = \sum_{(f_{x}, f_{y}, v)} {∣ H_{2, 1} (f_{x}, f_{y}, v) F_{1} (f_{x}, f_{y}) ψ (v) - G_{i} (f_{x}, f_{y}, v) ∣}^{2},

(13)

where the summation is over all spatial- and temporal-frequency samples, F ₁(f_x ,f_y ), given by

F_{1} (f_{x}, f_{y}) = \frac{\sum_{v} H_{2, 1}^{*} (f_{x}, f_{y}, v) ψ (v) G_{i} (f_{x}, f_{y}, v)}{ε + {\sum_{v} ∣ H_{2, 1} (f_{x}, f_{y}, v) ψ (v) ∣}^{2}},

(14)

which minimizes the value of E_ψ , analogous to Eq. (17) in Ref. [11

J. R. Fienup, “Invariant error metrics for image reconstruction,” Appl. Opt. 36, 8352–8357 (1997). [CrossRef]

], for any particular ψ(ν), and ε is a small number included in the denominator of Eq. (14) to prevent infinite values when numerically computing F ₁(f_x ,f_y ).

Additionally, the object spatial power spectrum Φ_o(f_x ,f_y ,ν), which is needed for constructing the Wiener filter, can be estimated from the unprocessed spectral imagery. Using the statistics of natural scenery [12-15

G. J. Burton and I. R. Moorhead, “Color and spatial structure in natural scenes,” Appl. Opt. 26, 157–170 (1987). [CrossRef] [PubMed]

] and the gray-world approximation, Φ_o(f_x ,f_y ,ν) is modeled as

Φ_{o} (f_{x}, f_{y}, v) = {\begin{matrix} \begin{matrix} A_{0}^{2} ψ^{2} (v) & for (f_{x}, f_{y}) = (0, 0) \end{matrix} \\ \begin{matrix} A^{2} ψ^{2} (v) {(f_{x}^{2} + f_{y}^{2})}^{- α} & for (f_{x}, f_{y}) \neq (0, 0) \end{matrix} \end{matrix} \begin{matrix}  \end{matrix}

(15)

where A ₀, A, and α are parameters of the model. The parameters A and α, as well as the noise power spectrum Φ_n (assumed to be white) are obtained by minimizing the cost function

E_{Φ} = \sum_{(f_{x}, f_{y}, v)} \frac{1}{\sqrt{f_{y}^{2} + f_{y}^{2}}} {\ln [∣ G_{i}^{(Re)} (f_{x}, f_{y}, v) ∣]

{- \ln [0.5 ∣ H_{1, 2} (f_{x}, f_{y}, v) + H_{2, 1} (f_{x}, f_{y}, v) ∣ \sqrt{Φ_{o} (f_{x}, f_{y}, v)} + \sqrt{Φ_{n}}]}}^{2},

(16)

where the dc spatial frequency samples, (f_x ,f_y )=(0,0), are not included in the summation. In practice, we find it necessary to blur |G_i ^(Re)(f_x ,f_y ,ν)| slightly (by convolution with a 5-8 pixel wide kernel) for use in Eq. (16) to prevent underestimating Φ_o(f_x ,f_y ,ν) caused by the logarithmic weighting of very small numbers. The parameter A ₀ is estimated as

A_{0} = \frac{G_{avg} (0, 0)}{\sum_{v} ψ (v)} .

(17)

The object and noise power spectra are then used to form a Wiener filter [10

C. W. Helstrom, “Image restoration by the method of least squares,” J. Opt. Soc. Am. 57, 297–303 (1967). [CrossRef]

], given by

W_{1} (f_{x}, f_{y}, v) = \frac{1}{0.5 [H_{1, 2} (f_{x}, f_{y}, v) + H_{2, 1} (f_{x}, f_{y}, v)]} \frac{{SNR}_{1} (f_{x}, f_{y}, v)}{1 + {SNR}_{1} (f_{x}, f_{y}, v)},

(18)

where SNR ₁(f_x ,f_y ,ν) is the estimated power signal-to-noise ratio (SNR) for G_i ^(Re)(f_x ,f_y ,ν), given by

{SNR}_{1} (f_{x}, f_{y}, v) = \frac{0.25 {∣ H_{1, 2} (f_{x}, f_{y}, v) + H_{2, 1} (f_{x}, f_{y}, v) ∣}^{2} Φ_{o} (f_{x}, f_{y}, v)}{Φ_{n}} .

(19)

The Wiener-filtered spectral image S ₁(x,y,ν) is calculated as the inverse spatial 2-D FFT of

G_{1} (f_{x}, f_{y}, v) = W_{1} (f_{x}, f_{y}, v) G_{i}^{(Re)} (f_{x}, f_{y}, v) .

(20)

Additionally, a gray-world object spectral density estimate S ₂(x,y,ν) is generated from G _avg(f_x ,f_y ) using the gray-world approximation. In the gray-world approximation, we can write

G_{avg} (f_{x}, f_{y}) \approx \sum_{v} [H_{1, 1} (f_{x}, f_{y}, v) + H_{2, 2} (f_{x}, f_{y}, v)] F_{2} (f_{x}, f_{y}) ψ (v)

(21)

where F ₂(f_x ,f_y ) is the 2D FFT of f(x,y) in Eq. (12). The power SNR of G _avg(f_x ,f_y ) can be estimated as

{SNR}_{2} (f_{x}, f_{y}) = \frac{{∣ \sum_{v} [H_{1, 2} (f_{x}, f_{y}, v) + H_{2, 1} (f_{x}, f_{y}, v)] \sqrt{Φ_{o} (f_{x}, f_{y}, v)} ∣}^{2}}{Φ_{n}} .

(22)

By inverting Eq. (21) and including a noise regularization term identical to that for the Wiener filter, F ₂(f_x ,f_y ) can be calculated as

F_{2} (f_{x}, f_{y}) = \frac{G_{avg} (f_{x}, f_{y})}{\sum_{v} [H_{1,1} (f_{x}, f_{y}, v) + H_{2,2} (f_{x}, f_{y}, v)] ψ (v)} \frac{{SNR}_{2} (f_{x}, f_{y})}{1 + {SNR}_{2} (f_{x}, f_{y})},

(23)

and the gray-world reconstruction is given by

G_{2} (f_{x}, f_{y}, v) = F_{2} (f_{x}, f_{y}) ψ (v) .

(24)

Finally, the reconstructions G ₁(f_x ,f_y ,ν) and G ₂(f_x ,f_y ,ν) are combined to yield a final spectral image G ₃(f_x ,f_y ,ν), given by

G_{3} (f_{x}, f_{y}, v) = \frac{1 + {SNR}_{1} (f_{x}, f_{y}, v)}{1 + {SNR}_{1} (f_{x}, f_{y}, v) + c^{- 1} {SNR}_{2} (f_{x}, f_{y})} G_{1} (f_{x}, f_{y}, v)

+ \frac{1 + {SNR}_{2} (f_{x}, f_{y})}{1 + c {SNR}_{1} (f_{x}, f_{y}, v) + {SNR}_{2} (f_{x}, f_{y})} G_{2} (f_{x}, f_{y}, v),

(25)

where c>1 is a weighting parameter used to emphasize the Wiener filtered spectral data G ₁(f_x ,f_y ,ν) over the gray-world reconstruction G ₂(f_x ,f_y ,ν). If G ₁(f_x ,f_y ,ν) and G ₂(f_x ,f_y ,ν) were independent Wiener filter reconstructions of G _o(f_x ,f_y ,ν), then a value of c=1 would yield a G ₃(f_x ,f_y ,ν) with minimum RMS error compared to G _o(f_x ,f_y ,ν) [16

L. P. Yaroslavsky and H. J. Caulfield, “Deconvolution of multiple images of the same object,” Appl. Opt. 33, 2157–2162 (1994). [CrossRef] [PubMed]

]. Since, however, for Fizeau FTIS G ₂(f_x ,f_y ,ν) is a gray-world reconstruction, we use a large value of c (=10³ in Section 4) to weight G ₁(f_x ,f_y ,ν) more heavily than G ₂(f_x ,f_y ,ν) in Eq. (25). With this weight, G ₃(f_x ,f_y ,ν) is essentially equal to G ₁(f_x ,f_y ,ν) throughout the support of the SOTF terms H _1,2(f_x ,f_y ,ν) and H _2,1(f_x ,f_y ,ν), except where SNR ₁(f_x ,f_y ,ν)≪SNR ₂(f_x ,f_y ,ν) and equal to G ₂(f_x ,f_y ,ν) at the remaining spatial frequencies passed by the conventional OTF H(f_x ,f_y ,ν).

4. Simulation results

Simulated FTIS data is used to illustrate the reconstruction algorithm. The optical system is an array of six telescopes in a ring configuration as shown in Fig. 1. Various optical system and detector parameters used in the simulation are given in Tables 1 and 2. The simulation is restricted to short wavelength infrared (SWIR) wavelengths λ=1.95–2.50 µm where the gray-world approximation is often quite good. To perform Fizeau FTIS, the subaperture telescopes are divided into two groups (q=1 or 2) as indicated by the shading in Fig. 1. During data collection the OPD between these groups of subapertures is modulated. Figure 2(a) shows an AVIRIS data set [17

Provided through the courtesy of Jet Propulsion Laboratory, California Institute of Technology, Pasadena, California, http://aviris.jpl.nasa.gov/.

] used for the object S_o(x,y,ν) in the simulation, while Figure 2(b) shows a narrowband reference image of the object S _ref(x,y,ν) used for comparison with reconstruction results later. The reference image was obtained by Wiener filtering the spectral image obtained from a noiseless Michelson FTIS simulation with the same optical system. Thus, S _ref(x,y,ν) is bandlimited to those spatial frequencies passed by the system OTF H(f_x ,f_y ,ν) and has the same spectral sampling and spectral resolution as the Fizeau

Fig. 1. Pupil of six-telescope array optical system, where the telescopes in group q=1 are white and those in group q=2 are gray.

Table 1. Optical System Parameters.

Parameter	Value	Units
Focal length, f	48	m
Encircled pupil diameter	1.7	m
Subaperture telescope diameter	0.50	m
Diameter of central obscurations	0.10	m
Optical efficiency/throughput	60	%

Table 2. Detector Parameters.

Parameter	Value	Units
Pixel pitch, Δ_x	27	µm
Fill factor	100	%
Quantum efficiency^a	0.50	photoelectrons/photon
Read noise	25	photoelectrons
Exposure time^b	0.108	sec

^a Uniform across the spectral bandwidth of simulation, λ=1.95–2.50 µm.

^b Chosen to yield an average of 10⁵ photoelectrons/pixel for each image frame.

Fig. 2. (a). Panchromatic view of object spectral density S _o(x,y,ν) in units of W/m²/sr, and (b) narrowband view of reference imagery S _ref(x,y,ν) in units of W/m²/sr/µm for ν=140 THz (λ=2.13 µm).

Fig. 3. Plot of the average radiance for each spectral band of S _o(x,y,ν). The errorbars represent +/− one standard deviation of the radiance for each spectral band.

FTIS imagery. Figure 3 gives an indication of the validity of the gray-world approximation for the AVIRIS data. As a side note, the AVIRIS data used for S _o(x,y,ν) has a ground sample distance (GSD) of 20m. The required object distance to match this GSD to the detector pixel pitch Δ_x is 35556 km, which is approximately equal to a geosynchronous orbital altitude of 35,786 km for a remote sensing system.

Table 3. FTIS parameters.

Parameter	Value	Units
Time delay sample spacing, Δτ	9.56	fsec
Corresponding OPD sample spacing, cΔτ	2.87	µm
Number of image frames, N_τ	64	-
Spectral resolution, Δν	1.63	THz
Corresponding wavenumber spectral resolution	54.5	cm⁻¹
Corresponding wavelength spectral resolution	20.3–35.4	nm

Fig. 4. Movie (668 KB) of Fizeau FTIS image intensity I(x,y,τ) [the still frame shows I(x,y,0)]. The visibility of the intensity modulation is very small for most of the movie, but is greatest near τ=0, which occurs halfway through the movie. [Media 1]

Table 3 lists the parameter values related to performing Fizeau FTIS. Since the spectral bandwidth of the simulation is limited to the SWIR wavelengths, I(x,y,τ) can be coarsely sampled with respect to τ [18

D. A. Naylor, B. G. Gom, M. K. Tahic, and G. R. Davis, “Astronomical spectroscopy using an aliased, step-and-integrate, Fourier transform spectrometer,” in Millimeter and Submillimeter Detectors for Astronomy II, J. Zmuidzinas, W. S. Holland, and S. Withington, eds., Proc. SPIE 5498-85 (2004). [CrossRef]

]. This reduces the size of the data set, and does not have any SNR penalty in the resulting spectral image for a fixed total data collection time in the shot noise limit [19

S. T. Thurman and J. R. Fienup, “Signal-to-noise ratio trade-offs associated with coarsely sampled Fourier transform spectroscopy,” J. Opt. Soc. Am. A 24, 2817–2821 (2007). [CrossRef]

]. Thus, note that the OPD sample spacing cΔτ=2.87 µm is almost a factor or 3 larger than the Nyquist sample spacing of min(λ)/2=0.975 µm. The total data collection time for making the FTIS measurements is given by t _exp N_τ =6.9 sec (neglecting time required for detector readout and OPD adjustments). The spectral resolution of the AVIRIS data used for S _o(x,y,ν) is Δλ=10 nm, and the FTIS parameters were chosen to yield a somewhat coarser spectral resolution than this to prevent artifacts due to the discrete nature of the simulation. Figure 4 shows a movie of simulated intensity measurements, I(x,y,τ), as a function of τ.

After computing the complex-valued spectral image S _i(x,y,ν), the first step in the reconstruction algorithm is to estimate the gray-world spectrum ψ(ν) from G _i(f_x ,f_y ,ν) by minimizing the objective function E_ψ given by Eq. (13). Figure 5(a) shows the resulting ψ(ν), as well as the average spectrum of the reference spectral image ψ _avg(ν), given by

ψ_{avg} (v) = \frac{1}{N_{x} N_{y}} \sum_{(x, y)} S_{ref} (x, y, v) .

(26)

Intuitively, one might expect that the best reconstruction results would be obtained when these spectra are equal, i.e., that the ideal gray world spectrum equals the average spectrum of the scene ψ(ν)=ψ _avg(ν). This, however, is not universally true. A normalized root mean square error (RMSE) metric can be used to evaluate the accuracy of the gray-world approximation. We define the normalized RMSE at each spatial frequency between the reference data G _ref(f_x ,f_y ,ν) [the 2D spatial Fourier transform of S _ref(x,y,ν)] and a gray-world approximation using ψ(ν) as

Fig. 5. (a). Gray-world spectrum ψ(ν) (solid line) estimated from the spectral image data G _i(f_x ,f_y ,ν) and the average spectrum of the reference spectral image ψ _avg(ν) (dotted line), and (b) azimuthal average of the spatial-frequency domain normalized RMSE E(f_x ,f_y ) corresponding to each gray world spectrum.

E (f_{x}, f_{y}) = \sqrt{\frac{\sum_{v} {∣ G_{ref} (f_{x}, f_{y}, v) - β (f_{x}, f_{y}) ψ (v) ∣}^{2}}{\sum_{v} {∣ G_{ref} (f_{x}, f_{y}, v) ∣}^{2}}},

(27)

where the coefficients β(f_x ,f_y ) are chosen to minimize E(f_x ,f_y ) [11

J. R. Fienup, “Invariant error metrics for image reconstruction,” Appl. Opt. 36, 8352–8357 (1997). [CrossRef]

] and are given by

β (f_{x}, f_{y}) = \frac{[\sum_{v} G_{ref} (f_{x}, f_{y}, v) ψ (v)]}{\sum_{v} {∣ ψ (v) ∣}^{2}} .

(28)

Figure 5(b) shows the azimuthal average of the normalized RMSE E(f_x ,f_y ) using the gray-world spectra shown in Fig. 5(a). While ψ _avg(ν) yields E(0,0)=0, i.e., ψ _avg(ν) matches the true object spectrum at the dc spatial frequency, the normalized RMSE of ψ _avg(ν) at the remaining spatial frequencies is ~15% worse than that obtained with ψ(ν), the gray-world spectrum estimated from the data. Note that the error shown in Fig. 5(b) increases dramatically for spatial frequencies beyond 0.3 cycles/pixel, which are outside the support of the system OTF.

The next step in the reconstruction algorithm is to estimate the object power spectrum model parameters A and α and the noise power spectrum Φ_n by minimizing E _Φ. Figure 6 shows comparisons of |G_i ^(Re)(f_x ,f_y ,ν)|² and |G _avg(f_x ,f_y )|² with the corresponding power spectrum model using the estimated model parameters. Notice that the smooth power spectrum model fits the spatial frequency data well on average. The power spectrum model is used only to estimate the Fourier-domain SNR quantities SNR₁(f_x ,f_y ,ν) and SNR₂(f_x ,f_y ), and an exact match between the power spectrum model and the Fourier data is not necessary. As a side note, the estimated value of α=1.07 is within the expected two-sigma range of α=1.2 ±0.3 for natural scenery [14

D. J. Tolhurst, Y. Tadmor, and T. Chao, “Amplitude spectra of natural images,” Ophthalic. Physiol. Opt. 12, 229–232 (1992). [CrossRef]

Next, the SNR estimates SNR ₁(f_x ,f_y ,ν) and SNR ₂(f_x ,f_y ) are used to compute the Wiener and gray-world reconstructions S ₁(x,y,ν) and S ₂(x,y,ν), respectively. These two estimates are then combined using Eq. (25) to form a final reconstruction S ₃(x,y,ν). Figure 7 shows each of these reconstructions in comparison to the unprocessed spectral image S _i ^(Re)(x,y,ν) at a particular spectral band. For comparison, the corresponding reference image S _ref(x,y,ν) is shown in Fig. 2(b). The missing low spatial-frequencies cause the amplitude of S _i ^(Re)(x,y,ν) to be appreciable only at edges within the scene, as is evident in Fig. 7(a). The Wiener filter sharpens the spectral image considerably, but cannot reconstruct the missing low spatial frequencies. Thus, S ₁(x,y,ν) is still zero mean and contains low spatial-frequency artifacts as seen in Fig. 7(b). The gray-world reconstruction S ₂(x,y,ν), shown in Fig. 7(c), has the general appearance of a conventional image [the small negative values in S ₂(x,y,ν) are Gibb’s ringing artifacts], but does not have much spectral utility since the spectrum at each pixel is identical. The reconstruction S ₃(x,y,ν) contains desirable properties from both S ₁(x,y,ν) and S ₂(x,y,ν), namely the features of conventional imagery and spectral utility.

Fig. 6. Power spectrum fitting results for: (a) |G _i ^(Re)(0,f_y ,ν)|² data (solid line) and the power spectrum fit (dotted line) for ν=140 THz and (b) |G _avg(f_x ,0)|² (solid line) and the corresponding gray-world power spectrum model.

Fig. 7. Spectral image data: (a) S _i ^(Re)(x,y,ν), (b) S ₁(x,y,ν), (c) S ₂(x,y,ν), and (d) S ₃(x,y,ν) for ν=140THz. The corresponding reference image S _ref(x,y,ν) is shown in Fig. 2(b).

Fig. 8. Normalized RMSE [calculated with respect to S _ref(x,y,ν)] for each image point for: (a) S _i ^(Re)(x,y,ν), (b) S ₁(x,y,ν), (c) S ₂(x,y,ν), and (d) S ₃(x,y,ν).

Figures. 8, 9, and 10(a) shows the normalized RMSE, between S _ref(x,y,ν) and each of S _i ^(Re)(x,y,ν), S ₁(x,y,ν), S ₂(x,y,n) and S ₃(x,y,ν) evaluated at each image point (x,y), at each spatial frequency (f_x ,f_y ), and for each spectral band ν, respectively. In each case, the normalized RMSE η is defined in a manner similar to Eq. (27). Additionally, Fig. 10(b) shows the standard RMSE for each spectral band of the various reconstructions and Table 4 lists the total reconstruction errors for each reconstruction. The quality of both S _i ^(Re)(x,y,ν) and the Wiener reconstruction S ₁(x,y,ν) appears to be rather poor in terms of η evaluated for each image point and each spectral band. Figure 9, however, indicates that the spectra at spatial frequencies passed by the SOTF terms H _1,2(f_x ,f_y ,ν) and H _2,1(f_x ,f_y ,ν) is fairly good. The quality of the gray-world reconstruction S ₂(x,y,ν) is fairly good in terms of η, but does not offer much spectral utility. It is interesting to note that better η values were obtained in the spatial frequency domain for G ₂(f_x ,f_y ,ν) than for G ₁(f_x ,f_y ,ν). We believe this results from the combination of three effects: (i) the error associated with the gray-world approximation is

Fig. 9. Normalized RMSE [calculated with respect to G _ref(f_x ,f_y ,ν)] for each spatial frequency for: (a) G _i ^(Re)(f_x ,f_y ,ν), (b) G ₁(f_x ,f_y ,ν), (c) G ₂(f_x ,f_y ,ν), and (d) G ₃(f_x ,f_y ,ν).

Fig. 10. (a). Normalized RMSE and (b) standard RMSE [calculated with respect to S _ref(x,y,ν)] for each spectral band for: S _i ^(Re)(x,y,ν), S ₁(x,y,ν), S ₂(x,y,ν), and S ₃(x,y,ν).

Table 4. Overall Error.

Reconstruction	RMSE [W m⁻² sr⁻¹ µm⁻¹]	Normalized RMSE
S _i ^(Re)(x,y,ν)	1.45	0.964
S ₁(x,y,ν)	1.38	0.950
S ₂(x,y,ν)	0.21	0.145
S ₃(x,y,ν)	0.15	0.102

Fig. 11. Comparison of Fizeau and Michelson FTIS reconstruction results: (a) RMSE and (b) fractional RMSE [calculated with respect to S _ref(x,y,ν)] for each spectral band.

smaller than the errors due to noise, (ii) when estimating the gray-world spectrum ψ(ν) via Eq. (13), the noise is reduced by effectively averaging over many spatial frequencies, and (iii) when using I _avg(x,y) to compute S ₂(x,y,ν), noise is also reduced by averaging over many time delays. Comparing the error metric values, the quality of the final reconstruction S ₃(x,y,ν) is better than that of S ₂(x,y,ν).

Finally, a Michelson FTIS simulation was performed for the same optical system and viewing conditions. The Michelson FTIS data was reconstructed using just the Wiener filter portion of the algorithm described in Section 3 since it does not suffer from missing low spatial frequencies. Figure 11 compares the RMSE for the final reconstructions obtained from the Fizeau and Michelson FTIS simulations. Note that the RMSEs are comparable at some wavelengths and the Michelson is up to two times better than the Fizeau at other wavelengths. For comparison, the Hyperion instrument, which is part of NASA’s Earth Observing-1 Mission [20

http://eo1.gsfc.nasa.gov/miscPages/home.html

], has a spectrally-averaged SNR of approximately 140 for the visible and 60 for the SWIR spectral regions [21

R. O. Green, B. E. Pavri, and T. G. Chrien, “On-orbit radiometric and spectral calibration characteristics of EO-1 Hyperion derived with an underflight of AVIRIS and in situ measurements at Salar de Arizaro, Argentina,” IEEE Trans. in Geosci. Remote Sens. , 41, 1194–1203 (2003). [CrossRef]

]. This implies that the fractional RMSE for Hyperion is approximately 2% in the SWIR region, compared with the fractional RMSEs of ~10% obtained from the Fizeau and Michelson FTIS simulations for the particular parameters we chose.

5. Discussion and summary

Fizeau FTIS poses a challenging image reconstruction problem because of missing low spatial frequencies. Unlike the optical superresolution problem [22

B. R. Hunt, “Super-resolution of images: algorithms, principles, performance,” Int. J. Imaging Syst. Technol. 6, 297–304 (1995). [CrossRef]

], in which one attempts to reconstruct high spatial-frequency information beyond the diffraction limit through extrapolation (which is very sensitive to noise), we wish to interpolate and reconstruct missing low spatial frequencies in the Fizeau FTIS problem. Furthermore, panchromatic information about the low spatial frequencies is available from the raw Fizeau FTIS measurements. The algorithm proposed here uses a gray-world approximation, which is typically applicable to natural scenery as long as a significant contrast reversal of the scene does not occur within the spectral bandwidth of the imagery. For example, while the gray-world approximation might be valid for imagery spanning the visible or near infrared wavelengths, the gray-world approximation is usually not valid for imagery containing vegetation that spans both the visible and near infrared wavelengths due to large contrast reversals associated with the chlorophyll reflectance spectrum. To overcome this problem, Fizeau FTIS could be performed twice for the same scene, using color filters to limit the spectral bandwidth of the imagery first to visible wavelengths and then to near infrared wavelengths.

Results with numerically simulated data indicate that moderate quality spectral imagery can be obtained through Fizeau FTIS using the proposed algorithm. The Fizeau FTIS results were found to be within approximately a factor of two of Michelson FTIS results (in terms of the estimated spectral SNR. The algorithm does assume knowledge of the system PSF, which may not be available in all cases, and the simulations did not explore the impact of incomplete knowledge or errors in the OPDs introduced for performing FTIS. While the performance of Fizeau FTIS is slightly less than Michelson FTIS, the ability of Fizeau FTIS to gather spectral information without the need for additional hardware makes it attractive for segmented-aperture or multiple telescope array systems.

Acknowledgement

This work was supported by Lockheed Martin Corporation.

References and links

1.	J. Kauppinen and J. Partanen, Fourier Transforms in Spectroscopy , (Wiley-VCH, Berlin, 2001). [CrossRef]
2.	N. J. E. Johnson, “Spectral imaging with the Michelson interferometer,” in Infrared Imaging Systems Technology, W. L Wolfe and J. Zimmerman, eds., Proc. SPIE 226, 2–9 (1980).
3.	C. L. Bennett, M. Carter, D. Fields, and J. Hernandez, “Imaging Fourier transform spectrometer,” in Imaging Spectrometry of the Terrestrial Environment, G. Vane, ed., Proc. SPIE 1937, 191–200 (1993). [CrossRef]
4.	M. R. Carter, C. L. Bennett, D. J. Fields, and F. D. Lee, “Livermore imaging Fourier transform infrared spectrometer,” in Imaging Spectrometry, M. R. Descour, J. M. Mooney, D. L. Perry, and L. R. Illing, eds., Proc. SPIE 2480, 380–386 (1995). [CrossRef]
5.	M. Frayman and J. A. Jamieson, “Scene imaging and spectroscopy using a spatial spectral interferometer,” in Amplitude and Intensity Spatial Interferometry, J. B. Breckingridge, ed., Proc. SPIE 1237, 585–603 (1990). [CrossRef]
6.	R. L. Kendrick, E. H. Smith, and A. L. Duncan, “Imaging Fourier transform spectrometry with a Fizeau interferometer,” in Interferometry in Space, M. Shao, ed., Proc. SPIE 4852, 657–662 (2003). [CrossRef]
7.	S. T. Thurman and J. R. Fienup, “Multi-aperture Fourier transform imaging spectroscopy: theory and imaging properties,” Opt. Express 13, 2160–2175 (2005). [CrossRef] [PubMed]
8.	J. Goodman, Introduction to Fourier Optics 2nd ed. , (McGraw-Hill, New York, 1996).
9.	Note that the additional factor of A/(λ ² f ² _i ) at the beginning of Eq. (9) was omitted from Eqs. (9), (16), and (17) of Ref. [6] in error.
10.	C. W. Helstrom, “Image restoration by the method of least squares,” J. Opt. Soc. Am. 57, 297–303 (1967). [CrossRef]
11.	J. R. Fienup, “Invariant error metrics for image reconstruction,” Appl. Opt. 36, 8352–8357 (1997). [CrossRef]
12.	G. J. Burton and I. R. Moorhead, “Color and spatial structure in natural scenes,” Appl. Opt. 26, 157–170 (1987). [CrossRef] [PubMed]
13.	D. J. Field, “Relations between the statistics of natural images and the response properties of cortical cells,” J. Opt. Soc. Am. A 4, 2379–2394 (1987). [CrossRef] [PubMed]
14.	D. J. Tolhurst, Y. Tadmor, and T. Chao, “Amplitude spectra of natural images,” Ophthalic. Physiol. Opt. 12, 229–232 (1992). [CrossRef]
15.	D. L. Ruderman and W. Bialek, “Statistics of natural images: scaling in the woods,” Phys. Rev. Lett. 73, 814–817 (1994). [CrossRef] [PubMed]
16.	L. P. Yaroslavsky and H. J. Caulfield, “Deconvolution of multiple images of the same object,” Appl. Opt. 33, 2157–2162 (1994). [CrossRef] [PubMed]
17.	Provided through the courtesy of Jet Propulsion Laboratory, California Institute of Technology, Pasadena, California, http://aviris.jpl.nasa.gov/.
18.	D. A. Naylor, B. G. Gom, M. K. Tahic, and G. R. Davis, “Astronomical spectroscopy using an aliased, step-and-integrate, Fourier transform spectrometer,” in Millimeter and Submillimeter Detectors for Astronomy II, J. Zmuidzinas, W. S. Holland, and S. Withington, eds., Proc. SPIE 5498-85 (2004). [CrossRef]
19.	S. T. Thurman and J. R. Fienup, “Signal-to-noise ratio trade-offs associated with coarsely sampled Fourier transform spectroscopy,” J. Opt. Soc. Am. A 24, 2817–2821 (2007). [CrossRef]
20.	http://eo1.gsfc.nasa.gov/miscPages/home.html
21.	R. O. Green, B. E. Pavri, and T. G. Chrien, “On-orbit radiometric and spectral calibration characteristics of EO-1 Hyperion derived with an underflight of AVIRIS and in situ measurements at Salar de Arizaro, Argentina,” IEEE Trans. in Geosci. Remote Sens. , 41, 1194–1203 (2003). [CrossRef]
22.	B. R. Hunt, “Super-resolution of images: algorithms, principles, performance,” Int. J. Imaging Syst. Technol. 6, 297–304 (1995). [CrossRef]

OCIS Codes
(070.2580) Fourier optics and signal processing : Paraxial wave optics
(100.3020) Image processing : Image reconstruction-restoration
(110.4850) Imaging systems : Optical transfer functions
(110.6770) Imaging systems : Telescopes
(120.3180) Instrumentation, measurement, and metrology : Interferometry
(120.6200) Instrumentation, measurement, and metrology : Spectrometers and spectroscopic instrumentation

ToC Category:
Fourier optics and signal processing

History
Original Manuscript: December 3, 2007
Revised Manuscript: March 24, 2008
Manuscript Accepted: March 25, 2008
Published: April 25, 2008

Citation
Samuel T. Thurman and James R. Fienup, "Fizeau Fourier transform imaging spectroscopy: missing data reconstruction," Opt. Express 16, 6631-6645 (2008)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-16-9-6631

Sort: Author | Year | Journal | Reset

References

J. Kauppinen and J. Partanen, Fourier Transforms in Spectroscopy (Wiley-VCH, Berlin, 2001). [CrossRef]
N. J. E. Johnson, "Spectral imaging with the Michelson interferometer," in Infrared Imaging Systems Technology, W. L Wolfe and J. Zimmerman, eds., Proc. SPIE 226, 2-9 (1980).
C. L. Bennett, M. Carter, D. Fields, and J. Hernandez, "Imaging Fourier transform spectrometer," in Imaging Spectrometry of the Terrestrial Environment, G. Vane, ed., Proc. SPIE 1937, 191-200 (1993). [CrossRef]
M. R. Carter, C. L. Bennett, D. J. Fields, and F. D. Lee, "Livermore imaging Fourier transform infrared spectrometer," in Imaging Spectrometry, M. R. Descour, J. M. Mooney, D. L. Perry, and L. R. Illing, eds., Proc. SPIE 2480, 380-386 (1995). [CrossRef]
M. Frayman and J. A. Jamieson, "Scene imaging and spectroscopy using a spatial spectral interferometer," in Amplitude and Intensity Spatial Interferometry, J. B. Breckingridge, ed., Proc. SPIE 1237, 585-603 (1990). [CrossRef]
R. L. Kendrick, E. H. Smith, and A. L. Duncan, "Imaging Fourier transform spectrometry with a Fizeau interferometer," in Interferometry in Space, M. Shao, ed., Proc. SPIE 4852, 657-662 (2003). [CrossRef]
S. T. Thurman and J. R. Fienup, "Multi-aperture Fourier transform imaging spectroscopy: theory and imaging properties," Opt. Express 13, 2160-2175 (2005). [CrossRef] [PubMed]
J. Goodman, Introduction to Fourier Optics 2nd ed. (McGraw-Hill, New York, 1996).
Note that the additional factor of A/(λ²fi²) at the beginning of Eq. (9) was omitted from Eqs. (9), (16), and (17) of Ref. [6] in error.
C. W. Helstrom, "Image restoration by the method of least squares," J. Opt. Soc. Am. 57, 297-303 (1967). [CrossRef]
J. R. Fienup, "Invariant error metrics for image reconstruction," Appl. Opt. 36, 8352-8357 (1997). [CrossRef]
G. J. Burton and I. R. Moorhead, "Color and spatial structure in natural scenes," Appl. Opt. 26, 157-170 (1987). [CrossRef] [PubMed]
D. J. Field, "Relations between the statistics of natural images and the response properties of cortical cells," J. Opt. Soc. Am. A 4, 2379-2394 (1987). [CrossRef] [PubMed]
D. J. Tolhurst, Y. Tadmor, and T. Chao, "Amplitude spectra of natural images," Ophthalmic Physiol. Opt. 12, 229-232 (1992). [CrossRef]
D. L. Ruderman and W. Bialek, "Statistics of natural images: scaling in the woods," Phys. Rev. Lett. 73, 814-817 (1994). [CrossRef] [PubMed]
L. P. Yaroslavsky and H. J. Caulfield, "Deconvolution of multiple images of the same object," Appl. Opt. 33, 2157-2162 (1994). [CrossRef] [PubMed]
Provided through the courtesy of Jet Propulsion Laboratory, California Institute of Technology, Pasadena, California, http://aviris.jpl.nasa.gov/.
D. A. Naylor, B. G. Gom, M. K. Tahic, and G. R. Davis, "Astronomical spectroscopy using an aliased, step-and-integrate, Fourier transform spectrometer," in Millimeter and Submillimeter Detectors for Astronomy II, J. Zmuidzinas, W. S. Holland, and S. Withington, eds., Proc. SPIE 5498-85 (2004). [CrossRef]
S. T. Thurman and J. R. Fienup, "Signal-to-noise ratio trade-offs associated with coarsely sampled Fourier transform spectroscopy," J. Opt. Soc. Am. A 24, 2817-2821 (2007). [CrossRef]
http://eo1.gsfc.nasa.gov/miscPages/home.html.
R. O. Green, B. E. Pavri, and T. G. Chrien, "On-orbit radiometric and spectral calibration characteristics of EO-1 Hyperion derived with an underflight of AVIRIS and in situ measurements at Salar de Arizaro, Argentina," IEEE Trans. in Geosci. Remote Sens., 41, 1194-1203 (2003). [CrossRef]
B. R. Hunt, "Super-resolution of images: algorithms, principles, performance," Int. J. Imaging Syst. Technol. 6, 297-304 (1995). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.

Click here to see a list of articles that cite this paper