Use DN instead of counts everywhere (#338)

In line with sunpy/sunpy#7585, this replaces the usage of `ct/count`
with `DN` everywhere. This may cause some friction in terms of
compatibility with sunpy as prior to v6, the units of any AIA map are
converted to `ct` rather than `DN`. I'm not sure the best way to resolve
that. We could of course just pin to >= v6, but that seems aggressive.

---------

Co-authored-by: Nabil Freij <nabil.freij@gmail.com>

aiapy/calibrate/tests/test_uncertainty.py

@@ -17,7 +17,7 @@
 
 @pytest.mark.parametrize("channel", CHANNELS)
 def test_error_all_channels(channel):
-    intensity = 10.0 * u.ct / u.pix
+    intensity = 10.0 * u.DN / u.pix
     error = estimate_error(intensity, channel, error_table=table_local)
     assert error.unit == intensity.unit
 
@@ -33,7 +33,7 @@ def test_error_all_channels(channel):
     ],
 )
 def test_counts_shapes(counts):
-    counts = np.abs(counts) * 1000 * u.ct / u.pix
+    counts = np.abs(counts) * 1000 * u.DN / u.pix
     errors = estimate_error(counts, 171 * u.angstrom, error_table=table_local)
     if counts.shape == ():
         assert errors.shape == (1,)
