From 525eec2e526d30e6947692dfa1fd9e8ff1ca2fbd Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Sun, 5 May 2024 20:08:02 -0400
Subject: [PATCH 1/6] Implemented test.

---
 romanisim/tests/test_image.py | 110 ++++++++++++++++++++++++++++++++++
 1 file changed, 110 insertions(+)

diff --git a/romanisim/tests/test_image.py b/romanisim/tests/test_image.py
index 8aff9e6a..9abdcfc2 100644
--- a/romanisim/tests/test_image.py
+++ b/romanisim/tests/test_image.py
@@ -746,3 +746,113 @@ def test_inject_source_into_image():
                    'flux': flux,
                    'tij': tij}
         af.write_to(os.path.join(artifactdir, 'dms231.asdf'))
+
+
+def test_inject_source_into_mosaic():
+    """Inject a source into a mosaic.
+    """
+
+    # Set constants and metadata
+    galsim.roman.n_pix = 100
+    rng_seed = 42
+    rng = galsim.UniformDeviate(rng_seed)
+    metadata = copy.deepcopy(parameters.default_parameters_dictionary)
+    metadata['instrument']['detector'] = 'WFI07'
+    metadata['instrument']['optical_element'] = 'F158'
+    metadata['exposure']['ma_table_number'] = 1
+    bandpass = [metadata['instrument']['optical_element']]
+    filter_name = metadata['instrument']['optical_element']
+    sca = 4
+
+    # Create PSF
+    twcs = wcs.get_wcs(metadata, usecrds=False)
+    psf = romanisim.psf.make_psf(sca, filter_name, wcs=twcs, chromatic=False, webbpsf=True)
+
+    # Create initial Level 2-like image
+
+    # Create Four-quadrant pattern of gaussian noise, centered around one
+    # Each quadrant's gaussin noise scales like total exposure time
+    # (total files contributed to each quadrant)
+
+    # Create gaussian noise generators
+    g1 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.01)
+    g2 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.02)
+    g3 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.05)
+    g4 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.1)
+
+    # Create cps Image data array
+    la_arr = np.zeros([200, 200])
+
+    # Populate the image array with gaussian noise from generators
+    g1.generate(la_arr[0:100, 0:100])
+    g2.generate(la_arr[0:100, 100:200])
+    g3.generate(la_arr[100:200, 0:100])
+    g4.generate(la_arr[100:200, 100:200])
+
+    # Poisson Noise of image array
+    var_poisson = (la_arr - 1)**2
+
+    # Set scaling factor for injected sources, one in the center of each quadrant
+    # Flux / sigma_p^2
+    Ct1 = la_arr[50][50] / var_poisson[50][50]
+    Ct2 = la_arr[50][150] / var_poisson[50][150]
+    Ct3 = la_arr[150][50] / var_poisson[150][50]
+    Ct4 = la_arr[150][150] / var_poisson[150][150]
+
+    # Create empty postage stamp galsim source images
+    sourcecounts1 = galsim.ImageF(100, 100, wcs=twcs, xmin=0, ymin=0)
+    sourcecounts2 = galsim.ImageF(100, 100, wcs=twcs, xmin=0, ymin=100)
+    sourcecounts3 = galsim.ImageF(100, 100, wcs=twcs, xmin=100, ymin=0)
+    sourcecounts4 = galsim.ImageF(100, 100, wcs=twcs, xmin=100, ymin=100)
+
+    # Create normalized psf source catalog (same source in each quadrant)
+    sc_dict = {"ra": [0.0], "dec": [0.0], "type": ["PSF"], "n": [-1.0],
+               "half_light_radius": [0.0], "pa": [0.0], "ba": [1.0], "F158": [1.0]}
+    source_cat = table.Table(sc_dict)
+    source_cat = catalog.table_to_catalog(source_cat, ["F158"])
+
+    # Simulate source postage stamps, one for each quadrant
+    image.add_objects_to_image(sourcecounts1, source_cat, xpos=[50], ypos=[50],
+                               psf=psf, flux_to_counts_factor=Ct1, bandpass="F158",
+                               filter_name="F158", rng=rng)
+
+    image.add_objects_to_image(sourcecounts2, source_cat, xpos=[50], ypos=[150],
+                               psf=psf, flux_to_counts_factor=Ct2, bandpass=bandpass,
+                               filter_name=filter_name, rng=rng)
+
+    image.add_objects_to_image(sourcecounts3, source_cat, xpos=[150], ypos=[50],
+                               psf=psf, flux_to_counts_factor=Ct3, bandpass=bandpass,
+                               filter_name=filter_name, rng=rng)
+
+    image.add_objects_to_image(sourcecounts4, source_cat, xpos=[150], ypos=[150],
+                               psf=psf, flux_to_counts_factor=Ct4, bandpass=bandpass,
+                               filter_name=filter_name, rng=rng)
+
+    # Scale the source images back by their quadrant's center flux ratios
+    sourcecounts1 /= Ct1
+    sourcecounts2 /= Ct2
+    sourcecounts3 /= Ct3
+    sourcecounts4 /= Ct4
+
+    # Add sources to the original image array -  one per quadrant
+    la_arr[0:100, 0:100] += sourcecounts1.array
+    la_arr[0:100, 100:200] += sourcecounts2.array
+    la_arr[100:200, 0:100] += sourcecounts3.array
+    la_arr[100:200, 100:200] += sourcecounts4.array
+
+    # Poisson Noise of the image with injected source
+    inject_var_poisson = (la_arr - 1)**2
+
+    # Ensure that the poisson variance of the source injected image in each quadrant
+    # is greater than variance of the original image
+    assert inject_var_poisson[150][150] > var_poisson[150][150]
+    assert inject_var_poisson[150][50] > var_poisson[150][50]
+    assert inject_var_poisson[50][150] > var_poisson[50][150]
+    assert inject_var_poisson[50][50] > var_poisson[50][50]
+
+    # Ensure that the poisson variance of the source injected image in each quadrant
+    # relatively scales with the exposure time
+    # Quadrants: 4 > 3 > 2 > 1
+    assert inject_var_poisson[150][150] > inject_var_poisson[150][50]
+    assert inject_var_poisson[150][50] > inject_var_poisson[50][150]
+    assert inject_var_poisson[50][150] > inject_var_poisson[50][50]

