RCAL-762: Simulate APT program

pyproject.toml

@@ -19,18 +19,17 @@ classifiers = [
 dependencies = [
     "asdf >=2.15.1",
     "astropy >=5.3",
-    "crds >=10.3.1",
+    "crds >=11.16.16",
     "defusedxml >=0.5.0",
     "galsim >=2.5.1",
-    #"rad >=0.18.0",
-    #"roman_datamodels >=0.18.0",
-    "rad @ git+https://github.com/spacetelescope/rad.git",
-    "roman_datamodels @ git+https://github.com/spacetelescope/roman_datamodels.git",
-    # dev versions of rad / roman_datamodels until released with DN/s units.
+    "rad >=0.19.2",
+    "roman_datamodels >=0.19.1",
+    # "rad @ git+https://github.com/spacetelescope/rad.git",
+    # "roman_datamodels @ git+https://github.com/spacetelescope/roman_datamodels.git",
     "gwcs >=0.18.1",
     "jsonschema >=4.8",
     "numpy >=1.22",
-    "webbpsf >=1.0.0",
+    "webbpsf >=1.2.1",
     "Cython >=0.29.21",
 ]
 dynamic = [

romanisim/ris_make_utils.py

@@ -4,6 +4,7 @@
 from copy import deepcopy
 import os
 import re
+import defusedxml.ElementTree
 import numpy as np
 import asdf
 from astropy import table
@@ -16,6 +17,7 @@
 from romanisim import catalog, image, wcs
 from romanisim import parameters, log
 
+NMAP = {'apt':'{http://www.stsci.edu/Roman/APT}'}
 
 def merge_nested_dicts(dict1, dict2):
     """
@@ -166,7 +168,7 @@ def create_catalog(metadata=None, catalog_name=None, bandpasses=['F087'],
 
 def parse_filename(filename):
     """
-    Try program / pass / visit / ... information out of the filename.
+    Pry program / pass / visit / ... information out of the filename.
 
     Parameters
     ----------
@@ -271,3 +273,81 @@ def simulate_image_file(args, metadata, cat, rng=None, persist=None):
     af.tree = {'roman': im, 'romanisim': romanisimdict}
 
     af.write_to(open(args.filename, 'wb'))
+
+
+def parse_apt_file(filename, csv_exposures):
+    """
+    Pry program / pass / visit / ... information out of the apt file.
+
+    Parameters
+    ----------
+    filename : str
+        filename of apt to parse
+
+    csv_exposures : int
+        number of exposures specified by the apt file
+
+    Returns
+    -------
+    list of name prefixes
+    """
+
+    # format is:
+    # r + PPPPPCCAAASSSOOOVVV_ggsaa_eeee_DET_suffix.asdf
+    # PPPPP = program
+    # CC = execution plan number
+    # AAA = pass number
+    # SSS = segment number
+    # OOO = observation number
+    # VVV = visit number
+    # gg = group identifier
+    # s = sequence identifier
+    # aa = activity identifier
+    # eeee = exposure number
+    # rPPPPPCCAAASSSOOOVVV_ggsaa_eeee
+
+    # Parse the xml
+    apt_tree = defusedxml.ElementTree.parse(filename)
+    program = apt_tree.find('.//{*}ProgramID', namespaces=NMAP).text \
+        if apt_tree.find('.//{*}ProgramID', namespaces=NMAP).text else 1
+
+    execution_plan = 1
+    pass_plan_tree = apt_tree.find('.//{*}PassPlans', namespaces=NMAP)
+    pass_plans = pass_plan_tree.findall('.//{*}PassPlan', namespaces=NMAP) \
+        if pass_plan_tree.findall('.//{*}PassPlan', namespaces=NMAP) else [-1]
+
+    total_apt_exposures = 0
+    for exp in apt_tree.findall('.//{*}NumberOfExposures', namespaces=NMAP):
+        total_apt_exposures += int(exp.text)
+
+    # Account for extra exposures due to dither pattern
+    if csv_exposures > total_apt_exposures:
+        avg_visits = int(csv_exposures / total_apt_exposures)
+    else:
+        avg_visits = 1
+
+    name_prefix_lst = []
+
+    # Set defaults
+    segment = 1
+    group = 1
+    sequence = 1
+    activity = 1
+
+    for pn in pass_plans:
+        pass_number = pn.get('Number')
+        obs_num = 0      
+
+        for on in pn.findall('.//{*}Observation', namespaces=NMAP):
+            obs_num += 1
+            exp_num = int(on.find('.//{*}NumberOfExposures', namespaces=NMAP).text)
+
+            for vn in range(avg_visits):
+                for en in range(exp_num):
+                    # Can make the Programmatic observation identifier list now
+                    prefix = "r" + str(program).zfill(5) + str(execution_plan).zfill(2) + str(pass_number).zfill(3) +\
+                             str(segment).zfill(3) + str(obs_num).zfill(3) + str(vn+1).zfill(3) + "_" +\
+                             str(group).zfill(2) + str(sequence).zfill(1) + str(activity).zfill(2) + "_" + str(en+1).zfill(4)
+                    name_prefix_lst.append(prefix)
+
+    return name_prefix_lst

scripts/romanisim-make-catalog

@@ -0,0 +1,65 @@
+#!/usr/bin/env python -u
+
+"""
+Gaia catalog generator for roman simulator.
+
+"""
+
+import argparse
+import astroquery
+from astroquery.gaia import Gaia
+from astropy.time import Time
+from romanisim import gaia
+import romanisim.parameters, romanisim.bandpass
+
+def main():
+    """Main function of catalog generator"""
+
+    parser = argparse.ArgumentParser(description='Command and options for estimating fit uncertainties')
+
+    # Command line Argument
+    parser.add_argument('radecrad',
+                        type=float,
+                        metavar='178.389741 -2.678655 1.0',
+                        nargs=3,
+                        help=f"The right ascension, declination, and radius of Gaia objects to catalog."
+                        )
+
+    # Optional output filename
+    parser.add_argument('-o', '--output',
+                        type=str,
+                        metavar='gaia-178-2-2027-06-01',
+                        help=f"Name prefix of gaia catalog file."
+                        )
+
+    # Optional Time
+    parser.add_argument('-t', '--time',
+                        type=str,
+                        metavar='2027-06-01',
+                        default='2027-06-01',
+                        help=f"Time of observations of sources in catalog."
+                        )
+
+    # Debug argument
+    parser.add_argument('-d', "--debug",
+                        action='store_true',
+                        help='Display verbose debug output')
+
+    # Collect the arguments
+    args = parser.parse_args()
+
+    if args.debug:
+        print("args = " + str(args))
+
+    q = f'select * from gaiadr3.gaia_source where distance({args.radecrad[0]}, {args.radecrad[1]}, ra, dec) < {args.radecrad[2]}'
+    job = Gaia.launch_job_async(q)
+    r = job.get_results()
+
+    outfile_name = f'gaia_{args.radecrad[0]:.2f}_{args.radecrad[1]:.2f}_{args.radecrad[2]:.2f}-{args.time}.ecsv'
+    if args.output:
+        outfile_name = args.output + ".ecsv"
+    gaia.gaia2romanisimcat(r, Time(args.time), fluxfields=set(romanisim.bandpass.galsim2roman_bandpass.values())).write(outfile_name, overwrite=True)
+
+# Call main if run (as opposed to imported as a module)
+if __name__ == "__main__":
+    main()

scripts/romanisim-make-stack

@@ -8,7 +8,7 @@ from astropy import units as u
 from astropy.table import Table
 from astropy.time import Time
 import galsim
-from romanisim import wcs, persistence, log, parameters
+from romanisim import wcs, persistence, log, parameters, apt
 from romanisim.parameters import default_parameters_dictionary
 from romanisim import ris_make_utils as ris
 
@@ -25,13 +25,19 @@ def main():
                         metavar='mosaic_list.csv',
                         help='Input (csv) file containing lists of observation parameters: '
                              'ra, dec, roll_angle, optical_element, date, overhead_time, '
-                             'ms_table_number')
+                             'ma_table_number')
 
     # WCS Object Catalog
     parser.add_argument('cat_file_name',
                         type=str,
                         metavar='object_table.csv',
                         help='Object catalog file for wcs matching (csv)')
