1.4. Calibration overview#

An image of the sky contains counts from several sources. The task of data reduction (another name for image calibration) is to remove all non-celestial counts from the image and to correct for non-uniform sensitivity.

At the end of the previous notebook we arrived at an expression for the counts in a science image in terms of the sources of counts:

\[ \text{raw image} = \text{bias} + \text{noise} + \text{dark current} + \text{flat} \times (\text{sky} + \text{stars}). \]

Solving for the counts just from the stars is as follows:

\[ \text{stars} + \text{noise} = \frac{\text{raw image} - \text{bias} - \text{dark current}}{\text{flat}} - \text{sky} \]

It is impossible to remove the noise from the raw image because the noise is random.

The dark current is typically calculated from a dark frame (aka dark image). Such an image has bias and read noise in it as well, so:

\[ \text{dark current} + \text{noise} = (\text{dark frame} - \text{bias})/(\text{dark exposure time}) \]

Once again, note that the noise cannot be removed.

1.4.1. This noise cannot be removed from CCD images#

To demonstrate that you cannot remove the noise from an image, let’s construct an image with just stars and noise and try to subtract a noise image created with the same parameters. The amount of noise here is exaggerated to make it clear in the images.

import numpy as np
%matplotlib inline
from matplotlib import pyplot as plt

from astropy.visualization import hist
from astropy.stats import histogram

import image_sim as imsim
from convenience_functions import show_image

---------------------------------------------------------------------------
ImportError                               Traceback (most recent call last)
Cell In[1], line 8
      5 from astropy.visualization import hist
      6 from astropy.stats import histogram
----> 8 import image_sim as imsim
      9 from convenience_functions import show_image

File ~/work/ccd-reduction-and-photometry-guide/ccd-reduction-and-photometry-guide/notebooks/image_sim.py:6
      3 import numpy as np
      5 from astropy.modeling.models import Gaussian2D, RickerWavelet2D, Const2D
----> 6 from photutils.datasets import (make_random_gaussians_table,
      7                                 make_gaussian_sources_image)
      8 from photutils.aperture import EllipticalAperture
     10 # To use a seed, set it in the environment. Useful for minimizing changes when
     11 # publishing the book.

ImportError: cannot import name 'make_random_gaussians_table' from 'photutils.datasets' (/usr/share/miniconda/envs/test/lib/python3.11/site-packages/photutils/datasets/__init__.py)

# Use custom style for larger fonts and figures
plt.style.use('guide.mplstyle')

1.4.1.1. First, some stars with noise#

The image below shows stars (the larger “blobs” in the image) but shows quite a bit of noise as well (the much smaller “dots”).

image = np.zeros([2000, 2000])
gain = 1.0
noise_amount = 1500 

stars_with_noise = imsim.stars(image, 50, max_counts=2000, fwhm=10) + imsim.read_noise(image, noise_amount, gain=gain)

show_image(stars_with_noise, cmap='gray', percu=99.9)
plt.title('Stars with noise')

1.4.1.2. Now an incorrect attempt at reducing noise#

Notice that the call to the noise function has exactly the same arguments as above, in much the same way your camera’s electronics will have the same noise properties every time you read out an image.

However, the amount of noise has increased, not decreased. It’s much harder to pick out the stars in this image.

incorrect_attempt_to_remove_noise = stars_with_noise - imsim.read_noise(image, noise_amount, gain=gain)

show_image(incorrect_attempt_to_remove_noise, cmap='gray', percu=99.9)

1.4.2. Every image has noise#

Every image, including calibration images like bias and dark frames, has noise. If we tried to calibrate images by taking a single bias image and a single dark image, the final result might well look worse than before the image is reduced.

For demonstration, we’ll see what happens below.

Note that here we construct realistic bias and dark, but leave read noise out of the flat; we’ll return to that point later.

1.4.2.1. First, set parameters for the CCD#

These are the same as in the previous notebook, except for the read noise, which is 700\(e-\), 100 times larger than in the previous notebook.

gain = 1.0
star_exposure = 30.0
dark_exposure = 60.0
dark = 0.1
sky_counts = 20
bias_level = 1100
read_noise_electrons = 700
max_star_counts = 2000

1.4.2.2. Generate the images, with noise#

bias_with_noise = (imsim.bias(image, bias_level, realistic=True) + 
                   imsim.read_noise(image, read_noise_electrons, gain=gain))

dark_frame_with_noise = (imsim.bias(image, bias_level, realistic=True) + 
                         imsim.dark_current(image, dark, dark_exposure, gain=gain, hot_pixels=True) +
                         imsim.read_noise(image, read_noise_electrons, gain=gain))

flat = imsim.sensitivity_variations(image)

realistic_stars = (imsim.stars(image, 50, max_counts=max_star_counts) +
                   imsim.dark_current(image, dark, star_exposure, gain=gain, hot_pixels=True) +
                   imsim.bias(image, bias_level, realistic=True) +
                   imsim.read_noise(image, read_noise_electrons, gain=gain)
                  )

1.4.2.3. Uncalibrated image#

Below we display the uncalibrated image; in a moment we’ll compare it to the calibrated version. Even though they don’t stand out there really are stars in it.

plt.figure(figsize=(12, 12))
show_image(realistic_stars, cmap='gray', percu=99.9, figsize=(9, 9))

1.4.2.4. Reduce (calibrate) the star image#

First we calculate the dark current, scaled to the exposure time of our light image.

scaled_dark_current = star_exposure * (dark_frame_with_noise - bias_with_noise) / dark_exposure

Next, we subtract the bias and dark current from the star image and then apply the flat correction.

calibrated_stars = (realistic_stars - bias_with_noise - scaled_dark_current) / flat

show_image(calibrated_stars, cmap='gray', percu=99.9)

1.4.2.5. Reducing the image cleans up the image a bit#

The stars stand more clearly than in the unreduced image.

This image does not look much better than the uncalibrated image, but remember that the read noise used in this simulated image, 700 \(e^-\) per pixel, is unrealistically high.

1.4.2.6. Reducing the image increases the noise in the image#

The histogram below shows pixel values before and after calibration. The width of the distribution is a measure of the read noise. As expected, reducing the image increases the read noise. One reason one takes several calibration images of each type is to reduce the amount of noise in the calibration image. That will, in turn, keep the noise in the final image as small as possible.

plt.figure(figsize=(9, 9))
hist(calibrated_stars.flatten(), bins='freedman', label='calibrated star image', alpha=0.5)
hist(stars_with_noise.flatten(), bins='freedman', label='raw star image', alpha=0.5)
plt.legend()
plt.grid()
plt.xlabel('Count level in image')
plt.ylabel('Number of pixels with that count');