From b4c15110dc4eecaea3c21ed72917ea1f800ede60 Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Wed, 15 May 2024 22:22:54 -0400
Subject: [PATCH 2/6] Addressed PR comments.

---
 romanisim/tests/test_image.py | 104 +++++++++++++---------------------
 romanisim/tests/test_wcs.py   |   8 +--
 2 files changed, 43 insertions(+), 69 deletions(-)

diff --git a/romanisim/tests/test_image.py b/romanisim/tests/test_image.py
index 9abdcfc2..7e42aa22 100644
--- a/romanisim/tests/test_image.py
+++ b/romanisim/tests/test_image.py
@@ -18,7 +18,7 @@
 import numpy as np
 import galsim
 from galsim import roman
-from romanisim import image, parameters, catalog, psf, util, wcs, persistence, ramp, l1
+from romanisim import image, parameters, catalog, psf, util, wcs, persistence, ramp, l1, l3
 from astropy.coordinates import SkyCoord
 from astropy import units as u
 from astropy.time import Time
@@ -30,6 +30,7 @@
 from metrics_logger.decorators import metrics_logger
 from romanisim import log
 from roman_datamodels.stnode import WfiScienceRaw, WfiImage
+import roman_datamodels.maker_utils as maker_utils
 import romanisim.bandpass
 
 
