Fringe Measurement and Control for the Keck Interferometer

The classical ABCD algorithm, as implemented on PTI and other instruments, is used on KI for real-time measurement of phase and group delay, as well as for estimation of fringe visibility for science. With an integrating detector, the flux in the 1/4-wave bins A, B, C, and D is computed from the differences of the detector reads z, a, b, c, and d following an array reset as the fringe is rapidly scanned over one wavelength, viz. A = a - z, etc. The fringe quadratures are computed as X = A - C and Y = B - D, and total flux as N = A + B + C + D.

The fringe phase is estimated as ϕ = arctan(Y/X), and the fringe power as NUM = X² + Y², from which the fringe visibility squared is computed, after bias correction of NUM, as V² = (π²/2)NUM_c/〈N〉², and the photon-noise-limited fringe signal-to-noise ratio squared as ( S / N ) pn 2 = 2     NUM c / 〈 N 〉 ; see F99 for more detail. Bias correction of NUM is discussed in § 2.2, and alternative estimators for NUM and an extended expression for (S/N)² are described in § 2.4.

When the length of the path length modulation stroke used to scan the fringe is different from the wavelength, the quantities X, Y, and N, above, differ from their true values. It is convenient, especially for real-time use, to define the "dewarping" coefficients (see P99)

where s = 2πd/λ, d is the length of the path length modulation stroke, and λ is the wavelength: s = 2π when the stroke length matches the wavelength. Using these coefficients, the corrected quadratures and flux are given by

where

Squared visibility estimates from the fringe quadratures require correction for biases attributable to photon noise and detector read noise (F99; Perrin & Ridgway 2005). However, after dewarping the phasors to correct for mismatched stroke versus wavelength, there remain cross terms that make the NUM bias phase dependent. We extend the analysis in F99 for this case, below.

The expectations on the squares of the corrected quadratures X_c and Y_c can be written

or

The photon-noise variances per bin are shown as , etc., and the read-noise variances per bin as , etc., each with value equal to the correlated-double-sampling (CDS) read variance . With an integrating detector, as assumed here, adjacent bins are correlated, viz. σ_AB ≠ 0, etc., as adjacent bins share a common raw detector read.

The only signal-dependent biases in equation (5) are from photon-noise terms. The read-noise biases calibrate simply from the dark values of and after correction for detector offsets: call this sum E. Thus the total signal-dependent bias on NUM is 2GN, where

so that an unbiased NUM estimator when using dewarped phasors is given by

where k is the pixel gain in dn per electron. In terms of dewarped quantities only, we substitute N = N_c + πγ(Y_c - X_c)/(2s), yielding

If s = 2π, then G = 1/2 and γ = 0, and this expression reduces to the standard result (see F99, eqs. [9]–[11])

For an interferometer observing through atmospheric turbulence with coherence diameter r₀ and wind speed W, the phase-difference coherence time (τ₀₂ = 0.207r₀/W) is 4 times smaller than the variance coherence time (T₀₂ = 0.815r₀/W) (see P99). In order to maintain adequate continuity between phase samples to allow for unwrapping, we typically require the spacing between samples δt to be ≲τ₀₂. Similarly, in order to avoid excess fringe blurring during the integration time, we typically require the frame time T to be ≲T₀₂. For a typical "full-frame" ABCD implementation, the period of the path length modulation stroke is equal to the frame time, and one phasor estimate is provided per frame. Thus, δt = T, and the frame time is typically limited by the smaller τ₀₂ constraint on continuity, which impacts instrument sensitivity.

The ABCD implementation also affects the performance of the system when the phase estimates are used in a feedforward system to stabilize a second optical path. As discussed in § 5, the error rejection for a feedforward system depends strongly on the total control latency Δ, which includes contributions from the fringe measurement, delay-line bandwidth, and data processing. The component of the control latency attributable to the fringe measurement has two components: (a) the data age, i.e., the average delay between when a photon is detected and when a phasor is computed; and (b) the average sampling delay, equal to one-half of the sampling period. Thus the fringe-measurement latency for a full-frame implementation is approximately the frame time: Δ_f ≃ T.

However, for the same frame time, it is possible to reduce the spacing between samples, as well as to reduce the control latency, by updating the phasor estimate as individual detector reads become available, rather than waiting for a full frame's worth of data. For example, when detector read c is available, a new value of the bin count C = c - b can be computed. Then using bin counts A and B from the current frame, and D from the previous frame, an updated phasor estimate can be computed using the new information. This approach, sliding-window estimation, is possible with a sawtooth path length modulation stroke, and is implemented on KI.

