1.1.

Wide Bandwidth Desired

Wide bandwidth high-resolution spectroscopy is one of the most valuable diagnostic tools for astronomy, from Doppler radial velocimetry in the search for exoplanets to studying the characteristics of stars through the size, breadth, and position of spectral features. A wide bandwidth is needed in radial velocimetry to collect as many photons as possible to reduce Doppler velocity noise. Wide bandwidth is desired for analysis of species population because the characteristic wavelengths are often distributed over a wide band.

However, it is difficult to achieve both wide bandwidth and high spectral resolution for a dispersive spectrograph, since this requires not only a large and expensive instrument but also a detector with many pixels. This is especially so for infrared wavelengths since the spectrograph size scales with resolving power times wavelength.

1.2.

Stability Against Point Spread Function Drifts Desired

Second, a serious challenge for high-resolution spectrographs is the unwanted shifting of the focal spot (point spread function, PSF) along the dispersion direction, such as due to air convention or thermomechanical drifts. While use of a calibrant spectrum measured in parallel or comingled with the stellar spectrum can remove some lower order components of the drift, the higher order components corresponding to gaps between strong features in the calibrant spectrum can remain, producing a wavenumber error. These environmental insults are frequently the limitation, rather than photon noise, to achieving the highest desired Doppler precision, or sufficiently accurate spectra to use as a template for subsequent Doppler calculations.

1.3.

Externally Dispersed Interferometry

A technique called externally dispersed interferometry¹ (EDI) can dramatically reduce the necessary size of spectrographs for Doppler radial velocimetry^1–9 and high-resolution spectroscopy.^10–15 Other workers have adopted the EDI method from our laboratories and demonstrated a Doppler planet detection.⁵ We have demonstrated EDI in both applications in the near-infrared at the Hale 5-m telescope at Mt. Palomar Observatory (Fig. 1) with a version called TEDI, (the T denotes the TripleSpec¹⁶ spectrograph to which it was coupled). An earlier report² describes the Doppler velocimetry, and this report describes the spectroscopy.

Fig. 1

(a) The T-EDI interferometer (black and silver) sits atop TripleSpec (native) spectrograph (blue cryogenic cylinder) which bolts to Cassegrain output (39 in. diameter cavity) of the Mt. Palomar 200 in. mirror. (b) TEDI 2.0¹ intercepts starlight beam (1) heading for TripleSpec and uses a dichroic split input tracking system (3) to mixing in ThAr lamp (5) to feed a symmetric, collimated interferometry cavity (beamsplitter 10). One interferometer output (13) is dichroic split (14) to a sensitive chopped IR diode (15) for flux maximization and a fringe tracking camera (16) for cavity nulling. The complementary interferometer output (17) is open-path relayed (18, 19) to the TripleSpec slit. (c) Externally dispersed interferometer (EDI) scheme—a fixed delay interferometer is crossed with a disperser. The sinusoidal transmission of the interferometer creates moiré patterns whose phase encodes the Doppler velocity and whose shape encodes the spectrum’s shape. (Schematic graphic from Ref. 1 reproduced with permission © The Astronomical Society of the Pacific. All rights reserved.)

This is the first time multiple delays have been used on starlight to recreate a high-resolution spectrum. Previous EDI use on stellar spectroscopy used a single delay.¹¹ Multiple delays had been used previously¹² only on laboratory sources.

The TEDI project was designed primarily to test Doppler measurements, and the selection of interferometer glass etalons and observing schedule were prioritized for this purpose. These were not optimal for testing high-resolution spectroscopy (having a gap in coverage of delay space when contiguous coverage is needed to avoid ringing in the lineshape). Nevertheless, the preliminary results for high-resolution spectroscopy using multiple delays are a resounding success in two important areas.

1.4.

Resolution, Robustness Demonstrated

We have demonstrated on starlight that the EDI technique can produce both high resolution, a factor of $4\times$ to $10\times$ beyond the native spectrograph resolution of $\sim 2700$, and have extremely wide simultaneous bandwidth, limited only by the bandwidth of the native spectrograph, in this case (0.95 to $2.45\mu \mathrm{m}$).

Second, the native spectrograph in our study suffered from extremely large and irregular PSF drifts that would normally have precluded the high resolution we obtained, even if by some other means the slit width and focal spot had been reduced by several times and the detector pixel density increased. We realize now that a great advantage of the EDI technique, especially when used with multiple delays rather than a single delay as used previously, is an order of magnitude or more improvement in the stability against PSF drifts.

We observe an approximate 20-fold decrease in the translational reaction of a ThAr lamp line to a translational insult along the dispersion direction. This robustness to PSF errors is explained theoretically and demonstrated with measurements. A recently realized method of further reducing the PSF shift sensitivity to zero is described theoretically and demonstrated in a simple simulation which produces a 350× times reduction.

Third, we demonstrate EDI’s dramatic robustness to fixed pattern (FP) noise such as from bad detector pixels. These pollute the native spectrum but not the EDI-derived spectrum measured simultaneously.

The TEDI data set is instrument limited rather than photon limited. A theoretical description of EDI photon-limited sensitivity compared to other spectroscopy techniques is presented in a companion paper.

2.1.

Externally Dispersed Interferometry Has Corrugated Point Spread Function

The EDI PSF¹⁰ (instrument response) is a corrugated peak or wavelet (bold curve of Fig. 3), having the envelope set by the native spectrograph (dashed curve), multiplied by a sinusoid set by the interferometer delay chosen by the user. The key point is that now the PSF slope can be much higher than the low slope of envelope (native spectrograph), being proportional to the interferometer delay. It is easy to implement the 1- to 5-cm delay needed to match the slope of a high-resolution spectrograph.

Fig. 3

The EDI PSF (instrument response) is a corrugated peak (bold curve), which has a much higher slope (and thus higher Doppler sensitivity) than the low resolution native spectrograph alone (dashes). The slope of the sinusoid is proportional to the interferometer delay. For use in high-resolution spectroscopy, the subject of this report, the central area (cross-hatched) is made unambiguous by combining several wavelets of different periodicities (delays), as shown in Figs. 4 and 5.

The PSF slope is now much more mathematically controlled, largely independent of the envelope shape. The interferometer has only three degrees of freedom (phase, period, and intensity offset), rather than the many degrees of freedom typical of a grating PSF (contributed by the aggregate of tens or hundreds of grating grooves). Second, since the interferometer component is small ($\sim 50\mathrm{cm}$), it is much easier to put in a stable environment than a stand-alone disperser, which for a 50,000 resolution can be several meters in size.

2.2.

Combining Several Delays Makes a Peaked Point Spread Function

For use in high-resolution spectroscopy, the subject of this report, the many redundant peaks in the corrugated PSF for a single delay would be a problem. The solution is to use several delays. The central area (cross-hatched) is made unambiguous by combining several wavelets of different periodicities (delays), as demonstrated in simulation in Fig. 4, and in practice in Fig. 5. The latter shows measurements of two ThAr lamp lines resolved by the EDI technique, but which otherwise are unresolved by the native spectrograph.

Fig. 4

Simulated summation of several simulated wavelets, weighted in Gaussian equalization scheme, produces a single Gaussian PSF peak with minimal ringing (bold red) and thinner than the native spectrograph (green curve) by factors of several.

Fig. 5

(a) Reconstruction of an otherwise unresolved ThAr doublet (7556, $7557.6{\mathrm{cm}}^{-1}$) from many wavelets measured at multiple delays. (b) The native spectrograph (green curve) cannot resolve the doublet. The EDI reconstructed spectrum (red curve), which is sum of wavelets, then equalized to $\mathrm{RES}=\mathrm{16,000}$, easily resolves the doublet.

3.1.

Data Analysis Overview

The data analysis pipeline is described in several sections:

4.1.

Two Stripes (A, B) for Each Echelle Order

After optional flat fielding to remove pixel-to-pixel gain variations, the 2-D data from echelle spectrograph detector, which are an FITS file, are converted to 16 1-D spectra ($4\text{fibers}\times 4\text{orders}$) by collapsing vertically across each of four stripes for each order. There are four useful orders, enumerated A to D starting from the top (high $Y$ pixel), and the fifth (E) is weak and mostly redundant in range with the D. These have wavenumber (wn or $\nu$) ranges of A: 4100 to $5300{\mathrm{cm}}^{-1}$, B: 5400 to 7100, C: 6700 to 8800, D: 8100 to 10,600, and E: 9400 to 10,700. These are wavelength ranges of A: 2.44 to $1.89\mu \mathrm{m}$, B: 1.85 to 1.41, C: 1.49 to 1.14, D: 1.23 to 0.94, and E: 1.06 to 0.93.

In the first version of the TEDI apparatus,⁶ both complementary interferometer outputs were sent through the spectrograph and detected. However, that interferometer design suffered from noncommon path errors. We rebuilt it in 2010 (version 2) using collimated beams, which eliminated the path errors. We also decided to use one of the two interferometer outputs to monitor the flux for maximizing starlight alignment at the input side (beam path described in Ref. 2). Although this prevents summing the two complementary outputs to create a nonfringing ordinary spectrum for each exposure, the latter is created by summing the 10 phase stepped exposures in the phase-stepping math process (Sec. 5) to create one nonfringing spectrum per set of exposures. The increase in flux due to better alignment of the starlight into the input fiber made up for the loss of half the outputs.