@@ -753,25 +754,20 @@ def test_inject_source_into_mosaic():
     """
 
     # Set constants and metadata
-    galsim.roman.n_pix = 100
+    galsim.roman.n_pix = 200
     rng_seed = 42
-    rng = galsim.UniformDeviate(rng_seed)
     metadata = copy.deepcopy(parameters.default_parameters_dictionary)
     metadata['instrument']['detector'] = 'WFI07'
     metadata['instrument']['optical_element'] = 'F158'
     metadata['exposure']['ma_table_number'] = 1
-    bandpass = [metadata['instrument']['optical_element']]
-    filter_name = metadata['instrument']['optical_element']
-    sca = 4
 
     # Create PSF
     twcs = wcs.get_wcs(metadata, usecrds=False)
-    psf = romanisim.psf.make_psf(sca, filter_name, wcs=twcs, chromatic=False, webbpsf=True)
 
-    # Create initial Level 2-like image
+    # Create initial Level 3-like image
 
     # Create Four-quadrant pattern of gaussian noise, centered around one
-    # Each quadrant's gaussin noise scales like total exposure time
+    # Each quadrant's gaussian noise scales like total exposure time
     # (total files contributed to each quadrant)
 
     # Create gaussian noise generators
@@ -780,75 +776,53 @@ def test_inject_source_into_mosaic():
     g3 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.05)
     g4 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.1)
 
-    # Create cps Image data array
-    la_arr = np.zeros([200, 200])
+    # Create level 3-like image model
+    l3_img = maker_utils.mk_level2_image(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
+
+    # Update metadata in the l3 model
+    for key in metadata.keys():
+        l3_img.meta[key].update(metadata[key])
 
     # Populate the image array with gaussian noise from generators
-    g1.generate(la_arr[0:100, 0:100])
-    g2.generate(la_arr[0:100, 100:200])
-    g3.generate(la_arr[100:200, 0:100])
-    g4.generate(la_arr[100:200, 100:200])
-
-    # Poisson Noise of image array
-    var_poisson = (la_arr - 1)**2
-
-    # Set scaling factor for injected sources, one in the center of each quadrant
-    # Flux / sigma_p^2
-    Ct1 = la_arr[50][50] / var_poisson[50][50]
-    Ct2 = la_arr[50][150] / var_poisson[50][150]
-    Ct3 = la_arr[150][50] / var_poisson[150][50]
-    Ct4 = la_arr[150][150] / var_poisson[150][150]
-
-    # Create empty postage stamp galsim source images
-    sourcecounts1 = galsim.ImageF(100, 100, wcs=twcs, xmin=0, ymin=0)
-    sourcecounts2 = galsim.ImageF(100, 100, wcs=twcs, xmin=0, ymin=100)
-    sourcecounts3 = galsim.ImageF(100, 100, wcs=twcs, xmin=100, ymin=0)
-    sourcecounts4 = galsim.ImageF(100, 100, wcs=twcs, xmin=100, ymin=100)
+    g1.generate(l3_img.data.value[0:100, 0:100])
+    g2.generate(l3_img.data.value[0:100, 100:200])
+    g3.generate(l3_img.data.value[100:200, 0:100])
+    g4.generate(l3_img.data.value[100:200, 100:200])
+
+    # Define Poisson Noise of image array
+    l3_img.var_poisson.value[0:100, 0:100] = 0.01**2
+    l3_img.var_poisson.value[0:100, 100:200] = 0.02**2
+    l3_img.var_poisson.value[100:200, 0:100] = 0.05**2
+    l3_img.var_poisson.value[100:200, 100:200] = 0.1**2
 
     # Create normalized psf source catalog (same source in each quadrant)
-    sc_dict = {"ra": [0.0], "dec": [0.0], "type": ["PSF"], "n": [-1.0],
-               "half_light_radius": [0.0], "pa": [0.0], "ba": [1.0], "F158": [1.0]}
+    sc_dict = {"ra": 4 * [0.0], "dec": 4 * [0.0], "type": 4 * ["PSF"], "n": 4 * [-1.0],
+               "half_light_radius": 4 * [0.0], "pa": 4 * [0.0], "ba": 4 * [1.0], "F158": 4 * [1.0]}
     source_cat = table.Table(sc_dict)
-    source_cat = catalog.table_to_catalog(source_cat, ["F158"])
-
-    # Simulate source postage stamps, one for each quadrant
-    image.add_objects_to_image(sourcecounts1, source_cat, xpos=[50], ypos=[50],
-                               psf=psf, flux_to_counts_factor=Ct1, bandpass="F158",
-                               filter_name="F158", rng=rng)
-
-    image.add_objects_to_image(sourcecounts2, source_cat, xpos=[50], ypos=[150],
-                               psf=psf, flux_to_counts_factor=Ct2, bandpass=bandpass,
-                               filter_name=filter_name, rng=rng)
 
-    image.add_objects_to_image(sourcecounts3, source_cat, xpos=[150], ypos=[50],
-                               psf=psf, flux_to_counts_factor=Ct3, bandpass=bandpass,
-                               filter_name=filter_name, rng=rng)
+    xpos, ypos = 50, 50
+    source_cat["ra"][0], source_cat["dec"][0] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    xpos, ypos = 50, 150
+    source_cat['ra'][1], source_cat['dec'][1] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    xpos, ypos = 150, 50
+    source_cat['ra'][2], source_cat['dec'][2] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    xpos, ypos = 150, 150
+    source_cat['ra'][3], source_cat['dec'][3] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
 
-    image.add_objects_to_image(sourcecounts4, source_cat, xpos=[150], ypos=[150],
-                               psf=psf, flux_to_counts_factor=Ct4, bandpass=bandpass,
-                               filter_name=filter_name, rng=rng)
-
-    # Scale the source images back by their quadrant's center flux ratios
-    sourcecounts1 /= Ct1
-    sourcecounts2 /= Ct2
-    sourcecounts3 /= Ct3
-    sourcecounts4 /= Ct4
+    source_cat = catalog.table_to_catalog(source_cat, ["F158"])
 
-    # Add sources to the original image array -  one per quadrant
-    la_arr[0:100, 0:100] += sourcecounts1.array
-    la_arr[0:100, 100:200] += sourcecounts2.array
-    la_arr[100:200, 0:100] += sourcecounts3.array
-    la_arr[100:200, 100:200] += sourcecounts4.array
+    # Add source_cat objects to l3_img
+    l3.add_objects_to_l3(l3_img, source_cat, seed=rng_seed)
 
     # Poisson Noise of the image with injected source
-    inject_var_poisson = (la_arr - 1)**2
+    inject_var_poisson = (l3_img.data.value - 1)**2
 
     # Ensure that the poisson variance of the source injected image in each quadrant
     # is greater than variance of the original image
-    assert inject_var_poisson[150][150] > var_poisson[150][150]
-    assert inject_var_poisson[150][50] > var_poisson[150][50]
-    assert inject_var_poisson[50][150] > var_poisson[50][150]
-    assert inject_var_poisson[50][50] > var_poisson[50][50]
+    assert inject_var_poisson[150][150] > l3_img.var_poisson.value[150][150]
+    assert inject_var_poisson[150][50] > l3_img.var_poisson.value[150][50]
+    assert inject_var_poisson[50][150] > l3_img.var_poisson.value[50][150]
+    assert inject_var_poisson[50][50] > l3_img.var_poisson.value[50][50]
 
     # Ensure that the poisson variance of the source injected image in each quadrant
     # relatively scales with the exposure time
diff --git a/romanisim/tests/test_wcs.py b/romanisim/tests/test_wcs.py
index 96236f44..d6f40a9c 100644
--- a/romanisim/tests/test_wcs.py
+++ b/romanisim/tests/test_wcs.py
@@ -13,10 +13,10 @@
 import pytest
 
 import roman_datamodels
-try:
-    import roman_datamodels.maker_utils as maker_utils
-except ImportError:
-    import roman_datamodels.testing.utils as maker_utils
+# try
+import roman_datamodels.maker_utils as maker_utils
+# except ImportError:
+#     import roman_datamodels.testing.utils as maker_utils
 
 
 def make_fake_distortion_function():

From 852e46b525a2ab6232647ec5fee674116f3b4b09 Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Fri, 17 May 2024 13:37:22 -0400
Subject: [PATCH 3/6] Added l3.py

---
 romanisim/l3.py | 73 +++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 73 insertions(+)
 create mode 100644 romanisim/l3.py

diff --git a/romanisim/l3.py b/romanisim/l3.py
new file mode 100644
index 00000000..650ad93d
--- /dev/null
+++ b/romanisim/l3.py
@@ -0,0 +1,73 @@
+"""Function for Level 3-like images.
+
+"""
+
+import numpy as np
+import galsim
+from . import log
+from . import image, wcs, psf
+
+
+def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
+    """Add objects to a Level 3-like image
+
+    Parameters
+    ----------
+    l3_img : ImageModel
+        Array containing a mosaic of images
+    source_cat : list
+        List of catalog objects to add to l3_img
+
+    Returns
+    -------
+    None
+        l3_img is updated in place
+    """
+
+    if rng is None and seed is None:
+        seed = 143
+        log.warning(
+            'No RNG set, constructing a new default RNG from default seed.')
+    if rng is None:
+        rng = galsim.UniformDeviate(seed)
+
+    # Obtain observation keywords
+    filter_name = l3_img.meta.instrument.optical_element
+    sca = int(l3_img.meta.instrument.detector[-2:])
+
+    # Generate wcs
+    twcs = wcs.get_wcs(l3_img.meta, usecrds=False)
+
+    # Create PSF
+    l3_psf = psf.make_psf(sca, filter_name, wcs=twcs, chromatic=False, webbpsf=True)
+
+    # Generate x,y positions for sources
+    coords = np.array([[o.sky_pos.ra.rad, o.sky_pos.dec.rad]
+                       for o in source_cat])
+    sourcecountsall = galsim.ImageF(l3_img.data.shape[0], l3_img.data.shape[1], wcs=twcs, xmin=0, ymin=0)
+    xpos, ypos = sourcecountsall.wcs._xy(coords[:, 0], coords[:, 1])
+    xpos = [round(x) for x in xpos]
+    ypos = [round(y) for y in ypos]
+
+    # Cycle over sources and add them to the image
+    for idx, (x, y) in enumerate(zip(xpos, ypos)):
+        # Set scaling factor for injected sources
+        # Flux / sigma_p^2
+        Ct = (l3_img.data[x][y].value / l3_img.var_poisson[x][y].value)
+
+        # Create empty postage stamp galsim source image
+        sourcecounts = galsim.ImageF(l3_img.data.shape[0], l3_img.data.shape[1], wcs=twcs, xmin=0, ymin=0)
+
+        # Simulate source postage stamp
+        image.add_objects_to_image(sourcecounts, [source_cat[idx]], xpos=[x], ypos=[y],
+                                   psf=l3_psf, flux_to_counts_factor=Ct, bandpass=[filter_name],
+                                   filter_name=filter_name, rng=rng)
+
+        # Scale the source image back by its flux ratios
+        sourcecounts /= Ct
+
+        # Add sources to the original image array
+        l3_img.data = (l3_img.data.value + sourcecounts.array) * l3_img.data.unit
+
+    # l3_img is updated in place, so no return
+    return None

From 700334a694adbfd57f72d216eadf5d91858e2907 Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Sun, 19 May 2024 17:38:18 -0400
Subject: [PATCH 4/6] Addressed PR comments.

---
 romanisim/l3.py               | 15 +++++-----
 romanisim/parameters.py       | 17 +++++++++++
 romanisim/tests/test_image.py | 53 +++++++++++++++++++++--------------
 romanisim/tests/test_wcs.py   |  3 --
 romanisim/wcs.py              | 41 +++++++++++++++++++++++++++
 5 files changed, 97 insertions(+), 32 deletions(-)

diff --git a/romanisim/l3.py b/romanisim/l3.py
index 650ad93d..0c3c0017 100644
--- a/romanisim/l3.py
+++ b/romanisim/l3.py
@@ -9,12 +9,12 @@
 
 
 def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
-    """Add objects to a Level 3-like image
+    """Add objects to a Level 3 mosaic
 
     Parameters
     ----------
-    l3_img : ImageModel
-        Array containing a mosaic of images
+    l3_img : MosaicModel
+        Mosaic of images
     source_cat : list
         List of catalog objects to add to l3_img
 
@@ -32,14 +32,13 @@ def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
         rng = galsim.UniformDeviate(seed)
 
     # Obtain observation keywords
-    filter_name = l3_img.meta.instrument.optical_element
-    sca = int(l3_img.meta.instrument.detector[-2:])
+    filter_name = l3_img.meta.basic.optical_element
 
-    # Generate wcs
-    twcs = wcs.get_wcs(l3_img.meta, usecrds=False)
+    # Generate WCS
+    twcs = wcs.get_mosaic_wcs(l3_img.meta)
 
     # Create PSF
-    l3_psf = psf.make_psf(sca, filter_name, wcs=twcs, chromatic=False, webbpsf=True)
+    l3_psf = psf.make_psf(filter_name=filter_name, sca=2, chromatic=False, webbpsf=True)
 
     # Generate x,y positions for sources
     coords = np.array([[o.sky_pos.ra.rad, o.sky_pos.dec.rad]
diff --git a/romanisim/parameters.py b/romanisim/parameters.py
index bf37d481..84930fac 100644
--- a/romanisim/parameters.py
+++ b/romanisim/parameters.py
@@ -58,6 +58,23 @@
                 },
 }
 
+# Default metadata for level 3 mosaics
+default_mosaic_parameters_dictionary = {
+    'basic':  {'time_mean_mjd': Time('2026-01-01T00:00:00').mjd,
+               'optical_element': 'F184',
+               },
+    'wcsinfo': {'ra_ref': 270.0,
+                'dec_ref': 66.0,
+                'v2_ref': 0,
+                'v3_ref': 0,
+                'roll_ref': 0,
+                'vparity': -1,
+                'v3yangle': -60.0,
+                # I don't know what vparity and v3yangle should really be,
+                # but they are always -1 and -60 in existing files.
+                },
+}
+
 reference_data = {
     "dark": 0.01 * u.electron / u.s,
     "distortion": None,
diff --git a/romanisim/tests/test_image.py b/romanisim/tests/test_image.py
index 6a8fb75f..d41468f2 100644
--- a/romanisim/tests/test_image.py
+++ b/romanisim/tests/test_image.py
@@ -757,6 +757,8 @@ def test_inject_source_into_image():
         af.write_to(os.path.join(artifactdir, 'dms231.asdf'))
 
 
+@metrics_logger("DMS232")
+@pytest.mark.soctests
 def test_inject_source_into_mosaic():
     """Inject a source into a mosaic.
     """
@@ -764,13 +766,11 @@ def test_inject_source_into_mosaic():
     # Set constants and metadata
     galsim.roman.n_pix = 200
     rng_seed = 42
-    metadata = copy.deepcopy(parameters.default_parameters_dictionary)
-    metadata['instrument']['detector'] = 'WFI07'
-    metadata['instrument']['optical_element'] = 'F158'
-    metadata['exposure']['ma_table_number'] = 1
+    metadata = copy.deepcopy(parameters.default_mosaic_parameters_dictionary)
+    metadata['basic']['optical_element'] = 'F158'
 
-    # Create PSF
-    twcs = wcs.get_wcs(metadata, usecrds=False)
+    # Create WCS
+    twcs = wcs.get_mosaic_wcs(metadata)
 
     # Create initial Level 3-like image
 
@@ -779,17 +779,18 @@ def test_inject_source_into_mosaic():
     # (total files contributed to each quadrant)
 
     # Create gaussian noise generators
-    g1 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.01)
-    g2 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.02)
-    g3 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.05)
-    g4 = galsim.GaussianDeviate(31415926, mean=1.0, sigma=0.1)
+    g1 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.01)
+    g2 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.02)
+    g3 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.05)
+    g4 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.1)
 
     # Create level 3-like image model
-    l3_img = maker_utils.mk_level2_image(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
+    l3_img = maker_utils.mk_level3_mosaic(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
 
     # Update metadata in the l3 model
     for key in metadata.keys():
-        l3_img.meta[key].update(metadata[key])
+        if key in l3_img.meta:
+            l3_img.meta[key].update(metadata[key])
 
     # Populate the image array with gaussian noise from generators
     g1.generate(l3_img.data.value[0:100, 0:100])
@@ -825,16 +826,26 @@ def test_inject_source_into_mosaic():
     # Poisson Noise of the image with injected source
     inject_var_poisson = (l3_img.data.value - 1)**2
 
-    # Ensure that the poisson variance of the source injected image in each quadrant
+    # Ensure that the total poisson variance of the source injected image in each quadrant
     # is greater than variance of the original image
-    assert inject_var_poisson[150][150] > l3_img.var_poisson.value[150][150]
-    assert inject_var_poisson[150][50] > l3_img.var_poisson.value[150][50]
-    assert inject_var_poisson[50][150] > l3_img.var_poisson.value[50][150]
-    assert inject_var_poisson[50][50] > l3_img.var_poisson.value[50][50]
+    assert np.sum(inject_var_poisson[0:100,0:100]) > np.sum(l3_img.var_poisson.value[0:100,0:100])
+    assert np.sum(inject_var_poisson[0:100,100:200]) > np.sum(l3_img.var_poisson.value[0:100,100:200])
+    assert np.sum(inject_var_poisson[100:200,0:100]) > np.sum(l3_img.var_poisson.value[100:200,0:100])
+    assert np.sum(inject_var_poisson[100:200,100:200]) > np.sum(l3_img.var_poisson.value[100:200,100:200])
 
-    # Ensure that the poisson variance of the source injected image in each quadrant
+    # Ensure that the total poisson variance of the source injected image in each quadrant
     # relatively scales with the exposure time
     # Quadrants: 4 > 3 > 2 > 1
-    assert inject_var_poisson[150][150] > inject_var_poisson[150][50]
-    assert inject_var_poisson[150][50] > inject_var_poisson[50][150]
-    assert inject_var_poisson[50][150] > inject_var_poisson[50][50]
+    assert np.sum(inject_var_poisson[100:200,100:200]) > np.sum(inject_var_poisson[100:200,0:100])
+    assert np.sum(inject_var_poisson[100:200,0:100]) > np.sum(inject_var_poisson[0:100,100:200])
+    assert np.sum(inject_var_poisson[0:100,100:200]) > np.sum(inject_var_poisson[0:100,0:100])
+
+    # Create log entry and artifacts
+    log.info('DMS232 successfully injected sources into a mosiac at points (50,50), (50,150), (150,50), (150,150).')
+
+    artifactdir = os.environ.get('TEST_ARTIFACT_DIR', None)
+    if artifactdir is not None:
+        af = asdf.AsdfFile()
+        af.tree = {'l3_img': l3_img,
+                   'inject_var_poisson': inject_var_poisson}
+        af.write_to(os.path.join(artifactdir, 'dms232.asdf'))
diff --git a/romanisim/tests/test_wcs.py b/romanisim/tests/test_wcs.py
index d6f40a9c..9b25067c 100644
--- a/romanisim/tests/test_wcs.py
+++ b/romanisim/tests/test_wcs.py
@@ -13,10 +13,7 @@
 import pytest
 
 import roman_datamodels
-# try
 import roman_datamodels.maker_utils as maker_utils
-# except ImportError:
-#     import roman_datamodels.testing.utils as maker_utils
 
 
 def make_fake_distortion_function():
diff --git a/romanisim/wcs.py b/romanisim/wcs.py
index f45a27cb..d6a4eaaf 100644
--- a/romanisim/wcs.py
+++ b/romanisim/wcs.py
@@ -404,3 +404,44 @@ def convert_wcs_to_gwcs(wcs):
     else:
         # make a gwcs WCS from a galsim.roman WCS
         return wcs_from_fits_header(wcs.header.header)
+
+
+def get_mosaic_wcs(mosaic):
+    """Get a WCS object for a given sca or set of CRDS parameters.
+
+    Parameters
+    ----------
+    mosaic : roman_datamodels.datamodels.MosaicModel or dict
+        Mosaic model or dictionary containing WCS parameters.
+
+    Returns
+    -------
+    galsim.CelestialWCS for the mosaic
+    """
+
+    # If sent a dictionary, create a temporary model for data interface
+    if (type(mosaic) is not roman_datamodels.datamodels.MosaicModel):
+        mosaic_node = maker_utils.mk_level3_mosaic()
+        for key in mosaic.keys():
+            if isinstance(mosaic[key], dict):
+                mosaic_node['meta'][key].update(mosaic[key])
+            else:
+                mosaic_node['meta'][key] = mosaic[key]
+        mosaic_node = roman_datamodels.datamodels.MosaicModel(mosaic_node)
+    else:
+        mosaic_node = mosaic
+
+    world_pos = astropy.coordinates.SkyCoord(
+        mosaic_node.meta.wcsinfo.ra_ref * u.deg,
+        mosaic_node.meta.wcsinfo.dec_ref * u.deg)
+
+    # Create a tangent plane WCS for the mosaic
+    # The affine parameters below should be reviewed and updated
+    affine = galsim.AffineTransform(
+        0.1, 0, 0, 0.1, origin=galsim.PositionI(mosaic_node.data.shape[0]/2,
+                                                mosaic_node.data.shape[1]/2,),
+        world_origin=galsim.PositionD(0, 0))
+    wcs = galsim.TanWCS(affine,
+                        util.celestialcoord(world_pos))
+
+    return wcs

From 3402ce41d4b506a95a105f168fca16b410341e72 Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Sun, 19 May 2024 17:53:36 -0400
Subject: [PATCH 5/6] Flake8 and comment updates.

---
 romanisim/l3.py               | 28 ++++++++++++++--------------
 romanisim/parameters.py       |  8 ++++----
 romanisim/tests/test_image.py | 30 +++++++++++++++---------------
 romanisim/wcs.py              |  8 ++++----
 4 files changed, 37 insertions(+), 37 deletions(-)

diff --git a/romanisim/l3.py b/romanisim/l3.py
index 0c3c0017..48654e99 100644
--- a/romanisim/l3.py
+++ b/romanisim/l3.py
@@ -8,20 +8,20 @@
 from . import image, wcs, psf
 
 
-def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
+def add_objects_to_l3(l3_mos, source_cat, rng=None, seed=None):
     """Add objects to a Level 3 mosaic
 
     Parameters
     ----------