+
+    # Create script making argument
+    parser.add_argument('-a', "--apt",
+                        type=str,
+                        metavar='small_dither_program.apt',
+                        help='APT file for metadata.')
 
     # Boresight
     parser.add_argument('-b', '--boresight',
@@ -64,7 +70,6 @@ def main():
                         help='Filename to output list of romanisim calls (sims.script) instead of '
                              'making simulation files (e.g. for cluster usage).')
 
-
     # Create script making argument
     parser.add_argument('-x', "--prefix",
                         type=str,
@@ -154,9 +159,13 @@ def main():
     if args.persistence:
         previous_file_name = {}
 
-    # Obtain pointing information from input file
+    # Open csv table pointing file    
     pointings = Table.read(args.pointing_file_name, comment="#", delimiter=" ")
 
+    if args.apt:
+        apt_idx = 0
+        apt_prefix = ris.parse_apt_file(args.apt, len(pointings))
+
     # Loop over pointings
     for entry_idx, entry in enumerate(pointings):
 
@@ -170,7 +179,10 @@ def main():
         # otherwise loop over all detectors
         while (sca <= parameters.NUMBER_OF_DETECTORS):
             # Create output file name
-            output_file_name = f"{args.prefix}_{entry_idx}_{entry['BANDPASS']}_WFI{sca:02}_{suffix}.asdf"
+            if args.apt:
+                output_file_name = f"{apt_prefix[apt_idx]}_WFI{sca:02}_{suffix}.asdf"
+            else:
+                output_file_name = f"{args.prefix}_{entry_idx}_{entry['BANDPASS']}_WFI{sca:02}_{suffix}.asdf"
 
             log.debug(f"output_file_name = {output_file_name}")
 
@@ -182,7 +194,7 @@ def main():
                 # Preserve relevant stack options for image lines
                 for option, value in options_dct.items():
                     if (type(value) != bool) and (option != 'sca') and (option != 'prefix') \
-                            and (option != 'cat_file_name'):
+                            and (option != 'cat_file_name') and (option != 'apt'):
                         line += f" --{option} {value}"
                     elif (option == 'sca'):
                         line += f" --{option} {sca}"
@@ -242,6 +254,8 @@ def main():
             else:
                 sca += 1
 
+        apt_idx += 1
+
         # Add time offset for the next exposure group
         time_offset += (float(entry['DURATION'])) * u.s