The KI fringe tracker servo runs at 5 times the frame rate, processing each of the 5 detector reads as they become available. Assume for simplicity an ideal 80% duty cycle path length modulation stroke as illustrated in Figure 1. In this example, at the time of the c read at t = 0.6T, a phasor estimate can be computed from the A, B, C, and D bins centered at 0.1T, 0.3T, 0.5T, and -0.3T, yielding an effective data age of 0.45T. A similar computation for all five reads {z,a,b,c,d} yields data ages of {0.6,0.55,0.5,0.45,0.4}T. Note that at the time of the z read, no new information is available, and the data age is just 0.2T longer than at the time of the d read. Adding the average 0.1T sampling delay of one half of the sample period, the fringe-measurement latency for the five reads is {0.7,0.65,0.6,0.55,0.5}T, with an average value Δ_f = 0.6T. By comparison, for a full frame readout (at the time of the d read) for this same 80% duty cycle example, the fringe-measurement latency is equal to the data age 0.4T plus the sampling delay of 0.5T, yielding Δ_f = 0.9T.

Zoom In Zoom Out Reset image size

In practice, both the sliding-window and full-frame values will be slightly larger when the actual duty cycle and camera readout time are included. We measured the latency at KI from the phase delay between time-tagged fringe-phase measurements and time-tagged delay-line position measurements for a 10 Hz sinusoidal disturbance. The measured value for the sliding-window implementation is 2.6 ms for a T = 4 ms frame time; this is ∼1.2 ms less than a comparable full-frame implementation would achieve.

We complete this section by returning to the issue of phase unwrapping. With the sliding-window implementation, samples are available at 5 times the frame rate, but we need to account for the different data ages between adjacent samples. Continuing the example in Figure 1, the effective delay from the previous sample for reads {z,a,b,c,d} is {0,0.25,0.25,0.25,0.25}T, with a maximum δt = 0.25T. This is a factor of ∼4 smaller than for the full-frame approach. This maximum value should also be closely achieved in practice, as additional overheads apply to each sample, affecting mostly the overall latency, and not the sample-to-sample spacing.

Thus, with a sliding-window implementation, longer frame times for improved sensitivity can be used while still maintaining adequate phase continuity for unwrapping; or for a fixed frame time, improved unwrapping and feedforward performance can be achieved. For cophasing, KI uses sliding-window estimation with a 4 ms frame time and an update rate of 1250 Hz. For noncophased K-band science measurements, typical frame times are 5 and 10 ms, with update rates of 1000 and 500 Hz, resp.

The fringe power NUM, as computed from the sum of the squares of the fringe quadratures (and in the absence of full-rate photometric monitoring), is biased by scintillation as the frequency of the scintillation approaches that of the path length modulation stroke. At KI, ordinary atmospheric scintillation is usually quite small, because of the near-infrared wavelength and the aperture averaging of the 10 m Keck telescopes. Most of the observed scintillation is caused by fluctuations of the coupling into the fringe tracker single-mode fibers caused by angular vibrations, guiding errors, and Strehl fluctuations.

KI implements a multistate fringe-acquisition algorithm, similar to that used at PTI, using filtered and unfiltered versions of the fringe S/N squared

to transition among states: thus, high levels of scintillation can cause false transitions. Some amelioration is provided by (S/N)² thresholds that include a fixed term along with a flux-dependent term that sets a lower limit on V². However, it is desirable to keep the coefficient of the flux-dependent term as small as possible to allow tracking of extended objects.

To reduce the effect of scintillation bias, the real-time system at KI uses an alternative phasor estimator to compute NUM to estimate (S/N)² using equation (10) for determining mode transitions. We now compare the scintillation properties of the standard phasor estimator and two alternatives.

2.4.1. Standard Phasor Estimator

With the standard phasor estimator, as described in § 2.1, the scintillation bias in the fringe power NUM = NUM₀ arises because the quadrature estimates X = A - C and Y = B - D are essentially derivatives of the incoming flux I(t), i.e., for a stroke period T, X(t) ∼ I(t) - I(t - T/2) and Y(t) ∼ I(t - T/4) - I(t - 3T/4). For this discussion, we assume a full-frame phasor estimator with a 100% duty cycle, and neglect for now the integration over the bin duration T/4. Thus, in the frequency domain, X(ω) ∼ [1 - exp(-jωT/2)]I(ω), and for unit scintillation at normalized frequency ω_n = 2πf_scint T, the bias on NUM₀ is given by

This has an intuitive form, with the bias from scintillation increasing quadratically as the scintillation frequency increases.

For comparison with next two estimators, the read-noise bias on NUM₀ is just the last term of equation (9) (assume pixel gain k = 1):