Hence we recorded four stripes per order, which we denoted sky A, fiber A, fiber B, and sky B. But since the companion sky stripes (which sensed blackbody background) were soon optionally subtracted from the data-carrying fiber stripes, we effectively had two stripes per order. (For the bluer orders, we could omit the subtraction.) The fibers A and B could contain either ThAr spectral lamp or starlight combined with ThAr light. We explored several data taking configurations. Typically, one fiber contained ThAr-only and the other ThAr + star. Often these roles were alternated during the night.

4.2.

Fringe Comb Visible for Small Delay

Figure 6 shows a portion of an FITS echelle spectrograph output around the C-order, and Fig. 7 shows three closeups showing the effect of interferometer delay stepping. In Figs. 6 and 7, note the periodic fringes imprinted on the stellar spectrum by interferometer acting on the continuum portion of the input spectrum. A small interferometer delay ($\sim 0.1\mathrm{cm}$, “E1”) was used here, which produces fringes wide enough to be resolvable by the native spectrograph (res $\sim 2700$). The sinusoidal fringes are not required to be resolvable, and for most of the higher delays ($>0.3\mathrm{cm}$), are not. The heterodyning occurring here creates moiré patterns that manifest high frequency information beaten down to lower, resolvable frequencies.

Fig. 6

Example C-order of the 2-D echelle fringing spectrum covering 6700 to $8800{\mathrm{cm}}^{-1}$ (left to right) over 2048 pixels. The three other orders (A, B, D) spanning 4100 to $\mathrm{11,000}{\mathrm{cm}}^{-1}$ (0.9 to $2.45\mu \mathrm{m}$) are shown in Fig. 2 of Ref. 2. Fiber A (white curves) is ThAr lamp alone, fiber B (red curves) is Star + ThAr. Note the sinusoidal comb from the continuum. Data from 2010Sep19, #0756, star HD962, etalon E1 ($\sim 0.1\mathrm{cm}$).

Fig. 7

Three closeups (a), (b), and (c) of the echelle C-order for ThAr lamp + star HD962, as interferometer delay changes for the E1 etalon ($\sim 0.1\mathrm{cm}$). Wavenumber is $\sim 7525{\mathrm{cm}}^{-1}$ or $1.33\mu \mathrm{m}$, so the $0.25\mu \mathrm{m}$ commanded delay change creates a 0.19 cycle phase change. Horizontal scaling $\sim 1{\mathrm{cm}}^{-1}$ per pixel. Each feature’s intensity varies sinusoidally with interferometer phase, but having various relative phases and magnitudes to neighboring features.

Figure 7 is a closeup in C-order of three exposures of the same set as Fig. 6 but having different phase. The top fiber (A) has only ThAr spectral lamp, and the bottom fiber (B) is a sum of stellar and ThAr. The change in fringe phase between the two exposures can be seen in the different intensities of the ThAr lines, comparing top to bottom panel, and in the apparent sideways shift of the periodic comb in fiber B stripe.

4.3.

Summing Across Stripes

For each of the four stripes of an order, the counts are summed across 15 vertical pixels to form a set of 1-D spectra outputted from stage 1. The pairs of white and red curves in the figures demark the 15-pixel vertical width of the stripes.

The stripe shape does not change significantly from night to night, and it merely needs to be offset vertically. The vertical summing does not correct for the effect of the slight slant in the spectrograph slit across the order, but the blurring this creates is insignificant since most of the counts are concentrated in 3 or 4 vertical pixels. However, the effect of the slant between the A and B stripes is more significant. Our routines compensate for this slight wavenumber shift between the two fibers in the pixel to wavenumber mapping done in stage 3, by having individual calibrations for each fiber. (However, our sky subtraction is currently unshifted.)

4.4.

Detector Artifact Cleaning and Flatfielding

An optional first step, before order extraction, is to run a routine that removes electronic ghosting artifacts specific to this detector. This routine was supplied by TripleSpec spectrograph group at Cornell.¹⁶ Then an optional flatfielding exposure divides out some FP noise due to pixel gain variation. Flatfielding is not critical for the fringing component since a built-in flatfielding-like process occurs already during phase-stepping mathematics that rejects stationary patterns. However, flatfielding benefits the ordinary spectrum.

4.5.

Stage 1: Output Is Two-Dimensional Array Stack Versus Delay and $\nu$

For each of the stripes, the 10 sequentially recorded phase-stepped 1-D spectra $B_n(\nu )$ are stacked into an array versus exposure number along $Y$-axis. Since exposure number is essentially delay, which is related to phase, this forms a pseudomultiphase format (Fig. 8). These arrays of 1-D uniphase data are thus in the same format as 2-D multiphase data of earlier EDIs (such as Fig. 5 of Ref. 1), where the phase spatially varies along slit length in each exposure. This is a visually intuitive and algorithmically convenient way of presenting the phase set data, and routines written earlier for multiphase data can be applied.

Fig. 8

The 10 sequentially recorded phase stepped 1-D spectra are stacked versus exposure # ($\sim \text{phase}$) along the vertical. (a) Small delay (E1, 0.1 cm), the comb period is wide enough to be resolved. (b) Larger delay (E3, 0.7 cm), the comb is unresolved. ThAr phase pattern constitutes a “fingerprint” that reveals the precise interferometer delay. Un-normalized intensity fluctuations due to atmospheric seeing cause horizontal artifacts in (b). Horizontal scaling $\sim 1{\mathrm{cm}}^{-1}$ per pixel.

Echelle spectrographs require uniphase mode because they usually do not have enough vertical pixels across an order to support varying phase spatially. Instead, the phase is varied temporally via interferometer mirror movement in 10 exposures. But there are some advantages to the uniphase mode that makes it attractive in any case. In uniphase mode, the same pixel is used repetitively, so FP noises are easily suppressed (because each pixel has a constant gain). Second, it allows using the spatial dimension to splay across the target, for imaging spectroscopy.

Comparing the two different delays shown in Figs. 8(a) and 8(b) note that the pattern of phase shifts for the ThAr lines is different. These constitute a “fingerprint” that reveals during analysis the precise delay that was used for each exposure, so we do not need to rely on knowledge of the commanded interferometer mirror position.

4.6.

Cosmic Ray Removal

For EDI, it is more reliable to do the cosmic ray (CR) removal when the data are in pseudomultiphase format, rather than earlier when individual exposures are in 2-D format. (In the latter, a bright ThAr line can be mistaken as a CR by a median filter.) We fit a sinusoid along each wavenumber column versus the phase. Unusually large residuals are deemed CRs and their values replaced by the sinusoidal fit value.

4.7.

Wavenumber Map

A pixel to wavenumber map is made from the ThAr line positions. This is used in stage 2 to convert from delay step to phase step through $\mathrm{\Delta }\varphi =\mathrm{\Delta }\tau /\lambda =\nu \mathrm{\Delta }\tau /{10}^4$, where phase $\varphi$ is in units of cycles or fringes, $\lambda$ is wavelength in microns, $\tau$ the delay in microns, and $\nu$ is wavenumber in ${\mathrm{cm}}^{-1}$.

5.1.

Brief Theory

Reference 11 describes single delay EDI theory regarding spectroscopy. Fringing $\mathbf{W}$ and nonfringing $B_{\mathrm{ord}}$ components can be computed from a series of phase stepped exposures $B_n$ of regular phase ${\varphi }_n$ (evenly distributed around the phase circle) by

5.2.

Overview

We now describe the phase-stepping algorithm we use to separate fringing ($\mathbf{W}$) and nonfringing ($B_{\mathrm{ord}}$) spectra from the spectral array data ($B_n$) of Fig. 8. This algorithm can handle unknown and irregular steps and intensities, and therefore it is useful to other interferometry applications generally (the first author uses it fruitfully for analyzing Doppler velocimetry interferograms having irregular visibilities and intensities in shock physics¹⁷). This distinguishes our algorithm from most others in the literature,¹⁸ which have trouble handling unknown and irregular steps.

5.2.1.

Alternate solving X and Y dependencies

The problem we face is to solve both $X$ and $Y$ behaviors. The $X$ axis is wavenumber $\nu$, and the $Y$ axis is delay step position $\delta {\tau }_n$. We alternate between analyzing data along a column using provisional knowledge of delay steps to produce “template” spectra ($\mathbf{W}$ and $B_{\mathrm{ord}}$), and then analyzing along rows using provisional knowledge of the template spectra to produce detailed delay step values (flowchart Fig. 9).

Fig. 9

Flowchart of iterative strategy for finding phases and amplitudes of each row (phase step, $Y$) of the interferogram, while also finding the ordinary spectrum (nonfringing intensity NF) and the fringing spectrum value $W(\nu )$ for each column (each $X$ or $\nu$). ${\text{Slit}}_y$ is the phase (or channel) dependent intensity variations, due to atmospheric seeing, which we normalize out.

The observed TEDI delay steps [Fig. 10(b)] are reasonably uniform so the initial calculation of the provisional templates by assuming constant $0.25\mu \mathrm{m}$ steps produces such an accurate result for W and B that the whole process converges rapidly in only 1 to 3 iterations, using only a single order. (A much tighter constraint on the delay step would result with the increased bandwidth of multiorder processing in parallel, but this has yet to be implemented.)

Fig. 10

