Information on how to build and run a solid-waffle script

A script can be run as follows:

from solid_waffle.correlation_run import run_ir_all, run_vis_all
run_ir_all("my_ir_config.txt") # this runs an IR script
run_vis_all("my_vis_config.txt") # this runs a visible script

The script file contains instructions on which darks and flats are to be analyzed, what format they are in, and what analysis options are desired.

A typical script might look something like the following:

DETECTOR: H4RG-18237
LIGHT:
  ../SCA18237/Set_002_Light_0001.fits
  ../SCA18237/Set_004_Light_0001.fits
  ../SCA18237/Set_006_Light_0001.fits
DARK:
  ../SCA18237/Set_001_Dark_0001.fits
  ../SCA18237/Set_001_Dark_0002.fits
  ../SCA18237/Set_001_Dark_0003.fits
FORMAT: 1

CHAR: Advanced 1 3 3 bfe

#NBIN: 128 16

TIME:    1 10 12 20
#TIME2A:  1  2  3 20
#TIME2B:  1  2  3 20
#TIME3:   1  2  3 20
OUTPUT:         out/adv2multi

Most commands in the script start with a command name (capital letters, followed by a colon, e.g., FORMAT:). Following the command is a list of arguments, which must be on the same line as the command except as specified below. Commands with no arguments do not include the colon. Blank lines and lines starting with a hash (#) are skipped.

Some script commands are required; others are optional, and come with defaults. There are additional commands required for a visible light run.

Required commands

LIGHT

This command contains a list of the flat field exposures. Input files in *.fits format are listed on subsequent lines as follows:

LIGHT:
  ../SCA18237/Set_002_Light_0001.fits
  ../SCA18237/Set_004_Light_0001.fits
  ../SCA18237/Set_006_Light_0001.fits

The list is terminated by the next command or by the end of the file.

DARK

This works the same way as the LIGHT command, but contains dark exposures. There must be the same number of dark exposures as flat field exposures, since internally in the script each flat is paired with a dark. The darks should be at least as long as the flats (in terms of number of frames), and need to be in the same format.

OUTPUT

This command tells the script where to save output files:

OUTPUT:         out/adv2multi

This is a stem, e.g. in this case the "summary.txt" output will be saved to out/adv2multi_summary.txt, etc.

Optional commands

FORMAT

The format of the input files:

FORMAT: 1

Note that aside from the array dimensions, the script does not use any information from the FITS headers, it only uses the image data.

The options fall into a few categories:

DCL and PPL test data (for use with ground test data from Roman and Euclid)
- 1 (default if not specified): H4RG-10, FITS cube, ascending ramps (signal level in DN increases during a flat). The data should be in HDU #0, with NAXIS = 3, NAXIS1 = NAXIS2 = 4096, and NAXIS3 equal to the number of frames in the flat.
- 2: same as 1, but for an H2RG (NAXIS1 = NAXIS2 = 2048).
- 3: H4RG-10, descending ramps (signal level in DN decreases during a flat). The images are in HDU #1 ... HDU #[last time slice], with NAXIS =2 and NAXIS1 = NAXIS2 = 4096.
- 4: H4RG-10, descending ramps (signal level in DN decreases during a flat). The images are in HDU #1, with NAXIS = 4, NAXIS1 = NAXIS2 = 4096, NAXIS3 = [number of time slices], and NAXIS4 = 1. (Used for some of the QE data.)
- 5: H4RG-10, descending ramps. The images are in HDU #1, with NAXIS1 = NAXIS2 = 4096, NAXIS3 = [number of time slices]. (Used by PPL.)
- 6: H4RG-10, ascending ramps. The images are in HDU #1, with NAXIS = 4, NAXIS1 = 4096 or 4224 (latter includes reference output as 33rd channel), NAXIS2 = 4096, NAXIS3 = [number of time slices], and NAXIS4 = 1. (Used for triplet testing.)
- 7: H2RG, similar to 3, but NAXIS1 = NAXIS2 = 2048. (Used for some of the Euclid detector data.)
Unit tests (for simulations of small-format arrays to enable fast tests)
- 1001: 512x512 array with 4 reference pixels, FITS cube, ascending ramps (signal level in DN increases during a flat). The data should be in HDU #0, with NAXIS = 3, NAXIS1 = NAXIS2 = 4096, and NAXIS3 equal to the number of frames in the flat.
Roman TVAC and flight-like data
- ( not yet enabled )

NBIN

The way the data are binned into super-pixels for statistical analysis:

NBIN: 128 16

This produces 128 super-pixels on the X axis, and 16 super-pixels on the Y axis. For an H4RG, this means that the super-pixels are 4096/128=32 pixels wide and 4096/16=256 pixels high. (Note that the fast read direction and the channel number direction are in X, and the slow read direction is in Y.)

The default is 32 super-pixels in each direction.

CHAR

The type of characterization to perform. Right now there are two options:

Basic, e.g., CHAR: Basic.

This is the default, and extracts basic information (gain, non-linearity, IPC) based on the statistics (medians, variances, and covariances) of two CDS images at different signal levels.
Advanced, e.g., CHAR: Advanced 1 3 3 bfe.

This takes more time, but fits gain, non-linearity, and IPC to a stack of CDS images extracted from different epochs in the ramp. The arguments of 'Advanced' are tchar1, tchar2, ncycle, and chartype. Here the CDS images have length tchar1 and tchar2 frames, and ncycle iterations of gain/non-linearity/IPC characterization and IPNL characterization are performed.

The type of IPNL characterization can be: 'none' or 'bfe'. If 'none' is selected, then IPNL is ignored during fits to gain, non-linearity, and IPC. If 'bfe' is selected, then the IPNL obtained from the non-overlapping correlation function is assumed to be due entirely to the BFE, and the resulting correction is applied to the statistical properties of the CDS images.

FULLNL

Controls the use of fully non-linear predictions for the flat correlation functions instead of the lowest order contributions in Hirata & Choi (2019). For example, in:

FULLNL: True True False

The 1st field will control use of the fully non-linear correlation function for the gain/IPC determination in advanced characterization (if on).
The 2nd field will control use of the fully non-linear correlation function for the IPNL kernel.
The 3rd field will control use of the full non-linearity curve instead of beta.

The default is all False.

TIME

Determines which time slices to use for initial characterization and IPNL correlations. For example, in

TIME:    1 10 12 20

there are four arguments, and they need to be in strictly ascending order. They represent frame numbers, with 1 being the first frame in the file. In the above example, the non-overlapping correlation function used for IPNL is the correlation of S_{1,10} * S_{12,20} (where recall S_{a,b} is the difference between frames a and b).

For the basic characterization (and the initial condition for advanced characterization), the non-linearity parameter would be based on frames 1, 10, and 20; the gain is determined from the variance and medians of S_{1,20} and S_{1,10}; and the IPC is determined from the autocorrelation of S_{1,20}.

TIME2A

Determines which time slices to use for the raw gain vs. interval duration test. If not present or commented out, the test is skipped. For example:

TIME2A:  1  2  3 20

There are four arguments; in this case, raw gain is computed for the frame triplets [1,2,3], [1,2,4], [1,2,5] ... [1,2,20]. Note that the first two frames are kept fixed, while the third frame is varied.

TIME2B

Determines which time slices to use for the raw gain vs. interval center test. If not present or commented out, the test is skipped. For example:

TIME2B:  1  2  3 20

There are four arguments; in this case, raw gain is computed for the frame triplets: [1,2,3], [2,3,4], [3,4,5] ... [18,19,20]. Note that the time differences used to compute raw gain are kept fixed, while the absolute time stamp is varied. If not present or commented out, the test is skipped.

TIME3

Determines which time slices to use for the IPC correlation vs. interval duration test. If not present or commented out, the test is skipped. For example:

TIME3:  1  2  3 20

There are four arguments; in this case, IPC is computed using basic characterization for the frame triplets: [1,2,3], [1,2,4], [1,2,5] ... [1,2,20]. Note that the first two frames are kept fixed, while the third frame is varied. The CDS autocorrelations of S_{1,3}, S_{1,4}, ... S_{1,20} are also computed.

TIMEREF

Determines which time slice is the "reference" (i.e., conversion gain is d[charge]/d[signal] at that signal level). The default is 1. If you want conversion gain reported at the reset level, and the 1st frame in the data cube is the frame after the reset, then include:

TIMEREF: 0

REF OFF

Turns off the left and right reference pixel subtraction. You may simply include the line REF OFF in the configuration.

NLPOLY

Allows a higher-order non-linear polynomial correction scheme for the classical non-linearity. We don't use this for any corrections yet, we just measure the curve and have the piping to put the results in the output. For example, NLPOLY: 4 2 32 allows a polynomial of 4th order, fit from time frames 2 .. 32 (inclusive).

DETECTOR

Name of the SCA (not used, except copied into the output data). The default is the null string ''. You may name your detector by including, e.g.:

DETECTOR: H4RG-18237

COLOR

Color scheme for output plots; for example:

COLOR: gnuplot

MASK

If provided, manually masks a super-pixel, e.g. MASK: 17 31 removes the (17,31) super-pixel. This is in X,Y order, with the lower-left being (0,0). The user may mask as many super-pixels as desired; if a super-pixel is repeated, this is equivalent to masking it once.

QUANTILE

For example, QUANTILE: 75. If provided, sets the quantile level used for estimating the variance in gain determination (the default is 75, which corresponds to using the inter-quartile range).

EPSILON

For example, EPSILON: .025. If provided, sets the clipping fraction for the IPC correlations (default: 0.01, meaning clip the top 1% and bottom 1% of the pixels before computing the covariance; >=0.50 will result in the IPC measurement failing).

IPCSUB

This feature can be turned on by setting IPCSUB: True. (The default is false.)

If True, then the correlation function used for IPC determination will subtract a 'baseline' region +/- 4 pixels to the left and right (fast-scan direction). This increases the noise but may be useful as a systematic check if you think your autocorrelation IPC measurement is affected by 1/f noise.

HOTPIX

For example, HOTPIX: 1500 3000 0.1 0.08. If provided, does a hot pixel-based IPC analysis. This uses only the dark frames.

Selected pixels should be:

in the range from 1500-3000 DN (in the example above)
'stable' (here meaning top to bottom of 0.1 or 10% of the median signal) and
'isolated' (here meaning no other pixels in a 5x5 box exceeding 0.08 or 8% of the selected pixel)
not within 2 pixels of a boundary or reference pixel

HOTPIX SLIDEMED

You may turn on this feature by including HOTPIX SLIDEMED.

This forces the median in the hot pixel IPC computation to be done using a "sliding median" instead of the traditional median. Recall that IPC is a ratio of the signal in nearest neighbors to total signal, i.e., we want y/x. The traditional median is simply the median of y/x. However in cases where x can fluctuate to negative values due to noise, one prefers the sliding median, which is the value of m where 50% of the y-mx values are positive and 50% are negative.

HOTPIX LOGTSPACE

You may turn on this feature by including HOTPIX LOGTSPACE.

In the sample number dependence of hot pixel IPC, uses binary logarithmic spacing of time stamps (i.e., 1, 3, 5, or 7 times a power of 2) instead of linear scalings (default; every 5th sample).

HOTREF AUTOCORR

You may turn on this feature by including HOTREF AUTOCORR.

In computing the sample number dependence of hot pixel IPC, uses the autocorrelation estimate of IPC as the 'reference' (for alpha_{this sample} - alpha_{reference}) instead of the last sample (default, in which case the last value of Delta alpha is trivially zero).

NARROWFIG

You may make narrow format output figures by including NARROWFIG (only some panels are displayed).

Outputs

The script outputs a set of summary files and figures.

`*_summary.txt`

This file contains the calibration data obtained on that SCA. It starts with a time stamp and version number:

# This summary created at Mon Jul  2 00:15:04 2018
# Uses pyirc v5

followed by a list of settings. This includes a list of the figures associated with the summary file. Then comes a table of the results in each super pixel; the number of columns may vary, but they always come with a description:

# Columns:
#   1, X (super pixel grid)
#   2, Y (super pixel grid)
...
#  37, BFE kernel K^2a (+NL-IPC) at ( 2, 2) (e^-1)

Following this is the actual list of data, with one line per super-pixel:

  0   0  1.5050000E+04  2.8682524E+00  ... -5.0412114E-07
...
 31  31  1.4963000E+04  2.8977701E+00  ...  7.2685133E-07

Failures or masked super-pixels have rows of zeros.

`*_multi.eps`

This is a 6-panel figure displaying the results of the characterization. In English reading order, the panels display:

The percentage of pixels that are good (max = 100%; note that reference pixels are rejected so you will see a lower percentage around the edges of the SCA).
The gain in e/DN, as corrected in the characterization chosen (this is always IPC and non-linearity corrected; the advanced characterization will be BFE-corrected if the 'bfe' option is chosen and the number of iterations is >=1).
The IPC in percent (nearest-neighbor, averaged over up, down, right, and left).
The quadratic non-linearity parameter beta (in ppm/e).
The charge per frame (in e).
The IPNL coefficient at zero lag (in ppm/e). See the paper for the exact definition, but roughly speaking this is the sum of the self-BFE coefficient (i.e. the fractional change in pixel area when 1 e of charge is deposited in it) and - 4 d[IPC alpha]/d[charge]. It may be interpreted as either BFE, NL-IPC, or both, depending on the outcome of other tests. It is usually negative, and with small numbers of flats may be very noisy.

`*_m23.eps`

This is a 4-panel figure displaying how raw gains, autocorrelations of synthetic CDS images (at 1-pixel lag), and basic characterization-inferred IPC coefficients vary with the intervals used. These panels all have slopes that depend on whether the IPNL is due to BFE or NL-IPC (or some combination). They can thus be used to distinguish the various effects.

The figures depend on the IPNL characterization methods used (2a, 2b, and 3), which can be controlled using the TIME2A, TIME2B, and TIME3 commands. If these commands are absent or invalid, then the corresponding figures are skipped, and a plot with fewer panels will be generated.

In each panel, the data points are the average over all super-pixels, with error bars computed by 'error on the mean' using the dispersion of the super-pixels. The predicted slopes of the relations are plotted, assuming that the IPNL is pure BFE (green) or pure NL-IPC (blue).

`*_hot.txt`

A table of output information on the hot pixels. This is a flat text file, with each row denoting a hot pixel. The columns are:

pixel x
pixel y
IPC alpha from auto-correlation in that pixel's super-pixel
S_{1,6} in that pixel
S_{1,6}, average in the 4 nearest neighbors
S_{1,6}, average in the 16-pixel "ring" region (5x5 minus 3x3 centered on that pixel)
S_{1,6}, asymmetry (left-right minus top-bottom, divided by 2)
S_{1,6}, average in the 4 diagonal neighbors

Then Columns 9-13 are for S_{1,11}, Columns 14-18 are for S_{1,16} ... all the way up to the last S_{1,1+5n} in the data cube.

`*_hotipc.eps`

A graphical display of the information on the hot pixels. It is a 4-panel plot. In English reading order, the panels display:

A map of the selected hot pixel locations on the SCA.
A plot of the alpha inferred from each hot pixel as a function of the signal level (i.e., how hot the pixel is).
A plot of Delta alpha: the IPC inferred from a partial ramp (so the hot pixel does not reach its final deposited charge) minus that inferred from the full ramp. For perfectly linear response (including linear IPC), this will be zero, with some scatter due to noise. The red points denote each sample; the black points are median stacks of the red points with binomial error bars. The last black point is zero by construction.
A scatter plot of the hot pixel and autocorrelation-based IPC measurements, binned into 16 'hyper-pixels', each N/4 x N/4 physical pixels (e.g. 1024x1024 for an H4RG). The red line indicates perfect agreement.

Commands for visible light runs

Some commands are used only with visible light runs (they have no effect on infrared runs):

VISLIGHT

This is the list of flats for the visible data, e.g.:

VISLIGHT:
../SCA20829/20191018_95K_1p1m0p1_ch21_500nm_gr3_filt6_shutter_open_20829_001.fits
../SCA20829/20191018_95K_1p1m0p1_ch21_500nm_gr3_filt6_shutter_open_20829_002.fits
../SCA20829/20191018_95K_1p1m0p1_ch21_500nm_gr3_filt6_shutter_open_20829_003.fits

VISDARK

This is the list of corresponding darks for the visible data, e.g.:

VISDARK:
../SCA20829/20191018_95K_1p1m0p1_ch0_1400nm_gr3_filt5_shutter_closed_20829_001.fits
../SCA20829/20191018_95K_1p1m0p1_ch0_1400nm_gr3_filt5_shutter_closed_20829_002.fits
../SCA20829/20191018_95K_1p1m0p1_ch0_1400nm_gr3_filt5_shutter_closed_20829_003.fits

The number of visible darks must equal the number of visible flats. It is allowed for IR darks to be repeated as visible darks (a dark is a dark, and has no illumination wavelength).

VISTIME

The timestamps to use in the visible data. There are 4 integer arguments, e.g.:

# visible times: start, end, dt1, dt2
VISTIME: 1 20 1 3

Here:

start = beginning of part of ramp to use (1st frame = 1)
end = ending of visible ramp to use
dt1, dt2 = time differences for visible autocorrelations (dt2>dt1; use dt2 for charge diffusion correlations, dt1 is for variance subtraction in determining the apparent gain and hence omega).

VISBFETIME

The timestamps to use for BFE determination in the visible data. (The code estimates BFE separately for visible and IR.) BFE measurement uses the non-overlapping correlation function C_{abcd}, where the timestamps a, b, c, and d are in this keyword:

VISBFETIME: 1 9 10 18

VISMEDCORR

If you add this keyword to the configuration file:

VISMEDCORR

then in the visible characterization, the variance C_{adad}(0,0)-C_{abab}(0,0) is re-scaled by the predicted / measured value of M_{bd} so that the slope of the photon transfer curve is corrected for errors in the current Ie. We recommend this option only as a test case.

COPYIRBFE