This document describes how to reduce
raw IVM magnetogram data using the IDL procedures in /solar/IVM/idl. The
program includes all the techniques that we have developed for handling
systematic errors; if all the procedures are used it may take a long time
to run. Omitting the corrections for scattered and parasitic light, and
for seeing, it is also the one we use for reducing the morning survey magnetograms
in near real-time. Capabilities are:
compute the flat-field corrections in the best way we know,
coalign successive images by rigid translation, using the geometry
images to determine the shifts,
correct for scattered light in the instrument,
correct for spectral leak,
reduce seeing-induced spurious polarization,
locate the line center spectral position pixel by pixel, and
compute the magnetic field parameters using a decent implementation
of the derivative method.
The software described here is the version of July 1996.
The magnetogram raw data file is named IMGRyymmdd.hhmm, with the lower-case
characters replaced by the date and time of the observation. The data file
is typically 62 MB, though sometimes smaller.
There are five ancillary files needed:
two dark frames, named IDK[01]yymmdd.hhmmss,
two flat-field data files named IFA[01]yymmdd.hhmmss, and
a calibration data file named ICARyymmdd.hhmm.
These files total about 20MB.
Output is to a FITS format file, named IMGMyymmdd.hhmm, containing
six images: continuum intensity, field line-of-sight component, field transverse
component, field transverse azimuth (referenced to celestial north), line
center offset in A, and line-of-sight velocity in km/sec.
Programs:
The driver is a function called rivm. It calls a number of other
routines, many of which are in the IVM directory, /solar/IVM/idl. Some
Yohkoh software routines are used. A precompiled, complete set of
needed routines is in rivm_run.exe. This file may be copied to another
machine and loaded to IDL with a restore command. Very handy for
running on other people's fast machines.
Calling syntax
gram = rivm(calfile,datafile[,fitmap])
The names of the calibration and data files are the only required arguments.
Default dark and flat files are the most recently observed dark and flat,
before the data observation. Their names are recorded in the data
file header and extracted by rivm. If different dark and flat files
are needed, they can be specified with keyword variables. The optional
argument fitmap is used to specify an output variable which will
hold a binary map of the success of the intensity profile fit. All other
arguments are keywords.
A useful way to remember the options in the call is the file runrivm.pro,
which may be edited and used as a batch file: @runrivm. An example
of this file contains the following (the keyword jump wouldn't
be used normally; this was for testing the deblur code):
Note that calf and datf must be defined as the appropriate filename
strings first.
Procedure
Rivm executes the following steps in the reduction process:
SETUP. Open a journal file in the current directory, named 'rivm_yymmdd_hhmm'.
The timestamp in the file name is the time of the data reduction. This
journal is useful later on, since it contains printouts of parameters in
the flat-field fits, the coregistration shifts, and sample results from
the magnetic-field derivation.
Read the header (Observation block or Oblk) from the data file. Extract
from it the image size, number of spectral steps, number of spectrum repeats,
Fabry-Perot settings for each spectral step, wavelength of the line. Determine
if 4-frame or 6-frame modulation sequence was used.
FLAT FIELD. Calculate flat-field images for each wavelength,
using procedure parfl2. The details are described in the header
of parfl2, but the intent is to separate the overall vignetting
from features which depend on the individual wavelength setting. This must
be done without distorting the line profile in the process.
Data flats: The data-camera flats file has one image for each wavelength
setting. First the continuum-wavelength image is read and dark-subtracted.
Next it is corrected for scattered light with a call to ivm_scatter,
if the keyword scat is set. A two-dimensional
parabolic fit to the image is computed, called conpar; it is used
later in computing the flats at other wavelengths. The continuum (inverse)
flat is simply the inverse of the scatter-corrected image, multiplied by
its mean value over its central quarter. Each of the other flats may have
a legitimate field angle variation if it's in the wing of the line where
brightness varies strongly with wavelength. To avoid removing this, we
divide the dark- and scatter-corrected flat image by the array conpar
mentioned above. Then fit a parabola to the remainder, divide it out
to leave only the high-frequency variation, then remultiply by conpar.
Finally, invert and multiply by the normalization constant.
Geometry flats: The geometry flats contain little spectral variation,
but there is some. We ignore it in the scattered light correction by computing
scattered light only for the mean of the images. Note:
check the normalization of scattered-light correction in this code in parfl2.
The flats are just the mean of the central area, divided by
the image.
CALIBRATION MATRIX. Call split_fil to separate the data-camera
images out of the calibration data file. This is currently a no-op, since
there are only data-camera images in the first place. Call ivm_calcal
to compute the Stokes calibration matrix, as follows: Apply the flat to
the calibration data, incorporating scattered light correction if keyword
scat is set. The continuum wavelength setting
is used for all the cal data, so only one of the flat frames is needed.
Demodulate the calibration data and construct the calibration matrix, which
consists of a 4x4 matrix at each pixel in the image plane. The matrix is
smoothed using a 16 x 16 median filter, since it is expected that pixel-to-pixel
variations will be only noise. Invert the 4x4 matrix at each spatial pixel
and write the result to a file named 'ICAR..._matinv'.
SPLIT FILES. Split the data file into two files, one for data-camera
images and one for geometry-camera images--unless there are only data-cam
images, as for some early datasets. The output files are named IMG0...
and IMG1..., with the remainder of the filenames copied from the data filename.
FLAT. Flat-field the data-camera images, subtracting a mean
dark level. Write the result to a file named 'IMG0..._df'. No scattered
light correction is done at this point.
SHIFTS. If there exist geometry images, flat-field them (without
scattered light correction) and then use them to compute the spatial shifts
of each image in the set, compared to a reference image (the one at the
middle of the set). The routine called is shifter_2d; it does an
FFT cross-correlation between the sample image and the reference, then
computes the shifts by finding the peak of the correlation function. Use
these shifts to coalign the data-camera images, in routine shifterd.
The assumption here is that the majority of the difference between successive
frames which is caused by 'seeing' is a rigid translation, the same for
the geometry and data cameras. And that the geometry camera sees no changes
due to polarization modulation, since it's looking at a wide wavelength
range. Any higher-order image distortion is ignored. We might be able to
reduce the seeing-noise by 'destretching' the images, but it would be time-consuming
and possibly subject to other noise in the data which would make destretching
ineffective. The files produced are called 'IMG1..._df' and 'IMG[01]..._df_re'.
The [x,y] pixel shifts are written to an IDL save file named 'shiftsyymmdd_hhmm',
identified by the time of processing. These data may be useful for debugging
if the results seem odd. The values are usually around a tenth of a pixel
if the guider was working well; a pixel or two if the guider was off. Many-pixel
values indicate pathology: occasionally the guider will jump out of the
spot for one frame, then return, but lots of large values probably means
the flat data are bad.
WARP. Compute the geometrical transformation matrix between
the data and geometry images. Routine ivm_matchcams uses setpts
to establish (with operator intervention) the coordinates of a set of features
in the average data image and the corresponding features in the average
geometry image. Calls caltrans to compute the best-fit linear transformation,
saves the matrix in a file with the suffix '_match'. Then uses poly_2d
and the transformation matrix to warp each of the geometry images.
SPECTRAL PARASITIC SUBTRACTION. Calls ivm_parasitic to
subtract a selected fraction of each geometry image from the corresponding
data image. This basically subtracts a constant from the spectrum at each
pixel, and is intended to compensate for leakage from the adjacent orders
of the Fabry-Perot which are not completely blocked by the prefilter. The
fraction used is 7 %, obtained from filter models and hard-coded into rivm.
CORRECT SCATTERED LIGHT. Optionally, subtract a scattered-light
image from each frame in the dataset. The amount is set by the keyword
SCAT. The purpose is to correct for light diffusely scattered in the IVM
instrument, as by the retarders, but it is a rather ad hoc correction based
on residual sky brightness in limb observations. Use with much caution.
If the correction is done, the output is written to 'IMG0..._df_re_av_sc',
and the '_sc' persists through successive stages. Routine ivm_scatter
is used to compute a correction image for each frame in the first spectral
repeat for the data camera, and one for each wavelength setting for the
geometry camera. Corrections are applied in routine descatter.
DEBLUR. This appears to be a powerful way of reducing spurious
polarization signals due to seeing. The idea is to use frame-to-frame changes
in the geometry data to characterize the seeing in such a way that we can
get a correction for the data-camera images taken at the same time. The
method is to difference a test image from a reference image, divide the
difference image into tiles, then for each tile model the difference signal
using the first and second spatial derivatives of one of the original images
as basis functions. Usually the X and Y gradients and the Laplacian are
used: the coefficients represent, in some sense, the vertical and horizontal
displacement of the tile between the two frames, and the difference in
blurring. The first step is to devise a reference geometry image. This
can be simply the average of all the geometry frames, but what is done
now is to compute the coefficient array representing the difference of
each frame from one arbitrary frame, then average the coefficients to compute
a correction which is then applied to the arbitrary test frame. This becomes
the standard. Then the difference of each geometry frame from the standard
is modeled, using the derivatives of the individual frame (not the standard)
as the basis set. Finally these model coefficients are applied to the derivatives
of the corresponding data-camera frame to compute a correction image, which
is then added to the data-camera frame.
AVERAGE. If the spectral scan is repeated in the observation,
the spectra are averaged. We have often used a setup with 20 spectral points,
repeated twice. Write the result to 'IMG0..._df_re_av'.
DEMODULATE. Demodulate the data, to obtain uncalibrated Stokes
images, using demod or demod2 depending on which modulation
sequence is appropriate. Write the result to 'IMG0..._df_re_av_st'.
CALIBRATE. Apply the inverted calibration matrix to the Stokes
data, to remove instrumental offsets and crosstalk, and correct the scaling
of the Stokes parameters. Write the result to 'IMG0..._df_re_av_st_ca'.
Optionally, subtract a background corresponding to leakage of the shutter
during CCD readout. This was a zero-order correction for the leakage of
the Ferroelectric liquid crystal shutters previously used. If applied,
the output file is 'IMG0..._df_re_av_st_ca_sh'.
MAGNETIC FIELD INFERENCE. Compute the magnetic field, using
the JLS method of comparing the polarizations at a chosen wing position
to the derivatives of the intensity profile at the same wavelength position.
Multiplies the normalized polarization parameters by intensity to recover
Stokes data. Then computes an average intensity spectral profile from a
region near the center of the image, and fits it with a Gaussian absorption
profile on a continuum which may have slope and curvature. This canonical
I fit is used as a starting point in the remainder of the processing. For
each spatial pixel, extract the Stokes profile, then call bjeff6,
which first fits the I profile. The fit at the previous point is used for
starting parameters in the least-squares fit, unless the previous fit was
invalid. The intensity derivative at the chosen offset DL is computed from
the fit. Smooth the Q,U,V profiles slightly using a convolution with a
symmetric kernel, then combine data from several wavelength points near
DL on both the red and blue side of the line, fitting points within 40
mA of DL with a quadratic curve. This little dance permits us to use several
spectral points to reduce the noise in the polarization data. Optionally
return an image of 1/0 points indicating the validity of the I profile
fit at each point.
Return the magnetic field data as a Magnetogram structure consistent
with Haleakala Polarimeter data and other reduction programs, first clipping
any aberrant data values which would invalidate the automatic scaling used
in creating the FITS file. Also write the data to a FITS file, with a header
intended to preserve much of the data in the original Oblk. The output
data includes maps of both line-center position and line-of-sight velocity,
and heliographic position (although the latter should not be taken too
seriously due to uncertainties in the basic IVM pointing methods).
Options
The options are many, but the default values are usually right so the
options are often not required. Fitmap is an optional argument; the remainder
are implemented as keywords. Listed alphabetically.
CALSMOOTH
Numeric entry which overrides the default 16-pixel smoothing which
is applied to the calibration matrix.
D0FIL, F0FIL
Use these keywords to override the dark and flat filenames in the header.
Useful for old data with these entries missing; also if the appropriate
dark/flat was obtained after the data file rather than just before it.
DETREND
Set to take limb-darkening out of correlation subimages in shifter_2d.
DL
Offset, in Angstroms, from line center at which the program is to measure
the polarizations and intensity derivative. The default value is 0.100
A, which appears to work the best for small fields, though a larger value
like 0.14 is less subject to saturation and magneto-optic effects.
FITMAP
Variable into which a yes/no map of the Intensity fit validity is returned.
FOURSTEP
Set to force demodulation using demod_2, assume the data and
calibration were obtained with four-step modulation sequence. This is the
default if the header item NFRAMES is 4.
GIF
Set to write gif-format images as well as FITS.
JUMP
Set to one of the strings 'shift', 'warp', 'para', 'scat', 'blur',
'avg', 'demod', 'calib' or 'bgram' to jump into the processing sequence,
skipping already executed stages. This may be useful if one wishes to do
the magnetic field computation at a different wing position, for example,
and the calibrated Stokes file has already been obtained.
KEEP_INT
If set, keep all the intermediate files. The default is to delete most
of them, retaining the inverse calibration matrix file, the inverse flats,
and the final calibrated Stokes file.
MATCHED
Set this to specify that the coalignment between the geometry and data
cameras has already been determined. The keyword should be a string containing
the name of the file which has the 4x2 transformation matrix.
NFUS
Numeric entry to specify the interval in Fabry-Perot units between
successive wavelength points. Normally extracted from the header.
NO_BGRAM
Abort processing after computing calibrated Stokes data; do not compute
magnetic fields.
NO_REG
Set to skip coalignment of frames. Usually coalignment is preferred.
NPTS
Numeric entry to specify the number of spectral positions, for files
where this number is not found in the header.
OLD_DATA
Set this keyword if the observations are from the very early days in
1992, before the header sizes were set. It causes RIVM to query the operator
for the correct header sizes.
PARASITIC
Set to correct for parasitic light: leakage from adjacent orders of
the FP into the spectrum.
PLOT
If set, the interval in sequential pixels between plots of the Stokes
profile fits in bjeff6. These plots are useful for estimating the
quality of the data and of the fits, though they increase processing time
somewhat. A value of 250 or 500 may be appropriate.
REDCAL
Set if the inverse calibration matrix is present in the working directory.
Skips processing of the calibration data file. Again, handy with flare-mode
data where one cal serves several magnetogram datasets.
REDFLAT
Set if the flat-field data have already been reduced and the inverse
flat files are present in the working directory. Skips processing of raw
flat data. This is helpful with flare-mode data, i.e. rapid sequences of
magnetograms depending on a common flat.
REGPARM
Array [bx,by,hx,hy] specifying subarray to use in coalignment. bx,
by: size of coalignment window (default: 128,128). hx, hy: centerpoint
of coalignment window in image frame (default: center of image).
SCAT
Set to apply the diffuse scattered light correction. Scat may be a
two-element array containing [level,power] or a scalar containing
level for ivm_scatter.
SHUT
Set to apply the shutter leakage correction using shuttersub.
Only needed for early data with LC shutter; doesn't work too well for that.
Set shut to a string containing the name of the correction file.
UNBLUR
Set to use geometry camera images to correct data for seeing distortions.
Set unblur to a value >3 to select a tile size.
WARP
Set to warp geometry images to coalign with data images.
Time and storage requirements
The intermediate files are by default deleted, but 130 MB of disk space
is required for short-term storage and the Stokes file. If the KEEP_INT
option is selected, about 400 MB of free space is required, in addition
to the 80 MB needed for the data files. The output FITS file is 660K for
a 256 x 256 image size, and can be compressed a little, to about 450K.
Processing takes about 35 minutes for a 256 x 256 file with 20 spectral
points and two spectral scans, on a Sparc-20/50. This time is for simple
processing without scattered light or unblur. They both add lots of additional
time.
Restrictions
Image is assumed to be square.
Assumes 6.0 mA per unit of F-P setting. This value has been determined
from the separation of the 630.15 and 630.25 nm lines. Note that the steps
in wavelength do not have to be uniform, but there is some code in xbgram6
(95Apr) that chooses which points in the I-profile to fit and that assumes
uniform wavelength steps.