-    l3_img : MosaicModel
+    l3_mos : MosaicModel
         Mosaic of images
     source_cat : list
-        List of catalog objects to add to l3_img
+        List of catalog objects to add to l3_mos
 
     Returns
     -------
     None
-        l3_img is updated in place
+        l3_mos is updated in place
     """
 
     if rng is None and seed is None:
@@ -31,11 +31,11 @@ def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
     if rng is None:
         rng = galsim.UniformDeviate(seed)
 
-    # Obtain observation keywords
-    filter_name = l3_img.meta.basic.optical_element
+    # Obtain optical element
+    filter_name = l3_mos.meta.basic.optical_element
 
     # Generate WCS
-    twcs = wcs.get_mosaic_wcs(l3_img.meta)
+    twcs = wcs.get_mosaic_wcs(l3_mos.meta)
 
     # Create PSF
     l3_psf = psf.make_psf(filter_name=filter_name, sca=2, chromatic=False, webbpsf=True)
@@ -43,19 +43,19 @@ def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
     # Generate x,y positions for sources
     coords = np.array([[o.sky_pos.ra.rad, o.sky_pos.dec.rad]
                        for o in source_cat])
-    sourcecountsall = galsim.ImageF(l3_img.data.shape[0], l3_img.data.shape[1], wcs=twcs, xmin=0, ymin=0)
+    sourcecountsall = galsim.ImageF(l3_mos.data.shape[0], l3_mos.data.shape[1], wcs=twcs, xmin=0, ymin=0)
     xpos, ypos = sourcecountsall.wcs._xy(coords[:, 0], coords[:, 1])
     xpos = [round(x) for x in xpos]
     ypos = [round(y) for y in ypos]
 
-    # Cycle over sources and add them to the image
+    # Cycle over sources and add them to the mosaic
     for idx, (x, y) in enumerate(zip(xpos, ypos)):
         # Set scaling factor for injected sources
         # Flux / sigma_p^2
-        Ct = (l3_img.data[x][y].value / l3_img.var_poisson[x][y].value)
+        Ct = (l3_mos.data[x][y].value / l3_mos.var_poisson[x][y].value)
 
         # Create empty postage stamp galsim source image
-        sourcecounts = galsim.ImageF(l3_img.data.shape[0], l3_img.data.shape[1], wcs=twcs, xmin=0, ymin=0)
+        sourcecounts = galsim.ImageF(l3_mos.data.shape[0], l3_mos.data.shape[1], wcs=twcs, xmin=0, ymin=0)
 
         # Simulate source postage stamp
         image.add_objects_to_image(sourcecounts, [source_cat[idx]], xpos=[x], ypos=[y],
@@ -65,8 +65,8 @@ def add_objects_to_l3(l3_img, source_cat, rng=None, seed=None):
         # Scale the source image back by its flux ratios
         sourcecounts /= Ct
 
-        # Add sources to the original image array
-        l3_img.data = (l3_img.data.value + sourcecounts.array) * l3_img.data.unit
+        # Add sources to the original mosaic data array
+        l3_mos.data = (l3_mos.data.value + sourcecounts.array) * l3_mos.data.unit
 
-    # l3_img is updated in place, so no return
+    # l3_mos is updated in place, so no return
     return None
diff --git a/romanisim/parameters.py b/romanisim/parameters.py
index 84930fac..221a14f2 100644
--- a/romanisim/parameters.py
+++ b/romanisim/parameters.py
@@ -55,14 +55,14 @@
                 },
     'aperture': {'name': 'WFI_CEN',
                  'position_angle': 0
-                },
+                 },
 }
 
 # Default metadata for level 3 mosaics
 default_mosaic_parameters_dictionary = {
-    'basic':  {'time_mean_mjd': Time('2026-01-01T00:00:00').mjd,
-               'optical_element': 'F184',
-               },
+    'basic': {'time_mean_mjd': Time('2026-01-01T00:00:00').mjd,
+              'optical_element': 'F184',
+              },
     'wcsinfo': {'ra_ref': 270.0,
                 'dec_ref': 66.0,
                 'v2_ref': 0,
diff --git a/romanisim/tests/test_image.py b/romanisim/tests/test_image.py
index d41468f2..eec08dd3 100644
--- a/romanisim/tests/test_image.py
+++ b/romanisim/tests/test_image.py
@@ -759,8 +759,8 @@ def test_inject_source_into_image():
 
 @metrics_logger("DMS232")
 @pytest.mark.soctests
-def test_inject_source_into_mosaic():
-    """Inject a source into a mosaic.
+def test_inject_sources_into_mosaic():
+    """Inject sources into a mosaic.
     """
 
     # Set constants and metadata
