Correcting for instrument degradation

This example demonstrates the degradation of the filters on AIA and how to correct it.

import astropy.time
import astropy.units as u
import matplotlib.pyplot as plt
from astropy.visualization import time_support
from sunpy.net import Fido
from sunpy.net import attrs as a

from aiapy.calibrate import degradation

# This lets you pass `astropy.time.Time` objects directly to matplotlib
time_support(format="jyear")

<astropy.visualization.time.time_support.<locals>.MplTimeConverter object at 0x7f32aa466b90>

First, let’s fetch the metadata for the 335 Å channel of AIA between 2021 and 2023 at a cadence of 7 days. We choose the 335 Å channel because it has experienced significant degradation compared to the other EUV channels.

results = Fido.search(
    a.Time("2021-01-01T00:00:00", "2023-01-01T00:00:00"),
    a.Sample(7 * u.day),
    a.jsoc.Series.aia_lev1_euv_12s,
    a.jsoc.Wavelength(335 * u.angstrom),
)

We only need the date and mean intensity columns from the metadata that was returned. We select those and nothing else.

table = results["jsoc"].show("DATE__OBS", "DATAMEAN")
table["DATAMEAN"].unit = u.DN
table["DATE_OBS"] = astropy.time.Time(table["DATE__OBS"], scale="utc")
del table["DATE__OBS"]

print(table)

DATAMEAN         DATE_OBS
   DN
-------- -----------------------
  0.8555 2021-01-01T00:00:00.620
  0.6709 2021-01-08T00:00:00.620
  0.6219 2021-01-15T00:00:00.630
  0.7758 2021-01-22T00:00:00.630
  0.7389 2021-01-29T00:00:00.630
  0.6407 2021-02-05T00:00:00.630
  0.6212 2021-02-12T00:00:00.630
  0.9368 2021-02-19T00:00:00.620
  0.9861 2021-02-26T00:00:00.630
     ...                     ...
  2.5536 2022-10-28T00:00:00.620
   2.682 2022-11-04T00:00:00.630
  2.5522 2022-11-11T00:00:00.630
  2.4727 2022-11-18T00:00:00.630
  2.1289 2022-11-25T00:00:00.630
  2.2118 2022-12-02T00:00:00.620
  2.6068 2022-12-09T00:00:00.620
  3.2372 2022-12-16T00:00:00.630
  2.5094 2022-12-23T00:00:00.630
  3.2047 2022-12-30T00:00:00.630
Length = 105 rows

Next, we pass the date column to the aiapy.calibrate.correct_degradation function. This function calculates the time-dependent correction factor based on the time and wavelength of the observation. We then divide the mean intensity by the correction factor to get the corrected intensity. For more details on how the correction factor is calculated, see the documentation for the aiapy.calibrate.degradation function.

correction_factor = degradation(335 * u.angstrom, table["DATE_OBS"])
table["DATAMEAN_DEG"] = table["DATAMEAN"] / correction_factor

To understand the effect of the degradation and the correction factor, we plot the corrected and uncorrected mean intensity as a function of time. Note that the uncorrected intensity decreases monotonically over time while the corrected intensity recovers to pre-2011 values in 2020.

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(table["DATE_OBS"], table["DATAMEAN"], label="mean", marker="o")
ax.plot(table["DATE_OBS"], table["DATAMEAN_DEG"], label="corrected mean", marker="+")
ax.set_title(f'{(335*u.angstrom).to_string(format="latex")} Channel Degradation')
ax.legend(frameon=False)

plt.show()