RCAL-942: Remove Units from Reference File Datamodels

changes/1474.general.rst

@@ -0,0 +1 @@
+Remove units from reference file datamodels.

pyproject.toml

@@ -24,9 +24,11 @@ dependencies = [
     "pyparsing >=2.4.7",
     "requests >=2.26",
     # "rad>=0.21.0,<0.22.0",
-    "rad @ git+https://github.com/spacetelescope/rad.git",
+    # "rad @ git+https://github.com/spacetelescope/rad.git",
+    "rad @ git+https://github.com/WilliamJamieson/rad.git@feature/remove_ref_units",
     # "roman_datamodels>=0.21.0,<0.22.0",
-    "roman_datamodels @ git+https://github.com/spacetelescope/roman_datamodels.git",
+    # "roman_datamodels @ git+https://github.com/spacetelescope/roman_datamodels.git",
+    "roman_datamodels @ git+https://github.com/WilliamJamieson/roman_datamodels.git@feature/remove_some_units",
     "scipy >=1.11",
 #    "stcal>=1.8.0,<1.9.0",
     "stcal @ git+https://github.com/spacetelescope/stcal.git@main",

romancal/dark_current/dark_current_step.py

@@ -54,7 +54,7 @@ def process(self, input):
             # Do the dark correction
             out_model = input_model
             nresultants = len(input_model.meta.exposure["read_pattern"])
-            out_model.data -= dark_model.data[:nresultants].value
+            out_model.data -= dark_model.data[:nresultants]
             out_model.pixeldq |= dark_model.dq
             out_model.meta.cal_step.dark = "COMPLETE"
 

romancal/dark_current/tests/test_dark.py

@@ -60,14 +60,14 @@ def test_dark_step_subtraction(instrument, exptype):
     # populate data array of science cube
     for i in range(0, 20):
         ramp_model.data[0, 0, i] = i
-        darkref_model.data[0, 0, i] = i * 0.1 * darkref_model.data.unit
+        darkref_model.data[0, 0, i] = i * 0.1
     orig_model = ramp_model.copy()
 
     # Perform Dark Current subtraction step
     result = DarkCurrentStep.call(ramp_model, override_dark=darkref_model)
 
     # check that the dark file is subtracted frame by frame from the science data
-    diff = orig_model.data - darkref_model.data.value
+    diff = orig_model.data - darkref_model.data
 
     # test that the output data file is equal to the difference found when subtracting
     # reffile from sci file

romancal/jump/jump_step.py

@@ -110,13 +110,13 @@ def process(self, input):
         gain_filename = self.get_reference_file(input_model, "gain")
         self.log.info("Using GAIN reference file: %s", gain_filename)
         gain_model = rdd.GainRefModel(gain_filename)
-        gain_2d = gain_model.data.value
+        gain_2d = gain_model.data
 
         readnoise_filename = self.get_reference_file(input_model, "readnoise")
         self.log.info("Using READNOISE reference file: %s", readnoise_filename)
         readnoise_model = rdd.ReadnoiseRefModel(readnoise_filename)
         # This is to clear the WRITEABLE=False flag?
-        readnoise_2d = np.copy(readnoise_model.data.value)
+        readnoise_2d = np.copy(readnoise_model.data)
 
         # DG 0810/21:  leave for now; make dqflags changes in a later,
         #              separate PR

romancal/jump/tests/test_jump_step.py

@@ -7,7 +7,6 @@
 import numpy as np
 import pytest
 import roman_datamodels.stnode as rds
-from astropy import units as u
 from astropy.time import Time
 from roman_datamodels import datamodels as rdm
 from roman_datamodels import maker_utils
@@ -42,15 +41,10 @@ def generate_wfi_reffiles(tmp_path_factory):
     meta["useafter"] = Time("2022-01-01T11:11:11.111")
 
     gain_ref["meta"] = meta
-    gain_ref["data"] = u.Quantity(
-        np.ones(shape, dtype=np.float32) * ingain, u.electron / u.DN, dtype=np.float32
-    )
+    gain_ref["data"] = np.ones(shape, dtype=np.float32) * ingain
+
     gain_ref["dq"] = np.zeros(shape, dtype=np.uint16)
-    gain_ref["err"] = u.Quantity(
-        (RNG.uniform(size=shape) * 0.05).astype(np.float64),
-        u.electron / u.DN,
-        dtype=np.float64,
-    )
+    gain_ref["err"] = (RNG.uniform(size=shape) * 0.05).astype(np.float64)
 
     gain_ref_model = GainRefModel(gain_ref)
     gain_ref_model.save(gainfile)
@@ -69,13 +63,10 @@ def generate_wfi_reffiles(tmp_path_factory):
     meta["exposure"]["p_exptype"] = "WFI_IMAGE|WFI_GRISM|WFI_PRISM|"
 
     rn_ref["meta"] = meta
-    rn_ref["data"] = u.Quantity(
-        np.ones(shape, dtype=np.float32), u.DN, dtype=np.float32
-    )
+    rn_ref["data"] = np.ones(shape, dtype=np.float32)
+
     rn_ref["dq"] = np.zeros(shape, dtype=np.uint16)
-    rn_ref["err"] = u.Quantity(
-        (RNG.uniform(size=shape) * 0.05).astype(np.float64), u.DN, dtype=np.float64
-    )
+    rn_ref["err"] = (RNG.uniform(size=shape) * 0.05).astype(np.float64)
 
     rn_ref_model = ReadnoiseRefModel(rn_ref)
     rn_ref_model.save(readnoisefile)

romancal/photom/photom.py

@@ -1,8 +1,6 @@
 import logging
 import warnings
 
-from astropy import units as u
-
 log = logging.getLogger(__name__)
 log.setLevel(logging.DEBUG)
 
@@ -28,22 +26,14 @@ def photom_io(input_model, photom_metadata):
 
     # Store the conversion factor in the meta data
     log.info(f"photmjsr value: {conversion:.6g}")
-    input_model.meta.photometry.conversion_megajanskys = conversion.value
-    input_model.meta.photometry.conversion_microjanskys = conversion.to(
-        u.microjansky / u.arcsecond**2
-    ).value
+    input_model.meta.photometry.conversion_megajanskys = conversion
 
     # Get the scalar conversion uncertainty factor
     uncertainty_conv = photom_metadata["uncertainty"]
 
     # Store the uncertainty conversion factor in the meta data
     log.info(f"uncertainty value: {uncertainty_conv:.6g}")
-    input_model.meta.photometry.conversion_megajanskys_uncertainty = (
-        uncertainty_conv.value
-    )
-    input_model.meta.photometry.conversion_microjanskys_uncertainty = (
-        uncertainty_conv.to(u.microjansky / u.arcsecond**2)
-    ).value
+    input_model.meta.photometry.conversion_megajanskys_uncertainty = uncertainty_conv
 
     # Return updated input model
     return input_model
@@ -64,12 +54,10 @@ def save_area_info(input_model, photom_parameters):
     """
 
     # Load the average pixel area values from the photom reference file header
-    area_ster = photom_parameters["pixelareasr"].value
-    area_a2 = photom_parameters["pixelareasr"].to(u.arcsecond**2).value
+    area_ster = photom_parameters["pixelareasr"]
 
     # Copy the pixel area values to the input model
-    log.debug(f"pixelarea_steradians = {area_ster}, pixelarea_arcsecsq = {area_a2}")
-    input_model.meta.photometry.pixelarea_arcsecsq = area_a2
+    log.debug(f"pixelarea_steradians = {area_ster}")
     input_model.meta.photometry.pixelarea_steradians = area_ster
 
     # Return updated input model

romancal/photom/tests/test_photom.py

@@ -46,19 +46,13 @@ def create_photom_wfi_image(min_r=3.1, delta=0.1):
     nrows = len(optical_element)
 
     # Create sample photometry keyword values
-    photmjsr = (
-        np.linspace(min_r, min_r + (nrows - 1.0) * delta, nrows)
-        * u.megajansky
-        / u.steradian
-    )
-    uncertainty = (
-        np.linspace(min_r / 20.0, min_r / 20.0 + (nrows - 1.0) * delta / 20.0, nrows)
-        * u.megajansky
-        / u.steradian
+    photmjsr = np.linspace(min_r, min_r + (nrows - 1.0) * delta, nrows)
+    uncertainty = np.linspace(
+        min_r / 20.0, min_r / 20.0 + (nrows - 1.0) * delta / 20.0, nrows
     )
 
     # Create sample area keyword values
-    area_ster = 2.31307642258977e-14 * u.steradian
+    area_ster = 2.31307642258977e-14
     pixelareasr = np.ones(nrows, dtype=np.float64) * area_ster
 
     # Bundle values into a list
@@ -138,63 +132,36 @@ def test_apply_photom1():
     output_model = photom.apply_photom(input_model, photom_model)
 
     # Set reference photometry
-    area_ster = 2.31307642258977e-14 * u.steradian
-    area_a2 = 0.000984102303070964 * u.arcsecond * u.arcsecond
+    area_ster = 2.31307642258977e-14
 
     # Tests for pixel areas
     assert np.isclose(
         output_model.meta.photometry.pixelarea_steradians,
-        area_ster.value,
-        atol=1.0e-7,
-    )
-    # assert output_model.meta.photometry.pixelarea_steradians.unit == area_ster.unit
-    assert np.isclose(
-        output_model.meta.photometry.pixelarea_arcsecsq,
-        area_a2.value,
+        area_ster,
         atol=1.0e-7,
     )
 
     # Set reference photometry
-    phot_ster = 3.5 * u.megajansky / u.steradian
-    phot_a2 = phot_ster.to(u.microjansky / u.arcsecond**2)
+    phot_ster = 3.5
 
     # Tests for photometry
     assert np.isclose(
         output_model.meta.photometry.conversion_megajanskys,
-        phot_ster.value,
-        atol=1.0e-7,
-    )
-    # assert output_model.meta.photometry.conversion_megajanskys.unit == phot_ster.unit
-    assert np.isclose(
-        output_model.meta.photometry.conversion_microjanskys,
-        phot_a2.value,
+        phot_ster,
         atol=1.0e-7,
     )
+
     # assert output_model.meta.photometry.conversion_microjanskys.unit == phot_a2.unit
 
     # Set reference photometric uncertainty
-    muphot_ster = 0.175 * u.megajansky / u.steradian
-    muphot_a2 = muphot_ster.to(u.microjansky / u.arcsecond**2)
+    muphot_ster = 0.175
 
     # Tests for photometric uncertainty
     assert np.isclose(
         output_model.meta.photometry.conversion_megajanskys_uncertainty,
-        muphot_ster.value,
-        atol=1.0e-7,
-    )
-    # assert (
-    #     output_model.meta.photometry.conversion_megajanskys_uncertainty.unit
-    #     == muphot_ster.unit
-    # )
-    assert np.isclose(
-        output_model.meta.photometry.conversion_microjanskys_uncertainty,
-        muphot_a2.value,
+        muphot_ster,
         atol=1.0e-7,
     )
-    # assert (
-    #     output_model.meta.photometry.conversion_microjanskys_uncertainty.unit
-    #     == muphot_a2.unit
-    # )
 
 
 def test_apply_photom2():

romancal/ramp_fitting/ramp_fit_step.py

@@ -3,7 +3,6 @@
 import logging
 
 import numpy as np
-from astropy import units as u
 from roman_datamodels import datamodels as rdd
 from roman_datamodels import maker_utils
 from roman_datamodels import stnode as rds
@@ -103,8 +102,8 @@ def ols_cas22(self, input_model, readnoise_model, gain_model, dark_model):
 
         resultants = input_model.data
         dq = input_model.groupdq
-        read_noise = readnoise_model.data.value
-        gain = gain_model.data.value
+        read_noise = readnoise_model.data
+        gain = gain_model.data
         read_time = input_model.meta.exposure.frame_time
 
         # Force read pattern to be pure lists not LNodes
@@ -115,9 +114,7 @@ def ols_cas22(self, input_model, readnoise_model, gain_model, dark_model):
         # add dark current back into resultants so that Poisson noise is
         # properly accounted for
         tbar = np.array([np.mean(reads) * read_time for reads in read_pattern])
-        resultants += (
-            dark_model.dark_slope.to(u.DN / u.s).value[None, ...] * tbar[:, None, None]
-        )
+        resultants += dark_model.dark_slope[None, ...] * tbar[:, None, None]
 
         # account for the gain
         resultants *= gain
@@ -142,7 +139,7 @@ def ols_cas22(self, input_model, readnoise_model, gain_model, dark_model):
         dq = output.dq.astype(np.uint32)
 
         # remove dark current contribution to slopes
-        slopes -= dark_model.dark_slope.to(u.DN / u.s).value * gain
+        slopes -= dark_model.dark_slope * gain
 
         # Propagate DQ flags forward.
         ramp_dq = get_pixeldq_flags(dq, input_model.pixeldq, slopes, err, gain)

romancal/ramp_fitting/tests/test_ramp_fit_cas22.py

@@ -2,7 +2,8 @@
 
 import numpy as np
 import pytest
-from astropy import units as u
+
+# from astropy import units as u
 from astropy.time import Time
 from roman_datamodels import maker_utils
 from roman_datamodels.datamodels import (
@@ -274,23 +275,11 @@ def generate_wfi_reffiles(
     gain_ref["meta"]["useafter"] = Time("2022-01-01T11:11:11.111")
 
     if randomize:
-        gain_ref["data"] = u.Quantity(
-            (rng.random(shape) * 0.5).astype(np.float32) * ingain,
-            u.electron / u.DN,
-            dtype=np.float32,
-        )
+        gain_ref["data"] = (rng.random(shape) * 0.5).astype(np.float32) * ingain
     else:
-        gain_ref["data"] = u.Quantity(
-            np.ones(shape).astype(np.float32) * ingain,
-            u.electron / u.DN,
-            dtype=np.float32,
-        )
+        gain_ref["data"] = np.ones(shape).astype(np.float32) * ingain
     gain_ref["dq"] = np.zeros(shape, dtype=np.uint16)
-    gain_ref["err"] = u.Quantity(
-        (rng.random(shape) * 0.05).astype(np.float32),
-        u.electron / u.DN,
-        dtype=np.float32,
-    )
+    gain_ref["err"] = (rng.random(shape) * 0.05).astype(np.float32)
 
     gain_ref_model = GainRefModel(gain_ref)
 
@@ -305,15 +294,9 @@ def generate_wfi_reffiles(
     rn_ref["meta"]["exposure"]["frame_time"] = 666
 
     if randomize:
-        rn_ref["data"] = u.Quantity(
-            (rng.random(shape) * rnoise).astype(np.float32),
-            u.DN,
-            dtype=np.float32,
-        )
+        rn_ref["data"] = ((rng.random(shape) * rnoise).astype(np.float32),)
     else:
-        rn_ref["data"] = u.Quantity(
-            np.ones(shape).astype(np.float32) * rnoise, u.DN, dtype=np.float32
-        )
+        rn_ref["data"] = np.ones(shape).astype(np.float32) * rnoise
 
     rn_ref_model = ReadnoiseRefModel(rn_ref)
 
@@ -329,9 +312,7 @@ def generate_wfi_reffiles(
     dark_ref["meta"]["exposure"]["type"] = "WFI_IMAGE"
     dark_ref["meta"]["exposure"]["frame_time"] = 666
 
-    dark_ref["dark_slope"] = u.Quantity(
-        np.zeros(shape).astype(np.float32) * 0.01, u.DN / u.s, dtype=np.float32
-    )
+    dark_ref["dark_slope"] = np.zeros(shape).astype(np.float32) * 0.01
 
     dark_ref_model = DarkRefModel(dark_ref)
 

romancal/saturation/saturation.py

@@ -41,7 +41,7 @@ def flag_saturation(input_model, ref_model):
     pdq = input_model.pixeldq[np.newaxis, :]
 
     # Copy information from saturation reference file
-    sat_thresh = ref_model.data.value.copy()
+    sat_thresh = ref_model.data.copy()
     sat_dq = ref_model.dq.copy()
 
     # Obtain dq arrays updated for saturation