Allow SCA -1 to indicate that multiple images should be made.

.github/workflows/ci.yml

@@ -71,8 +71,7 @@ jobs:
       cache-path: ${{ needs.data.outputs.path }}
       cache-key: data-${{ needs.data.outputs.hash }}
       envs: |
-        - linux: py39-oldestdeps-cov-xdist
-        - linux: py39-xdist
+        - linux: py310-oldestdeps-cov-xdist
         - linux: py310-xdist
         - linux: py3-xdist
         - macos: py3-xdist

.github/workflows/ci_cron.yml

@@ -60,6 +60,5 @@ jobs:
       cache-path: ${{ needs.data.outputs.path }}
       cache-key: data-${{ needs.data.outputs.hash }}
       envs: |
-        - macos: py39-xdist
         - macos: py310-xdist
         - linux: py3-devdeps-xdist

pyproject.toml

@@ -2,7 +2,7 @@
 name = 'romanisim'
 description = 'Nancy Grace Roman Space Telescope WFI Simulator'
 readme = 'README.md'
-requires-python = '>=3.9'
+requires-python = '>=3.10'
 license = { file = 'LICENSE' }
 authors = [{ name = 'STScI', email = '[email protected]' }]
 classifiers = [
@@ -13,7 +13,7 @@ classifiers = [
     'Operating System :: POSIX',
     'Programming Language :: Python :: 3',
     'Programming Language :: Python :: 3 :: Only',
-    'Programming Language :: Python :: 3.9',
+    'Programming Language :: Python :: 3.10',
 ]
 dependencies = [
     'asdf >=2.15.1',

romanisim/image.py

@@ -458,7 +458,9 @@ def simulate_counts_generic(image, exptime, objlist=None, psf=None,
 
     if not np.all(flat == 1):
         image.quantize()
-        image.array[:, :] = np.random.binomial(
+        rng_numpy_seed = rng.raw()
+        rng_numpy = np.random.default_rng(rng_numpy_seed)
+        image.array[:, :] = rng_numpy.binomial(
             image.array.astype('i4'), flat / maxflat)
 
     if dark is not None:
@@ -800,6 +802,7 @@ def simulate(metadata, objlist,
             extras["simulate_reffiles"][key] = value
 
     extras['simcatobj'] = simcatobj
+    extras['wcs'] = wcs.convert_wcs_to_gwcs(counts.wcs)
     log.info('Simulation complete.')
     return im, extras
 
@@ -878,12 +881,7 @@ def make_asdf(slope, slopevar_rn, slopevar_poisson, metadata=None,
         out['meta'].update(metadata)
 
     if imwcs is not None:  # add a WCS
-        if isinstance(imwcs, wcs.GWCS):
-            out['meta'].update(wcs=imwcs.wcs)
-        else:
-            # make a gwcs WCS from a galsim.roman WCS
-            out['meta'].update(
-                wcs=wcs.wcs_from_fits_header(imwcs.header.header))
+        out['meta'].update(wcs=wcs.convert_wcs_to_gwcs(imwcs))
 
     out['data'] = slope
     out['dq'] = np.zeros(slope.shape, dtype='u4')

romanisim/psf.py

@@ -152,9 +152,10 @@ def make_psf(sca, filter_name, wcs=None, webbpsf=True, pix=None,
                             pix=pix, chromatic=chromatic, **kw)
     elif pix is not None:
         raise ValueError('cannot set both pix and variable')
-    buf = 40
+    buf = 49
     # WebbPSF complains if we get too close to (0, 0) for some reason.
     # For other corners one can go to within a fraction of a pixel.
+    # if we go larger than 49 we have to change some of the tests, which use a 100x100 image.
     corners = dict(
         ll=[buf, buf], lr=[roman.n_pix - buf, buf],
         ul=[buf, roman.n_pix - buf], ur=[roman.n_pix - buf, roman.n_pix - buf])

romanisim/ris_make_utils.py

@@ -4,6 +4,7 @@
 from copy import deepcopy
 import os
 import re
+import numpy as np
 import asdf
 from astropy import table
 from astropy import time
@@ -13,7 +14,7 @@
 import roman_datamodels
 from roman_datamodels import stnode
 from romanisim import catalog, image, wcs
-from romanisim import parameters
+from romanisim import parameters, log
 
 
 def merge_nested_dicts(dict1, dict2):
@@ -148,6 +149,17 @@
             coord, bandpasses=bandpasses, nobj=nobj, rng=rng)
     else:
         cat = table.Table.read(catalog_name, comment="#", delimiter=" ")
+        bandpass = [f for f in cat.dtype.names if f[0] == 'F']
+        bad = np.zeros(len(cat), dtype='bool')
+        for b in bandpass:
+            bad |= ~np.isfinite(cat[b])
+            if hasattr(cat[b], 'mask'):
+                bad |= cat[b].mask
+        cat = cat[~bad]
+        nbad = np.sum(bad)
+        if nbad > 0:
+            log.info(f'Removing {nbad} catalog entries with non-finite or '
+                     'masked fluxes.')
 
     return cat
 
@@ -251,7 +263,7 @@
             args.pretend_spectral.upper())
 
     drop_extra_dq = getattr(args, 'drop_extra_dq', False)
-    if drop_extra_dq:
+    if drop_extra_dq and ('dq' in romanisimdict):
         romanisimdict.pop('dq')
 
     # Write file

romanisim/wcs.py

@@ -384,3 +384,22 @@
     gw.bounding_box = ((-0.5, nx - 0.5), (-0.5, ny - 0.5))
 
     return gw
+
+
+def convert_wcs_to_gwcs(wcs):
+    """Convert a GalSim WCS object into a GWCS object.
+
+    Parameters
+    ----------
+    wcs : gwcs.wcs.WCS or wcs.GWCS
+        input WCS to convert
+
+    Returns
+    -------
+    wcs.GWCS corresponding to wcs.
+    """
+    if isinstance(wcs, GWCS):
+        return wcs.wcs
+    else:
+        # make a gwcs WCS from a galsim.roman WCS
+        return wcs_from_fits_header(wcs.header.header)

scripts/romancal_make_regtest_l1s.sh

@@ -28,3 +28,4 @@ romanisim-make-image --radec 270.00 66.00 --level 1 --sca 1 --bandpass F158 --ca
 # truncated spectroscopic image
 romanisim-make-image --radec 270.00 66.00 --level 1 --sca 1 --bandpass F158 --catalog gaia-270-66-2027-06-01.ecsv --webbpsf --usecrds --ma_table_number 110 --date 2027-06-01T00:30:00 --rng_seed 7 --drop-extra-dq r0000201001001001001_01101_0004_WFI01_uncal.asdf --pretend-spectral GRISM &
 
+wait

scripts/romanisim-make-image

@@ -8,6 +8,60 @@ from astropy import units as u
 import galsim
 from romanisim import log, wcs, persistence, parameters
 from romanisim import ris_make_utils as ris
+from copy import copy
+
+
+def go(args):
+    if args.config is not None:
+        # Open and parse overrides file
+        with open(args.config, "r") as config_file:
+            config = yaml.safe_load(config_file)
+            combo_dict = parameters.__dict__
+            ris.merge_nested_dicts(combo_dict, config)
+            parameters.__dict__.update(combo_dict)
+    elif args.usecrds:
+        # don't use default values
+        for k in parameters.reference_data:
+            parameters.reference_data[k] = None
+
+    if args.sca == -1:
+        # simulate all 18 SCAs sequentially
+        origfilename = copy(args.filename)
+        origprevious = copy(args.previous)
+        for i in range(1, 19):
+            args.sca = i
+            args.filename = origfilename.format(f'WFI{args.sca:02d}')
+            if origprevious is not None:
+                args.previous = origprevious.format(f'WFI{args.sca:02d}')
+            go(args)
+        return
+
+    metadata = ris.set_metadata(
+        date=args.date, bandpass=args.bandpass,
+        sca=args.sca, ma_table_number=args.ma_table_number,
+        truncate=args.truncate)
+
+    if args.radec is not None:
+        coord = SkyCoord(ra=args.radec[0] * u.deg, dec=args.radec[1] * u.deg,
+                         frame='icrs')
+        wcs.fill_in_parameters(metadata, coord, boresight=args.boresight, pa_aper=args.roll)
+
+    rng = galsim.UniformDeviate(args.rng_seed)
+
+    # Create catalog
+    cat = ris.create_catalog(metadata=metadata, catalog_name=args.catalog,
+                             bandpasses=[args.bandpass], rng=rng,
+                             nobj=args.nobj, usecrds=args.usecrds)
+
+    # Create persistence object
+    if args.previous is not None:
+        persist = persistence.Persistence.read(args.previous)
+    else:
+        persist = persistence.Persistence()
+
+    # Simulate image and write to file
+    ris.simulate_image_file(args, metadata, cat, rng, persist)
+
 
 if __name__ == '__main__':
     parser = argparse.ArgumentParser(
@@ -61,40 +115,4 @@ if __name__ == '__main__':
                 "packages like galsim's roman package or STIPS may better "
                 "serve such purposes.")
 
-    if args.config is not None:
-        # Open and parse overrides file
-        with open(args.config, "r") as config_file:
-            config = yaml.safe_load(config_file)
-            combo_dict = parameters.__dict__ 
-            ris.merge_nested_dicts(combo_dict, config)
-            parameters.__dict__.update(combo_dict)
-    elif args.usecrds:
-        # don't use default values
-        for k in parameters.reference_data:
-            parameters.reference_data[k] = None
-
-    metadata = ris.set_metadata(
-        date=args.date, bandpass=args.bandpass,
-        sca=args.sca, ma_table_number=args.ma_table_number,
-        truncate=args.truncate)
-
-    if args.radec is not None:
-        coord = SkyCoord(ra=args.radec[0] * u.deg, dec=args.radec[1] * u.deg,
-                         frame='icrs')
-        wcs.fill_in_parameters(metadata, coord, boresight=args.boresight, pa_aper=args.roll)
-
-    rng = galsim.UniformDeviate(args.rng_seed)
-
-    # Create catalog
-    cat = ris.create_catalog(metadata=metadata, catalog_name=args.catalog,
-                             bandpasses=[args.bandpass], rng=rng,
-                             nobj=args.nobj, usecrds=args.usecrds)
-
-    # Create persistence object
-    if args.previous is not None:
-        persist = persistence.Persistence.read(args.previous)
-    else:
-        persist = persistence.Persistence()
-
-    # Simulate image and write to file
-    ris.simulate_image_file(args, metadata, cat, rng, persist)
+    go(args)