Example diagnostics during phase-stepping math, etalon E4 ($\sim 1\mathrm{cm}$) A-order data. (a) Delay values (blue dots) $\delta \tau$ found from ThAr data. Gray line is initial guess having $0.25\mu \mathrm{m}$ steps. (b) Interval in delay $\mathrm{\Delta }{\tau }_n=\delta {\tau }_{n+1}-\delta {\tau }_n$. Deviations typically $\sim 0.01\mu \mathrm{m}$. Our algorithm does not require uniform steps, but larger deviations flag presence of CRs or other data failure. (c) Phase at $x=1000$, showing that after sorting (solid squares to open circles), even a constant step in delay causes an irregularly spaced phase step, due to wrapping around circle.

Figures 10 and 11 show example results for $B_{\mathrm{ord}}(\nu )$ and $\mathbf{W}(\nu )$ for star *11UMi plus weak ThAr lamp, for delay E4 ($\sim 1\mathrm{cm}$). The magnitude of the whirl is typically just a few percent of the ordinary nonfringing spectrum for the higher delays where the sinusoidal comb is unresolved. For the lowest delay (E1, $\sim 0.1\mathrm{cm}$) where the sinusoidal comb is fully resolved the magnitude can approach 100%.

Fig. 11

(a) An ordinary (nonfringing) spectrum and (b) complex spectrum (fringing) called a “whirl” are both produced from phase-stepping math for each multiphase input array. Red and blue designate real and imaginary parts. Whirl magnitude is typically a few percent of ordinary spectrum. The target is star *11UMi plus weak ThAr lamp, and etalon E4 (delay $\sim 1\mathrm{cm}$) is used. Vertical axis is intensity, in thousands of counts, and for bottom panel is change in intensity. Note a ThAr line at pixel $X=1402$ ($4849{\mathrm{cm}}^{-1}$).

Pathologies in the iterative process can be detected by monitoring the plot of $\delta {\tau }_n$ in Fig. 10(b). Since the typical TEDI delay step is fairly uniform, if the observed $\delta {\tau }_n$ is erratic this indicates a problem in the data. Often this is due to an uncorrected CR.

5.3.

Fitting Along a Column Assuming Y-Knowledge

5.3.1.

Sine fitting versus vector fitting

For each (wavenumber) column of Fig. 8, the intensity versus phase is essentially fitted to a sinusoid. The sine and cosine amplitudes represent the imaginary (blue) and real (red) parts of the complex spectrum $\mathbf{W}$ [Fig. 11(b)], called “whirl,” which represents the moiré pattern. The vertical offset to the sine fit represents the contribution of the nonfringing ordinary spectrum $B_{\mathrm{ord}}$, also called the native spectrum [Fig. 11(a)].

In the case of uniform visibility fringes and knowledge of phase step positions, a standard sine fitting routine can be used. But if the visibility is nonuniform [as in Fig. 12(a)], such as under atmospheric seeing induced intensity variations, then a standard sine fitting routine can give an erroneous result. (Magnifying a weak signal also magnifies the photon noise, so fixing severe drops in fringe visibility only through normalization is not recommended.)

Fig. 12

(a) Sine fitting is equivalent to a (b) vector construction where step phase and visibility are designated by vector length and angle, which can both be irregular. To isolate the nonfringing from fringing components the vector sum must be zero. The values of weights to achieve this are calculated from phase of ThAr features. Whatever arithmetic is performed on vectors is also performed on data.

Rather than sine fitting, we use a vector combination method instead [Fig. 12(b)], which has greater flexibility in handling irregular fringe visibility and irregular phases, and has the ability for massively parallel computations over many columns simultaneously—useful for rapidly processing the wide bandwidths of our data. To separate the fringing and nonfringing component, the vector sum must be zero, and weightings are applied to achieve this.

The method of choosing the optimal weightings is discussed in Sec. 5.3.3 on “healing.” Notably, a “clumping” technique is used so that if some of the vectors have near zero magnitude, they do not spoil the value of the clump. So dropouts in signal due to a passing cloud, e.g., are tolerated well in our scheme.

5.3.2.

Designating vectors

It is useful to introduce a set of designating vectors (Fig. 13), whose phase and magnitude represent the interferometer comb phase and visibility for a given pixel or group of pixels. Whatever operation is performed on the designating vectors is also performed on the underlying data it represent.

Fig. 13

Regularizing an irregular phase set of fringing data (a) into a regular set (b) involves, for each input vector, adding small amounts of other data from the same set. In the presence of noise, it is beneficial to add averages (d), “clumps,” of those data (c) facing similar directions.

We visualize Eq. (4) for the ordinary spectrum $B_{\mathrm{ord}}$ as a sum of vectors [Fig. 13(d)] such that the fringing component is cancelled because the vector sum is zero. We are allowed to weight the vectors. We find a combination of weights (there are many) that produces vector cancellation, and then apply the same weights and arithmetic to the data to produce the ordinary nonfringing spectrum.

Analogously, we visualize Eq. (3) for the fringing component $\mathbf{W}$ as a series of rotations $e^{-i2\pi {\varphi }_n}$, tailored to each exposure, to align all the vectors to point in the same direction. These software rotations are ideally exactly counter-rotating to the delay-induced rotations. Now, it is the nonfringing terms which appear to rotate and vector cancel, leaving the purely fringing component.

5.4.

Fitting Along a Row Assuming X-Knowledge

5.4.4.

Normalization problems lead to pollution of whirl with ordinary

We found that if the normalization was not working properly, this harms the phase-stepping math, and pollutes the whirl with a small fraction of ordinary component, causing the whirl to have a pedestal appearance [Fig. 14(c)], seen between the real and imaginary parts and having the same shape as the ordinary spectrum. This is not critical for the ThAr because it is large. But for the Star whirl which is only a few percent, a small fractional pollution of the $\sim 50\times$ larger ordinary component produces a significant absolute amount of pedestal.

Fig. 14

(a) Earlier spectral reconstruction results had excessive ringing (yellow curve), compared to telluric theory¹⁹ (dashed). Bad normalization polluted the fringing signal (whirl) with nonfringing component, which created a pedestal-like offset in the whirl (c), which became ringing during heterodyning reversal frequency shift. (b) Ringing problem was solved (red curve) using a different normalization method.

This pedestal is a DC-like (zero frequency) component that gets upshifted during heterodyning reversal to contribute a ringing artifact [Fig. 14(a)]. This was most easily detected in places of the telluric spectrum which are heavily absorbing and thus should be relatively flat (zero intensity). When we eliminated the pollution of the whirl by changing normalization schemes, we got a beautiful ring-free result [Fig. 14(b)]. (We could also solve this by manually subtracting an adjustable amount of ordinary spectrum from the whirl, but that did not address the root problem.)

The pollution and ringing problem disappeared when we changed the method of computing normalization for the StarThar to the simple sum method (after masking away ThAr lines).

6.1.

Rebinning

The horizontal pixels are rebinned more finely $\sim 10\times$, to a uniform grid of $0.05{\mathrm{cm}}^{-1}$. The criteria are that half the reciprocal of the final grid spacing (called the Nyquist frequency) is 10 cm, and this needs to be larger than the largest delay we intend to use (4.6 cm of E8).

We found that a standard cubic spline does a better job than our custom-written FFT method, regarding avoiding any overshoot artifacts produced by the narrowest features.

6.2.

Pixel to Wavenumber Map Calibration

During rebinning, we map the data from pixels to wave numbers using an improved calibration curve (map) which we fine-tune here. We use individual calibrations for fibers A and B that accounts for the slit slant between those fibers. We add to the original polynomial fit a line segmented “tweak map” which forces the map to exactly pass through each important ThAr line or other feature. We primarily use the magnitude of the whirl for the comparison because it is much cleaner than the ordinary spectrum (since the former is insensitive to FP noise). In sparse ThAr regions we use telluric features.

6.3.

Bulk Shift of Map versus (Time or Etalon)

We find that the map has an etalon-dependent offset, which we call shear. This is due to a time-dependent drift, probably thermomechanical from the changed gravity vector after repointing the telescope. The magnitude of the shift is $\sim 0.3{\mathrm{cm}}^{-1}$ over 1 or 2 h.

After removing the monolithic shift in each ordinary spectrum, we notice that there remains a significant misalignment which changes its amount rapidly along the band, and even changes polarity. Examples are shown in the instrument noise section.

6.4.

Lowpass Filtering to Prevent Nyquist Artifacts

Figures 15(a) and 15(b) show aliasing artifacts in data from etalon E5 (1.27 cm), that occurs when the detector pixel density is too sparse. Energy from the continuum fringe comb spike leaks beyond the Nyquist frequency ($0.5{\text{pixel}}^{-1}$, or $\rho \sim 1\mathrm{cm}$), wrapping to the other frequency branch on the negative frequency side. This causes an erroneous signal.

Fig. 15

Aliasing problem because detector pixels too sparse. (a) For delay E5 the FFT has spatial frequencies above $0.5{\mathrm{pix}}^{-1}$, so it wraps to the negative frequency axis. (b) FFT of data after horizontal linearization—aliased spike easier to see at $-0.65\mathrm{cm}$ (which belongs at 1.27 cm). (c) Lowpass filtering removes most aliased spikes.

Our remedy [Fig. 15(c)] is to restrict frequencies to $<0.5\mathrm{cm}$, then apply a standard cubic spline interpolation routine to rebin. The latter has good behavior similar to a Hanning-shaped lowpass filter.

6.5.

Detwisting the Whirl, Removing Glass Dispersion

