what’s your bias? (bias correction)

/ November 20, 2020/ prepare to astronomize, research log

In this series, I’ll walk through each of the steps necessary to go from taking data at an (optical) telescope to doing science with that data. This article touches on why astronomical images need to be corrected before use, and goes into detail about the bias correction process. If you’d like to follow along and create your own image of the NGC 663 cluster you can clone the github repository. I’ll be learning something too, because I haven’t studied the NGC 663 cluster in detail (we’ll get to discuss it when we begin to analyze our data later). The data I use in this explainer series was taken from the Amherst College Observatory in October 2019, from a 11in Cassegrain telescope and a SBIG CCD, and the sorted raw data is available to download from this website as a tarball. Much of the code I present here is modeled after code handed down to me by Sarah Betti and Kim Ward-Duong, who taught me most of what I know (thanks guys!).

Some basics

Astronomical images, when they first descend from the heavens and grace the camera at the bottom of an astronomer’s telescope, don’t look pretty. They pass through the atmosphere, which distorts and smudges them, before entering a telescope aperture. The telescope aperture causes the light from the sky to diffract, which causes each star in the image to take on a particular shape, called a Point Spread Function (PSF). The camera introduces its own problems into the image, which have to be corrected. In this series we’ll discuss the ways that astronomers correct for these problems in order to collect images which enable them to collect the information needed to answer specific science questions.

The PSF of a telescope depends on the shape of its aperture and any other optics that the light passes through before it is collected by a camera. The atmosphere smudges the PSFs of smaller telescopes, turning them into boring, near-normal distributions. Large, complicated telescopes which mitigate the atmosphere (either by flying above it, or by using a technique called adaptive optics, which corrects for the atmospheric distortion) have strangely shaped PSFs. Below are three PSFs from three telescopes: WIYN 0.9m, Keck 10m, and the Roman Space Telescope.

WIYN 0.9m PSF (author’s image)

Keck Adaptive Optics PSF (source)

The real @NASARoman coronagraph is more complicated than this but the basic principal is shaping where interference happens. See if you can spot the planet in this version!

Originally tweeted by Bruce Macintosh (@bmac_astro) on November 18, 2020.

Images from optical telescopes are captured by a Charged Coupled Device (CCD) camera, which acts like a grid of “photon buckets” which collect photons from the telescope image plane in potential wells. This is how digital cameras like the one in your smartphone operate. The important thing to know is that this camera, while incredibly powerful for being able to measure astronophysical signal, introduces some sources of error into the resultant images.

For optical images, there are three dominant sources of noise (and therefore error, sometimes called uncertainty) in our observations. There are therefore three correction steps we take before analyzing our images. The CCD introduces “bias” when it reads out the image to a computer; we discuss bias correction in this next post. The CCD also introduces “dark current”: excess photons which come from the heat generated by the CCD itself, and are therefore not part of the original light from the sky; we discuss dark subtraction in the next post. After correcting the CCD’s noise, we attempt to set every pixel in our images on equal footing, correcting for some of the effects of the optics in the telescope and any variations or defects in the sensitivity of the pixels in the CCD; we discuss this “flat fielding” in a future post. After that, we’ll calibrate the brightness of our images using a standard star, align them to a reference coordinate system, and finally analyze the data to answer some questions we have about the cluster NGC 663.

Bias correction

Biasing is an electronics term that describes the application of a voltage/current to a circuit or piece of circuit to create more optimal operating conditions. A CCD used for astronomy applies a bias to its pixels in order to avoid negative values. The CCD wants to avoid negative values because a set of negative and positive pixel values takes twice as much space to store as a set of purely positive values.

In a perfect world, this would result in every pixel in any image having some offset from zero of one value; we could simply subtract this value from our images to get the “true” pixel values for a given image. A bias image (or “frame” as images are occasionally referred to by astronomers) essentially attempts to capture this value. By taking a zero second exposure, letting no light or heat build up on the detector, the given bias frame should have pixel values equal to only the bias voltage applied; you can then subtract this uniform image from a science image in order to remove the bias voltage. In reality, bias voltage isn’t uniformly applied, and a CCD’s amplifier cannot convert voltages to digital units with perfect accuracy. In this way, a bias frame captures more information to help us mitigate these sources of error.

Readout refers to the CCD’s process of recording the voltage induced in a given pixel and converting that voltage into some unit of ‘counts’ (which attempt to measure how many photons a certain pixel has captured). A signal is induced in the pixel by a number of photons from an astronomical source striking the pixel and (via the photoelectric effect) generating some voltage across the pixel. There is some inherent uncertainty in the measurement of said voltage. Depending on how long it takes for a given pixel to be readout (in a camera with thousands of pixels) the voltage might change. These are some sources of read noise, which is mitigated by measuring and subtracting bias.

A bias frame is then a way to capture any spurious readout patterns or affects and remove them from science images, in addition to correcting for the applied bias voltage. Bias frames are taken by reading out the CCD when it is empty, that is, when there has been no exposure to outside light, and the exposure time of every bias frame is zero.

