Enable slitless spectroscopic features (prism and grism) with Roman optical models

.gitignore

@@ -0,0 +1,3 @@
+*.pyc
+*checkpoint*
+*egg-info*

config/prism.yaml

@@ -0,0 +1,172 @@
+# Default settings for roman simulation
+# Includes creation of noisless oversampled images (including PSF)
+#  -- processing of other detector and instrument effects are still handled in the
+#     python postprocessing layer to enable things not currently in galsim.roman
+
+modules:
+
+    # Including galsim.roman in the list of modules to import will add a number of Roman-specific
+    # functions and classes that we will use here.
+    - roman_imsim
+    - galsim.roman
+
+    # We need this for one of our Eval items.  GalSim does not by default import datetime into
+    # the globals dict it uses when evaluating Eval items, so we can tell it to import it here.
+    - datetime
+
+# Define some other information about the images
+image:
+
+    # A special Image type that knows all the Roman SCA geometry, WCS, gain, etc.
+    # It also by default applies a number of detector effects, but these can be turned
+    # off if desired by setting some parameters (given below) to False.
+    type: roman_sca
+
+    wcs:
+        type: RomanWCS
+        SCA: '@image.SCA'
+        ra: { type: ObSeqData, field: ra }
+        dec: { type: ObSeqData, field: dec }
+        pa: { type: ObSeqData, field: pa }
+        # mjd: { type: ObSeqData, field: mjd }
+        mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+        max_sun_angle: 50
+        force_cvz: True
+
+    bandpass:
+        type: RomanBandpass
+        name: { type: ObSeqData, field: filter }
+
+    # When you want to have multiple images generate the same random galaxies, then
+    # you can set up multiple random number generators with different update cadences
+    # by making random_seed a list.
+    # The default behavior is just to have the random seeds for each object go in order by
+    # object number across all images, but this shows how to set it up so we use two separate
+    # cadences.
+    # The first one behaves normally, which will be used for things like noise on the image.
+    # The second one sets the initial seed for each object to repeat to the same starting value
+    # at the start of each filter.  If we were doing more than 3 total files, it would then
+    # move on to another sequence for the next 3 and so on.
+    random_seed:
+        # Used for noise and nobjects.
+        - { type: ObSeqData, field: visit }
+
+        # Used for objects.  Repeats sequence for each filter
+        # Note: Don't use $ shorthand here, since that will implicitly be evaluated once and then
+        # treated the same way as an integer (i.e. making a regular sequence starting from that
+        # value).  Using an explicit dict with an Eval type means GalSim will leave it alone and
+        # evaluate it as is for each object.
+
+
+    # We're just doing one SCA here.
+    # If you wanted to do all of them in each of three filters (given below), you could use:
+    #
+    # SCA:
+    #     type: Sequence
+    #     first: 1
+    #     last: 18
+    #     repeat: 3  # repeat each SCA num 3 times before moving on, for the 3 filters.
+    #
+    SCA: 16
+    # mjd: { type: ObSeqData, field: mjd }
+    mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+    filter: { type: ObSeqData, field: filter }
+    exptime: { type: ObSeqData, field: exptime }
+
+    # Photon shooting is way faster for chromatic objects than fft, especially when most of them
+    # are fairly faint.  The cross-over point for achromatic objects is generally of order
+    # flux=1.e6 or so (depending on the profile).  Most of our objects here are much fainter than
+    # that.  The fft rendering for chromatic is a factor of 10 or so slower still, whereas
+    # chromatic photon shooting is only slighly slower than achromatic, so the difference
+    # is even more pronounced in this case.
+    draw_method: 'phot'
+
+    # These are all by default turned on, but you can turn any of them off if desired:
+    ignore_noise: True
+    stray_light: False
+    thermal_background: False
+    reciprocity_failure: False
+    dark_current: False
+    nonlinearity: False
+    ipc: False
+    read_noise: False
+    sky_subtract: False
+
+stamp:
+    type: Roman_stamp
+    world_pos:
+        type: SkyCatWorldPos
+    exptime: { type: ObSeqData, field: exptime }
+    min_size: 1024
+    flux_cap: 1000000
+
+    photon_ops:
+        -
+            type: SlitlessSpec
+        #-
+           # type: ChargeDiff
+
+# psf:
+#     type: roman_psf
+#     # If omitted, it would figure this out automatically, because we are using the RomanSCA image
+#     # type.  But if we weren't, you'd have to tell it which SCA to build the PSF for.
+#     SCA: '@image.SCA'
+#     # n_waves defines how finely to sample the PSF profile over the bandpass.
+#     # Using 10 wavelengths usually gives decent accuracy.
+#     n_waves: 10
+
+# Define the galaxy type and positions to use
+gal:
+    type: SkyCatObj
+
+input:
+    obseq_data:
+        #file_name: /hpc/group/cosmology/repos/prism_sim/Roman_TDS_obseq_prism.fits
+        file_name: /hpc/group/cosmology/OpenUniverse2024/RomanTDS/Roman_TDS_obseq_11_6_23.fits
+        visit: 14234
+        SCA: '@image.SCA'
+    roman_psf:
+        SCA: '@image.SCA'
+        n_waves: 5
+    sky_catalog:
+        file_name: /hpc/group/cosmology/repos/roman_imsim/roman_rubin_cats_v1.1.2_faint/skyCatalog.yaml 
+        edge_pix: 512
+        # mjd: { type: ObSeqData, field: mjd }
+        mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+        exptime: { type: ObSeqData, field: exptime }          
+        obj_types: ["snana"]  # diffsky_galaxy, star
+
+output:
+
+    nfiles: 1
+    dir: RomanTDS_prism/images/truth
+    file_name:
+        type: FormattedStr
+        format: "Roman_TDS_truth_test_%s_%i_%i.fits.gz"
+        items:
+            - { type: ObSeqData, field: filter }
+            - { type: ObSeqData, field: visit }
+            - '@image.SCA'
+
+    truth:
+        dir: RomanTDS_prism/truth
+        file_name:
+            type: FormattedStr
+            format: "Roman_TDS_index_test_%s_%i_%i.txt"
+            items:
+                - { type: ObSeqData, field: filter }
+                - { type: ObSeqData, field: visit }
+                - '@image.SCA'
+        columns:
+            object_id: "@object_id"
+            ra: "$sky_pos.ra.deg"
+            dec: "$sky_pos.dec.deg"
+            x: "$image_pos.x"
+            y: "$image_pos.y"
+            realized_flux: "@realized_flux"
+            flux: "@flux"
+            mag: "@mag"
+            obj_type: "@object_type"