The interferometer delay is a function of wavenumber due to the refractive index dispersion of the glass etalon. This creates a phase twist versus wavenumber, ${\theta }_{\text{glass}}(\nu )$, associated with each etalon. Data whirls are untwisted via $\mathbf{W}(\nu )e^{-i2\pi {\theta }_{\text{glass}}(\nu )}$. Otherwise, the FFT will not produce a narrow spike [Fig. 16(b)] for the fringe comb.

Fig. 16

The effect of linearizing the horizontal axis from pixels to wavenumbers, and compensating for glass dispersion, is to narrow the comb spike. (a) The Fourier transform of whirl for delay E2 (0.33 cm) in $1/\text{pixel}$ frequency space, showing a broad continuum comb spike. (b) Linearized data, the spike is extremely narrow.

For each etalon, we find ${\theta }_{\text{glass}}(\nu )$ [Fig. 17(a)] by comparing measured to calculated phases of ThAr lines. Fitting a polynomial to ${\theta }_{\text{glass}}(\nu )$ produces a baseline fit good to $\sim 0.1$ cycle in chunks of 20 to 200 pixels. The precision is further improved later to better than 0.01 cycle using the ThAr or telluric features of each target’s data (Fig. 18).

Fig. 17

The dispersion characteristics of set of glass etalons were calibrated over a 4000 to $\mathrm{11,000}{\mathrm{cm}}^{-1}$ bandwidth by phase shift comparison with ThAr emission line model.²⁰ (a) The resulting phase versus wavenumber curve is then used to “untwist” stellar fringing spectra so that a constant base delay ${\tau }_0$ can be used during reversal of heterodyning effect. (b) The derivative of the phase data provides delay versus wavenumber curves. The dispersion characteristic of the beamsplitter substrate (most obvious in the smallest delay E1) is a common component of all. Delay values for E1-8 are $\sim 0.1$, 0.3, 0.7, 1, 1.3, 1.7, 2.9, 4.6 cm, and 2.4 cm for E6.5 (not shown).

Fig. 18

Aligning a data wavelet (light green) to theory wavelet (black) in wavenumber and phase is easy to do precisely by using Lissajous presentation (bottom row), which plots data versus theory, rather than data and theory versus wavenumber (top row). One restricts bandwidth to a small ($\sim 10{\mathrm{cm}}^{-1}$) region otherwise the loops obscure each other. (b) Effect of a $0.3{\mathrm{cm}}^{-1}$ change in wavenumber. (c) Effect of a 0.02 cycle phase change. (a) Perfect alignment. Data of ThAr lamp from September 19, 2010 etalon E4, A-order.

For the spectral reconstruction, we only need this phase curve. However, it is useful to produce a delay versus $\nu$ curve [Fig. 17(b)] by taking the $\nu$-derivative. This is useful for detecting skipped integer fringes of phase, lost spanning band regions that lack features.

6.6.

Fine Delay Determination

Due to our wide bandwidth, by observing the twist in phase along the band, we could determine the delay unique to the nearest fringe (nearest $2\mu \mathrm{m}$), i.e., with no integer ambiguity. Then within that fringe we further determined the delay to a fraction of a fringe and removed the last fractional phase shift.

Due to mirror repositioning errors, the delays would change a few microns from night to night; we could also detect a delay difference of a fraction of a micron between A and B fibers. Since these were so small, we kept the original delay value and instead slightly adjusted the whirl phase.

The values of the eight delays in our quiver for 2010 September were for E1-E8 0.083, 0.34, 0.66, 0.96, 1.27, 1.75, 2.92, 4.63 cm, and 2.4 cm for E6.5, added in June 2011. (The latter was used to fill the gap between 1.7 and 3 cm, and swapped for the E1 delay on the eight-position filter wheel.)

6.7.

Calculating Phase Difference Between Theory and Data

We have two methods for calculating the phase difference between theory and data. (1) The dot product method is automated and thus used initially. However, it has the drawback that it blindly uses the entire theory spectra without regard to any incorrect components to it. (2) Then we fine-tune the result using a manual Lissajous method which has the ability for the user to ignore components suspected to be incorrect.

7.1.

Summary of Spectral Reconstruction

In stage 4, “spectral reconstruction,” the whirls are numerically shifted in frequency space (heterodyning reversal, Fig. 19) to undo the analog frequency shift imparted by the interferometer. This produces a wavelet signal ${\mathrm{\Omega }}_n(\nu )$ for each etalon, covering a different spatial frequency range centered on its delay value. These wavelets are summed to produce a provisional output spectrum ${\mathrm{SR}}_{\mathrm{prov}}(\nu )$. This is equalized to remove the bumps in the net frequency response, so that the latter is Gaussian with a user-chosen final resolution, to form a final spectral output ${\mathrm{SR}}_{\text{final}}(\nu )$ having minimal ringing. Comparison of ThAr data to theory for an isolated peak provides the equalization curve used for both ThAr and Star. The choice of final resolution is limited by the maximum delay which produces a contiguous range of delays, as $\sim {\tau }_{\max }/\lambda$.

Fig. 19

Stages of heterodyning reversal processing. (a) FFT of whirl, here for a small delay E2 (0.34 cm). Continuum comb spike at 0.34 cm is within range of native spectrograph ($-0.5<\rho <0.5{\mathrm{cm}}^{-1}$). Horizontal axis is frequency in features per ${\mathrm{cm}}^{-1}$, having units of cm. Vertical axis limited to $\pm 2400$ to more easily see details, and the spike is about $\pm \mathrm{10,000}$ high. (b) The FFT is shifted left by the same amount as the delay, to move the spike to zero frequency and other signal energy to left branch. (For larger delay etalons there may be no comb spike.) The real part of the inverse FFT of (b) produces the spectrum (d), which may look like wavelets for higher delays, but for lower delays will contain the continuum shape. Equivalently, one could copy the complex conjugate to the opposite branches to form (c) then apply the inverse FFT to form (d).

7.2.

Bell-Like Weighting to Reduce Net Photon Noise

Prior to heterodyning reversal, frequencies outside a set value (0.52 cm) were deleted, since they contain little signal. In addition, a bell-shaped weighting $h_0(\rho )$ was sometimes applied $\mathbf{w}\to h_0(\rho )\mathbf{w}(\rho )$.

For photon-limited noise, the ideal shape is that of the native sensitivity peak, but for instrument limited noise (the TEDI situation) there is little benefit. If we used a weighting, it was a Hanning function since it goes to zero absolutely at its edges (set to 0.52 cm).

7.3.

Reversal of Heterodyning

Figure 19(b) shows the effect of reversal of the heterodyning, where the signal has been shifted in the direction that moves the continuum comb to zero. The effect is illustrated in frequency space, but is more accurately computed in the software in $\nu$ space as

Our convention has the continuum fringe comb starting at a positive frequency $\rho$. Hence during numerical heterodyning reversal, Figs. 19(a) and 19(b), the signals shift to left toward negative $\rho$. This convention is arbitrary since after the real part is taken (d) both frequency branches (c) have energy.

7.4.

Wavelets Produced from Real Part of Shifted Whirl

The output spectrum we are reconstructing must be purely real. We take the real part of the heterodyne reversed signal (in $\nu$ space) to form what wavelet signal

In practice, the above equation might actually require using the imaginary part of $\mathbf{W}$, depending on how the phase is defined. The test is that the continuum comb portion of $\mathbf{W}$ (produced by a small delay) should be, after heterodyning reversal, mostly real and positive.

7.5.

Summation of Wavelets to Form Provisional Spectral Reconstruction

The wavelets are summed to form a provisional SR spectral output called ${\mathrm{SR}}_{\text{raw}}$

7.5.2.

Ordinaries not included in TEDI output spectra

Initially, we included the ordinary native spectrum in the output spectrum, to provide the lowest frequencies, but we discovered for TEDI data we get cleaner results NOT including the ordinary spectrum. Instead, we rely on fringing from the low delay etalons E1 and E2 for the lowest frequencies information. (However, we still include it in plots of the SR output, for comparison.)

Figure 20 demonstrates in stellar data that the fringing-derived spectrum (red curve) agrees extremely well with the telluric model (black curve). In some cases, it is more accurate than the native spectrum (green curve), such as when FP noise is present. The EDI is extremely robust to FP noise because only changes between exposures contribute to the whirl.

Fig. 20

The EDI result (red, $1\times$) can provide an accurate measurement of the low-resolution information usually provided by the native spectrograph (green). Here we deliberately blur the fringing EDI result to the same resolution as native spectrograph. It accurately follows the details of telluric model (black). Green and black curves offset vertically. Target is star HD219134 + ThAr lamp, September 19, 2010.

7.6.

Modulation Transfer Function Displays Net Frequency Response

The net frequency response of the EDI system is probed by a modulation transfer function (MTF), which is the ratio out versus in of complex amplitude, for each frequency $\rho$

7.6.2.

Measured MTF of externally dispersed interferometry data

Red curve in Fig. 21 shows an MTF measured for a star. Note that the positions of the peaks (associated with each delay) are not ideal. We desire shoulder-to-shoulder coverage to produce a ring-free net lineshape. Unfortunately, the set of TEDI delays had a gap between E7 and E8 (3 and 4.6 cm), since TEDI delays were selected for Doppler radial velocimetry purposes, not spectroscopy, and there were only eight positions in the rotary etalon holder.

Fig. 21