For the standard phasor estimator, the standard deviation of NUM in the read-noise limit is just E₀ (see F99). Thus, we use the read-noise bias on NUM as a proxy for the noise on the (S/N)² estimator.

2.4.2. Modified Phasor Estimator 1

We can create quadrature estimates that are less sensitive to scintillation by subtracting out a local trend from adjacent bins, i.e., we can compute X₁ = 2[(B + D)/2 - C] and Y₁ = 2[B - (A + C)/2]. With respect to scintillation, these look like centered derivatives, and thus reject low-frequency terms. The bias for unit scintillation on the fringe power estimate computed from these quadratures is given by

Compared with the standard NUM estimator, the scintillation bias dependence with frequency is now to the fourth power, and the bias drops rapidly with frequency below the corner frequency ω_n = 8.

However, this improved rejection comes at the expense of read-noise bias. With nondestructive reads, adjacent bins are correlated. Writing NUM₁ in terms of detector reads yields NUM₁ = (-a + 3b - 3c + d)² + (z - 3a + 3b - c)². The reset noise drops out in this expression, so that the effective noise variance per detector read is 1/2 of the CDS value, and the bias on NUM₁ is given by

Because of its ease of computation—it just needs one frame's worth of bin data—and its scintillation-rejection properties, this algorithm was used on KI for the first several years in the computation of (S/N)² for mode transitions.

2.4.3. Modified Phasor Estimator 2

An alternative NUM estimator can be created which retains the scintillation rejection of previous version, but with a slightly larger leading coefficient on the bias expression. This estimator uses data from adjacent frames to compute the necessary fringe quadratures. Thus for each new bin count A, B, C, and D, we compute a flux trend by comparison with the prior value of that bin, i.e., on a new A, we compute I^' = A - A_last, etc. Thus, in the sliding-window computations, we update the total flux N using the sum of the four most recent bin counts, and update the quadrature estimates for each bin as

On a new A, X₂ = A - C - α(A - A_last)/2;

On a new B, Y₂ = B - D - α(B - B_last)/2;

On a new C, X₂ = A - C - α(C - C_last)/2;

On a new D, Y₂ = B - D - α(D - D_last)/2.

The factor of 1/2 converts the trend measured over 4 bins to the value for a 2 bin difference, and the factor α is for duty cycle correction; with 100% duty cycle (integration time/frame time), α = 1; more typically, α = 0.8–0.9. The bias, for unit scintillation, for a fringe power estimate computed using these phasor estimates (assuming α = 1) is

While not quite as effective at rejecting scintillation as the previous estimator, this NUM estimator does retain the fourth power dependence with frequency, and has much lower noise. On a new D, NUM₂ = (A - (C + C_last)/2)² + (B - (D + D_last)/2)², and thus

which is even smaller than for the standard estimator because of the effective averaging of C and D. The tradeoff for low noise is that the effective frame time is slightly longer than the standard estimator. We find that this NUM estimator works well, and after bias correction using equation (8), is our default in the real-time system for estimating (S/N)² using equation (10). However, all other phasor computations use the standard estimator.

The scintillation bias, for unit scintillation, of the three NUM estimators discussed in this section is plotted in Figure 2. The plot also includes the common factor of sinc²(πf_scint T/4) to account for integration over the bin duration T/4.

Zoom In Zoom Out Reset image size

The KI fringe tracker uses two post-combination single-mode fibers on the two outputs of the Michelson beamsplitter to feed an infrared camera,⁷ and the camera is typically configured to provide a white-light channel on one output and a spectrometer channel on the other (Vasisht et al. 2003). The white-light channel is typically used for phase estimation for high bandwidth fringe tracking, while the spectrometer channel is used for group delay estimation for central fringe identification. This approach for phase and group delay measurement is similar to that used at PTI.

At KI, we use the sliding-window estimator described in § 2.3 to provide sets of dark-corrected bin counts A, B, C, and D for the white-light pixel and each spectrometer pixel. With the sliding-window estimator, the computations are done at 5 times the frame rate.

For the low spectral resolution modes, we can improve the S/N for phase estimation by combining the white-light and spectral channels to create a synthetic white-light channel. We can create a set of synthetic white-light bin counts by summing the spectrometer bin counts and combining with the white-light bin counts as⁸

From either the white-light or synthetic white-light bin counts, we compute the complex phasor Γ_WL(n) = X_WL(n) + jY_WL(n), dewarp, and compute the phase as ϕ(n) = arg[Γ_WL(n)], as in § 2.1. The raw phase is unwrapped using a first-difference continuity unwrapper with a small amount of memory, viz.