@@ -59,7 +59,7 @@ def test_counts_shapes(counts):
 def test_flags(include_preflight, include_eve, include_chianti, expectation):
     with expectation:
         errors = estimate_error(
-            1 * u.ct / u.pix,
+            1 * u.DN / u.pix,
             94 * u.angstrom,
             error_table=table_local,
             include_eve=include_eve,

aiapy/calibrate/uncertainty.py

@@ -16,7 +16,7 @@
 @u.quantity_input
 @validate_channel("channel")
 def estimate_error(
-    counts: u.ct / u.pix,
+    counts: u.DN / u.pix,
     channel: u.angstrom,
     *,
     n_sample=1,
@@ -25,7 +25,7 @@ def estimate_error(
     include_chianti=False,
     error_table=None,
     **kwargs,
-) -> u.ct / u.pix:
+) -> u.DN / u.pix:
     """
     Given an observed number of counts estimate the associated errors.
 
@@ -81,24 +81,24 @@ def estimate_error(
     # Dark noise
     # NOTE: The dark error of 0.18 is from an analysis of long-term trends in the residual
     # dark error from 2015.
-    dark = 0.18 * u.ct / u.pix
+    dark = 0.18 * u.DN / u.pix
 
     # Read noise
     if kwargs.get("compare_idl", False):
         # The IDL version hardcodes the read noise as 1.15 DN / pixel so we
         # want to use this when comparing against that code and not at any
         # other time. This is why this option is not documented.
-        read_noise = 1.15 * u.ct / u.pix
+        read_noise = 1.15 * u.DN / u.pix
     else:
         # NOTE: This lookup table comes from Table 6 of Boerner et al. (2012)
         # (http://adsabs.harvard.edu/abs/2012SoPh..275...41B). Each
         # channel is attached to one of the four cameras such that the
         # read noise for a given channel depends on what camera it is on.
         read_noise = {
-            1: 1.18 * u.ct / u.pix,
-            2: 1.20 * u.ct / u.pix,
-            3: 1.15 * u.ct / u.pix,
-            4: 1.14 * u.ct / u.pix,
+            1: 1.18 * u.DN / u.pix,
+            2: 1.20 * u.DN / u.pix,
+            3: 1.15 * u.DN / u.pix,
+            4: 1.14 * u.DN / u.pix,
         }[telescope_number(channel)]
     read = read_noise / np.sqrt(n_sample)
 
@@ -108,7 +108,7 @@ def estimate_error(
     # with the signal and has an approximately uniform distribution. The RMS value is thus the
     # standard deviation of this distribution, approximately 1/sqrt(12) LSB.
     # See https://en.wikipedia.org/wiki/Quantization_(signal_processing)
-    quant_rms = (1 / np.sqrt(12)) * u.ct / u.pix
+    quant_rms = (1 / np.sqrt(12)) * u.DN / u.pix
     quant = quant_rms / np.sqrt(n_sample)
 
     # Onboard compression

aiapy/calibrate/util.py

@@ -244,7 +244,7 @@ def get_error_table(error_table=None):
         warnings.simplefilter("ignore", category=ErfaWarning)
         table["T_STOP"] = Time(table["T_STOP"], scale="utc")
     table["WAVELNTH"] = u.Quantity(table["WAVELNTH"], "Angstrom")
-    table["DNPERPHT"] = u.Quantity(table["DNPERPHT"], "ct photon-1")
+    table["DNPERPHT"] = u.Quantity(table["DNPERPHT"], "DN photon-1")
     return table
 
 

aiapy/response/channel.py

@@ -337,7 +337,7 @@ def eve_correction(self, obstime, **kwargs) -> u.dimensionless_unscaled:
     @u.quantity_input
     def gain(
         self,
-    ) -> u.count / u.ph:
+    ) -> u.DN / u.ph:
         r"""
         Gain of the CCD camera system.
 
@@ -359,7 +359,7 @@ def gain(
         electron_per_ev = self._data["elecperev"] * _e / u.eV
         energy_per_photon = const.h * const.c / self.wavelength / u.ph
         electron_per_photon = (electron_per_ev * energy_per_photon).to(_e / u.ph)
-        return electron_per_photon / (self._data["elecperdn"] * _e / u.count)
+        return electron_per_photon / (self._data["elecperdn"] * _e / u.DN)
 
     @u.quantity_input
     def wavelength_response(
@@ -369,7 +369,7 @@ def wavelength_response(
         include_eve_correction=False,
         include_crosstalk=True,
         **kwargs,
-    ) -> u.count / u.ph * u.cm**2:
+    ) -> u.DN / u.ph * u.cm**2:
         r"""
         The wavelength response function is the product of the gain and the
         effective area.

aiapy/response/tests/test_channel.py

@@ -177,8 +177,8 @@ def test_wavelength_response_no_idl(channel):
 def test_wavelength_response_uncorrected(channel, idl_environment):
     r = channel.wavelength_response()
     ssw = idl_environment.run("r = aia_get_response(/area,/dn,evenorm=0)", save_vars=["r"], verbose=False)
-    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.count / u.ph
-    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.count / u.ph)
+    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.DN / u.ph
+    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.DN / u.ph)
 
 
 def test_wavelength_response_no_crosstalk(channel, idl_environment):
@@ -188,8 +188,8 @@ def test_wavelength_response_no_crosstalk(channel, idl_environment):
         save_vars=["r"],
         verbose=False,
     )
-    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.count / u.ph
-    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.count / u.ph)
+    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.DN / u.ph
+    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.DN / u.ph)
 
 
 @pytest.mark.parametrize("include_eve_correction", [False, True])
@@ -217,8 +217,8 @@ def test_wavelength_response_time(channel, idl_environment, include_eve_correcti
         },
         verbose=False,
     )
-    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.count / u.ph
-    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.count / u.ph)
+    r_ssw = ssw["r"][f"A{channel.name}"][0]["ea"][0] * u.cm**2 * u.DN / u.ph
+    assert u.allclose(r, r_ssw, rtol=1e-4, atol=0.0 * u.cm**2 * u.DN / u.ph)
 
 
 @pytest.mark.remote_data()

aiapy/tests/test_idl.py

@@ -14,10 +14,10 @@
 
 @pytest.mark.parametrize(
     ("channel", "counts", "include_eve", "include_preflight", "include_chianti"),
-    [[c, 10 * u.ct / u.pixel] + 3 * [False] for c in CHANNELS]
+    [[c, 10 * u.DN / u.pixel] + 3 * [False] for c in CHANNELS]
     + [
-        [171 * u.angstrom, 1000 * u.ct / u.pix, True, False, False],
-        [171 * u.angstrom, 1000 * u.ct / u.pix, False, False, True],
+        [171 * u.angstrom, 1000 * u.DN / u.pix, True, False, False],
+        [171 * u.angstrom, 1000 * u.DN / u.pix, False, False, True],
     ],
 )
 def test_error_consistent(idl_environment, channel, counts, include_eve, include_preflight, include_chianti):

changelog/338.breaking.rst

@@ -0,0 +1 @@
+All inputs and outputs that previously used "counts" or "ct" now use "DN".

examples/correct_degradation.py

@@ -35,7 +35,7 @@
 # metadata that was returned. We select those and nothing else.
 
 table = results["jsoc"].show("DATE__OBS", "DATAMEAN")
-table["DATAMEAN"].unit = u.ct
+table["DATAMEAN"].unit = u.DN
 table["DATE_OBS"] = astropy.time.Time(table["DATE__OBS"], scale="utc")
 del table["DATE__OBS"]