(a) Net MTF versus frequency $\rho$ for the TEDI system near $7540{\mathrm{cm}}^{-1}$ ($1.34\mu \mathrm{m}$) shows peaks for each delay (red), having same shape as native peak (green) but shifted by ${\tau }_n$. Some of peaks have blended together (the ideal situation). (b) Equalization curve (purple) is Gaussian shape of $\mathrm{res}=\mathrm{27,000}$ (thick black curve) divided by red curve MTF. If additional etalons filled the gap between 3 and 4.6 cm, resolution of 45,000 at $7500{\mathrm{cm}}^{-1}$ is possible. Data of for star ${\kappa }_{\mathrm{CrB}}$ on June 19, 2011. The E6.5 etalon (2.4 cm) swapped with E1 (0.1 cm) in the rotary holder. Most data were taken without E6.5, thereby producing a gap at 2.4 cm which limited ring-free resolution.

The presence of the 3 to 4.6 cm gap limits the effective resolution, because attempting to use a wide Gaussian MTF that extends over this region (such as the $\mathrm{Res}=\mathrm{36,000}$ curve in Fig. 21) would produce ringing in the final lineshape in $\nu$-space. If additional etalons filled the gap between 3 and 4.6 cm, a resolution of 45,000 at $7500{\mathrm{cm}}^{-1}$ is possible.

For data taken prior to June 11, 2011, there was a gap at 2.4 cm between E6 and E7. This was remedied by a purchase of a new etalon, called E6.5, of 2.4 cm delay which occupied the position in the eight-position “filter wheel” previously held by the smallest etalon E1. (E2 can cover the range of E1.) Only a few exposures have been taken with the E6.5, during what turned out to be the last run of TEDI, June 11 to 19, 2011.

7.7.

Equalization Applies Taper to Modulation Transfer Function

The purple curve in Fig. 21(b) is the equalization curve used to reshape the MTF into a smoother Gaussian. It is the ratio between the desired MTF (black curve) and the actual MTF (red curve)

The EQ is applied to the FFT of the data, ${\mathrm{sr}}_{\text{raw}}(\rho )$ via

The choice of goal resolution for the Gaussian in Eq. (33) is a tradeoff between resolution and ringing. In the presence of gaps in the delay coverage, increasing resolution increases the ringing.

8.1.

Comparing Two Fibers on Same Star

Figures 26 and 27 show EDI measurements comparing the two fibers A (blue curve) and B (red curve) measuring the same source, their difference curve (gold), and compared to the telluric model (gray). The two fibers used different detector pixels in the $Y$ direction and were exposed at times $\sim 10\min$ different because the role of the fibers was being alternated between Star and Star + lamp. Hence the fibers experienced different PSF distortions to the native spectrograph (see Sec. 9.4). Yet the resulting absorption lines yield reasonably accurate lines because it sums over several wavelets.

Fig. 26

Comparison of SR result from two fibers A (thin blue) and B (thick red), thus on different pixels, looking at same stellar source HD219134 in A-order at 4x boost, resolution of 11,000. The difference between fibers is the tan curve. The ordinary spectrum in native pixels is green curve. That and the telluric + lamp model (gray) are offset vertically for clarity. SR results detects the telluric and lamp features in the model, as well as many stellar features not in the model. The two fibers were measured a few minutes apart, on different exposures, since the roles (ThAr-only, ThAr + star) were alternated between the fibers. This figure shows 4715 to $4775{\mathrm{cm}}^{-1}$. (b) and (c) show the wavelet arrays for the 6 significant etalons used. Fiber B appears to have better phase alignment than fiber A.

Fig. 27

Additional regions of the A-order when comparing two fibers A (thin blue) and B (thick red) and their difference (tan); same data as Fig. 26. (a) 4240 to $4354{\mathrm{cm}}^{-1}$, (b) 4916 to $4989{\mathrm{cm}}^{-1}$, (c) 5200 to $5272{\mathrm{cm}}^{-1}$. Telluric model (gray) offset vertically.

8.2.

Measurement Comparing Star to Literature Spectrum

Figure 28 compares our EDI measurement of HD219134 with a spectrum obtained by another group with a conventional spectrograph, the NIRSPEC at Keck Observatory.²¹ All the stellar NIRSPEC features are seen in the EDI output, but many are unresolved in the TripleSpec native spectrum (green curve).

Fig. 28

Comparison of SR result (red) at 11,000 resolution on star HD219134 on September 19, 2010, to data measured independently by the NIRSPEC conventional spectrograph²¹ at Keck Observatory on June 21, 2003. Black curve is the (star*telluric + lamp) model, where the star curve is NIRSPEC data blurred from its original 18,000 resolution to 11,000 resolution (purple), and is multiplied by the telluric and ThAr lamp model (gray). Native TripleSpec spectrum (green) cannot resolve stellar features (vertical lines) now resolved by the EDI. Wavelength shifts due to Earth and solar system motion not corrected.

8.3.

Measurement Comparing Two Stars: Unexpectedly Deep Lines and Doppler Shift Value

Figures 29 and 30 compare two stars measured on same night, star HD219134 (blue curve) and HD962 (red curve). We observed the latter star because the SIMBAD database described it as a featureless A-type, and thus could indicate the telluric component of the spectrum. However, our observations disagree with SIMBAD in two ways for HD962, while agreeing for HD219134.

Fig. 29

Science result comparing spectra of two stars, in A-order, measured nearly simultaneously star HD219134 (blue curve) and HD962 (red curve). The panel shows only a very small $53{\mathrm{cm}}^{-1}$ region (4715 to $4768{\mathrm{cm}}^{-1}$) of the $7200{\mathrm{cm}}^{-1}$ wide spectrum measured simultaneously. Telluric + lamp model¹⁹ is black curve. This region is interesting because it contains a ThAr spectral emission line ($4764{\mathrm{cm}}^{-1}$) which appears in all three curves, telluric absorption features (e.g., $4758{\mathrm{cm}}^{-1}$) which appear in all curves, and three stellar absorption features (4724, 4740, and $4747{\mathrm{cm}}^{-1}$) which appear only in the stellar data not the model and which appear Doppler shifted one star relative to the other. Data were taken on both fibers and averaged. Effective resolution here was set to only $4\times$ boost or 11,000 to reduce ringing to an extreme minimum, which aids identifying weak lines. Higher resolution can be used if some ringing is tolerated.

Fig. 30

D-order data version of the two stars of Fig. 29, comparing two stars measured nearly simultaneously, HD962 (red curve), and HD219134 (blue curve) on September 19, 2015, around UT of 10 h. They have resolution of 19,000, a $6\times$ boost over native spectra (thin and thick green curves are of HD962 and HD219134). The SIMBAD database classifies HD962 as a featureless A-type star. However, we find stellar absorption lines of similar depth as to HD219134, which is K-type. Second, we measure for HD962 a Doppler velocity that differs significantly from its SIMBAD value, ($29\mathrm{km}/\mathrm{s}$ too red—the HD962 lines should be bluer than HD219134, not redder), whereas for HD219134 it agrees perfectly. The ThAr and telluric features of HD962 agree nicely with theory (see 9188 and $9206{\mathrm{cm}}^{-1}$), and the native spectra are consistent with the EDI result, so we believe the measurement is accurate. The difference in Doppler velocity observed here in the D-order is the same as observed in the A-order in Fig. 29.

First, the depth of narrow stellar absorption features of HD962 is similar to that of HD219134, which is a K-type. So the designation of A-type for HD962 is suspect. This is seen clearly in Fig. 31 which compares the native spectra of HD962 to several A- and B-type stars (HR1178, Gam Oph, HD7891) and the K-type star HD219134.

Fig. 31

Comparing the native spectra (star + ThAr lamp) of the disputed A-type star HD962 (bold black) with three other “featureless” (B and A type) stars (thin blue, red, green for HR1178, gam Oph, HD7891), and with K-type star HD219134 (thin black) shows that the disputed star is more similar to the K-type star regarding the narrowest stellar absorption lines (such as at $9450{\mathrm{cm}}^{-1}$). However, there is a wider feature at $9140{\mathrm{cm}}^{-1}$ in the undisputed A stars that is weakly evident in the disputed star and not in the K-star. Curves are vertically offset for clarity, and the original continuum intensities for a single phase stepped exposure are 13.5k, 75k, 26.5k, 7.4k, 28.5k counts, for HR1178, Gam Oph, HD7891, HD962, HD219134.

All were measured on the same night September 19, 2010. Strikingly, HD962 has more similarity to the K-type star than the A- or B-type stars, in that it shares the narrow stellar absorption lines at $9450{\mathrm{cm}}^{-1}$, which are absent for the three telluric calibrator stars. Interestingly, the three telluric calibrator stars all show (Fig. 31) a broader absorption feature at $9140{\mathrm{cm}}^{-1}$ not seen in the K-type, and the dispute star HD962 also shows this, although weakly.

Second, we observe a relative Doppler shift between the stellar lines of the two stars, but it is of the wrong polarity. According to SIMBAD, we should be seeing HD962 blue shifted relative to HD219134. Instead we see HD962 redder. The Doppler velocity of HD219134 that we measure agrees perfectly with the SIMBAD database when the barycentric motion is accounted for. The two stars were measured within an hour of each other, using the same data analysis procedure, and the telluric and ThAr features match up perfectly (see dips at 9188 and $9206{\mathrm{cm}}^{-1}$). We checked the RA and Dec values of the FITS file header to confirm the telescope was pointed at the right star. Thus we believe the SIMBAD database for HD962 is in error.