@@ -772,7 +772,7 @@ def test_inject_source_into_mosaic():
     # Create WCS
     twcs = wcs.get_mosaic_wcs(metadata)
 
-    # Create initial Level 3-like image
+    # Create initial Level 3 mosaic
 
     # Create Four-quadrant pattern of gaussian noise, centered around one
     # Each quadrant's gaussian noise scales like total exposure time
@@ -784,7 +784,7 @@ def test_inject_source_into_mosaic():
     g3 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.05)
     g4 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.1)
 
-    # Create level 3-like image model
+    # Create level 3 mosaic model
     l3_img = maker_utils.mk_level3_mosaic(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
 
     # Update metadata in the l3 model
@@ -792,13 +792,13 @@ def test_inject_source_into_mosaic():
         if key in l3_img.meta:
             l3_img.meta[key].update(metadata[key])
 
-    # Populate the image array with gaussian noise from generators
+    # Populate the mosaic data array with gaussian noise from generators
     g1.generate(l3_img.data.value[0:100, 0:100])
     g2.generate(l3_img.data.value[0:100, 100:200])
     g3.generate(l3_img.data.value[100:200, 0:100])
     g4.generate(l3_img.data.value[100:200, 100:200])
 
-    # Define Poisson Noise of image array
+    # Define Poisson Noise of mosaic
     l3_img.var_poisson.value[0:100, 0:100] = 0.01**2
     l3_img.var_poisson.value[0:100, 100:200] = 0.02**2
     l3_img.var_poisson.value[100:200, 0:100] = 0.05**2
@@ -820,25 +820,25 @@ def test_inject_source_into_mosaic():
 
     source_cat = catalog.table_to_catalog(source_cat, ["F158"])
 
-    # Add source_cat objects to l3_img
+    # Add source_cat objects to mosaic
     l3.add_objects_to_l3(l3_img, source_cat, seed=rng_seed)
 
-    # Poisson Noise of the image with injected source
+    # Poisson Noise of the mosaic with injected source
     inject_var_poisson = (l3_img.data.value - 1)**2
 
     # Ensure that the total poisson variance of the source injected image in each quadrant
     # is greater than variance of the original image
-    assert np.sum(inject_var_poisson[0:100,0:100]) > np.sum(l3_img.var_poisson.value[0:100,0:100])
-    assert np.sum(inject_var_poisson[0:100,100:200]) > np.sum(l3_img.var_poisson.value[0:100,100:200])
-    assert np.sum(inject_var_poisson[100:200,0:100]) > np.sum(l3_img.var_poisson.value[100:200,0:100])
-    assert np.sum(inject_var_poisson[100:200,100:200]) > np.sum(l3_img.var_poisson.value[100:200,100:200])
+    assert np.sum(inject_var_poisson[0:100, 0:100]) > np.sum(l3_img.var_poisson.value[0:100, 0:100])
+    assert np.sum(inject_var_poisson[0:100, 100:200]) > np.sum(l3_img.var_poisson.value[0:100, 100:200])
+    assert np.sum(inject_var_poisson[100:200, 0:100]) > np.sum(l3_img.var_poisson.value[100:200, 0:100])
+    assert np.sum(inject_var_poisson[100:200, 100:200]) > np.sum(l3_img.var_poisson.value[100:200, 100:200])
 
     # Ensure that the total poisson variance of the source injected image in each quadrant
     # relatively scales with the exposure time
     # Quadrants: 4 > 3 > 2 > 1
-    assert np.sum(inject_var_poisson[100:200,100:200]) > np.sum(inject_var_poisson[100:200,0:100])
-    assert np.sum(inject_var_poisson[100:200,0:100]) > np.sum(inject_var_poisson[0:100,100:200])
-    assert np.sum(inject_var_poisson[0:100,100:200]) > np.sum(inject_var_poisson[0:100,0:100])
+    assert np.sum(inject_var_poisson[100:200, 100:200]) > np.sum(inject_var_poisson[100:200, 0:100])
+    assert np.sum(inject_var_poisson[100:200, 0:100]) > np.sum(inject_var_poisson[0:100, 100:200])
+    assert np.sum(inject_var_poisson[0:100, 100:200]) > np.sum(inject_var_poisson[0:100, 0:100])
 
     # Create log entry and artifacts
     log.info('DMS232 successfully injected sources into a mosiac at points (50,50), (50,150), (150,50), (150,150).')
diff --git a/romanisim/wcs.py b/romanisim/wcs.py
index d6a4eaaf..ec8f5055 100644
--- a/romanisim/wcs.py
+++ b/romanisim/wcs.py
@@ -77,8 +77,8 @@ def fill_in_parameters(parameters, coord, pa_aper=0, boresight=True):
         parameters['pointing']['dec_v1'])
 
     # Romanisim uses ROLL_REF = PA_APER - V3IdlYAngle
