The newest version of the defringeflat package is available at http://www.das.uchile.cl/~pato/sw. For a more detailed description of the algorithm, please refer to Rojo & Harrington 2006 (ApJ submitted).
This routine removes a fringe pattern with repetitive characteristics from a flat field.
A typical flat field array from a VLT observation using the ISAAC instrument is shown in Figure 1. In the upper frame a periodic fringe pattern is clearly visible. It can also be seen that the period is not constant throughout the array. Instead, it changes visibly, for example, in the upper right quadrant.
This change in periodicity gets confirmed in the lower frame of Figure 1 where only shortwards of pixel ~400 averaging brings up the signal-to-noise. Longwards of that same pixel the average gets flatten because of the fringe period's changes. region.
Also, the smooth change in period inmediately explains while you
cannot use methods like Fast Fourier Transform (FFT) where only one
frequency per row can be found.
Figure 1: Upper frame is ISAAC's flat field image obtained after
subtraction of dark frame with same exposure time. Lower frame is
the average of every row.
Template id 'ISAACSW_spec_cal_NightCalib'. Integration time: '4
sec'. FITS file 'flat.fits' is supplied with defringeflat package.
A more detailed, quantitative introduction to the wavelet transform is given by C. Torrence and G. Compo, only a brief introduction is presented below.
In short, the wavelet transform calculates the frequency spectrum in a region around each point of a one-dimensional array. The resulting array shows the spectrum's variation with location in the array. The main and critical difference from the windowed Fourier transform.
Figure 2 shows an example of a periodic function (a sine function whose period changes linearly from 70 to 40 pixels) and its wavelet transform. The lower frame shows the amplitude of the wavelet transform. The region in the shaded area is inside the "cone of influence" of the edges of the data, values in this region cannot be trusted.
A trace in the transform that peaks at the expected frequency is
clearly visible. However, this trace has a shape with some
non-zero width. As can be seen in cross-section (Figure 3), the shape is fairly approximated by a
Gaussian function.
Figure
2. Upper frame: sinusoidal function with changing period and
amplitude. Lower frame: amplitude of a Morlet wavelet transform of the
spectrum.
Figure 3. Vertical cross section of the wavelet transform's
amplitude (solid line) and the best-fit Gaussian (dotted line).
The routine defringeflat uses the following procedure to compute the fringe pattern of a flat field. A priory, the user must identify the boundaries of the sub-array containing the significant values. Also, the user must orient the frame such that, only in a per-row basis, the period of the fringes is at least several pixels but may not exceed about a quarter of the row width.
To improve the signal-to-noise ratio and remove "wild" pixels, a 1-D median average of a specified number of rows is computed for each of the rows (Figure 4). For each of the filtered rows, a polynomial fit is subtracted to obtain an enhanced row.
A wavelet transform is performed according to specified wavelet parameters, and then its amplitude possibly interpolated in period (Figure 5).
Then, the trace of the wavelet amplitudes that corresponds to the fringe system is identified. A specified function (see next section: Trace fitting methods) is fit to the trace profile using a user-identified fringe period (at some reference image coordinates) as the initial guess (Figure 6). Then we follow the local maxima across the columns. Repeating this procedure for each row, we obtain a fitted set of parameters (e.g., center, width, and amplitude of a Gaussian) for each pixel in the array.
If specified, each parameter array is compared to a median-filtered version of itself to identify and replace outliers, and then it is mean-smoothed (Figure 7).
Then, the wavelet transform amplitude is reconstructed from the (smoothed) parameters obtaining the fringe's wavelet transform (Figure 8). When reconstructing, the background of the fitted functions (whenever available) can be kept or not. Not keeping it is the default because it allows for some non-fringe component in the fringe's period to stay in the frame.
Finally, an inverse wavelet transform is performed to obtain the
computed fringe pattern (Figure 9), which is then
removed from the flat (Figure 10).
Figure 4. Flat field in Fig. 1 after being median
averaged in the vertical direction. In this particular example, each row
is the median average of the 40 rows surrounding it.
Figure 5. Upper frame shows middle row of Fig. 4. Middle frame shows the
amplitude of the wavelet transform in the relevant periods. Lower
frame: Same as middle frame but interpolated in period to 10 times the
original period resolution.
Figure 6. Cross section of wavelet transform's middle column, original
points are always represented by small crosses. Upper frame shows
original points and the interpolated profile (solid line). Middle and
lower frames show the profiles fit (solid line) when using methods
'Gaussian' and 'true shape', respectively.
Figure 7. Parameter smoothing. From left to right the parameters
are amplitude, center, width, and background of a Gaussian fit of
the trace. From top to bottom, they are shown as raw
fitted, patched (bad pixels replaced), and smoothed parameters.
Figure 8. Fringe's wavelet. Amplitude of the reconstructed fringe
trace.
Figure 9. Fringe obtained from inverse transform of Figure 8.
Figure 10. Flat field from Figure 1 after being cleaned by
defringeflat.
To run this example you would need the wavelet routines from C. Torrence and G. Compo and the astrolib routines from NASA.
IDL> flat = readfits('example/flat.fits', h)
IDL> func = 'gauss'
IDL> clean = defringeflat(flat, 40, binw=40, trim=[10,145,950,950], $
osamp=10, /smooth, fitname=func, fringe=fringe)
IDL> tvscl, flat
IDL> tvscl, fringe
IDL> tvscl, clean
We tested the algorithm's performance at different noise levels
by using a synthetic image consisting of a fringe pattern, a
background intensity, and random noise with a Gaussian distribution that
mimics pixel-to-pixel flat-field variations and photon noise.
We define noise strength as the standard deviation of the
Gaussian noise divided by the standard deviation of the noiseless fringe
pattern. Figure 11 shows the fraction of remnant
fringe after running the algorithm on simulated data with different
fitting functions and varying noise strength. The remnant fringe level
is not strongly dependent on noise strength.
Remaining fringe for
varying noise strengths. The two nokeep methods give the best
results. In general, all methods give similar results for low noise
strength, but smoothed nokeep methods are better for high noise
strength. Top: Fraction of fringe remaining. Bottom: Absolute fringe
remaining when noise level is scaled to 1.