where mod _-w/2,w/2(x) is the modulo of x over [-w/2,w/2) and 〈f(n)〉_N is an N-point boxcar low-pass filter. Typically N = 2 for this unwrapper. The phase delay x(n) is just Φ(n)λ/2π.⁹

We compute the group delay g(n) using the spectrometer phasors Γ(n; k) as a function of wavenumber k. To improve S/N, the spectrometer phasors are phase referenced to the WL phase and averaged as

Typically M = 60, or 50 ms for the nuller configuration. The group delay g(n) is then computed from the Fourier transform vs. wavenumber of the average, phase-referenced spectrometer phasor array (see P99). The DFT is zero padded by typically 8 times to avoid scalloping losses, and a quadratic interpolation about the peak is done to improve the group delay estimate.

KI implements different fringe tracker servos for V² and nulling. The purpose of the tracker is not only to center the phase track point near the peak of the coherence envelope, but also to correct for fringe hops in the continuity unwrapper.

The servo implementation for V² observations uses a fast inner loop and a slow outer loop, as illustrated in Figure 3, to generate servo targets c_K for the fast delay line. The inner loop uses the phase delay x with an integral gain a₁; a limiter preceding the integrator constrains the maximum slew rate to prevent fringe blurring. Additional vibration rejection is provided by control block H_n, which will be discussed in § 3.4.¹⁰ We typically use an integral gain a₁ = 0.1, which provides, for a 5 ms frame time, a 16 Hz closed-loop bandwidth. It is worth noting that for V² measurements, only the motion during a frame time reduces the fringe visibility, and thus high bandwidth is not required: bandwidth is only required to stay within the coherence length of the channel used for the V² measurement.

Zoom In Zoom Out Reset image size

The outer loop uses the group delay g with an integral gain a₀; typically a₀ ∼ a₁/10, providing an outer loop bandwidth of 1/10 of the inner loop, typically 1.6 Hz for a 5 ms frame time at KI. The outer loop implements the fringe centering, by providing the target for the phase-tracking inner loop. Thus the outer loop corrects for the arbitrary zero point of the unwrapped phase, as well as for fringe hops which may occur in the continuity unwrapper. The correction occurs slowly, at the time constant of the outer loop, so there is insignificant fringe blurring from the changing phase target. By centering continuously, rather than stepwise, no synchronization or coordination of centering commands is needed, and this simplification is one of the advantages of this implementation.

The outer loop also accepts an external target. For V² measurements, we typically introduce a triangular group delay target modulation, m, with a one-wave peak-to-peak amplitude and a period of 5 s; the period is chosen to be short compared with the length of a science scan. The effect of the target modulation is that the average V² for the scan includes measurements over a near-uniform distribution of phases. This greatly reduces the effect of errors that are periodic with a wavelength, such as those that are attributable to errors in the assumed wavelength, or residual nonlinearities in the path length modulation stroke.

It is important to note that while phase delay and group delay are shown as logically equivalent in Figure 3, they differ in practice due to systematic and random dispersion. The systematic dispersion is attributable to average air and material dispersion, with a sidereal component from changing air paths as the delay lines move. The random dispersion is primarily attributable to water vapor seeing (see Colavita et al. 2004a). In the inner-loop/outer-loop implementation, the dispersion gets pushed into the phase target. As this is a low-bandwidth effect, it is not especially important for V². In fact, for V², it is the proper behavior, as we wish to track to the center of the coherence envelope, including the effects of dispersion, in order to maximize V².

For phase referencing, as implemented for the KI nuller, we wish to track to constant phase, rather than to constant group delay, as in the previous V² implementation. The implementation we use for phase referencing is shown in Figure 4. In the figure, the inner loop, using the phase delay, remains unchanged from Figure 3. What does change is the source of the phase target, which was previously a continuous target from an outer-loop group delay servo, and included dispersion as well as integer fringe hops. In Figure 4, the phase target is now always an integer number of wavelengths, so we always track to "zero" phase. In principle, we could implement this simply, adding or subtracting wavelengths to the phase target to keep the group delay within ± 0.5λ. However, with varying dispersion, especially from water vapor seeing, the phase target would continually toggle. For phase referencing, at least for reasonably short tracks where the change in the systematic dispersion due to sidereal motion is less than a wavelength, it is better to follow a single phase point as it wanders within the coherence envelope. For longer tracks, it is necessary to periodically reset the track point in order to not drift too far from the center of the envelope.¹¹

Zoom In Zoom Out Reset image size