-    V3IdlYAngle = -60 # this value should eventually be taken from the SIAF
-    parameters['wcsinfo']['roll_ref'] = pa_aper - V3IdlYAngle 
+    V3IdlYAngle = -60  # this value should eventually be taken from the SIAF
+    parameters['wcsinfo']['roll_ref'] = pa_aper - V3IdlYAngle
 
     if boresight:
         parameters['wcsinfo']['v2_ref'] = 0
@@ -438,8 +438,8 @@ def get_mosaic_wcs(mosaic):
     # Create a tangent plane WCS for the mosaic
     # The affine parameters below should be reviewed and updated
     affine = galsim.AffineTransform(
-        0.1, 0, 0, 0.1, origin=galsim.PositionI(mosaic_node.data.shape[0]/2,
-                                                mosaic_node.data.shape[1]/2,),
+        0.1, 0, 0, 0.1, origin=galsim.PositionI(mosaic_node.data.shape[0] / 2.0,
+                                                mosaic_node.data.shape[1] / 2.0,),
         world_origin=galsim.PositionD(0, 0))
     wcs = galsim.TanWCS(affine,
                         util.celestialcoord(world_pos))

From 0a1aac517b844a4e94d7a4636ecdce2b20455388 Mon Sep 17 00:00:00 2001
From: PaulHuwe <phuwe@stsci.edu>
Date: Tue, 21 May 2024 19:15:12 -0400
Subject: [PATCH 6/6] Adressed PR comments.