Some details of the Doppler calculation–Solar spectrum from FTS measurements at NIST are shown shifted by velocities 5.7 and $26.4\mathrm{km}/\mathrm{s}$ to match the position of absorption lines. We observe geocentric blueshifts of 5.7 and $26.4\mathrm{km}/\mathrm{s}$ for HD962 and HD219134, respectively. Given Earth’s motion of 11.679 and $7.566\mathrm{km}/\mathrm{s}$ toward HD962 and HD219134 at UT times 10:49 and 9:50, respectively, this implies heliocentric blueshifts of 5.7 to $11.679=-5.98\mathrm{km}/\mathrm{s}$ for HD962 and 26.4 to $7.566=18.83\mathrm{km}/\mathrm{s}$ for HD219134. SIMBAD blueshift values are 23.0 and $18.83\mathrm{km}/\mathrm{s}$ for HD962 and HD219134. Thus they agree nicely for HD219134, but are off by $23+5.98=29\mathrm{km}/\mathrm{s}$ for HD962.

8.4.

High Dynamic Range and No Fixed Pattern Noise

Figure 32 shows the extremely high dynamic range and robustness to FP noise of which EDI is capable. These are of our ThAr lamp to a $6\times$ boosted resolution of 19,000 in the D-order. We resolve weak peaks of heights $\sim 0.04\%$ of the brightest line in that order. (Note the weak but distinct peaks at 9171 and $9240{\mathrm{cm}}^{-1}$, for example.) The data are normalized so that the brightest line of the order, at $9548{\mathrm{cm}}^{-1}$, defines unity height. The $\sim 1{\mathrm{e}}^{-4}$ noise at the base of the 9188 line (used to calibrate the equalization) is a dynamic range of $\sim \mathrm{10,000}$, which is excellent. There is very little ringing ($\sim 0.01\%$).

Fig. 32

Measurement (red curve) of our ThAr lamp to a $6\times$ boosted resolution of 19,000, compared to NIST measurement²⁰ of a different ThAr lamp blurred to same resolution (black dashes), D-order. Purple text is their assignment of species. Native spectrograph (green curve) has resolution $\sim 3300$. Note disagreement in heights of Th lines, whereas Ar lines agree well, due to individuality of each ThAr lamp. Note extremely high dynamic range of measurement—$\sim 0.1\%$ lines are easily observed (heights are fraction of $9548{\mathrm{cm}}^{-1}$ line). The EDI curve is robust to FP noise such as false peak in the native spectrum at 9174 and $9214{\mathrm{cm}}^{-1}$ due to bad pixels at $X=1033$ and 1066 (inset shows detector there). Bad pixels are constant and thus do not affect whirls which look at changes between exposures.

This data also illustrates EDI’s robustness to FP noise. The inset is a small section of the 2-D raw detector image at the D-order associated with this data, showing a couple of bad pixels. These produce erroneous bumps in the native spectrum (green curve) at 9174 and $9214{\mathrm{cm}}^{-1}$. The red EDI output spectrum does not react to these bad pixels, because they are constant and thus do not respond to the sinusoidal interferometer phase dithering.

9.1.

Fixed Pattern Noise Rejection Advantage

9.2.

Severe Distortions of Native Spectrograph

Figure 34(a) shows that high-resolution spectra can still be produced in spite of native PSF drifts that are extremely large and irregular—locally as much as $\sim 0.6{\mathrm{cm}}^{-1}$ as shown in Fig. 35, and averaged over $10{\mathrm{cm}}^{-1}$ in Fig. 37. Figure 35 is a stack of the eight native spectra showing drifts along the dispersion direction. There is a bulk drift versus time (Fig. 36) which we call shear, which is removed prior to heterodyning reversal. Note the significant residual irregular component, which even causes the shear polarity to flip over small wavenumber distances. Figure 38 shows how this horizontal noise is equivalent to an intensity noise of $\sim 3\%$.

Fig. 35

Large PSF irregularity in the native spectrograph to which EDI is robust. Row is native spectra of A-order versus etalon # ($\sim 10\min$ apart). The PSF shift $\sim 0.6{\mathrm{cm}}^{-1}$ versus etalon # is called “shear.” Bulk shear has already been removed here, and an irregularity remains. Note irregular change in shear polarity. Fortunately, for the EDI the native lineshape only affects the envelope of the wavelets, but not the internal shape.

Fig. 36

Bulk translation of the native spectrum along the dispersion axis due to thermomechanical or changing gravity vector, and due to slit slant acting on vertical separation of fibers A and B, for star HD219134 on September 19, 2010.

Fig. 37

Irregular PSF drift of native spectrograph, versus wavenumber, correlated in $10{\mathrm{cm}}^{-1}$ chunks between E6 and E1 (red dots), and E2 and E1 (blue open squares), for the A-order and September 19, 2010, data of HD219134. The time-dependent bulk offset (Fig. 36) has already been removed.

Fig. 38

(c) Measured amount of (vertical) instrumental noise for A-order, as a fraction of continuum, due to (horizontal) PSF drift or width change, by subtracting the native spectrum averaged over all eight etalons (d) from each spectrum (a) measured individually for each etalon to produce about a 3% residual (b).

9.3.

Robustness to Point Spread Function Widening

The EDI is very robust to both changes in the PSF position and width, such as due to atmospheric seeing disturbances or thermomechanical drifts. Let us first discuss the effect of an unexpected change in the width of the native lineshape. For concreteness, let us assume the net focal spot size is a sum in quadrature of the native slit width and some fixed additional blur, fixed in detector space, and that we have a $10\times$ boost so that the slit width of the EDI native spectrograph is $10\times$ wider than the conventional spectrograph it is trying to match.

This is numerically simulated in Fig. 39. The result is that EDI lineshape (red) is dramatically more robust than the conventional spectrograph (black) because EDI data are sensed at spatial frequencies defined by a physical delay, not controlled by focal spot size or slit width. The conventional spectrograph lineshape (black solid to dash) doubles and its MTF is halved in width. Whereas the effect on the EDI MTF is just to have slightly narrower bumps [solid to dashed red curves of (b)]. Since the original bumps are equalized away, the EDI final resolution is mostly unchanged by the blurring.

Fig. 39

Simulated effect of adding focal blur to the native instrument lineshape for conventional spectrograph (black), and an EDI (red) having a $10\times$ boost—the EDI is much more robust. The conventional spectrograph lineshape (black solid to dash) doubles (a) and its MTF (b) is halved in width. For EDI native spectrograph, already $10\times$ wider slit width, relative broadening is very small ($\sqrt{1^2+{0.173}^2}=1.015$), so change in bump shapes (b) also small. Note, the EDI still retains significant high frequency sensitivity. The 1.5% change in bump width creates a 5% error in equalization in removing the bumps, creating two $\sim 1\%$ height satellite peaks at $3.5{\mathrm{cm}}^{-1}$ from central lineshape (off graph).

9.4.

Robustness to Point Spread Function Drift of Native Spectrograph

9.4.1.

Translational reaction coefficient for single delay externally dispersed interferometry

Now we will discuss, at length, the EDI robustness to an unwanted translation or “insult” $\mathrm{\Delta }x$ (in wavenumber units) of the native spectrograph PSF along the dispersion direction on the detector, producing a perceived shift $\mathrm{\Delta }\nu$ in the final output spectra. We define a translational reaction coefficient (TRC) $\equiv \mathrm{\Delta }\nu /\mathrm{\Delta }x$. For a conventional dispersive spectrograph, $\mathrm{TRC}=1$.

For a single delay EDI, the TRC is smaller than unity by a factor of several, and can be zero, depending on the shape of the spectrum where TRC is being evaluated. For example, at an isolated narrow feature such as a ThAr line, it is zero. This is because the moiré pattern has a phase which has no slope along the dispersion direction, so translating the moiré pattern does not change the perceived phase, and hence perceived $\nu$.

However, for a pair or more of lines that are partially blended together by the native blur, there can be a slope to the moiré pattern. This slope, acting against a translation $\mathrm{\Delta }x$, can produce a phase change which leads to a nonzero TRC (but still much smaller than unity by a factor of several). How TRC behaves on average for an arbitrary spectrum is best seen in Fourier space, by multiplying the peak shape at a delay $\tau$ of $\mathrm{psf}(\rho -\tau )$ by a line $(\rho -\tau )$ running through the center of the peak [Fig. 40(b)]. A similar line $(\rho -0)$ runs through the native spectrograph $\mathrm{psf}(\rho )$ in (a) to calculate its reaction to $\mathrm{\Delta }x$. We have colored green similar regions of $\rho$ in the panels to suggest how the net reaction for the native spectrograph (a) is several times that of the EDI peak in (b), because the location of multiplier line is at much higher $\rho$ for the EDI, because the heterodyning effect moves it, along with the peak, to be located at $\rho =\tau$.

Fig. 40

(a) The reaction of a conventional spectrograph to an unwanted displacement $\mathrm{\Delta }x$ along the dispersion direction at the detector, in Fourier space versus frequency $\rho$, is to multiply the $\mathrm{psf}(\rho )$ or $\mathrm{MTF}(\rho )$ by $\rho$ (red line). Science frequencies typically are high (green band), producing a large reaction. (b) In contrast, for the EDI both the sensitivity peak and red line are shifted to high frequency by the interferometer delay, and the reaction flips polarity as $\rho$ crosses the center of the peak. This can reduce the integrated reaction by an order of magnitude. (c) With a multiple delay EDI the net reaction is further reduced when the positive and negative portions of the peaks overlap.