Bias frames are relatively easy to take, and can be taken at any time during the observing run. The only time constraint on taking these calibration frames is that they take as long to create as your CCD takes to read out all the pixels (usually a few seconds, depending on how big your detector is).

A single bias frame.

Once we’ve taken a sizable number of bias frames (usually at least 11), we take the “median combination” of the bias frames to create a “master bias” which we subtract from each of our images. This median combination means that the master bias captures the most common bias signal for each pixel, and is therefore more accurate than just taking one bias frame which may have some variation from the average.

To take a median combination of astronomical images, I’ve written a python function:

mediancombine.py
# import statements
import numpy as np
from astropy.io import fits


# function
def mediancombine(filelist):
    '''
    median combine frames
    '''
    # gather some information about the images
    n = len(filelist) 
    first_frame_data = fits.getdata(filelist[0])
    imsize_y, imsize_x = first_frame_data.shape
    # construct an empty cube with proper dimensions 
    fits_stack = np.zeros((imsize_y, imsize_x, n))
    # fill cube with images in filelist
    for ii in range(0, n):
        im = fits.getdata(filelist[ii])
        fits_stack[:, :, ii] = im
    # take the median combination of the images along the 
    # correct axis
    med_frame = np.median(fits_stack, axis=2)
    return med_frame
An eleven bias median “master bias” calibration frame.

Looking at the average pixel values in this master bias, it appears the bias that is applied to this CCD was about 900 counts. To subtract the master bias from each of our science images, I’ve written the following functions into a python script to automate the task. By running “python bias.py” from the command line, the code is executed.

bias.py
# import statements
import os
import glob
from astropy.io import fits
from mediancombine import mediancombine


# functions
def bias_subtract(filename, path_to_bias, outpath):
    '''
    bias subtract frames
    '''
    # gather image, header, and bias image
    targetdata = fits.getdata(filename)
    target_header = fits.getheader(filename)
    biasdata = fits.getdata(path_to_bias)

    # perform subtraction
    b_data = targetdata - biasdata

    # check if files exist in directory other than notebook
    fitsname = filename.split('\\')[-1]
    # save image to disk
    fits.writeto(outpath + "/" + 'b_' + fitsname, b_data,
                 target_header, overwrite=True)
    return


def run_master_bias(fil):
    '''
    create master bias
    '''
    # check that a folder for the master calibration frames
    # exists or not, if not, create it
    if not os.path.exists(f'{fil}/Masters'):
        os.makedirs(f'{fil}/Masters')
        print(f'Masters folder created at: {fil}/Masters')

    # store path for final image as string
    masterbiaspath = fil + '/Masters/MasterBias.fit'
    # if the master bias has not been created, create it
    # and write to disk
    if not os.path.exists(f'{fil}/Masters/MasterBias.fit'):
        print('Making Master Bias')
        # create master bias
        bias_fits = glob.glob(fil + '/bias/*.fit')

        median_bias = mediancombine(bias_fits)
        fits.writeto(masterbiaspath, median_bias,
                     header=fits.getheader(bias_fits[0]), overwrite=True)
    # if the master bias already exists, point user to it
    else:
        print(f'Master Bias in: {masterbiaspath}')
    return


def bias_sub_images(fil, targs):
    '''
    bias correct science targets
    '''
    # for each target we have, run the reduction
    for target in targs:
        print()
        print('------------o------------')
        print('target: ', target)
        print()
        # store master bias path as string
        masterbiaspath = fil + '/Masters/MasterBias.fit'

        # get lists of target images
        bias_images = glob.glob(f'{fil}\\{target}\\b_*.fit')
        scidata = glob.glob(fil + '\\' + target + '\\*.fit')
        # if images have not been bias subtracted, subtract them
        if len(bias_images) != len(scidata):
            [os.remove(im) for im in bias_images]
            print('Bias subtracting ')
            for sci_image in scidata:
                # read in data
                filtername = fits.getheader(sci_image)['FILTER'][0]
                # store filepath as string
                sci_outpath = fil + '\\' + target + '\\' + filtername + 'band'
                # if the folder does not exist already, make it
                if not os.path.exists(sci_outpath):
                    os.makedirs(sci_outpath)
                # run bias subtraction function
                bias_subtract(sci_image, masterbiaspath, sci_outpath)

    return('Bias subtracted target images')

if __name__ == '__main__':
    # Running as a script
    root = os.getcwd()
    fil = root+'\\Data\\20191025'
    targs = ['NGC 663', 'Entered Coordinates']
    run_master_bias(fil)
    bias_sub_images(fil, targs)

From here, our images are ready to be calibrated further.

The next post details “dark current” and the dark subtraction process. Until then, you can star the github repo or follow this blog. Clear skies, and happy astronomizing!

Share this Post

About willb

I'm Will, an undergraduate astronomer studying transition disks, direct imaging, and planet accretion and formation at the Follette Lab at Amherst College. I use they/them/theirs pronouns.

Leave a Comment

Your email address will not be published. Required fields are marked *

You may use these HTML tags and attributes: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>
*
*