---
 romanisim/l3.py               | 16 ++++++---
 romanisim/tests/test_image.py | 68 +++++++++++++++++------------------
 2 files changed, 45 insertions(+), 39 deletions(-)

diff --git a/romanisim/l3.py b/romanisim/l3.py
index 48654e99..0b47022c 100644
--- a/romanisim/l3.py
+++ b/romanisim/l3.py
@@ -45,11 +45,14 @@ def add_objects_to_l3(l3_mos, source_cat, rng=None, seed=None):
                        for o in source_cat])
     sourcecountsall = galsim.ImageF(l3_mos.data.shape[0], l3_mos.data.shape[1], wcs=twcs, xmin=0, ymin=0)
     xpos, ypos = sourcecountsall.wcs._xy(coords[:, 0], coords[:, 1])
-    xpos = [round(x) for x in xpos]
-    ypos = [round(y) for y in ypos]
+    xpos_idx = [round(x) for x in xpos]
+    ypos_idx = [round(y) for y in ypos]
+
+    # Create overall scaling factor map
+    Ct_all = (l3_mos.data.value / l3_mos.var_poisson)
 
     # Cycle over sources and add them to the mosaic
-    for idx, (x, y) in enumerate(zip(xpos, ypos)):
+    for idx, (x, y) in enumerate(zip(xpos_idx, ypos_idx)):
         # Set scaling factor for injected sources
         # Flux / sigma_p^2
         Ct = (l3_mos.data[x][y].value / l3_mos.var_poisson[x][y].value)
@@ -58,7 +61,7 @@ def add_objects_to_l3(l3_mos, source_cat, rng=None, seed=None):
         sourcecounts = galsim.ImageF(l3_mos.data.shape[0], l3_mos.data.shape[1], wcs=twcs, xmin=0, ymin=0)
 
         # Simulate source postage stamp
