Basejoint

Procedure to make a monthyl bias for the STIS CCD.

Python translation of a python translation of an original IRAF cl script to create superbias reference file from a list of bias frames.

The input image (with multiple extensions) is overscan-subtracted and cosmic-ray-rejected using the CALSTIS algorithms within STSDAS. The cosmic-ray-rejected image is divided by the number of imsets present in the input image (since ocrreject adds up the individual imsets). After that, the superbias is median filtered using a window of 15 x 1 pixels. The median-filtered is subtracted from the superbias to produce a “residual” image containing hot columns and such. This “residual” image is averaged along rows and replicated back to the original image size, so that hot columns clearly show up. After that, the image values in hot columns and pixels of the original superbias image (defined as those pixels having values greater than (mean + 5 sigma of Poisson noise) are replaced by those in the median-filtered bias image. Plots are made of the row- and column-averaged superbias, with plotting scales appropriate to the gain and binning settings of the superbias.

refstis.basejoint.average_biases(bias_list)[source]

Create a weighted sum of the individual input files.

First make sure all individual input files have been ocrrejected.

Parameters

bias_list (list) – list of the input biases

Returns

  • mean_file (str) – name of the averaged filename

  • totalweight (float) – sum of the NCOMBINE header keywords from the data

refstis.basejoint.calibrate(input_file)[source]

calibrate input file

refstis.basejoint.make_basebias(input_list, refbias_name='basebias.fits')[source]

Make the basebias for an anneal month

1- Calbrate each bias in the list 2- Average together the biases 3- Replace pixels and colums with median values 4- Set header keywords

Parameters
  • input_list (list) – list of input bias files.

  • refbias_name (str) – name of the output reference file.

refstis.basejoint.replace_hot_cols(mean_bias, median_image, residual_image, yfrac=1)[source]

Replace hot columns in the mean_bias as identified from the residual image with values from the bias_median

‘hot’ is 3* sigma

mean_bias will be updated in place

refstis.basejoint.replace_hot_pix(mean_bias, median_image)[source]

Replace image values in residual single hot pixels

defined as those having values greater than (mean + 5 sigma of Poisson noise) by those in median-filtered bias image. This represents the science extension of the final output reference superbias.

mean_bias will be updated in place.

Parameters
  • mean_bias (str) – name of the mean bias bias

  • median_image (np.ndarray) – 2d median image of the bias