Figure 40 is notional and does not yet multiply the peak shape against the line. That is done in Fig. 41(a). This latter graph actually depicts the net reaction (bold black curve) for multiple delays, but if you look at the reaction curve under one of the two highest delays, E7 or E8, they are isolated, and thus depict what happens for a single delay. This net reaction curve (bold black) in Fig. 41 is a phase change in $\rho$ space, which manifests translations in $\nu$ space. The black and gray dashed slopes depict how that phase would change for $\mathrm{TRC}=1$ and $\mathrm{TRC}=0.05$, respectively. Note that the reaction curve for E7 at 3 cm is roughly below the $\mathrm{TRC}=0.05$ slope, and would be under $\mathrm{TRC}=0.02$ for the highest delay E8 at 4.5 cm. The TRC scales as $1/\tau$ for a single delay, and would be of order 0.15 for the smaller delay E4 at 1 cm.

Fig. 41

Net reaction (bold black curves) to a PSF $\mathrm{\Delta }x$ shift for set of TEDI peaks, for (a) original lineshape and (b) peaks reshaped into a sinc function crossfading lineshape (triangle shape also works). Reaction is peak $n$ shape (colored) multiplied by line $(\rho -{\tau }_n)$ passing through peak center, then summed over $n$. (a) For original peak shapes, net reaction is consistent with same slope (gray dashes) of $\mathrm{TRC}=0.05$ ($\mathrm{\Delta }\nu /\mathrm{\Delta }x$) seen in Fig. 42. (b) Reshaped peaks cancel each other’s reaction so $\mathrm{TRC}=0$. Remarkably, PSF can now drift in wavenumber without shifting the output spectrum—not possible in conventional spectroscopy where $\mathrm{TRC}=1$ (black dashes). Missing delays at 2.3 and 3.8 cm prevent crossfading above 2 cm (but TRC still less than 0.05). The $\mathrm{TRC}=0$ region 0 to 1.8 cm allows $6\times$ boost (16,000 resolution) in the A-order at $4700{\mathrm{cm}}^{-1}$. Crossfaded lineshapes are asymmetrical to achieve foot to center positioning for irregular delays ${\tau }_n$.

9.4.2.

Translational reaction coefficient for multiple delays

We find that with the use of multiple delays, and without any special lineshape reshaping, we can reduce the $\mathrm{\Delta }\nu$ per $\mathrm{\Delta }x$ by roughly a factor of 20. The idea is that having multiple overlapping peaks will cause the negative reaction on the left side of one peak to cancel with the positive reaction on the right side of a neighboring peak. This is notionally shown in Fig. 40(c), and calculated for actual TEDI peak positions and original shapes in Fig. 41(a). The bold black curve is consistent with $\mathrm{TRC}=0.05$ in the region of the overlapped peaks E1 to E6, 0 to 2 cm. This agrees with the simulation of Fig. 42, even though this was done in a different manner—it calculates movement of a specific feature rather than estimate for an arbitrary spectrum, and does it in wavenumber space rather than Fourier space. It follows a ThAr line in TEDI data and computes $\mathrm{\Delta }\nu$ when all etalon rows of the moiré data are artificially shifted by $\mathrm{\Delta }x=0.5{\mathrm{cm}}^{-1}$ to the left. The result is a $\mathrm{\Delta }\nu =0.025{\mathrm{cm}}^{-1}$ shift, for a $\mathrm{TRC}=1/20$. This means the Doppler noise would be $20\times$ smaller than if using the native spectrograph alone (if $\mathrm{\Delta }x$ is the only noise insult type).

Fig. 42

Simulated EDI reaction to a native spectrograph PSF translation using ThAr data and TEDI data processing. Moiré input data become wavelets—a stack of seven is shown in (b), which sum to produce a narrow peak in (a) (thin red curve) at $\sim \mathrm{16,000}$ resolution ($6\times$ boost). (Unequalized gap between E6 and E7 delays produces slight ringing.) The effect of a PSF shift is simulated by translating every moiré input by $\mathrm{\Delta }x=0.5{\mathrm{cm}}^{-1}$ to the left (c). Note shift in the envelope of the wavelet stack yet internal peaks are stationary. The summed peak (thin black) moves only $\mathrm{\Delta }\nu =0.025{\mathrm{cm}}^{-1}$ relative to the unshifted version (red peak), hence the TRC, $\mathrm{\Delta }\nu /\mathrm{\Delta }x=\mathrm{TRC}=0.05$, a $20\times$ reduction.

10.1.

General Lineshape for Two Crossfading Peaks

To engineer perfect cancellation of the net reaction of overlapped peaks, the shape and width of the peaks are modified into a so-called “crossfading” shape. The original psf is reshaped to a new psf by multiplying the data in Fourier space with a multiplier $[{\mathrm{psf}}_{\text{new}}(\rho )/{\mathrm{psf}}_{\mathrm{orig}}(\rho )]$, as in Fig. 43(a) from a typical peak with gradual wings into a triangle crossfading shape having a well-defined width. For the dual peak crossfade discussed here (mixing more peaks is possible), the most important consideration is that the edge of the ${\text{peak}}_1$ (called a “foot”) goes to zero magnitude at the center position ${\tau }_2$ of the neighboring ${\text{peak}}_2$, and stays zero beyond that.

Fig. 43

Achieving perfect PSF drift robustness between two delays ${\tau }_1$ and ${\tau }_2$ by crossfading using a triangular lineshape. This achieves zero net translational reaction ($\mathrm{TRC}=0$) for frequency $\rho$ between the peaks (frequency and delay have same units). (a) The $\mathrm{psf}(\rho )$ or MTF for each peak is reshaped to a crossfading shape, here a triangle, by multiplication of data in $\rho$-space by $[{\mathrm{MTF}}_{\text{new}}(\rho )/{\mathrm{MTF}}_{\mathrm{orig}}(\rho )]$. (b) For dual-peak crossfading, foot of ${\text{peak}}_1$ (black dot) must be 0% at ${\tau }_2$, which is maximum of neighboring ${\text{peak}}_2$ (green dashed triangle), and vice versa for ${\text{peak}}_2$ foot (green dotted) at ${\tau }_1$. Reaction to a wavenumber shift $\mathrm{\Delta }x$ is estimated by multiplying MTF by line through center of each peak. (c) For triangular peak this creates a sine-like odd-symmetry shape consisting of two opposite parabolas. (For sinc function, it makes an exact sinusoid.) When overlaid with the identical but shifted sine-like shape of the neighboring peak (dashed), the net reaction (bold blue segment) is zero between ${\tau }_1$ and ${\tau }_2$. (The $\mathrm{\Delta }x$ presumed same for both peaks.) Scheme is extended to all pairs of peaks ${\tau }_n$ to ${\tau }_{n+1}$ to cover all delay space with $\mathrm{TRC}=0$. For unevenly spaced peaks, reshaped peaks are asymmetrical.

From the isolated peaks E7 and E8 of Fig. 41(a), one can see that each peak produces an odd symmetry sine-like contribution to the reaction curve, with the right side being positive and the left is negative. The new idea is to reshape the peaks to overlap so that the positive contribution of one peak perfectly cancels the negative contribution of the neighboring peak. This is shown in Figs. 43(b) and 43(c) for the simple case of a triangle lineshape. Between the two peaks, i.e., between ${\tau }_1$ and ${\tau }_2$, a crossfading is said to occur.

There are a family of curves that satisfy this condition of producing $\mathrm{TRC}=0$ in this crossfade region between the peaks. Let us derive the equation that describes them. The net reaction to a horizontal PSF drift is the product of the line $(\rho -{\tau }_n)$ through the center of the peak, times each peak. For $n$ peaks this is

10.2.

Case of Triangle Lineshape

The simplest case is a triangle lineshape, separated by distance 1 and height 1, which thus has linear crossfading of unity slope between the peaks. Then the right side of the left peak is $f(r)=1-r$, so the first term of Eq. (36) is $r(1-r)$, and the left side of the right peak is $r$ so the second term is $(r-1)r$. Hence Eq. (36) is satisfied since $r(1-r)+(r-1)r=0$, and in the region between the two peaks $\mathrm{TRC}=0$.

We apply this crossfade reshaping to every pair of overlapping peaks to fill a contiguously overlapped $\rho$ region entirely with $\mathrm{TRC}=0$. For the case of irregular ${\tau }_n$, the peaks will be asymmetrical, treated as half-peaks, in order to always place the foot on the center of the neighboring peak and vice versa. That is, the left half of peak at ${\tau }_{n+1}$ is the reflection of right half of peak at ${\tau }_1$, and these can be a different shape than neighboring pairs of half-peaks. So one should view a manifold of peaks as instead a manifold of crossfading pairs of half-peaks.

10.3.

Other Crossfading Lineshapes

Colleague Eric Linder pointed out²² that a sinc function, $f(r)=\sin (\pi r)/r$, for just the region between its first zeros, also satisfies Eq. (36). This produces $F(r)=\sin (\pi r)$, which has the required property that shifting by 1 flips its polarity. This has the advantage of a rounded peak that avoids some possible ringing produced by the sharp tip of a triangle. Interestingly, the sinc function (but with its many satellite peaks) is the Fourier transform of a rectangle and hence is naturally the psf of a slit opened much wider than its spectrograph focal blur. We are studying the utility of other crossfading shapes that satisfy Eq. (36).