config/roman_cutouts.yaml

@@ -0,0 +1,181 @@
+# Default settings for roman simulation
+# Includes creation of noisless oversampled images (including PSF)
+#  -- processing of other detector and instrument effects are still handled in the
+#     python postprocessing layer to enable things not currently in galsim.roman
+
+modules:
+
+    # Including galsim.roman in the list of modules to import will add a number of Roman-specific
+    # functions and classes that we will use here.
+- roman_imsim
+- galsim.roman
+
+    # We need this for one of our Eval items.  GalSim does not by default import datetime into
+    # the globals dict it uses when evaluating Eval items, so we can tell it to import it here.
+- datetime
+
+# Define some other information about the images
+image:
+
+    # A special Image type that knows all the Roman SCA geometry, WCS, gain, etc.
+    # It also by default applies a number of detector effects, but these can be turned
+    # off if desired by setting some parameters (given below) to False.
+  type: roman_sca
+
+  wcs:
+    type: RomanWCS
+    SCA: '@image.SCA'
+    ra: {type: ObSeqData, field: ra}
+    dec: {type: ObSeqData, field: dec}
+    pa: {type: ObSeqData, field: pa}
+    mjd: &id001
+      type: Sum_float
+      items:
+      - type: ObSeqData_float
+        field: mjd
+      - 105.76
+        # mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+    max_sun_angle: 50
+    force_cvz: true
+
+  bandpass:
+    type: RomanBandpass
+    name: {type: ObSeqData, field: filter}
+
+    # When you want to have multiple images generate the same random galaxies, then
+    # you can set up multiple random number generators with different update cadences
+    # by making random_seed a list.
+    # The default behavior is just to have the random seeds for each object go in order by
+    # object number across all images, but this shows how to set it up so we use two separate
+    # cadences.
+    # The first one behaves normally, which will be used for things like noise on the image.
+    # The second one sets the initial seed for each object to repeat to the same starting value
+    # at the start of each filter.  If we were doing more than 3 total files, it would then
+    # move on to another sequence for the next 3 and so on.
+  random_seed:
+        # Used for noise and nobjects.
+  - {type: ObSeqData, field: visit}
+
+        # Used for objects.  Repeats sequence for each filter
+        # Note: Don't use $ shorthand here, since that will implicitly be evaluated once and then
+        # treated the same way as an integer (i.e. making a regular sequence starting from that
+        # value).  Using an explicit dict with an Eval type means GalSim will leave it alone and
+        # evaluate it as is for each object.
+
+
+    # We're just doing one SCA here.
+    # If you wanted to do all of them in each of three filters (given below), you could use:
+    #
+    # SCA:
+    #     type: Sequence
+    #     first: 1
+    #     last: 18
+    #     repeat: 3  # repeat each SCA num 3 times before moving on, for the 3 filters.
+    #
+  SCA: 5
+  mjd: *id001
+    # mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+  filter: {type: ObSeqData, field: filter}
+  exptime: {type: ObSeqData, field: exptime}
+
+    # Photon shooting is way faster for chromatic objects than fft, especially when most of them
+    # are fairly faint.  The cross-over point for achromatic objects is generally of order
+    # flux=1.e6 or so (depending on the profile).  Most of our objects here are much fainter than
+    # that.  The fft rendering for chromatic is a factor of 10 or so slower still, whereas
+    # chromatic photon shooting is only slighly slower than achromatic, so the difference
+    # is even more pronounced in this case.
+  draw_method: 'phot'
+
+    # These are all by default turned on, but you can turn any of them off if desired:
+  ignore_noise: true
+  stray_light: false
+  thermal_background: false
+  reciprocity_failure: false
+  dark_current: false
+  nonlinearity: false
+  ipc: false
+  read_noise: false
+  sky_subtract: false
+
+stamp:
+  type: Roman_stamp
+  world_pos:
+    type: SkyCatWorldPos
+  exptime: {type: ObSeqData, field: exptime}
+  min_size: 1024
+  flux_cap: 10000000000.0
+
+  photon_ops:
+  - type: SlitlessSpec
+  #- type: ChargeDiff
+
+# psf:
+#     type: roman_psf
+#     # If omitted, it would figure this out automatically, because we are using the RomanSCA image
+#     # type.  But if we weren't, you'd have to tell it which SCA to build the PSF for.
+#     SCA: '@image.SCA'
+#     # n_waves defines how finely to sample the PSF profile over the bandpass.
+#     # Using 10 wavelengths usually gives decent accuracy.
+#     n_waves: 10
+
+# Define the galaxy type and positions to use
+gal:
+  type: SkyCatObj
+
+input:
+  obseq_data:
+        #file_name: /hpc/group/cosmology/repos/prism_sim/Roman_TDS_obseq_prism.fits
+    file_name: /hpc/group/cosmology/OpenUniverse2024/RomanTDS/Roman_TDS_obseq_11_6_23.fits
+    visit: 55726
+    SCA: '@image.SCA'
+  roman_psf:
+    SCA: '@image.SCA'
+    n_waves: 5
+  sky_catalog:
+    file_name: /hpc/group/cosmology/repos/roman_imsim/roman_rubin_cats_v1.1.2_faint/skyCatalog.yaml
+    edge_pix: 512
+    mjd: *id001
+        # mjd: { type: Sum_float, items: [ { type: ObSeqData_float, field: mjd }, 170 ] }
+
+    exptime: {type: ObSeqData, field: exptime}
+    obj_types: ["snana"]      # diffsky_galaxy, star
+
+label:
+  lab: "1*exp_p+15d"
+
+
+output:
+
+  nfiles: 1
+  dir: /hpc/home/jr542/SNPIT_OU24_extract1d/exptime/truth
+  file_name:
+    type: FormattedStr
+    format: "Roman_TDS_truth_%s_%s_%i_%i.fits.gz"
+    items:
+    - '@label.lab'
+    - {type: ObSeqData, field: filter}
+    - {type: ObSeqData, field: visit}
+    - '@image.SCA'
+
+  truth:
+    dir: /hpc/home/jr542/SNPIT_OU24_extract1d/exptime/index
+    file_name:
+      type: FormattedStr
+      format: "Roman_TDS_index_%s_%s_%i_%i.txt"
+      items:
+      - '@label.lab'
+      - {type: ObSeqData, field: filter}
+      - {type: ObSeqData, field: visit}
+      - '@image.SCA'
+    columns:
+      object_id: "@object_id"
+      ra: "$sky_pos.ra.deg"
+      dec: "$sky_pos.dec.deg"
+      x: "$image_pos.x"
+      y: "$image_pos.y"
+      realized_flux: "@realized_flux"
+      flux: "@flux"
+      mag: "@mag"
+      obj_type: "@object_type"

det.py

@@ -0,0 +1,21 @@
+import sys
+
+import roman_imsim
+
+params = sys.argv[1]
+visit = int(sys.argv[2])
+sca = int(sys.argv[3])
+if len(sys.argv) > 4:
+    dither_from_file = sys.argv[4]
+else:
+    import warnings
+
+    warnings.warn("Ditherlist not provided. Persistence effect will not simulated.")
+    dither_from_file = None
+
+if len(sys.argv) > 5:
+    sca_path = sys.argv[5]
+else:
+    sca_path = None
+
+roman_imsim.detector_physics.modify_image(params, visit, sca, dither_from_file, sca_filepath=sca_path)