-        image.add_objects_to_image(sourcecounts, [source_cat[idx]], xpos=[x], ypos=[y],
+        image.add_objects_to_image(sourcecounts, [source_cat[idx]], xpos=[xpos[idx]], ypos=[ypos[idx]],
                                    psf=l3_psf, flux_to_counts_factor=Ct, bandpass=[filter_name],
                                    filter_name=filter_name, rng=rng)
 
@@ -68,5 +71,10 @@ def add_objects_to_l3(l3_mos, source_cat, rng=None, seed=None):
         # Add sources to the original mosaic data array
         l3_mos.data = (l3_mos.data.value + sourcecounts.array) * l3_mos.data.unit
 
+        # Note for the future - other noise sources (read and flat) need to be implemented
+
+    # Set new poisson variance
+    l3_mos.var_poisson = l3_mos.data.value / Ct_all
+
     # l3_mos is updated in place, so no return
     return None
diff --git a/romanisim/tests/test_image.py b/romanisim/tests/test_image.py
index eec08dd3..554b2227 100644
--- a/romanisim/tests/test_image.py
+++ b/romanisim/tests/test_image.py
@@ -785,60 +785,56 @@ def test_inject_sources_into_mosaic():
     g4 = galsim.GaussianDeviate(rng_seed, mean=1.0, sigma=0.1)
 
     # Create level 3 mosaic model
-    l3_img = maker_utils.mk_level3_mosaic(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
+    l3_mos = maker_utils.mk_level3_mosaic(shape=(galsim.roman.n_pix, galsim.roman.n_pix))
 
     # Update metadata in the l3 model
     for key in metadata.keys():
-        if key in l3_img.meta:
-            l3_img.meta[key].update(metadata[key])
+        if key in l3_mos.meta:
+            l3_mos.meta[key].update(metadata[key])
 
     # Populate the mosaic data array with gaussian noise from generators
-    g1.generate(l3_img.data.value[0:100, 0:100])
-    g2.generate(l3_img.data.value[0:100, 100:200])
-    g3.generate(l3_img.data.value[100:200, 0:100])
-    g4.generate(l3_img.data.value[100:200, 100:200])
+    g1.generate(l3_mos.data.value[0:100, 0:100])
+    g2.generate(l3_mos.data.value[0:100, 100:200])
+    g3.generate(l3_mos.data.value[100:200, 0:100])
+    g4.generate(l3_mos.data.value[100:200, 100:200])
 
     # Define Poisson Noise of mosaic
-    l3_img.var_poisson.value[0:100, 0:100] = 0.01**2
-    l3_img.var_poisson.value[0:100, 100:200] = 0.02**2
-    l3_img.var_poisson.value[100:200, 0:100] = 0.05**2
-    l3_img.var_poisson.value[100:200, 100:200] = 0.1**2
+    l3_mos.var_poisson.value[0:100, 0:100] = 0.01**2
+    l3_mos.var_poisson.value[0:100, 100:200] = 0.02**2
+    l3_mos.var_poisson.value[100:200, 0:100] = 0.05**2
+    l3_mos.var_poisson.value[100:200, 100:200] = 0.1**2
 
     # Create normalized psf source catalog (same source in each quadrant)
     sc_dict = {"ra": 4 * [0.0], "dec": 4 * [0.0], "type": 4 * ["PSF"], "n": 4 * [-1.0],
                "half_light_radius": 4 * [0.0], "pa": 4 * [0.0], "ba": 4 * [1.0], "F158": 4 * [1.0]}
-    source_cat = table.Table(sc_dict)
+    sc_table = table.Table(sc_dict)
 
     xpos, ypos = 50, 50
-    source_cat["ra"][0], source_cat["dec"][0] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    sc_table["ra"][0], sc_table["dec"][0] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
     xpos, ypos = 50, 150
-    source_cat['ra'][1], source_cat['dec'][1] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    sc_table['ra'][1], sc_table['dec'][1] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
     xpos, ypos = 150, 50
-    source_cat['ra'][2], source_cat['dec'][2] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    sc_table['ra'][2], sc_table['dec'][2] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
     xpos, ypos = 150, 150
-    source_cat['ra'][3], source_cat['dec'][3] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
+    sc_table['ra'][3], sc_table['dec'][3] = (twcs._radec(xpos, ypos) * u.rad).to(u.deg).value
 
-    source_cat = catalog.table_to_catalog(source_cat, ["F158"])
+    source_cat = catalog.table_to_catalog(sc_table, ["F158"])
 
-    # Add source_cat objects to mosaic
-    l3.add_objects_to_l3(l3_img, source_cat, seed=rng_seed)
+    # Copy original Mosaic before adding sources
+    l3_mos_orig = l3_mos.copy()
 
-    # Poisson Noise of the mosaic with injected source
-    inject_var_poisson = (l3_img.data.value - 1)**2
+    # Add source_cat objects to mosaic
+    l3.add_objects_to_l3(l3_mos, source_cat, seed=rng_seed)
 
-    # Ensure that the total poisson variance of the source injected image in each quadrant
-    # is greater than variance of the original image
-    assert np.sum(inject_var_poisson[0:100, 0:100]) > np.sum(l3_img.var_poisson.value[0:100, 0:100])
-    assert np.sum(inject_var_poisson[0:100, 100:200]) > np.sum(l3_img.var_poisson.value[0:100, 100:200])
-    assert np.sum(inject_var_poisson[100:200, 0:100]) > np.sum(l3_img.var_poisson.value[100:200, 0:100])
-    assert np.sum(inject_var_poisson[100:200, 100:200]) > np.sum(l3_img.var_poisson.value[100:200, 100:200])
+    # Ensure that every data pixel value has increased or
+    # remained the same with the new sources injected
+    assert np.all(l3_mos.data.value >= l3_mos_orig.data.value)
 
-    # Ensure that the total poisson variance of the source injected image in each quadrant
-    # relatively scales with the exposure time
-    # Quadrants: 4 > 3 > 2 > 1
-    assert np.sum(inject_var_poisson[100:200, 100:200]) > np.sum(inject_var_poisson[100:200, 0:100])
-    assert np.sum(inject_var_poisson[100:200, 0:100]) > np.sum(inject_var_poisson[0:100, 100:200])
-    assert np.sum(inject_var_poisson[0:100, 100:200]) > np.sum(inject_var_poisson[0:100, 0:100])
+    # Ensure that every pixel's poisson variance has increased or
+    # remained the same with the new sources injected
+    # Numpy isclose is needed to determine equality, due to float precision issues
+    close_mask = np.isclose(l3_mos.var_poisson.value, l3_mos_orig.var_poisson.value, rtol=1e-06)
+    assert np.all(l3_mos.var_poisson.value[~close_mask] > l3_mos_orig.var_poisson.value[~close_mask])
 
     # Create log entry and artifacts
     log.info('DMS232 successfully injected sources into a mosiac at points (50,50), (50,150), (150,50), (150,150).')
@@ -846,6 +842,8 @@ def test_inject_sources_into_mosaic():
     artifactdir = os.environ.get('TEST_ARTIFACT_DIR', None)
     if artifactdir is not None:
         af = asdf.AsdfFile()
-        af.tree = {'l3_img': l3_img,
-                   'inject_var_poisson': inject_var_poisson}
+        af.tree = {'l3_mos': l3_mos,
+                   'l3_mos_orig': l3_mos_orig,
+                   'source_cat_table': sc_table,
+                   }
         af.write_to(os.path.join(artifactdir, 'dms232.asdf'))