Thus, in Figure 4, we compare the group delay with a slow, nonwrapping average of the dispersion 〈g - x〉; we describe how this is computed in § 4. A comparator, with thresholds of typically ± 0.75λ, increments or decrements a counter in units of λ. The output of this counter passes though a short boxcar filter, B, which limits the slew rate of the step, and then to the input of the phase loop; typically B is 20 samples long. The rate at which the wavelength counter can change is limited to allow fringe hops to be fully reflected in the group delay, including the time it takes to fill the group delay coherent average filter, before another hop can occur. Typically, the limit is 100 samples between hops.

With typical closed-loop bandwidths of 16 Hz, the basic fringe tracker has limited ability to attenuate high frequency disturbances. However, at KI there are residual disturbances at high frequencies that are not observed by the laser metrology or accelerometers. Most of these disturbances are narrowband vibrations, and typical sources are induction motors which introduce vibrations at ∼29 and ∼58 Hz (for US 60 Hz ac power), as well as "mystery" sources at other frequencies. For such narrowband disturbances, it is possible to actively suppress them to improve fringe tracker performance, even though their frequencies exceed the ordinary fringe servo closed-loop bandwidth.

The control block used for vibration suppression is shown in Figure 5, and implements higher harmonic control (HHC; Shaw & Albion 1981). This block is inserted in the fringe tracker forward path as shown in Figures 3 and 4, and provides loop gain at the frequency of the vibrations. As shown in Figure 5, the overall control block has a pass-through path on the bottom, and multiple HHC blocks for each vibration frequency under control. While each control block typically controls five frequencies at KI, we will generally assume a single controlled frequency for the rest of this discussion.

Zoom In Zoom Out Reset image size

In Figure 5, quadrature oscillators configured to the nominal vibration frequency mix the input signal to baseband. The phase shift ϕ₀ on the input oscillators compensates for phase shifts in the forward path of the fringe tracker at the frequency of the vibration, including the phase of the ordinary fringe tracker integral controller. The controllers on the baseband disturbance signals are implemented as leaky integrators with gain k₀. The small leak δ is so that the oscillators decay to zero when the fringe tracker enters a hold state on loss of lock, although it does limit the ultimate vibration rejection. The integrator outputs are then mixed back up to the original frequency and summed with the pass-through signal.

The response of the HHC controller is given by Hall & Wereley (1989; summarized in Sievers & Von Flotow 1992). Including the pass-through, and the finite leak in the HHC integrator, the transfer function of the vibration control block is

where a = cos ϕ and b = sin ϕ. Thus, H_n(ω) provides high gain at the configured frequency, and is approximately unity elsewhere.

At KI, H_n(ω) is used in series with an ordinary servo controller G(ω). While the approach is generalizable, we are normally interested in providing suppression above the unity-gain frequency of the ordinary servo, i.e., where G = |G(ω₀)| ≪ 1; we assume G incorporates any dynamics of the plant. In this case, near the configured frequency, the second term of equation (20) is ≫1, and so we can approximate H_n(ω) near ω₀ as

In our application, we typically set ϕ ≃ arg[(G(ω₀)]. The 3 dB half width of these active notches Δω, for Δω≫δ, depends on the product of the HHC gain and the ordinary loop gain as Δω ≃ Gk₀ rad/s, with peak suppression ≃Δω/δ.

To first order, the maximum notch width is limited by the requirement that the ordinary controller and plant be approximately constant in phase and amplitude over the notch width in order to control the phase at the loop crossovers on either side of the notch. This constraint is similar to that which limits the bandwidth of the ordinary controller, and thus the notch widths are typically limited to a fraction of the ordinary closed-loop bandwidth.

We find these modules to be robust and effective in suppressing vibrations with well-known frequencies, which is the case for the largest-amplitude vibrations at KI.¹² The HHC approach also allows for accurate single-parameter tuning of the notch frequency. For a different approach to narrowband suppression, see Di Lieto et al. (2008), which describes the adaptive algorithm used on VLTI's FINITO fringe tracker.

As the ordinary controller gain, while less than unity at the frequencies where these notches are employed, is still finite, and there are multiple notches which can interact, we usually tune the notches empirically. In our implementation, we include a test oscillator in the HHC block to measure arg[(G(ω₀)], and we set the notch gain based on frequency-response measurements in order to balance notch depth with band-edge peaking. Typically bandwidths of order 1 Hz are used, which allows for some uncertainty in the vibration frequency, e.g., to account for shifts in the frequency of induction motors with changing loads. For example, with a 1 Hz notch width, and δ ≪ 1 Hz, a factor of ∼4 suppression is possible for uncertainties in the vibration frequency of ± 0.125 Hz.