10.4.

Crossfaded Multiple Delays Produce TRC$=$0

Figure 41(b) is the same calculation of net reaction as (a), but using sinc function peaks instead of original peak shapes, and having asymmetrical shapes to satisfying the foot to center condition for the irregular delay spacing. Remarkably, this produces $\mathrm{TRC}=0$ over the contiguously overlapped region from E1 to E6 (up to 1.8 cm). This means that for resolutions that only need 0 to 2 cm of frequency ($6\times$ boost, or $\sim \mathrm{16,000}$ in the A-order), the interferometer is completely decoupled from the native spectrograph and not disturbed by changes in the native PSF location. (In principle, three new etalons at 2.3, 3.4, and 4.1 cm to fill gaps between E6, E7, and E8 could be added to extend the contiguous delay region to 4.6 cm, greatly increasing resolution for $\mathrm{TRC}=0$ operation to $\sim \mathrm{45,000}$).

Crossfading can include the native peak, but that is not done here since the smallest delay E1 provides a cleaner signal of the lowest $\rho$ than native spectrograph peak (robust to the FP noise seen in the native).

Crossfading has not yet been applied to bulk TEDI data. How well one can model the actual peak shape, and how much it varies across the band is under study.

10.5.

Simple Crossfading Simulation Using Two Delays

Figure 44 shows a simple demonstration using a synthetic source spectrum consisting a group of evenly spaced lines (a), so that a reasonably narrow frequency range around 1.4 cm is produced, so that only two peaks that straddle this frequency, at 1.2 and 1.6 cm, are needed. The source spectrum and moiré pattern are calculated twice, one in the original position, and again after artificially shifted by $\mathrm{\Delta }x=0.75{\mathrm{cm}}^{-1}$ to the left.

Fig. 44

Demonstration of the crossfading weighting technique on a simulated spectrum (a) consisting of several rectangular lines spaced $0.7{\mathrm{cm}}^{-1}$, so that Fourier energy (at $\sim 1.4\mathrm{cm}$) is between overlapping delay peaks at 1.2 and 1.6 cm. Moirés for the (a) smaller and (c) larger delays have opposite polarity slopes, and thus react oppositely in phase to a simulated PSF drift of $\mathrm{\Delta }x=0.75{\mathrm{cm}}^{-1}$ to left. The net reaction can be cancelled by strategic shaping of lineshape prior to combining (sinc function shape was used, triangle produces slightly more ringing). (d) Peaks in output spectra for 1.2 cm (green dash) and 1.6 cm (black) delays react oppositely (and differently than the envelope which shifts to left). Average of the 1.2 and 1.6 curves forms the final outputs (e), for original (blue dash) and artificially shifted (black) moiré positions. Peaks appear in their original positions, robust to the moiré shift. Change in peak position between blue dashed and black curves is only $0.0022{\mathrm{cm}}^{-1}$, 350 times smaller than the simulated PSF drift of $0.75{\mathrm{cm}}^{-1}$. That is, TRC $\sim 0.003$.

A key point is that the moiré pattern changes its slope (phase versus wavenumber) when the delay passes through 1.4 cm, from 1.2 cm (b) compared to 1.6 cm (c). Since the moiré flips its slope polarity, the idea is to mix the two results to produce a net cancellation to the insult of a $\mathrm{\Delta }x$ translation. The mixing happens not to the data in moiré form, but after it has been processed into an output spectrum (d). Both a triangle and sinc function crossfading shape were tested (sinc results are shown), and the triangle produces slightly more ringing.

The average of the 1.2 and 1.6 outputs form the final outputs shown in (e), for original (blue dash) and forcibly shifted (black) cases. These two cases overlap in position extremely well, producing a perceived shift of only $\mathrm{\Delta }\nu =0.0022{\mathrm{cm}}^{-1}$ for a $\mathrm{TRC}=1/350$, which is remarkably small.

10.6.

Correlated Point Spread Function Drifts Needed

The above crossfading method presumes that the same $\mathrm{\Delta }x$ PSF drift applies to the two overlapped delay peaks. Otherwise the cancellation will be incomplete. This is naturally satisfied by long-term drifts in a TEDI style apparatus which measures delays sequentially on the same pixels using different glass etalons that are rotated into the interferometer beam, where “long” means longer than the time to measure each delay (typically several minutes).

We propose that shorter term drifts would be more correlated between peaks using a stepped mirror (Fig. 45) which can be designed to manifest several delays simultaneously from the same source, although on different detector pixel rows. Then comparison of original (having shift) to crossfaded data (having no shift) yields information on the PSF shift $\mathrm{\Delta }x$, which thus could be at short times scales as well as long. However, this would be susceptible to any row dependence to the $\mathrm{\Delta }x$ drifts (such as thermally induced on the detector, and thus likely to be long time scale).

Fig. 45

Notional diagram of stepped mirror scheme⁴ for measuring multiple delays simultaneously with the same disperser, so that short-term PSF spectrograph drifts can be shared among at least two delays, to cancel drift using crossfaded data mixing. Staircase stepped shape of mirror and glass etalon of interferometer (dashed square) images different delays to different detector pixels, but passes through same disperser. This is undesirable if drifts are pixel row dependent. However, the stepped mirror could be combined with TEDI’s swapped delay scheme that uses same pixel rows, to discover both long- and short-term drifts on the same pixels for at least two delays.

By combining in the same apparatus a stepped mirror with a sequentially changed etalon similar to TEDI, information on short timescale variations from the stepped mirrors but on different pixel rows could be joined with that sequentially measured on the same rows. We speculate that this could produce a near zero TRC apparatus for Doppler velocimetry or high-resolution spectroscopy that was robust both on long and short timescales. The details of how to combine the data to best tease out the desired information are still under investigation.

10.7.

Implications for Precision Radial Velocimetry

In light of the great advantage of having at least two overlapping delays to cancel the effect of PSF drift, we encourage users of single delay EDI for velocimetry to consider modifying their apparatus to measure two or more delays quasi-simultaneously, such as in Fig. 45. Or if their existing single delay overlaps their native peak, they can apply crossfading to mix native and single delay signals, for frequencies up to but not beyond the delay. The decrease in TRC should directly reduce a horizontal (wavelength) noise, translating into reduced Doppler velocity noise if PSF translation is the dominant type of insult, and if not dominated by photon or read noise. (Example noise sources not addressed by the TRC could be a telluric contribution that is time varying and blended into the spectrum, and spatial drifts in the reference ThAr beam relative to the stellar beam, causing slight differential phase shifts in the interferometer.)

11.1.

Disadvantages of Externally Dispersed Interferometry Implementation

Disadvantages of implementing EDI include increase number of optics that the beam must pass through or reflect, and thus increases the parasitic loss, which decreases the flux. This is only relevant for photon-limited applications. We note that researchers have passed starlight through fiber scramblers²⁴ to reduce focal spot mode variation, while accepting some flux loss such scrambling entails. The increase in PSF stability is worth some flux loss, because it is critical to attaining precision high-resolution spectroscopy or Doppler velocimetry.

Another EDI disadvantage is that there are two interferometer outputs, which (unless polarization is used) necessarily travel on two separate beam directions. So additional optical engineering is required to detect both outputs, otherwise a 50% flux hit is suffered. This could be with optics that redirects the two outputs to different portions of the same spectrograph slit, or to two different spectrographs. For TEDI version 1 both outputs were detected,⁶ but for version 2 the complementary output was found² more useful as a monitoring port to maximize the amount of flux entering the input fiber.

Third, additional read noise is generated because multiple exposures (at least three) are needed. However, this disadvantage may often be immediately outweighed by the advantage of FP noise rejection that such phase stepping creates (where two exposures recently taken are compared for differences), since FP noises can often be more significant than photon noise. (Read noise behaves like detector noise, which is described in the companion paper.²³) Also, a recent paper⁹ proposes using a readout noise-free detector (EMCCD) with fast scanning of interferometer phase to improve Doppler precision.

11.2.

Advantages of Externally Dispersed Interferometry

Summarizing, the EDI method offers three significant advantages: (1) resolution increases beyond classical limits (Fig. 46), and (2) dramatic robustness to PSF irregularities (Figs. 41 and 42), and (3) superb rejection of FP noises such as from bad pixels (Fig. 33). While the resolution increase tends to be discussed most often, what may prove to be the greatest advantage is the robustness to PSF irregularities, since this is the most serious practical limitation to many precision Doppler velocimeters or high-resolution spectrographs in astronomy today.

Fig. 46

Notional diagram: the EDI surmounts classical limits to resolution from native spectrograph lens aberrations or detector pixel Nyquist sampling requirements. EDI resolution $\mathrm{Res}\sim {\tau }_{\max }/\lambda$ is independent of slit width.

11.3.

Summary

In summary, the EDI allows one to make high-resolution measurements not possible or practical with your current apparatus due to physical laws of optics, such as a nonzero focal spot size, sparse pixel spacing, bad pixels, and severe drifts and irregularities of the focal position and diameter. The bandwidth is set by the native spectrograph, which tends to be wide because it has lower resolution. Hence the bandwidth resolution product can be several times larger than in a conventional spectrograph, which is limited by the number of pixels. Spacecraft and aircraft have severe weight restrictions that make large spectrographs impractical. Since the cost and weight of a dispersive spectrograph grows in a power law with resolution, any technology that reduces the native spectrograph size by a factor of several to achieve a given resolution goal, offers a significant opportunity.