From 46eb1eea1c8ad7d1c3e2d97ad28717cb43359947 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 7 Aug 2026 18:25:47 -0400
Subject: [PATCH 1/9] Rename OpenUniverse2024 tutorials to consistent filenames

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 ...eview_Firefly.md => openuniverse2024_firefly_visualization.md} | 0
 ...imedomainsurvey.md => openuniverse2024_roman_tds_supernova.md} | 0
 2 files changed, 0 insertions(+), 0 deletions(-)
 rename tutorials/simulated-data/OpenUniverse2024/{OpenUniverse2024Preview_Firefly.md => openuniverse2024_firefly_visualization.md} (100%)
 rename tutorials/simulated-data/OpenUniverse2024/{openuniverse2024_roman_simulated_timedomainsurvey.md => openuniverse2024_roman_tds_supernova.md} (100%)

diff --git a/tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
similarity index 100%
rename from tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md
rename to tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
similarity index 100%
rename from tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md
rename to tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md

From 1df412adb04ee07b680f6b36fee53709a28eecc6 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 7 Aug 2026 18:28:02 -0400
Subject: [PATCH 2/9] Update OpenUniverse2024 tutorials to use the full dataset

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .gitignore                                    |   5 +
 deploy_to_fornax_manifest.in                  |   6 +-
 ignore_tutorials/ignore_circleci_testing      |   3 +-
 notebook_metadata.yml                         |  26 +-
 toc.yml                                       |  24 +-
 .../openuniverse2024_SED_fit.md               |  11 +-
 .../openuniverse2024_firefly_visualization.md | 392 +++++++++++++-----
 .../openuniverse2024_quickstart.md            |   2 +-
 ...2024_roman_demodata_20000808_instances.csv |  40 --
 ...erse2024_roman_simulated_wideareasurvey.md | 174 --------
 .../openuniverse2024_roman_tds_supernova.md   | 168 +++++---
 tutorials/simulated-data/simulated.md         |  19 +-
 12 files changed, 439 insertions(+), 431 deletions(-)
 delete mode 100644 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_demodata_20000808_instances.csv
 delete mode 100644 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md

diff --git a/.gitignore b/.gitignore
index 5fb8486d..d9bb91cf 100644
--- a/.gitignore
+++ b/.gitignore
@@ -17,6 +17,7 @@ downloaded_hdf5_files/
 *.ecsv
 *.fits
 *.gzip
+*.h5
 *.jpg
 *.jpeg
 *.log
@@ -53,3 +54,7 @@ _build/
 
 #pycache
 __pycache__
+
+# Local working notes, not part of the tutorials
+CLAUDE.md
+OU_PR_REVIEW.md
diff --git a/deploy_to_fornax_manifest.in b/deploy_to_fornax_manifest.in
index 547bfad9..57db50d2 100644
--- a/deploy_to_fornax_manifest.in
+++ b/deploy_to_fornax_manifest.in
@@ -8,11 +8,9 @@ tutorials/euclid/Euclid_ERO.md
 tutorials/euclid/euclid-cloud-access.md
 tutorials/euclid/merged-objects-hats-catalog/1-euclid-q1-hats-intro.md
 tutorials/euclid/merged-objects-hats-catalog/4-euclid-q1-hats-magnitudes.md
-tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md
+tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
 tutorials/simulated-data/cosmoDC2_TAP_access.md
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_demodata_20000808_instances.csv
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
+tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
diff --git a/ignore_tutorials/ignore_circleci_testing b/ignore_tutorials/ignore_circleci_testing
index 54a81f21..036d747f 100644
--- a/ignore_tutorials/ignore_circleci_testing
+++ b/ignore_tutorials/ignore_circleci_testing
@@ -1,8 +1,7 @@
 # ignore these as they run out of resources on CircleCI.
 # Note that we don't duplicate the listing of notebooks that are ignored for all types of testing.
 #
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey
+tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova
 tutorials/simulated-data/roman_hlss_number_density
 tutorials/euclid/euclid-cloud-access
 tutorials/techniques-and-tools/cloud-access-intro
diff --git a/notebook_metadata.yml b/notebook_metadata.yml
index bbf963af..12a8f26f 100644
--- a/notebook_metadata.yml
+++ b/notebook_metadata.yml
@@ -42,34 +42,30 @@
   file: tutorials/euclid/merged-objects-hats-catalog/4-euclid-q1-hats-magnitudes.md
   section: Euclid
   description: Review the types of flux measurements available, load template-fit and aperture magnitudes, and plot distributions and comparisons for different object types.
-- title: OpenUniverse2024 Visualization
-  file: tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md
+- title: Firefly Visualization
+  file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
   section: Simulated Data
-  description: Use Firefly to get an overview of survey structure and visualize content.
+  description: Use Firefly to explore the structure of the OpenUniverse2024 Roman and Rubin simulations and visualize their content.
 - title: CosmoDC2 Data Access
   file: tutorials/simulated-data/cosmoDC2_TAP_access.md
   section: Simulated Data
   description: Access, query, and visualize the CosmoDC2 catalog.
-- title: OpenUniverse2024 Quickstart
+- title: Quickstart
   file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
   section: Simulated Data
   description: Access OpenUniverse2024 Roman and Rubin images and catalogs.
-- title: OpenUniverse2024 TDE Light Curve
+- title: TDE Light Curve
   file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
   section: Simulated Data
-  description: Locate a simulated TDE, retrieve Roman images, and build a multi-epoch light curve.
-- title: OpenUniverse2024 SED Fitting
+  description: Locate a simulated TDE in OpenUniverse2024, retrieve Roman images, and build a multi-epoch light curve.
+- title: SED Fitting
   file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
   section: Simulated Data
-  description: Fit spectral energy distributions for supernova host galaxies using Prospector.
-- title: OpenUniverse2024 Time Domain
-  file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md
+  description: Fit spectral energy distributions for OpenUniverse2024 supernova host galaxies using Prospector.
+- title: Roman TDS Supernova
+  file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
   section: Simulated Data
-  description: Access and analyze the simulated time-domain OpenUniverse2024 survey.
-- title: OpenUniverse2024 Roman Coadds
-  file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
-  section: Simulated Data
-  description: Access OpenUniverse2024 wide-area simulated survey data.
+  description: Explore the OpenUniverse2024 Roman Time Domain Survey observation sequence and animate a simulated SN Ia across epochs.
 - title: Roman HLSS Number Density
   file: tutorials/simulated-data/roman_hlss_number_density.md
   section: Simulated Data
diff --git a/toc.yml b/toc.yml
index c214df18..969e3887 100644
--- a/toc.yml
+++ b/toc.yml
@@ -67,22 +67,22 @@ project:
     - title: Simulated Data
       file: tutorials/simulated-data/simulated.md
       children:
+        - title: OpenUniverse2024
+          children:
+            - title: Quickstart
+              file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
+            - title: Firefly Visualization
+              file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+            - title: Roman TDS Supernova
+              file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
+            - title: TDE Light Curve
+              file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
+            - title: SED Fitting
+              file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
         - title: Roman HLSS Number Density
           file: tutorials/simulated-data/roman_hlss_number_density.md
-        - title: OpenUniverse2024 Roman Coadds
-          file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
-        - title: OpenUniverse2024 Visualization
-          file: tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md
-        - title: OpenUniverse2024 Time Domain
-          file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md
         - title: CosmoDC2 Data Access
           file: tutorials/simulated-data/cosmoDC2_TAP_access.md
-        - title: OpenUniverse2024 Quickstart
-          file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
-        - title: OpenUniverse2024 TDE Light Curve
-          file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
-        - title: OpenUniverse2024 SED Fitting
-          file: tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
     - title: Techniques & Tools
       file: tutorials/techniques-and-tools/techniques.md
       children:
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
index 9bf9704e..0ac19228 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
@@ -209,9 +209,16 @@ def assemble_SN_data(sn_flux_file, galaxy_flux_file, galaxy_info_file, *, region
     sn_df = pq.read_table(sn_flux_file, filesystem=fs).to_pandas()
     print(f"  Loaded {len(sn_df)} SN entries")
 
-    # Load galaxy flux and info tables
+    # Load galaxy flux and info tables. The flux table is all fluxes and we want
+    # every band, but the info table carries a lot of shape and shear columns we
+    # never touch, so we ask Parquet for just the three we do use. That leaves
+    # about a sixth of the data to pull across and hold in memory.
     gal_flux = pq.read_table(galaxy_flux_file, filesystem=fs).to_pandas()
-    gal_info = pq.read_table(galaxy_info_file, filesystem=fs).to_pandas()
+    gal_info = pq.read_table(
+        galaxy_info_file,
+        filesystem=fs,
+        columns=["galaxy_id", "redshift", "um_source_galaxy_obs_sm"],
+    ).to_pandas()
 
     # Join host-galaxy flux and info tables on galaxy_id
     gal_joined = gal_flux.merge(gal_info, on="galaxy_id", how="inner")
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
index 52742549..22e87587 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
@@ -15,7 +15,7 @@ authors:
   - name: IRSA Team
 ---
 
-# Using Firefly to Explore OpenUniverse2024 Data Preview Simulated Roman and Rubin Images
+# Using Firefly to Explore OpenUniverse2024 Simulated Roman and Rubin Images
 
 +++
 
@@ -25,6 +25,8 @@ By the end of this tutorial, you will:
 
 - Learn how to access cloud-hosted Roman and Rubin simulated images.
 
+- Learn how to locate a Roman coadd block from a sky position, and how to find Rubin images covering that position with the IRSA Simple Image Access (SIA) service.
+
 - Learn how to launch an interactive Firefly instance inside JupyterLab.
 
 - Learn how to use the Firefly Jupyterlab extension to visualize cloud-hosted simulated images, overplot ds9 regions, overplot catalogs in Parquet format, and create 3 color images.
@@ -33,7 +35,7 @@ By the end of this tutorial, you will:
 
 ## Introduction
 
-The purpose of this tutorial is to become familiar with simulated Roman and Rubin simulated data published through the OpenUniverse 2024 data preview, and to become familiar with the Firefly JupyterLab Extension for visualizing astronomical data products.
+The purpose of this tutorial is to become familiar with the simulated Roman and Rubin data published through OpenUniverse2024, and to become familiar with the Firefly JupyterLab Extension for visualizing astronomical data products.
 
 OpenUniverse2024 is a project to simulate spatially overlapping imaging surveys to be carried out by the Nancy Grace Roman Telescope and the Vera C. Rubin Observatory. The simulations were carried out on Argonne's Theta cluster and consist of:
 
@@ -43,7 +45,7 @@ OpenUniverse2024 is a project to simulate spatially overlapping imaging surveys
 - Overlapping Roman Wide-Area Survey (WAS) in the same region
 - A deep-field calibration region of the Roman WAS in the same region
 
-This data preview release consists of a subset of data from each of the five categories above. More information about data preview can be found at [IRSA holding of this dataset](https://irsa.ipac.caltech.edu/data/theory/openuniverse2024/overview.html). 
+This tutorial works with the full simulation, which covers the entire survey footprint. More information about the dataset can be found at [IRSA's holding of this dataset](https://irsa.ipac.caltech.edu/data/theory/openuniverse2024/overview.html), and the [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) describes how the full simulation differs from the smaller preview subset that is also available.
 
 Firefly is an open-source web-based UI library for astronomical data archive access and visualization developed at Caltech and used by multiple space- and ground-based astrophysics archives. More information on Firefly can be found [here](https://github.com/Caltech-IPAC/firefly/blob/dev/README.md).
 
@@ -51,95 +53,162 @@ In addition to being used to make web applications, Firefly can be used from Pyt
 
 The Firefly JupyterLab Extension makes it particularly easy to use Firefly to efficiently visualize cloud-hosted astronomical data using JupyterLab instances running locally or on cloud. More information on Firefly JupyterLab Extension can be found [here](https://github.com/Caltech-IPAC/jupyter_firefly_extensions/blob/master/README.md).
 
+If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart) tutorial introduces the directory layout, the parquet catalogs, and the SIA image search that this notebook builds on.
+
 +++
 
 ## Imports
 
 - astropy.io.fits for accessing FITS files
 - numpy for numerical computing
-- s3fs for browsing cloud buckets
 - matplotlib.pyplot for creating static visualizations of FITS images
 - matplotlib.patches for annotating visualizations of FITS images
 - astropy.wcs for dealing with astronomical world coordinate systems
 - astropy.units for dealing with astronomical units
 - astropy.coordinates.SkyCoord for dealing with astronomical coordinates
+- astroquery.ipac.irsa.Irsa for asking IRSA which images cover a position
+- hpgeom for converting a sky position into a HEALPix region index
+- json for reading the cloud location returned by the image search
 - firefly_client.FireflyClient for using the Firefly python client
 - astropy.nddata.Cutout2D for making image cutouts
-- itertools.product to support looping over Roman blocks
 - reproject.reproject_interp to convert Roman coadds from STG to TAN projection
 - io.BytesIO for writing a fits file to an in-memory stream
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-# !pip install numpy astropy s3fs matplotlib firefly_client reproject
+# !pip install numpy astropy matplotlib firefly_client reproject astroquery hpgeom
 ```
 
 ```{code-cell} ipython3
 from astropy.io import fits
 import numpy as np
-import s3fs
 from matplotlib import pyplot as plt
 from matplotlib import patches
 from astropy import wcs
 from astropy import units as u
 from astropy.coordinates import SkyCoord
+from astroquery.ipac.irsa import Irsa
+import hpgeom
+import json
 from firefly_client import FireflyClient
 from astropy.nddata import Cutout2D
-from itertools import product
 from reproject import reproject_interp
 from io import BytesIO
 ```
 
 ## 1. Learn where the OpenUniverse2024 data are hosted in the cloud
 
-The OpenUniverse2024 data preview is hosted in the cloud via Amazon Web Services (AWS). To access these data, you need to create a client to read data from Amazon's Simple Storage Service (s3) buckets, and you need to know some information about those buckets. The OpenUniverse2024 data preview contains simulations of the Roman Wide-Area Survey (WAS) and the Roman Time Domain Survey (TDS). In this tutorial, we will focus on the WAS.
+The OpenUniverse2024 data are hosted in the cloud via Amazon Web Services (AWS). To access these data, you need to create a client to read data from Amazon's Simple Storage Service (s3) buckets, and you need to know some information about those buckets. OpenUniverse2024 contains simulations of the Roman Wide-Area Survey (WAS) and the Roman Time Domain Survey (TDS). In this tutorial, we will focus on the WAS.
+
+The two telescopes need two different approaches. The Roman coadds sit on a fixed grid of sky positions, so once we know the grid we can work out a file path ourselves. The Rubin images are individual visits scattered across the sky, so instead we ask IRSA's Simple Image Access (SIA) service which ones cover the position we care about. We set up both routes here.
 
 ```{code-cell} ipython3
 BUCKET_NAME = "nasa-irsa-simulations"
-ROMAN_PREFIX = "openuniverse2024/roman/preview"
-ROMAN_COADD_PATH = f"{ROMAN_PREFIX}/RomanWAS/images/coadds"
-TRUTH_FILES_PATH = f"{ROMAN_PREFIX}/roman_rubin_cats_v1.1.2_faint"
+OU_PREFIX = "openuniverse2024"
+ROMAN_COADD_PATH = f"{OU_PREFIX}/roman/full/RomanWAS/images/coadd"
+TRUTH_FILES_PATH = f"{OU_PREFIX}/roman/full/roman_rubin_cats_v1.1.2_faint"
+```
+
+```{code-cell} ipython3
+# Point the astroquery IRSA client at the simulated-data services, which are
+# separate from the ones serving IRSA's observed data.
+Irsa.sia_url = "https://irsa.ipac.caltech.edu/simulated/SIA"
+Irsa.tap_url = "https://irsa.ipac.caltech.edu/simulated/TAP"
+
+OU_RUBIN_SIA_COLLECTION = "simulated_rubin_openuniverse2024"
 
-RUBIN_PREFIX = "openuniverse2024/rubin/preview"
-RUBIN_COADD_PATH = f"{RUBIN_PREFIX}/u/descdm/preview_data_step3_2877_19_w_2024_12/20240403T150003Z/deepCoadd_calexp/2877/19"
+# A small radius is all we need, since we only want images containing a given point.
+SEARCH_RADIUS = 1 * u.arcsec
 ```
 
 ## 2. Roman Coadds
 
-The Nancy Grace Roman Space Telescope will carry out a wide-area survey (WAS) in the near infrared. The OpenUniverse2024 data preview includes coadded mosaics of simulated WAS data, created with the IMCOM algorithm (Rowe et al. 2011). Bands include F184, H158, J129, K213, Y106. In this section, we define some functions that make it convenient to retrieve a given cloud-hosted simulated Roman coadd based on position and filter.
+The Nancy Grace Roman Space Telescope will carry out a wide-area survey (WAS) in the near infrared. OpenUniverse2024 includes coadded mosaics of simulated WAS data, created with the IMCOM algorithm (Rowe et al. 2011). Bands include F184, H158, J129, K213, Y106. In this section, we define some functions that make it convenient to retrieve a given cloud-hosted simulated Roman coadd based on position and filter.
 
 +++
 
-### Define Roman Simulated "Blocks"
+### Describe the grid of Roman simulated "blocks"
 
-The simulated Roman coadds are arranged in 100 arcsecond blocks, as described in Hirata et al. 2024. Below we define the RA and Dec of the block centers. This cell should not be altered.
+The simulated Roman coadds are arranged in blocks, as described in Hirata et al. 2024. Rather than listing the position of every block, we can describe the entire grid with a handful of numbers, because all 1296 blocks share one projection centered on the survey and differ only in where that center falls within each block. These values are recorded in the header of every coadd file, and we check them against a real header once we have opened one.
 
 ```{code-cell} ipython3
-#Centers of roman data preview blocks. Do not alter.
-ra_block_centers = np.array([9.76330352298415, 9.724522605135252, 9.68574158906671,
-                        9.646960496603766, 9.608179349571955, 9.56939816979703,
-                        9.530616979104877, 9.491835799321422, 9.453054652272561,
-                        9.414273559784032, 9.375492543681393, 9.336711625789874]) * u.deg
-dec_block_centers = np.array([-44.252584927082495, -44.22480733304182, -44.197029724175756,
-                            -44.16925210374898, -44.14147447502621, -44.11369684127218,
-                            -44.08591920575162, -44.05814157172923, -44.03036394246976,
-                            -44.0025863212379, -43.974808711298394, -43.94703111591591]) * u.deg
-block_size = 100 * u.arcsec # each block is 100 arcsec across
+SURVEY_CRVAL = (9.55, -44.1)         # deg, the projection center shared by every block
+PIXEL_SCALE = 1.0850694444444e-05    # deg/pixel
+BLOCK_NPIX = 2688                    # pixels along each side of a block
+BLOCK_STEP = 2560                    # pixel offset between the centers of adjacent blocks
+N_BLOCKS = 36                        # the survey is a 36 x 36 grid of blocks
+REF_BLOCK = 18                       # the block whose CRPIX sits at the survey center
+REF_CRPIX = 64.5                     # that block's CRPIX, in the 1-based FITS convention
 ```
 
-### Define a function that returns the simulated Roman block column/row for a given RA/Dec.
-Blocks are arranged in a grid with axes parallel to RA and Dec. A given RA represents a column of blocks and a given Dec represents a row of blocks. This function returns the block column or row given the block centers from the cell above and a given RA or Dec, respectively.
+Each block is wider than the spacing between blocks, so neighbors overlap slightly rather than butting up against each other.
 
 ```{code-cell} ipython3
-def get_block_axis(block_centers, coord, ra_or_dec):
-    ra_or_dec_coord = getattr(coord, ra_or_dec)
-    block_dist_array = np.absolute(block_centers - ra_or_dec_coord)
-    closest_block_idx = block_dist_array.argmin()
-    if (ra_or_dec_coord < block_centers.min()-block_size/2 
-        or ra_or_dec_coord > block_centers.max()+block_size/2):
-        raise ValueError(f"Chosen {ra_or_dec}: {ra_or_dec_coord} not covered by OpenUniverse 2024 data preview simulated Roman coadds")
-    else:
-        return closest_block_idx + 12 # preview covers central 12 rows 12 columns, in a grid of 36x36 blocks
+block_size = (BLOCK_NPIX * PIXEL_SCALE * u.deg).to(u.arcsec)
+block_spacing = (BLOCK_STEP * PIXEL_SCALE * u.deg).to(u.arcsec)
+
+print(f"Each block is {block_size:.1f} across, laid down every {block_spacing:.1f}")
+```
+
+### Define a function that returns the WCS of any block
+
+Because the blocks differ only in their CRPIX, we can write down the WCS of any block without opening its file.
+
+```{code-cell} ipython3
+def block_wcs(col, row):
+    """
+    Build the WCS of the Roman WAS coadd block in a given column and row.
+
+    Parameters
+    ----------
+    col, row : int
+        Block indices, each running from 0 to 35.
+
+    Returns
+    -------
+    astropy.wcs.WCS
+        The WCS of that block.
+    """
+    w = wcs.WCS(naxis=2)
+    w.wcs.ctype = ["RA---STG", "DEC--STG"]
+    w.wcs.crval = list(SURVEY_CRVAL)
+    w.wcs.cdelt = [-PIXEL_SCALE, PIXEL_SCALE]
+    w.wcs.crpix = [REF_CRPIX - (col - REF_BLOCK) * BLOCK_STEP,
+                   REF_CRPIX - (row - REF_BLOCK) * BLOCK_STEP]
+    w.array_shape = (BLOCK_NPIX, BLOCK_NPIX)
+    return w
+```
+
+### Define a function that returns the block containing a given sky position
+
+To find the block holding a position, we express that position in the pixel grid of the reference block, then count how many block widths it lands away from that block's center.
+
+```{code-cell} ipython3
+def find_block(coord):
+    """
+    Find the Roman WAS coadd block that contains a sky position.
+
+    Parameters
+    ----------
+    coord : astropy.coordinates.SkyCoord
+        The position of interest.
+
+    Returns
+    -------
+    tuple of int
+        The (column, row) of the block containing that position.
+    """
+    x, y = block_wcs(REF_BLOCK, REF_BLOCK).world_to_pixel(coord)
+    block_center = (BLOCK_NPIX - 1) / 2
+
+    col = REF_BLOCK + int(np.round((x - block_center) / BLOCK_STEP))
+    row = REF_BLOCK + int(np.round((y - block_center) / BLOCK_STEP))
+
+    if not (0 <= col < N_BLOCKS and 0 <= row < N_BLOCKS):
+        raise ValueError(f"{coord.to_string('hmsdms')} is not covered by the "
+                         "OpenUniverse2024 simulated Roman coadds")
+
+    return col, row
 ```
 
 ### Define a function that retrieves a Roman simulated coadd given a sky position and filter.
@@ -150,18 +219,19 @@ Each of the cloud-hosted simulated Roman coadds can be accessed via a S3 filepat
 
 ```{code-cell} ipython3
 def get_roman_coadd_fpath(coord, filter):
-    col = get_block_axis(ra_block_centers, coord, 'ra')
-    row = get_block_axis(dec_block_centers, coord, 'dec')
-    
-    # Construct the coadd filename from the chosen filter, row, and column.
-    coadd_fname_root = f"prod_{filter[0]}_{col}_{row}_map.fits"
-    coadd_fpath = f"{BUCKET_NAME}/{ROMAN_COADD_PATH}/{filter}/Row{row}/{coadd_fname_root}"
+    col, row = find_block(coord)
+
+    # Construct the coadd filename from the chosen filter, column, and row.
+    coadd_fname_root = f"prod_{filter[0]}_{col:02d}_{row:02d}_map.fits.gz"
+    coadd_fpath = f"{BUCKET_NAME}/{ROMAN_COADD_PATH}/{filter}/{coadd_fname_root}"
     return coadd_fpath
 ```
 
-Now we use this access path and prefix it with `s3://` and use astropy.fits to extract a subset of it (for more info see [this section](https://docs.astropy.org/en/stable/io/fits/usage/cloud.html#subsetting-fits-files-hosted-in-amazon-s3-cloud-storage) of astropy docs).
+Now we prefix that path with `s3://` and open it with astropy.
+
+We use [`.section`](https://docs.astropy.org/en/stable/io/fits/usage/cloud.html#subsetting-fits-files-hosted-in-amazon-s3-cloud-storage) to pull out just the science plane as a 2D `numpy.array`, so the full 15-plane cube is never assembled in memory, and we take the WCS from the fits header as an `astropy.wcs.WCS` object. The function below returns a dictionary of both, along with the header itself so we can look at it later.
 
-We use `.section` to retrieve just the science image data from fits HDU as a 2D `numpy.array`, and extract WCS information from fits header as `astropy.wcs.WCS` object. The following function returns a dictionary of both.
+A note on speed: `.section` can also cut down how much of a cloud-hosted file travels across the network, but only when the file can be read out of order. These coadds are gzipped, and a compressed stream has to be read from the beginning, so astropy works through the whole file to reach the plane we want. Expect each of these reads to take twenty to thirty seconds.
 
 ```{code-cell} ipython3
 def get_roman_coadd(coord, filter):
@@ -171,23 +241,27 @@ def get_roman_coadd(coord, filter):
 
     with fits.open(coadd_s3_uri, fsspec_kwargs={"anon": True}) as hdul:
         # retrieve science data from coadd fits
-        coadd_data = hdul[0].section[0,0, :, :]  # has (2688, 2688, 15, 1) shape, with 0th layer in the cube as science image
+        coadd_data = hdul[0].section[0, 0, :, :]  # has (1, 15, 2688, 2688) shape, with 0th layer in the cube as science image
 
         # make wcs using header
         coadd_wcs = wcs.WCS(hdul[0].header, naxis=2)
 
-        return {'data': coadd_data, 'wcs': coadd_wcs}
+        return {'data': coadd_data, 'wcs': coadd_wcs, 'header': hdul[0].header}
 ```
 
 ### Inspect a simulated Roman Coadd
 
 +++
 
-Choose a filter and position that lies within the data preview region
+Choose a filter and a position within the survey footprint
 
 ```{code-cell} ipython3
 coord = SkyCoord(ra=9.6055383, dec=-44.1895542, unit="deg")
-filter_roman = 'H158' #F184, H158, J129, K213, and Y106 are available in the data preview
+filter_roman = 'H158' #F184, H158, J129, K213, and Y106 are available
+```
+
+```{code-cell} ipython3
+print("This position falls in block (column, row):", find_block(coord))
 ```
 
 Retrieve the data and header information from the simulated Roman coadd corresponding to the chosen position and filter.
@@ -196,6 +270,15 @@ Retrieve the data and header information from the simulated Roman coadd correspo
 coadd_roman = get_roman_coadd(coord, filter_roman)
 ```
 
+With a real header in hand, we can confirm that the grid description above reproduces it exactly.
+
+```{code-cell} ipython3
+predicted = block_wcs(*find_block(coord))
+
+print("CRPIX from the file:", coadd_roman['header']['CRPIX1'], coadd_roman['header']['CRPIX2'])
+print("CRPIX from our grid:", *predicted.wcs.crpix)
+```
+
 ### Understand the size of a simulated Roman coadd.
 
 ```{code-cell} ipython3
@@ -218,6 +301,54 @@ The field of view of Roman coadd is ~100 arcsec.
 
 +++
 
+### Look at everything else the coadd file contains
+
+So far we have read a single plane out of the primary HDU, but each coadd file carries a good deal more than the science image. Listing the extensions shows the full picture. Finding each extension means walking past the one before it, so this cell has to work through the entire compressed file and takes a few seconds.
+
+```{code-cell} ipython3
+with fits.open(f"s3://{get_roman_coadd_fpath(coord, filter_roman)}",
+               fsspec_kwargs={"anon": True}) as hdul:
+    hdul.info()
+```
+
+The Primary HDU for the coadded image is a cube with 15 layers, i.e., its shape is 1x15x2688x2688. The layers are as follows:
+
+0 = simulated "Science" image (Roman+Rubin simulation, units of e/(0.11 arcsec)^2/exposure)
+
+1 = lab noise: based on dark frames from the April 2023 test, masked at 3 e/p/s. Units: e/(0.11 arcsec)^2/s
+
+2 = GalSim stars, on HEALPix resolution 14 grid, normalized to total flux of 1
+
+3 = noisy stars, on HEALPix resolution 14 grid, normalized to total flux of 2.4e5 e with self-Poisson noise, including 86 e^2/input pixel background variance
+
+4 = stars, on HEALPix resolution 14 grid, total flux 1, but on in only one of the passes (to test transient response)
+
+5 = stars, on HEALPix resolution 14 grid, with total flux that varies by 5% from center to edge of the focal plane (to test what happens when the filter bandpass varies; 5% is highly exaggerated)
+
+6 = GalSim extended objects, on HEALPix resolution 14 grid, right now exponential profiles. The scale radius is log-distributed between 0.125 and 0.500 arcsec, and the ellipticity (g1,g2) is uniformly distributed in the disc of radius 0.5, i.e., g1^2+g2^2<0.5^2.
+
+7,8,9 = same objects as layer 6, but with applied shear of 0.02. The shear orientations are spaced by 60° in tangent vector space, so that in the (g1,g2)-space they are spaced by 120° and can be used for finite differences. Specifically, the directions are: layer 7 -> in East-West direction (shear PA = 270°). (g1,g2) = (0.02,0) layer 8 -> in NNW-SSE direction (shear PA = 330°). (g1,g2) = (-0.02/2,0.02√3/2) layer 9 -> in NNE-SSW direction (shear PA = 30°). (g1,g2) = (-0.02/2,-0.02√3/2)
+
+10 = coadded 1/f noise map, normalized to variance per ln f of 1
+
+11,12,13,14 = coadded white noise maps, different seeds, normalized to variance of 1 in each input pixel
+
+The remaining HDUs contain additional information:
+
+CONFIG = the configuration file
+
+INDATA = the input images used, as a binary table. The columns are: obsid (int32) -> observation ID sca (int16) -> SCA (1 through 18, inclusive) ra (float64) -> right ascension of pointing center in degrees dec (float64) -> declination of pointing center in degrees pa (float64) -> position angle of pointing in degrees valid (logical) -> input science data file is valid (should be True)
+
+INWEIGHT = the mean input weights for how much each 1.25x1.25 arcsec postage stamp depends on each input exposure. The shape is 1 x Nin x 84 x 84, where Nin is the number of input images listed in INDATA. Note that each postage stamp is 32 output pixels, so 84x32=2688.
+
+If summed on axis 1, this will normally be something close to 1. Deviations of ~10% are common, due to plate scale variations and the normalization issues introduced by diffraction spikes.
+
+INWTFLAT = a reshape of INWEIGHT suitable for display as a single image in a viewer such as DS9. The contributions from the Nin input exposures are rearranged into a 1 x 84 x (N_in*84) array.
+
+FIDELITY, SIGMA, INWTSUM, EFFCOVER = maps of U_alpha/C, S_alpha, sum_i T_{alpha i}, and the effective coverage as rescaled int16 maps. See Rowe et al. (2011) for details on the definitions of these quantities. The comment in the 'UNIT' keyword indicates how to rescale these.
+
++++
+
 ### Use the WCS from the Roman simulated header to convert the specified coordinate into a pixel position.
 
 ```{code-cell} ipython3
@@ -241,78 +372,88 @@ plt.imshow(coadd_roman['data'], origin='lower',
 plt.plot(*coord_arr_idx, 'r+', markersize=15)
 ```
 
-## 3. Rubin Coadds
+## 3. Rubin Images
 
-The OpenUniverse2024 data preview includes coadded mosaics in the following filters: u, g, r, i, z, y. In this section, we define some functions that make it convenient to retrieve a given cloud-hosted simulated Roman coadd based on position and filter.
+OpenUniverse2024 includes simulated Rubin images in the following filters: u, g, r, i, z, y. These are individual visits rather than a fixed grid of mosaics, so instead of building a path we let the SIA service tell us which images cover our position. In this section, we define functions that retrieve a Rubin image for a chosen position and filter, returning data in the same structure as the functions we defined above for Roman.
 
 +++
 
-### Retrieve Rubin Coadds
+### Retrieve Rubin images
 
-The OpenUniverse2024 data preview includes only one simulated Rubin coadd per filter. Below we define functions that make it convenient to retrieve the simulated Rubin coadd corresponding to the desired filter. We return data in same structure as the functions we defined above for Roman.
+The image search returns one row per image, with the cloud location of each tucked inside a `cloud_access` JSON string. Many visits cover any given position, so we sort by observation time and take the earliest; that keeps the notebook reproducible from one run to the next.
 
 ```{code-cell} ipython3
-def get_rubin_coadd_fpath(filter): 
-    coadd_fname_root = f"deepCoadd_calexp_2877_19_{filter}_DC2_u_descdm_preview_data_step3_2877_19_w_2024_12_20240403T150003Z.fits"
-    coadd_fpath = f"{BUCKET_NAME}/{RUBIN_COADD_PATH}/{filter}/{coadd_fname_root}"
-    return coadd_fpath
+def get_rubin_image_fpath(coord, filter):
+    results = Irsa.query_sia(pos=(coord, SEARCH_RADIUS),
+                             collection=OU_RUBIN_SIA_COLLECTION)
+
+    # Band names come back like "r_57", so compare only the part before the underscore.
+    in_band = results[[str(r['energy_bandpassname']).split('_')[0] == filter
+                       for r in results]]
+    if len(in_band) == 0:
+        raise ValueError(f"No Rubin {filter}-band images cover "
+                         f"{coord.to_string('hmsdms')}")
+
+    in_band.sort('t_min')
+    cloud_info = json.loads(in_band['cloud_access'][0])['aws']
+    return f"{cloud_info['bucket_name']}/{cloud_info['key']}"
 ```
 
 ```{code-cell} ipython3
-def get_rubin_coadd(filter):
-    coadd_s3_fpath = get_rubin_coadd_fpath(filter)
+def get_rubin_image(coord, filter):
+    image_s3_fpath = get_rubin_image_fpath(coord, filter)
 
-    with fits.open(f"s3://{coadd_s3_fpath}", fsspec_kwargs={"anon": True}) as hdul:
-        # retrieve science data from coadd fits
-        coadd_data = hdul[1].section[:,:]
+    with fits.open(f"s3://{image_s3_fpath}", fsspec_kwargs={"anon": True}) as hdul:
+        # retrieve science data, which sits in the first extension
+        image_data = hdul[1].data
 
         # make wcs using header
-        coadd_wcs = wcs.WCS(hdul[1].header)
+        image_wcs = wcs.WCS(hdul[1].header)
 
-        return {'data': coadd_data, 'wcs': coadd_wcs}
+        return {'data': image_data, 'wcs': image_wcs}
 ```
 
-### Inspect a simulated Rubin Coadd
+### Inspect a simulated Rubin image
 
 +++
 
-Choose a filter and retrieve the data and header information from the simulated Rubin coadd corresponding to that filter.
+Choose a filter and retrieve the data and header information from the simulated Rubin image covering our position.
 
 ```{code-cell} ipython3
 filter_rubin = 'r'
-coadd_rubin = get_rubin_coadd(filter_rubin)
+image_rubin = get_rubin_image(coord, filter_rubin)
 ```
 
-### Understand the size of a simulated Rubin coadd.
+### Understand the size of a simulated Rubin image.
 
 ```{code-cell} ipython3
 # Number of pixels (Y, X)
-coadd_rubin['data'].shape
+image_rubin['data'].shape
 ```
 
 ```{code-cell} ipython3
 # Pixel size (scale Y, scale X) [degrees/pixel]
-coadd_rubin['wcs'].proj_plane_pixel_scales()
+image_rubin['wcs'].proj_plane_pixel_scales()
 ```
 
 ```{code-cell} ipython3
-# Coadd size (FOV Y, FOV X)
+# Image size (FOV Y, FOV X)
 [(num * size).to('arcsec') for num, size in zip(
-    coadd_rubin['data'].shape, coadd_rubin['wcs'].proj_plane_pixel_scales())]
+    image_rubin['data'].shape, image_rubin['wcs'].proj_plane_pixel_scales())]
 ```
 
-The field of view of Rubin coadd is 840 arcsec.
+A single Rubin detector covers a far larger patch of sky than one Roman coadd block, at a coarser pixel scale.
 
 +++
 
-### Use matplotlib imshow to create a static visualization of the Rubin simulated coadd and overplot the selected position.
+### Use matplotlib imshow to create a static visualization of the Rubin simulated image and overplot the selected position.
 
 ```{code-cell} ipython3
-plt.imshow(coadd_rubin['data'], origin='lower', 
-           clim=stretch_color(coadd_rubin['data'], 1)
+plt.imshow(image_rubin['data'], origin='lower', 
+           clim=stretch_color(image_rubin['data'], 1)
            )
 
-plt.plot(*coord_to_xy(coadd_rubin['wcs'], coord), 'r+', markersize=15)
+plt.plot(*coord_to_xy(image_rubin['wcs'], coord), 'r+', markersize=15)
 ```
 
 ### Define a function that returns the URL for a given S3 filepath
@@ -324,11 +465,11 @@ def https_url(s3_fpath):
     return f"https://{BUCKET_NAME}.s3.amazonaws.com/{s3_fpath_without_bucket}"
 ```
 
-Let's generate URL for the Rubin coadd we plotted above. Clicking on the returned URL will allow you to download this image locally.
+Let's generate URL for the Rubin image we plotted above. Clicking on the returned URL will allow you to download this image locally.
 
 ```{code-cell} ipython3
-coadd_s3_fpath_rubin = get_rubin_coadd_fpath(filter_rubin)
-https_url(coadd_s3_fpath_rubin)
+image_s3_fpath_rubin = get_rubin_image_fpath(coord, filter_rubin)
+https_url(image_s3_fpath_rubin)
 ```
 
 ## 4. Compare simulated Roman and Rubin cutouts for a selected position
@@ -345,7 +486,7 @@ cutout_size = 50*u.arcsec
 
 ```{code-cell} ipython3
 cutout_roman = Cutout2D(coadd_roman['data'], coord, size=cutout_size, wcs=coadd_roman['wcs'])
-cutout_rubin = Cutout2D(coadd_rubin['data'], coord, size=cutout_size, wcs=coadd_rubin['wcs'])
+cutout_rubin = Cutout2D(image_rubin['data'], coord, size=cutout_size, wcs=image_rubin['wcs'])
 ```
 
 ### Use matplotlib imshow to plot static side-by-side comparisons of the cutouts
@@ -372,7 +513,7 @@ plt.tight_layout(rect=[0, 0, 1, 0.97])
 
 ## 5. Use Firefly to interactively identify a blended source
 
-Clearly, the simulated Roman coadd has higher spatial resolution than the Rubin simulated coadd. Let's try to locate blended objects to compare in the simulated Rubin and Roman images. We will use Firefly's interactive visualization to make this task easier.
+Clearly, the simulated Roman coadd has higher spatial resolution than the simulated Rubin image. Let's try to locate blended objects to compare in the simulated Rubin and Roman images. We will use Firefly's interactive visualization to make this task easier.
 
 +++
 
@@ -391,15 +532,15 @@ fc = FireflyClient.make_client(url="https://irsa.ipac.caltech.edu/irsaviewer")
 fc.reinit_viewer() # to clean the state, if this cell ran earlier
 ```
 
-### Send the simulated Rubin coadd to Firefly using show_fits.
+### Send the simulated Rubin image to Firefly using show_fits.
 
-For displaying the FITS image of Rubin coadd in Firefly, we use [`show_fits`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits):
+For displaying the FITS image of the Rubin visit in Firefly, we use [`show_fits`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits):
 
 ```{code-cell} ipython3
-coadd_ff_id_rubin = 'rubin-coadd-filter-r'
-fc.show_fits(url=https_url(coadd_s3_fpath_rubin),
-             plot_id=coadd_ff_id_rubin,
-             Title="Rubin Coadd"
+image_ff_id_rubin = 'rubin-image-filter-r'
+fc.show_fits(url=https_url(image_s3_fpath_rubin),
+             plot_id=image_ff_id_rubin,
+             Title="Rubin Image"
              )
 ```
 
@@ -410,11 +551,25 @@ overlay a region layer on the loaded FITS images, we can use [`overlay_region_la
 
 Region data is defined in ds9 region syntax that can be found [here](https://ds9.si.edu/doc/ref/region.html).
 
+Drawing all 1296 blocks would clutter the display, so we outline only the ones in the neighborhood of our target. Their centers come straight from the WCS of each block.
+
+```{code-cell} ipython3
+def block_center(col, row):
+    """Sky position of the center of a given coadd block."""
+    return block_wcs(col, row).pixel_to_world((BLOCK_NPIX - 1) / 2,
+                                              (BLOCK_NPIX - 1) / 2)
+
+col, row = find_block(coord)
+nearby_blocks = [(c, r)
+                 for c in range(max(col - 3, 0), min(col + 4, N_BLOCKS))
+                 for r in range(max(row - 3, 0), min(row + 4, N_BLOCKS))]
+```
+
 ```{code-cell} ipython3
 # mark the roman coadd blocks as boxes
 roman_regions = [
-    f'icrs;box {ra_block_center.value}d {dec_block_center.value}d {block_size.value}" {block_size.value}" 0d'
-    for (ra_block_center, dec_block_center) in product(ra_block_centers, dec_block_centers)
+    f'icrs;box {center.ra.deg}d {center.dec.deg}d {block_size.value}" {block_size.value}" 0d'
+    for center in (block_center(c, r) for c, r in nearby_blocks)
 ]
 
 roman_regions_id = 'roman_regions'
@@ -443,7 +598,7 @@ coords_of_interest = SkyCoord('0h38m25.35s -44d00m10.1s', frame='icrs') # locate
 coords_of_interest
 ```
 
-We can now use this [astropy `SkyCoord` object](https://docs.astropy.org/en/stable/api/astropy.coordinates.SkyCoord.html#astropy.coordinates.SkyCoord) to compare our coadds.
+We can now use this [astropy `SkyCoord` object](https://docs.astropy.org/en/stable/api/astropy.coordinates.SkyCoord.html#astropy.coordinates.SkyCoord) to compare our images.
 
 +++
 
@@ -460,6 +615,7 @@ fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id)
 
 ```{code-cell} ipython3
 coadd_roman = get_roman_coadd(coords_of_interest, filter_roman)
+image_rubin = get_rubin_image(coords_of_interest, filter_rubin)
 ```
 
 ```{code-cell} ipython3
@@ -468,7 +624,7 @@ cutout_size = 20*u.arcsec
 
 ```{code-cell} ipython3
 cutout_roman = Cutout2D(coadd_roman['data'], coords_of_interest, size=cutout_size, wcs=coadd_roman['wcs'])
-cutout_rubin = Cutout2D(coadd_rubin['data'], coords_of_interest, size=cutout_size, wcs=coadd_rubin['wcs'])
+cutout_rubin = Cutout2D(image_rubin['data'], coords_of_interest, size=cutout_size, wcs=image_rubin['wcs'])
 ```
 
 ```{code-cell} ipython3
@@ -505,34 +661,40 @@ plt.tight_layout(rect=[0, 0, 1, 0.97])
 # plt.savefig("plot.pdf", bbox_inches='tight', pad_inches=0.2)
 ```
 
-## 7. Use Firefly to visualize the OpenUniverse2024 data preview catalogs
-Let's inspect the properties of sources in the Rubin coadd image. For this we will use the input truth files present in S3 bucket.
+## 7. Use Firefly to visualize the OpenUniverse2024 catalogs
+Let's inspect the properties of sources in the Rubin image. For this we will use the input truth files present in S3 bucket.
 
-The OpenUniverse2024 data preview includes the input truth files that were used to create the simulated images. These files are in Parquet and HDF5 format, and include information about the properties of galaxies, stars, and transients.
+OpenUniverse2024 includes the input truth files that were used to create the simulated images. These files are in Parquet and HDF5 format, and include information about the properties of galaxies, stars, and transients.
 
-To list the available files, we use [s3fs](https://s3fs.readthedocs.io/en/latest/index.html). We open a file-system like connection 
-to AWS S3 and then simply list (`ls`) truth files directory:
+The catalogs are split by HEALPix sky region (nside=32, RING ordering) and the region index appears in each filename, so we can convert our sky position into that index with [`hpgeom`](https://hpgeom.readthedocs.io/en/latest/) rather than hunting through the directory listing. This is the same approach the [Quickstart](openuniverse2024_quickstart) tutorial takes.
 
 ```{code-cell} ipython3
-s3 = s3fs.S3FileSystem(anon=True) # to browse s3 bucket
-s3.ls(f"{BUCKET_NAME}/{TRUTH_FILES_PATH}")
+region = hpgeom.angle_to_pixel(32, coord.ra.deg, coord.dec.deg, lonlat=True, nest=False)
+print(f"HEALPix region for our position: {region}")
 ```
 
 ```{code-cell} ipython3
 # Catalog table of star properties (in parquet format)
-pointsource_cat_path = f"{BUCKET_NAME}/{TRUTH_FILES_PATH}/pointsource_10307.parquet"
+pointsource_cat_path = f"{BUCKET_NAME}/{TRUTH_FILES_PATH}/pointsource_{region}.parquet"
 pointsource_cat_path
 ```
 
 ```{code-cell} ipython3
 # Catalog table of galaxy properties (in parquet format)
-galaxy_cat_path = f"{BUCKET_NAME}/{TRUTH_FILES_PATH}/galaxy_10307.parquet"
+galaxy_cat_path = f"{BUCKET_NAME}/{TRUTH_FILES_PATH}/galaxy_{region}.parquet"
 galaxy_cat_path
 ```
 
-### Use Firefly's show_table to overlay the catalogs on interactive image of coadd
+### Use Firefly's show_table to overlay the catalogs on the interactive image
+
+Each catalog spans a whole HEALPix region, which is much larger than the image we are looking at, so we filter the table down to the sky area the Rubin image actually covers. (Note: you can remove filters through the table UI if you wish to see the entire data)
 
-The input truth files cover a region much larger than the data preview, so we define filters on table to constraint the sources in catalog within the ra and dec bounds of preview data. (Note: you can remove filters through the table UI if you wish to see the entire data)
+```{code-cell} ipython3
+# Bound the catalogs by the footprint of the Rubin image we sent to Firefly.
+footprint = image_rubin['wcs'].calc_footprint()
+ra_bounds = (footprint[:, 0].min(), footprint[:, 0].max())
+dec_bounds = (footprint[:, 1].min(), footprint[:, 1].max())
+```
 
 +++
 
@@ -540,8 +702,8 @@ You can visualize catalogs interactively with Firefly using [`show_table`](https
 
 ```{code-cell} ipython3
 cat_filters = [
-    f'("ra" >= {ra_block_centers.min().value} AND "ra" <= {ra_block_centers.max().value})', 
-    f'("dec" >= {dec_block_centers.min().value} AND "dec" <= {dec_block_centers.max().value})'
+    f'("ra" >= {ra_bounds[0]} AND "ra" <= {ra_bounds[1]})', 
+    f'("dec" >= {dec_bounds[0]} AND "dec" <= {dec_bounds[1]})'
 ]
 cat_filters
 ```
@@ -588,18 +750,18 @@ ROMAN_RGB_FILTERS = ['H158', 'J129', 'Y106']
 RUBIN_RGB_FILTERS = ['r', 'g', 'u']
 ```
 
-We already have Rubin coadd with catalog overlaid, let's make a 3 color image to see colors of marked objects more clearly. For this we will use [`show_fits_3color`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits_3color) method:
+We already have a Rubin image with the catalog overlaid, let's make a 3 color image to see colors of marked objects more clearly. For this we will use [`show_fits_3color`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits_3color) method:
 
 ```{code-cell} ipython3
-coadd_ff_id_rubin_3color = 'rubin-coadd-3color'
+image_ff_id_rubin_3color = 'rubin-image-3color'
 threeC = [
-    dict(url=https_url(get_rubin_coadd_fpath(filter_name)),
-         Title="Rubin Coadd 3 color")
+    dict(url=https_url(get_rubin_image_fpath(coord, filter_name)),
+         Title="Rubin 3 color")
     for filter_name in RUBIN_RGB_FILTERS
 ]
 
 fc.show_fits_3color(three_color_params=threeC,
-                    plot_id=coadd_ff_id_rubin_3color)
+                    plot_id=image_ff_id_rubin_3color)
 ```
 
 ### Use Firefly's interactivity to identify a region of interest 
@@ -627,7 +789,7 @@ coadd_roman['wcs']
 
 ### Prepare Roman coadds for displaying in Firefly
 
-Roman coadds have STG projection which cannot be read by Firefly yet. Firefly can display any FITS image but it needs to read the WCS for overlaying catalogs and other interactive features. So unlike Rubin 3 color image where we directly passed URL of coadd files to Firefly, we will read Roman coadd files in Python, reproject them from STG to TAN, and write them back to FITS to pass them to Firefly.
+Roman coadds have STG projection which cannot be read by Firefly yet. Firefly can display any FITS image but it needs to read the WCS for overlaying catalogs and other interactive features. So unlike the Rubin 3 color image where we directly passed URLs to Firefly, we will read Roman coadd files in Python, reproject them from STG to TAN, and write them back to FITS to pass them to Firefly.
 
 Let's first define functions to do so:
 
@@ -660,7 +822,7 @@ def get_fits_stream(coadd_roman):
     return fits_stream
 ```
 
-Then we perform all 3 operations we mentioned above for the RGB filters of Roman:
+Then we perform all 3 operations we mentioned above for the RGB filters of Roman. This reads another coadd for each band, so the cell takes a couple of minutes:
 
 ```{code-cell} ipython3
 coadds_rgb = []
@@ -700,7 +862,7 @@ fc.show_fits_3color(three_color_params=three_color_params,
                     plot_id=coadd_ff_id_roman_3color)
 ```
 
-We can see 3 color image of Roman coadd containing the high-redshift galaxy sources. Try panning and zomming out, you can notice it spans over one block compared to the Rubin coadd which is much larger.
+We can see 3 color image of Roman coadd containing the high-redshift galaxy sources. Try panning and zomming out, you can notice it spans over one block compared to the Rubin image which is much larger.
 
 ### Use Firefly's pan/zoom/align methods to locate high redshift sources
 Now, let's pan & zoom to the region where we located high-redshift galaxy sources. Also align & lock all images being displayed by WCS. For these operations we use these 3 methods (click on them to see their documentation):
@@ -727,6 +889,6 @@ fc.set_stretch(plot_id=coadd_ff_id_roman_3color, stype='sigma', algorithm='squar
 
 ## About This Notebook
 
-**Updated:** 2024-12-19
+**Updated:** 2026-08-05
 
 **Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems.
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
index 36c30c72..327eed4b 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
@@ -36,7 +36,7 @@ This tutorial is a focused introduction to data access only. It covers the three
 2. **Parquet catalogs** — transient (SNANA), galaxy, and galaxy-flux tables, indexed by HEALPix sky region
 3. **Image search via SIA** — querying which images cover a given sky position using astroquery and the IRSA Simple Image Access service
 
-No astrophysical analysis is performed here. For science workflows that build on these access patterns, see the [TDE Light Curve](TDE_light_curve) and [SED Fitting](SED_fit) tutorials in this repository.
+No astrophysical analysis is performed here. For science workflows that build on these access patterns, see the [TDE Light Curve](openuniverse2024_TDE_light_curve) and [SED Fitting](openuniverse2024_SED_fit) tutorials in this repository.
 
 ### Instructions
 
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_demodata_20000808_instances.csv b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_demodata_20000808_instances.csv
deleted file mode 100644
index a6a3610b..00000000
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_demodata_20000808_instances.csv
+++ /dev/null
@@ -1,40 +0,0 @@
-,filter,pointing,sca
-14,R062,407,7
-19,R062,4257,10
-8,R062,4647,6
-6,R062,6182,5
-27,R062,8112,15
-31,R062,9652,17
-0,R062,11197,1
-1,R062,11582,1
-17,R062,15037,8
-23,R062,16962,11
-20,R062,18502,10
-9,R062,18507,6
-10,R062,20812,6
-37,R062,21592,18
-32,R062,21977,17
-28,R062,22357,15
-33,R062,23512,17
-2,R062,25442,1
-25,R062,27367,13
-15,R062,28512,7
-21,R062,32362,10
-11,R062,32752,6
-7,R062,34287,5
-29,R062,36217,15
-34,R062,37757,17
-3,R062,39302,1
-4,R062,39687,1
-18,R062,43142,8
-24,R062,45067,11
-22,R062,46607,10
-12,R062,46612,6
-13,R062,48917,6
-38,R062,49697,18
-35,R062,50082,17
-30,R062,50462,15
-36,R062,51617,17
-5,R062,53547,1
-26,R062,55472,13
-16,R062,56617,7
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
deleted file mode 100644
index 6e0544e2..00000000
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
+++ /dev/null
@@ -1,174 +0,0 @@
----
-jupytext:
-  text_representation:
-    extension: .md
-    format_name: myst
-    format_version: 0.13
-    jupytext_version: 1.16.2
-kernelspec:
-  display_name: science_demo
-  language: python
-  name: python3
-authors:
-  - name: Vandana Desai
-  - name: IPAC Science Platform team
----
-
-# Analyzing cloud-hosted simulated Roman coadded images
-
-+++
-
-## Learning Goals
-
-By the end of this tutorial, you will learn how to do the following:
-
-- Identify the row and column that contains a particular ra, dec coordinate.
-- Browse the bucket containing the simulated Roman coadds.
-- Identify a simulated Roman coadd by filter, row, and column.
-- Take a closer look at a simulated coadd file.
-
-## Install and Import required modules
-
-```{code-cell} ipython3
-# Uncomment the next line to install dependencies if needed.
-# !pip install astropy numpy s3fs
-```
-
-```{code-cell} ipython3
-#Import modules
-from astropy.io import fits
-import numpy as np
-import s3fs  # browse buckets
-```
-
-## Identify the row and column that contains a particular ra, dec coordinate.
-
-The full simulation covers a 1 deg x 1 deg area centered around RA = 9.5 deg and
-Dec = -44.1 deg. This region is divided into 1296 blocks (36 rows and 36 columns),
-each 100 arcsec across.
-
-The data preview presented here covers the central 20x20 arcmin, corresponding to 144 blocks (12 rows and 12 columns).
-
-```{code-cell} ipython3
-#Choose an RA, DEC of interest.
-ra = 9.5981595
-dec = -44.2026950
-
-#Centers of data preview blocks. Do not alter.
-ra_block_centers = np.array([9.76330352298415, 9.724522605135252, 9.68574158906671,
-                    9.646960496603766, 9.608179349571955, 9.56939816979703,
-                    9.530616979104877, 9.491835799321422, 9.453054652272561,
-                    9.414273559784032, 9.375492543681393, 9.336711625789874])
-dec_block_centers = np.array([-44.252584927082495, -44.22480733304182, -44.197029724175756,
-                              -44.16925210374898, -44.14147447502621, -44.11369684127218,
-                              44.08591920575162, -44.05814157172923, -44.03036394246976,
-                              -44.0025863212379, -43.974808711298394, -43.94703111591591])
-
-ra_difference_array = np.absolute(ra_block_centers-ra)
-ra_block_centers_index = ra_difference_array.argmin()
-closest_ra_center = ra_block_centers[ra_block_centers_index]
-ra_dist = 3600. * ra_difference_array[ra_block_centers_index]
-if ra_dist > 50:
-    print("Chosen ra not covered by OpenUniverse 2024 data preview simulated Roman coadds")
-else:
-    COLUMN = ra_block_centers_index + 12
-    print("COLUMN:", COLUMN)
-    print("")
-
-dec_difference_array = np.absolute(dec_block_centers-dec)
-dec_block_centers_index = dec_difference_array.argmin()
-closest_dec_center = dec_block_centers[dec_block_centers_index]
-dec_dist = 3600. * dec_difference_array[dec_block_centers_index]
-
-if dec_dist > 50:
-    print("Chosen dec not covered by OpenUniverse 2024 data preview simulated Roman coadds")
-else:
-    ROW = dec_block_centers_index + 12
-    print("ROW:", ROW)
-```
-
-## Browse the bucket containing the simulated Roman coadds.
-
-```{code-cell} ipython3
-s3 = s3fs.S3FileSystem(anon=True) # create an S3 client
-
-BUCKET_NAME = "nasa-irsa-simulations"
-ROMAN_PREFIX = "openuniverse2024/roman/preview"
-COADD_PATH = f"{ROMAN_PREFIX}/RomanWAS/images/coadds"
-
-s3.ls(f"{BUCKET_NAME}/{COADD_PATH}")
-```
-
-## Identify a Roman simulated coadd by filter, row, and column.
-
-A simulated coadd can be uniquely identified by filter, row, and column.
-
-```{code-cell} ipython3
-#Choose a filter, row, and column
-FILTER = 'H158' #Filters F184, H158, J129, K213, and Y106 are available in the data preview.
-#ROW = 12 #Rows 12-23 are available in the data preview.
-#COLUMN = 12 #Columns 12-23 are available in the data preview.
-
-#Construct the coadd filename from the chosen filter, row, and column.
-filename_root = f"prod_{FILTER[0]}_{COLUMN}_{ROW}_map.fits"
-
-#Construct the full coadd path from the chosen filter, row, and column.
-s3_uri = f"s3://{BUCKET_NAME}/{COADD_PATH}/{FILTER}/Row{ROW}/{filename_root}"
-
-#List this filename to make sure it is found.
-s3.ls(s3_uri)
-```
-
-## Take a closer look at the simulated coadd file you identified.
-
-```{code-cell} ipython3
-#Show a summary of extensions for this file.
-
-with fits.open(s3_uri, fsspec_kwargs={"anon": True}) as hdul:
-    hdul.info()
-```
-
-The Primary HDU for the coadded image is a cube with 15 layers, i.e., its shape is 1x15x2688x2688. The layers are as follows:
-
-0 = simulated "Science" image (Roman+Rubin simulation, units of e/(0.11 arcsec)^2/exposure)
-
-1 = lab noise: based on dark frames from the April 2023 test, masked at 3 e/p/s. Units: e/(0.11 arcsec)^2/s
-
-2 = GalSim stars, on HEALPix resolution 14 grid, normalized to total flux of 1
-
-3 = noisy stars, on HEALPix resolution 14 grid, normalized to total flux of 2.4e5 e with self-Poisson noise, including 86 e^2/input pixel background variance
-
-4 = stars, on HEALPix resolution 14 grid, total flux 1, but on in only one of the passes (to test transient response)
-
-5 = stars, on HEALPix resolution 14 grid, with total flux that varies by 5% from center to edge of the focal plane (to test what happens when the filter bandpass varies; 5% is highly exaggerated)
-
-6 = GalSim extended objects, on HEALPix resolution 14 grid, right now exponential profiles. The scale radius is log-distributed between 0.125 and 0.500 arcsec, and the ellipticity (g1,g2) is uniformly distributed in the disc of radius 0.5, i.e., g1^2+g2^2<0.5^2.
-
-7,8,9 = same objects as layer 6, but with applied shear of 0.02. The shear orientations are spaced by 60° in tangent vector space, so that in the (g1,g2)-space they are spaced by 120° and can be used for finite differences. Specifically, the directions are: layer 7 -> in East-West direction (shear PA = 270°). (g1,g2) = (0.02,0) layer 8 -> in NNW-SSE direction (shear PA = 330°). (g1,g2) = (-0.02/2,0.02√3/2) layer 9 -> in NNE-SSW direction (shear PA = 30°). (g1,g2) = (-0.02/2,-0.02√3/2)
-
-10 = coadded 1/f noise map, normalized to variance per ln f of 1
-
-11,12,13,14 = coadded white noise maps, different seeds, normalized to variance of 1 in each input pixel
-
-The following HDUs contain additional information:
-
-CONFIG = the configuration file
-
-INDATA = the input images used, as a binary table. The columns are: obsid (int32) -> observation ID sca (int16) -> SCA (1 through 18, inclusive) ra (float64) -> right ascension of pointing center in degrees dec (float64) -> declination of pointing center in degrees pa (float64) -> position angle of pointing in degrees valid (logical) -> input science data file is valid (should be True)
-
-INWEIGHT = the mean input weights for how much each 1.25x1.25 arcsec postage stamp depends on each input exposure. The shape is 1 x Nin x 84 x 84, where Nin is the number of input images listed in INDATA. Note that each postage stamp is 32 output pixels, so 84x32=2688.
-
-If summed on axis 1, this will normally be something close to 1. Deviations of ~10% are common, due to plate scale variations and the normalization issues introduced by diffraction spikes.
-
-INWEIGHTFLAT = a reshape of INWEIGHT suitable for display as a single image in a viewer such as DS9. The contributions from the Nin input exposures are rearranged into a 1 x 84 x (N_in*84) array.
-
-FIDELITY, SIGMA, INWTSUM, EFFCOVER = maps of U_alpha/C, S_alpha, sum_i T_{alpha i}, and the effective coverage as rescaled int16 maps. See Rowe et al. (2011) for details on the definitions of these quantities. The comment in the 'UNIT' keyword indicates how to rescale these.
-
-+++
-***
-
-## About this notebook
-
-**Updated:** 2024-06-10
-
-**Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems.
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
index 9c104976..a91adf09 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
@@ -13,7 +13,7 @@ authors:
   - name: Lauren Aldoroty
 ---
 
-# Analyzing cloud-hosted simulated Roman Time Domain Survey images
+# Following a simulated supernova through the Roman Time Domain Survey
 
 +++
 
@@ -21,10 +21,21 @@ authors:
 
 By the end of this tutorial, you will:
 
-1. learn more about the "observations" that make up the simulated Roman TDS preview.
-2. learn how to find the locations of simulated supernovae in the preview data.
-3. learn how to create aligned cutouts of simulated Roman images.
-4. learn how to make an animated gif from these cutouts.
+1. learn more about the "observations" that make up the simulated Roman Time Domain Survey (TDS).
+2. learn how to find the locations of simulated supernovae in the transient catalog.
+3. learn how to ask IRSA which simulated Roman images cover a given position and time.
+4. learn how to create aligned cutouts of simulated Roman images.
+5. learn how to make an animated gif from these cutouts.
+
++++
+
+## Introduction
+
+The Roman Time Domain Survey revisits the same patch of sky over and over, which is what makes it possible to watch a transient appear and fade. This notebook picks one simulated Type Ia supernova out of the OpenUniverse2024 transient catalog, collects every Roman image that covers it while it is bright, and stacks those images into a short movie.
+
+The survey stores its images by pointing and detector rather than by sky position, so which files contain a particular supernova is not something you can work out from the file paths. IRSA's Simple Image Access (SIA) service answers exactly that question, and we use it here to assemble the list of images to stack.
+
+If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart) tutorial introduces the directory layout, the parquet catalogs, and the image search used below.
 
 +++
 
@@ -32,23 +43,25 @@ By the end of this tutorial, you will:
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-# !pip install astropy matplotlib numpy pandas pyarrow s3fs scipy
+# !pip install astropy matplotlib numpy pandas pyarrow s3fs scipy astroquery hpgeom
 ```
 
 ```{code-cell} ipython3
 # Import modules
 import warnings
+import json
 
 import astropy.units as u
 import matplotlib.pyplot as plt
 import numpy as np
 import pandas as pd
+import hpgeom
 from astropy.coordinates import SkyCoord
 from astropy.io import fits
 from astropy.nddata import Cutout2D
 from astropy.nddata.utils import NoOverlapError
-from astropy.table import Table
 from astropy.wcs import WCS, FITSFixedWarning
+from astroquery.ipac.irsa import Irsa
 from matplotlib import animation
 from scipy.ndimage import rotate
 
@@ -61,28 +74,13 @@ s3 = s3fs.S3FileSystem(anon=True)  # create an S3 client
 warnings.simplefilter('ignore', category=FITSFixedWarning)
 ```
 
-## Define a module to get the date (mjd) of a particular pointing.
-
 ```{code-cell} ipython3
-def get_mjd(pointing,
-            obseq_path=f's3://nasa-irsa-simulations/openuniverse2024/roman/preview/RomanTDS/Roman_TDS_obseq_11_6_23.fits'):
-
-    """
-    Retrieve MJD of a given pointing.
-
-    :param pointing: Pointing ID.
-    :type pointing: int
-    :param obseq_path: Path to obseq file Roman_TDS_obseq_11_6_23.fits.
-    :type obseq_path: str, optional
-    :return: MJD of specified pointing.
-    :rtype: float
-    """
-
-    with fits.open(obseq_path, fsspec_kwargs={"anon": True}) as obs:
-        obseq = Table(obs[1].data)
-    mjd = float(obseq['date'][int(pointing)])
+# Point the astroquery IRSA client at the simulated-data services, which are
+# separate from the ones serving IRSA's observed data.
+Irsa.sia_url = "https://irsa.ipac.caltech.edu/simulated/SIA"
+Irsa.tap_url = "https://irsa.ipac.caltech.edu/simulated/TAP"
 
-    return mjd
+OU_ROMAN_SIA_COLLECTION = 'simulated_roman_openuniverse2024'
 ```
 
 ## Define a module to create an animated gif from a collection of cutouts.
@@ -137,7 +135,7 @@ def animate_stamps(stamps, savepath, no_whitespace=True,
 # Read in the (simulated) Observation Sequence File.
 
 BUCKET_NAME = 'nasa-irsa-simulations'
-ROMAN_PREFIX = 'openuniverse2024/roman/preview'
+ROMAN_PREFIX = 'openuniverse2024/roman/full'
 
 ROMAN_TDS_PATH = f'{ROMAN_PREFIX}/RomanTDS'
 FILENAME = 'Roman_TDS_obseq_11_6_23.fits'
@@ -149,7 +147,7 @@ obseq = pd.DataFrame(obseq_hdu[1].data)
 print(obseq)
 ```
 
-## What is the spatial and temporal coverage of the openuniverse2024 Roman TDS data preview?
+## What is the spatial and temporal coverage of the openuniverse2024 Roman TDS?
 
 ```{code-cell} ipython3
 # Find the ranges of RA, Dec, and date listed in the observation sequence file.
@@ -165,12 +163,17 @@ print("mjd_min, mjd_max:", mjd_min, mjd_max)
 
 ## Read in the Supernova Analysis (SNANA) file.
 
+The transient catalogs are split by HEALPix sky region (nside=32, RING ordering), with the region index in the filename. We convert the center of the Roman TDS into that index rather than guessing at the file name.
+
 ```{code-cell} ipython3
-parquet_file = f's3://{BUCKET_NAME}/{ROMAN_PREFIX}/roman_rubin_cats_v1.1.2_faint/snana_10307.parquet'
+# The Roman Time-Domain Survey is centered near the LSST ELAIS-S1 Deep Drilling Field.
+region = hpgeom.angle_to_pixel(32, 9.45, -44.02, lonlat=True, nest=False)
+
+parquet_file = f's3://{BUCKET_NAME}/{ROMAN_PREFIX}/roman_rubin_cats_v1.1.2_faint/snana_{region}.parquet'
 transients = pd.read_parquet(parquet_file, filesystem=s3)
 ```
 
-## Let's find a relatively nearby SN Ia that lies within the region of the data preview.
+## Let's find a relatively nearby SN Ia that the survey actually watched go off.
 
 ```{code-cell} ipython3
 #List the unique models in the SNANA file.
@@ -186,55 +189,95 @@ print('Number of SN1a in SNANA file: ', len(sn1a))
 ```
 
 ```{code-cell} ipython3
-# Choose the SNIa that overlap with the spatial extent of the OpenUniverse2024 Roman TDS data preview.
+# Choose the SNIa that overlap with the spatial and temporal extent of the survey.
 ra_mask = np.logical_and(sn1a['ra'] > ra_min, sn1a['ra'] < ra_max)
 dec_mask = np.logical_and(sn1a['dec'] > dec_min, sn1a['dec'] < dec_max)
 mjd_mask = np.logical_and(sn1a['start_mjd'] > mjd_min, sn1a['end_mjd'] < mjd_max)
 all_mask = np.logical_and.reduce((ra_mask,dec_mask,mjd_mask))
-preview_sn1a = sn1a[all_mask]
-print('Number of SNIa in OpenUniverse2024 data preview:', len(preview_sn1a))
+covered_sn1a = sn1a[all_mask]
+print('Number of SNIa within the survey:', len(covered_sn1a))
 ```
 
 ```{code-cell} ipython3
-# Choose the SNIa in the data preview that are nearby, at redshifts less than 0.7.
-nearby_preview_sn1a = preview_sn1a[preview_sn1a['z_CMB'] < 0.7]
-print('Number of nearby SNIa in OpenUniverse2024 data preview:', len(nearby_preview_sn1a))
+# Choose the SNIa that are nearby, at redshifts less than 0.5, so they are bright enough
+# to stand out clearly against their host galaxy.
+nearby_sn1a = covered_sn1a[covered_sn1a['z_CMB'] < 0.5]
+print('Number of nearby SNIa:', len(nearby_sn1a))
 ```
 
+A supernova is only worth animating if Roman happened to be looking at that patch of sky while it was bright. The catalog records the date each one peaks in `peak_mjd`, and the survey visits any given field in bursts rather than continuously, so we pick an object whose peak falls inside a well-visited stretch of the survey.
+
 ```{code-cell} ipython3
-# Let's choose SN 20000808.
-oid = 20000808
-chosen_object = nearby_preview_sn1a[nearby_preview_sn1a['id'] == oid]
-ra = chosen_object.get('ra')
-dec = chosen_object.get('dec')
-ra, dec = 9.619282, -44.313894
+# Let's choose SN 20131477, which peaks while its field is being visited regularly.
+oid = 20131477
+chosen_object = nearby_sn1a[nearby_sn1a['id'] == oid].iloc[0]
+
+ra, dec = chosen_object['ra'], chosen_object['dec']
+peak_mjd = chosen_object['peak_mjd']
 coord = SkyCoord(ra*u.deg, dec*u.deg)
+
+print(f"SN {oid}: RA={ra:.6f}, Dec={dec:.6f}, z={chosen_object['z_CMB']:.3f}, peaks at MJD {peak_mjd:.1f}")
 ```
 
-## Read in the auxiliary file that lists the simulated Roman TDS images covering the chosen SNIa.
+## Ask IRSA which simulated Roman images cover the chosen SNIa.
+
+We hand the position to IRSA's image search, which returns one row per image along with the time it was taken and where the file lives in the cloud.
 
 ```{code-cell} ipython3
-# The auxiliary file contains all the images this thing is in.
-# If you need to download this file, see https://irsa.ipac.caltech.edu/docs/notebooks/.
-csvfile = './openuniverse2024_roman_demodata_20000808_instances.csv'
-instances = pd.read_csv(csvfile, usecols=['filter','pointing','sca'])
-instances
+sia_results = Irsa.query_sia(pos=(coord, 1 * u.arcsec),
+                             collection=OU_ROMAN_SIA_COLLECTION)
+
+print(f"Images covering this position: {len(sia_results)}")
 ```
 
-## Create cutouts of the chosen SNIa in the band of your choice.
+That covers every band and both Roman surveys, so we narrow it down to the Time Domain Survey images in a single band.
 
 ```{code-cell} ipython3
 band = 'R062'
-instances = instances[instances['filter'] == band]
+
+is_tds = np.array(['TDS_simple_model' in str(obs_id) for obs_id in sia_results['obs_id']])
+is_band = np.char.strip(np.array(sia_results['energy_bandpassname'], dtype=str)) == band
+
+instances = sia_results[is_tds & is_band]
+instances.sort('t_min')
+
+print(f"{band} images covering SN {oid}: {len(instances)}")
+```
+
+Finally we keep only the epochs around the peak. Frames from years before the explosion would add nothing to the movie beyond download time, so we take a window that starts shortly before the supernova appears and runs until it has faded.
+
+```{code-cell} ipython3
+epoch_mjd = np.asarray(instances['t_min'], dtype=float)
+in_window = (epoch_mjd >= peak_mjd - 40) & (epoch_mjd <= peak_mjd + 80)
+
+instances = instances[in_window]
+epoch_mjd = epoch_mjd[in_window]
+
+print(f"Epochs to animate: {len(instances)}, "
+      f"spanning MJD {epoch_mjd.min():.1f} to {epoch_mjd.max():.1f}")
 ```
 
+The cloud location of each image arrives as a JSON string, which we unpack into an S3 path.
+
 ```{code-cell} ipython3
-#Make the cutouts; this will take a couple of minutes.
+def get_s3_fpath(cloud_access):
+    """Extract the S3 URI from the cloud_access JSON string in an image search result."""
+    cloud_info = json.loads(cloud_access)['aws']
+    return f"s3://{cloud_info['bucket_name']}/{cloud_info['key']}"
+
+image_paths = [get_s3_fpath(row['cloud_access']) for row in instances]
+image_paths[:3]
+```
+
+## Create cutouts of the chosen SNIa.
+
+Roman's focal plane sits at a different angle on the sky at each visit, so the same patch of sky arrives rotated differently in every image. Before the frames can be stacked into a movie we rotate each one so that north points the same way throughout. The rotation angle is recorded in the header, and detectors in every third slot are mounted flipped, which we correct for as well.
+
+```{code-cell} ipython3
+#Make the cutouts; this will take a few minutes.
 stamps = []
 mjd = []
-for i, row in enumerate(instances.itertuples()):
-    band, pointing, sca = row[1], row[2], row[3]
-    imgpath = f's3://{BUCKET_NAME}/{ROMAN_TDS_PATH}/images/simple_model/{band}/{pointing}/Roman_TDS_simple_model_{band}_{pointing}_{sca}.fits.gz'
+for imgpath, epoch in zip(image_paths, epoch_mjd):
     print(imgpath)
     with fits.open(imgpath, fsspec_kwargs={"anon": True}) as hdu:
         img = hdu[1].data
@@ -246,10 +289,10 @@ for i, row in enumerate(instances.itertuples()):
         CDmat = np.array([header['CD1_1'], header['CD1_2'],
                           header['CD2_1'], header['CD2_2']]).reshape(2,2)
 
-        orientation = hdu[0].header['ORIENTAT']
+        orientation = header['ORIENTAT']
 
         # These chips are "flipped".
-        if sca % 3 == 0:
+        if header['SCA_NUM'] % 3 == 0:
             orientation += 180
 
         # Build rotation matrix.
@@ -274,7 +317,6 @@ for i, row in enumerate(instances.itertuples()):
 
         # Rotate the image.
         rot_img = rotate(img,angle=orientation,reshape=False,cval=np.nan)
-        hdu[1].data = rot_img
 
         rot_wcs = WCS(header)
 
@@ -282,23 +324,27 @@ for i, row in enumerate(instances.itertuples()):
             # Make cutout around SN Ia location.
             cutout = Cutout2D(rot_img,coord,100,wcs=rot_wcs,mode='partial')
             stamps.append(cutout.data)
-            mjd.append(get_mjd(pointing))
+            mjd.append(epoch)
         except NoOverlapError:
             pass
+
+print(f"Collected {len(stamps)} cutouts")
 ```
 
 ## Make an animated gif out of the cutouts.
 
+Watch the supernova brighten near the middle of the sequence and fade away again.
+
 ```{code-cell} ipython3
 savepath = f'SN{oid}.gif'
 savepath
-animate_stamps(stamps, savepath, labels=mjd)
+animate_stamps(stamps, savepath, labels=[f'MJD {m:.1f}' for m in mjd])
 ```
 
 ***
 
 ## About this notebook
 
-**Updated:** 2024-06-10
+**Updated:** 2026-08-05
 
 **Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems.
diff --git a/tutorials/simulated-data/simulated.md b/tutorials/simulated-data/simulated.md
index c5ac243f..c3acc46b 100644
--- a/tutorials/simulated-data/simulated.md
+++ b/tutorials/simulated-data/simulated.md
@@ -5,13 +5,22 @@ Because this collection is heterogeneous in coverage, structure, and intended us
 Access methods are tailored to the structure and scale of each product.
 These tutorials are designed to help users get started with accessing, visualizing, and analyzing simulated datasets hosted at IRSA.
 
+## OpenUniverse2024
+
+[OpenUniverse2024](https://arxiv.org/abs/2501.05632) simulates the overlapping imaging surveys to be carried out by the Nancy Grace Roman Space Telescope and the Vera C. Rubin Observatory, covering roughly 70 square degrees of matched optical and infrared sky.
+Start with the Quickstart tutorial, which introduces the directory layout, the parquet catalogs, and the image search that the other tutorials build on.
+
 ```{notebook-gallery} notebook_metadata.yml
-tutorials/simulated-data/roman_hlss_number_density.md
-tutorials/simulated-data/cosmoDC2_TAP_access.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
-tutorials/simulated-data/OpenUniverse2024/OpenUniverse2024Preview_Firefly.md
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_wideareasurvey.md
+tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
 tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
-tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_simulated_timedomainsurvey.md
+```
+
+## Other Simulated Datasets
+
+```{notebook-gallery} notebook_metadata.yml
+tutorials/simulated-data/roman_hlss_number_density.md
+tutorials/simulated-data/cosmoDC2_TAP_access.md
 ```

From baaddeef254a908fbb5e96e7a9a066da08fa46dd Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 7 Aug 2026 17:54:45 -0400
Subject: [PATCH 3/9] Add playback controls and visit-coverage guidance to the
 supernova notebook

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../openuniverse2024_roman_tds_supernova.md      | 16 +++++++++++++++-
 1 file changed, 15 insertions(+), 1 deletion(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
index a91adf09..e2d00bc6 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
@@ -64,6 +64,7 @@ from astropy.wcs import WCS, FITSFixedWarning
 from astroquery.ipac.irsa import Irsa
 from matplotlib import animation
 from scipy.ndimage import rotate
+from IPython.display import HTML
 
 # Needed to access data in the cloud
 import s3fs
@@ -127,6 +128,11 @@ def animate_stamps(stamps, savepath, no_whitespace=True,
     writer = animation.PillowWriter()
     anim = animation.FuncAnimation(fig, animate, interval=600, frames=len(stamps))
     anim.save(savepath, writer=writer)
+
+    # Close the figure so the notebook does not also print a static copy of the first frame,
+    # and hand the animation back so it can be displayed with playback controls.
+    plt.close(fig)
+    return anim
 ```
 
 ## Read in the Observation Sequence File to learn more about the "observations" that make up the simulated Roman Time Domain Survey.
@@ -207,6 +213,8 @@ print('Number of nearby SNIa:', len(nearby_sn1a))
 
 A supernova is only worth animating if Roman happened to be looking at that patch of sky while it was bright. The catalog records the date each one peaks in `peak_mjd`, and the survey visits any given field in bursts rather than continuously, so we pick an object whose peak falls inside a well-visited stretch of the survey.
 
+The cuts above cannot check this for us: they confirm a supernova went off somewhere inside the survey's footprint and date range, not that the telescope was pointed at it at the time. If you change `oid` to explore a different object, compare the dates it was bright, `start_mjd` through `end_mjd`, against the dates Roman actually visited its position, which arrive as the `t_min` column of the image search a few cells below. Objects that pass every cut here but were only observed long before or after they erupted are common; for those, the epoch window a few cells below comes back empty and the notebook cannot build an animation at all.
+
 ```{code-cell} ipython3
 # Let's choose SN 20131477, which peaks while its field is being visited regularly.
 oid = 20131477
@@ -338,7 +346,13 @@ Watch the supernova brighten near the middle of the sequence and fade away again
 ```{code-cell} ipython3
 savepath = f'SN{oid}.gif'
 savepath
-animate_stamps(stamps, savepath, labels=[f'MJD {m:.1f}' for m in mjd])
+anim = animate_stamps(stamps, savepath, labels=[f'MJD {m:.1f}' for m in mjd])
+```
+
+The saved gif loops without stopping, which makes a brief brightening hard to follow. Displaying the animation instead gives playback controls: pause it, step one epoch at a time with the arrows, drag the slider to any frame, and read the date in the corner as you go.
+
+```{code-cell} ipython3
+HTML(anim.to_jshtml(default_mode='once'))
 ```
 
 ***

From 50c034e8485a8367a42cd887ef9f6481335e7a9a Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Tue, 11 Aug 2026 17:03:27 -0400
Subject: [PATCH 4/9] Fold helper function definitions into hidden cells

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../openuniverse2024_SED_fit.md               |  2 +-
 .../openuniverse2024_TDE_light_curve.md       |  4 ++--
 .../openuniverse2024_quickstart.md            | 19 +++++++++++++++++--
 .../openuniverse2024_roman_tds_supernova.md   | 12 ++++++++++++
 4 files changed, 32 insertions(+), 5 deletions(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
index 0ac19228..07a369fb 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_SED_fit.md
@@ -713,7 +713,7 @@ def plot_many_seds(df, rubin_bands, roman_bands, n_galaxies=10, loglog=False):
 plot_many_seds(df_scaled, rubin_bands, roman_bands, n_galaxies=50, loglog=True)
 ```
 
-+++ {"jupyter": {"source_hidden": true}}
++++
 
 ## 3. SED fitting
 In this section we will use [Prospector](https://prospect.readthedocs.io/en/stable/index.html), a Bayesian SED fitting code built on [FSPS](https://dfm.io/python-fsps/current/), to infer stellar population parameters for our galaxies by fitting their SEDs.  We choose Prospector because it is an open-source, powerful, yet flexible package to infer stellar population properties that is widely used in the community.
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
index 8a37e265..1b0e4919 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_TDE_light_curve.md
@@ -82,11 +82,11 @@ import itertools
 
 ## 1. Explore the OpenUniverse2024 data directories
 
-This section of the tutorial demonstrates how to explore the OpenUniverse2024 data directories directly on S3 and inspect simulated Roman and Rubin images without downloading large datasets locally. It establishes a connection to the public NASA IRSA simulations bucket using s3fs, defines key directory paths for the full Roman and Rubin simulations (not the preview subsets), and illustrates how to browse image files for a selected band and pointing. The accompanying functions — summarize_fits_files() and show_gallery() — provide tools for quickly summarizing FITS file metadata (e.g., number of extensions, pointing information, pixel scale) and for visualizing a small gallery of example images from the chosen directory.
+This section of the tutorial demonstrates how to explore the OpenUniverse2024 data directories directly on S3 and inspect simulated Roman and Rubin images without downloading large datasets locally. It establishes a connection to the public NASA IRSA simulations bucket using s3fs, defines key directory paths for the full Roman and Rubin simulations, and illustrates how to browse image files for a selected band and pointing. The accompanying functions — summarize_fits_files() and show_gallery() — provide tools for quickly summarizing FITS file metadata (e.g., number of extensions, pointing information, pixel scale) and for visualizing a small gallery of example images from the chosen directory.
 
 In the prefix you will see that we choose "simple_model" simulations and not "truth" simulations because the simple_model images are the ones with noise and real effects, while "Truth" are noise free, perfect images.
 
-Also in the prefix you will see that we choose the full simulation, not the preview simulation for both Roman and Rubin. Differences between the "full" and "preview" simulations are clarified in the [this](https://arxiv.org/abs/2501.05632) publication
+The prefix also selects the full simulation, which covers the entire survey footprint, for both Roman and Rubin. The [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) describes how the simulation was produced.
 
 ```{code-cell} ipython3
 # Setup
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
index 327eed4b..51f2877d 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_quickstart.md
@@ -77,7 +77,7 @@ from astropy.io import fits
 
 The OpenUniverse2024 data live on the cloud in a public AWS S3 bucket and can be accessed anonymously using `s3fs`. This section shows how to establish that connection, navigate the directory tree, and inspect the contents of a FITS image file.
 
-In the path below, `simple_model` refers to the simulated images with noise and realistic instrument effects, as opposed to `truth` images which are noise-free. The `full` simulation covers the complete survey footprint; a smaller `preview` subset is also available. See the [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) for details on the differences. A `pointing` is a unique Roman observation visit — each pointing corresponds to one placement of the 18-detector focal plane on the sky, producing up to 18 individual FITS files (one per detector).
+In the path below, `simple_model` refers to the simulated images with noise and realistic instrument effects, as opposed to `truth` images which are noise-free. The `full` simulation covers the complete survey footprint. See the [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) for details on how it was produced. A `pointing` is a unique Roman observation visit — each pointing corresponds to one placement of the 18-detector focal plane on the sky, producing up to 18 individual FITS files (one per detector).
 
 ```{code-cell} ipython3
 # Create an anonymous (public read-only) connection to the NASA IRSA S3 bucket.
@@ -130,6 +130,11 @@ Let's display a gallery of example images to get a sense of the data.
 Note this gallery can take about a minute to build.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def show_gallery(files, max_images=9):
     """
     Display a gallery of FITS images.
@@ -205,10 +210,15 @@ print("Galaxy flux file: ", gal_flux_path)
 
 ### 2.1 Inspect the SNANA Transient Catalog
 
-`inspect_parquet_columns()` reads only the Parquet metadata footer to print the row count and column names — no data is loaded into memory.
+`inspect_parquet_files()` reads only the Parquet metadata footer to print the row count and column names — no data is loaded into memory.
 We use it here for the SNANA catalog and repeat it for the galaxy info and flux catalogs below.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def inspect_parquet_files(s3_path, *, region='us-east-1'):
     """
     Print the structure of a Parquet file on S3 without reading its data.
@@ -305,6 +315,11 @@ OU_RUBIN_SIA_COLLECTION = 'simulated_rubin_openuniverse2024'
 ```
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def get_s3_fpath(cloud_access):
     """Extract the S3 URI from the cloud_access JSON string in an SIA result."""
     cloud_info = json.loads(cloud_access)
diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
index e2d00bc6..deadbb43 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_roman_tds_supernova.md
@@ -87,6 +87,11 @@ OU_ROMAN_SIA_COLLECTION = 'simulated_roman_openuniverse2024'
 ## Define a module to create an animated gif from a collection of cutouts.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def animate_stamps(stamps, savepath, no_whitespace=True,
                    labels=[],labelxy=(0.05,0.95),
                    **kwargs):
@@ -268,11 +273,18 @@ print(f"Epochs to animate: {len(instances)}, "
 The cloud location of each image arrives as a JSON string, which we unpack into an S3 path.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def get_s3_fpath(cloud_access):
     """Extract the S3 URI from the cloud_access JSON string in an image search result."""
     cloud_info = json.loads(cloud_access)['aws']
     return f"s3://{cloud_info['bucket_name']}/{cloud_info['key']}"
+```
 
+```{code-cell} ipython3
 image_paths = [get_s3_fpath(row['cloud_access']) for row in instances]
 image_paths[:3]
 ```

From 32a6893671e1fbf25bd0457bad55673548b19e81 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Tue, 11 Aug 2026 18:04:50 -0400
Subject: [PATCH 5/9] Rework the firefly notebook to use single exposures

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../openuniverse2024_firefly_visualization.md | 465 +++++++-----------
 1 file changed, 184 insertions(+), 281 deletions(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
index 22e87587..a7429f74 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
@@ -1,4 +1,8 @@
 ---
+authors:
+- name: Jaladh Singhal
+- name: Vandana Desai
+- name: IRSA Team
 jupytext:
   text_representation:
     extension: .md
@@ -6,13 +10,9 @@ jupytext:
     format_version: 0.13
     jupytext_version: 1.16.3
 kernelspec:
-  display_name: Python 3 (ipykernel)
-  language: python
   name: python3
-authors:
-  - name: Jaladh Singhal
-  - name: Vandana Desai
-  - name: IRSA Team
+  display_name: python3
+  language: python
 ---
 
 # Using Firefly to Explore OpenUniverse2024 Simulated Roman and Rubin Images
@@ -25,7 +25,7 @@ By the end of this tutorial, you will:
 
 - Learn how to access cloud-hosted Roman and Rubin simulated images.
 
-- Learn how to locate a Roman coadd block from a sky position, and how to find Rubin images covering that position with the IRSA Simple Image Access (SIA) service.
+- Learn how to find the Roman and Rubin images covering a sky position with the IRSA Simple Image Access (SIA) service.
 
 - Learn how to launch an interactive Firefly instance inside JupyterLab.
 
@@ -45,7 +45,7 @@ OpenUniverse2024 is a project to simulate spatially overlapping imaging surveys
 - Overlapping Roman Wide-Area Survey (WAS) in the same region
 - A deep-field calibration region of the Roman WAS in the same region
 
-This tutorial works with the full simulation, which covers the entire survey footprint. More information about the dataset can be found at [IRSA's holding of this dataset](https://irsa.ipac.caltech.edu/data/theory/openuniverse2024/overview.html), and the [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) describes how the full simulation differs from the smaller preview subset that is also available.
+This tutorial works with the full simulation, which covers the entire survey footprint. More information about the dataset can be found at [IRSA's holding of this dataset](https://irsa.ipac.caltech.edu/data/theory/openuniverse2024/overview.html), and the [OpenUniverse2024 paper](https://arxiv.org/abs/2501.05632) describes how the simulation was produced.
 
 Firefly is an open-source web-based UI library for astronomical data archive access and visualization developed at Caltech and used by multiple space- and ground-based astrophysics archives. More information on Firefly can be found [here](https://github.com/Caltech-IPAC/firefly/blob/dev/README.md).
 
@@ -71,12 +71,11 @@ If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart
 - json for reading the cloud location returned by the image search
 - firefly_client.FireflyClient for using the Firefly python client
 - astropy.nddata.Cutout2D for making image cutouts
-- reproject.reproject_interp to convert Roman coadds from STG to TAN projection
 - io.BytesIO for writing a fits file to an in-memory stream
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-# !pip install numpy astropy matplotlib firefly_client reproject astroquery hpgeom
+# !pip install numpy astropy matplotlib firefly_client astroquery hpgeom
 ```
 
 ```{code-cell} ipython3
@@ -92,7 +91,6 @@ import hpgeom
 import json
 from firefly_client import FireflyClient
 from astropy.nddata import Cutout2D
-from reproject import reproject_interp
 from io import BytesIO
 ```
 
@@ -100,12 +98,11 @@ from io import BytesIO
 
 The OpenUniverse2024 data are hosted in the cloud via Amazon Web Services (AWS). To access these data, you need to create a client to read data from Amazon's Simple Storage Service (s3) buckets, and you need to know some information about those buckets. OpenUniverse2024 contains simulations of the Roman Wide-Area Survey (WAS) and the Roman Time Domain Survey (TDS). In this tutorial, we will focus on the WAS.
 
-The two telescopes need two different approaches. The Roman coadds sit on a fixed grid of sky positions, so once we know the grid we can work out a file path ourselves. The Rubin images are individual visits scattered across the sky, so instead we ask IRSA's Simple Image Access (SIA) service which ones cover the position we care about. We set up both routes here.
+Both telescopes are reached the same way. Neither survey stores its images by sky position, so which file covers a given patch of sky is not something you can work out from a file path. Instead we ask IRSA's Simple Image Access (SIA) service, which answers exactly that question for Roman and Rubin alike.
 
 ```{code-cell} ipython3
 BUCKET_NAME = "nasa-irsa-simulations"
 OU_PREFIX = "openuniverse2024"
-ROMAN_COADD_PATH = f"{OU_PREFIX}/roman/full/RomanWAS/images/coadd"
 TRUTH_FILES_PATH = f"{OU_PREFIX}/roman/full/roman_rubin_cats_v1.1.2_faint"
 ```
 
@@ -115,141 +112,74 @@ TRUTH_FILES_PATH = f"{OU_PREFIX}/roman/full/roman_rubin_cats_v1.1.2_faint"
 Irsa.sia_url = "https://irsa.ipac.caltech.edu/simulated/SIA"
 Irsa.tap_url = "https://irsa.ipac.caltech.edu/simulated/TAP"
 
+OU_ROMAN_SIA_COLLECTION = "simulated_roman_openuniverse2024"
 OU_RUBIN_SIA_COLLECTION = "simulated_rubin_openuniverse2024"
 
 # A small radius is all we need, since we only want images containing a given point.
 SEARCH_RADIUS = 1 * u.arcsec
 ```
 
-## 2. Roman Coadds
-
-The Nancy Grace Roman Space Telescope will carry out a wide-area survey (WAS) in the near infrared. OpenUniverse2024 includes coadded mosaics of simulated WAS data, created with the IMCOM algorithm (Rowe et al. 2011). Bands include F184, H158, J129, K213, Y106. In this section, we define some functions that make it convenient to retrieve a given cloud-hosted simulated Roman coadd based on position and filter.
-
-+++
-
-### Describe the grid of Roman simulated "blocks"
-
-The simulated Roman coadds are arranged in blocks, as described in Hirata et al. 2024. Rather than listing the position of every block, we can describe the entire grid with a handful of numbers, because all 1296 blocks share one projection centered on the survey and differ only in where that center falls within each block. These values are recorded in the header of every coadd file, and we check them against a real header once we have opened one.
-
-```{code-cell} ipython3
-SURVEY_CRVAL = (9.55, -44.1)         # deg, the projection center shared by every block
-PIXEL_SCALE = 1.0850694444444e-05    # deg/pixel
-BLOCK_NPIX = 2688                    # pixels along each side of a block
-BLOCK_STEP = 2560                    # pixel offset between the centers of adjacent blocks
-N_BLOCKS = 36                        # the survey is a 36 x 36 grid of blocks
-REF_BLOCK = 18                       # the block whose CRPIX sits at the survey center
-REF_CRPIX = 64.5                     # that block's CRPIX, in the 1-based FITS convention
-```
-
-Each block is wider than the spacing between blocks, so neighbors overlap slightly rather than butting up against each other.
-
-```{code-cell} ipython3
-block_size = (BLOCK_NPIX * PIXEL_SCALE * u.deg).to(u.arcsec)
-block_spacing = (BLOCK_STEP * PIXEL_SCALE * u.deg).to(u.arcsec)
-
-print(f"Each block is {block_size:.1f} across, laid down every {block_spacing:.1f}")
-```
-
-### Define a function that returns the WCS of any block
+## 2. Roman Images
 
-Because the blocks differ only in their CRPIX, we can write down the WCS of any block without opening its file.
+The Nancy Grace Roman Space Telescope will carry out a wide-area survey (WAS) in the near infrared. OpenUniverse2024 simulates that survey as individual exposures, each covering one detector's worth of sky. Bands include F184, H158, J129, K213, W146, Y106.
 
-```{code-cell} ipython3
-def block_wcs(col, row):
-    """
-    Build the WCS of the Roman WAS coadd block in a given column and row.
-
-    Parameters
-    ----------
-    col, row : int
-        Block indices, each running from 0 to 35.
+The survey stores those exposures by pointing and detector rather than by sky position, so in this section we define functions that ask the image search which exposure covers a position, and then read it.
 
-    Returns
-    -------
-    astropy.wcs.WCS
-        The WCS of that block.
-    """
-    w = wcs.WCS(naxis=2)
-    w.wcs.ctype = ["RA---STG", "DEC--STG"]
-    w.wcs.crval = list(SURVEY_CRVAL)
-    w.wcs.cdelt = [-PIXEL_SCALE, PIXEL_SCALE]
-    w.wcs.crpix = [REF_CRPIX - (col - REF_BLOCK) * BLOCK_STEP,
-                   REF_CRPIX - (row - REF_BLOCK) * BLOCK_STEP]
-    w.array_shape = (BLOCK_NPIX, BLOCK_NPIX)
-    return w
-```
++++
 
-### Define a function that returns the block containing a given sky position
+### Define a function that returns the file path of a Roman image given a sky position and filter.
 
-To find the block holding a position, we express that position in the pixel grid of the reference block, then count how many block widths it lands away from that block's center.
+The Roman collection holds both the Wide-Area Survey and the Time Domain Survey, so we keep only the WAS exposures and then only the band we want. Many exposures cover any given position, so we sort by observation time and take the earliest; that keeps the notebook reproducible from one run to the next.
 
 ```{code-cell} ipython3
-def find_block(coord):
-    """
-    Find the Roman WAS coadd block that contains a sky position.
-
-    Parameters
-    ----------
-    coord : astropy.coordinates.SkyCoord
-        The position of interest.
-
-    Returns
-    -------
-    tuple of int
-        The (column, row) of the block containing that position.
-    """
-    x, y = block_wcs(REF_BLOCK, REF_BLOCK).world_to_pixel(coord)
-    block_center = (BLOCK_NPIX - 1) / 2
-
-    col = REF_BLOCK + int(np.round((x - block_center) / BLOCK_STEP))
-    row = REF_BLOCK + int(np.round((y - block_center) / BLOCK_STEP))
-
-    if not (0 <= col < N_BLOCKS and 0 <= row < N_BLOCKS):
-        raise ValueError(f"{coord.to_string('hmsdms')} is not covered by the "
-                         "OpenUniverse2024 simulated Roman coadds")
-
-    return col, row
-```
-
-### Define a function that retrieves a Roman simulated coadd given a sky position and filter.
-
-+++
-
-Each of the cloud-hosted simulated Roman coadds can be accessed via a S3 filepath. This function returns the access path for the simulated Roman coadd that includes a specified position on the sky and desired filter.
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
+def get_roman_image_fpath(coord, filter):
+    results = Irsa.query_sia(pos=(coord, SEARCH_RADIUS),
+                             collection=OU_ROMAN_SIA_COLLECTION)
 
-```{code-cell} ipython3
-def get_roman_coadd_fpath(coord, filter):
-    col, row = find_block(coord)
+    in_was = results[['WAS_simple_model' in str(r['obs_id']) for r in results]]
+    in_band = in_was[[str(r['energy_bandpassname']).strip() == filter
+                      for r in in_was]]
+    if len(in_band) == 0:
+        raise ValueError(f"No Roman WAS {filter}-band images cover "
+                         f"{coord.to_string('hmsdms')}")
 
-    # Construct the coadd filename from the chosen filter, column, and row.
-    coadd_fname_root = f"prod_{filter[0]}_{col:02d}_{row:02d}_map.fits.gz"
-    coadd_fpath = f"{BUCKET_NAME}/{ROMAN_COADD_PATH}/{filter}/{coadd_fname_root}"
-    return coadd_fpath
+    in_band.sort('t_min')
+    cloud_info = json.loads(in_band['cloud_access'][0])['aws']
+    return f"{cloud_info['bucket_name']}/{cloud_info['key']}"
 ```
 
-Now we prefix that path with `s3://` and open it with astropy.
-
-We use [`.section`](https://docs.astropy.org/en/stable/io/fits/usage/cloud.html#subsetting-fits-files-hosted-in-amazon-s3-cloud-storage) to pull out just the science plane as a 2D `numpy.array`, so the full 15-plane cube is never assembled in memory, and we take the WCS from the fits header as an `astropy.wcs.WCS` object. The function below returns a dictionary of both, along with the header itself so we can look at it later.
+### Define a function that retrieves a Roman simulated image given a sky position and filter.
 
-A note on speed: `.section` can also cut down how much of a cloud-hosted file travels across the network, but only when the file can be read out of order. These coadds are gzipped, and a compressed stream has to be read from the beginning, so astropy works through the whole file to reach the plane we want. Expect each of these reads to take twenty to thirty seconds.
+Now we prefix that path with `s3://` and open it with astropy. The science image sits in the first extension, and we take the WCS from the same header as an `astropy.wcs.WCS` object. The function below returns a dictionary of both, along with the header itself so we can look at it later.
 
 ```{code-cell} ipython3
-def get_roman_coadd(coord, filter):
-    # retrive fits file of block/tile from the coadd mosiac
-    coadd_s3_fpath = get_roman_coadd_fpath(coord, filter)
-    coadd_s3_uri = f"s3://{coadd_s3_fpath}"
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
+def get_roman_image(coord, filter):
+    image_s3_fpath = get_roman_image_fpath(coord, filter)
 
-    with fits.open(coadd_s3_uri, fsspec_kwargs={"anon": True}) as hdul:
-        # retrieve science data from coadd fits
-        coadd_data = hdul[0].section[0, 0, :, :]  # has (1, 15, 2688, 2688) shape, with 0th layer in the cube as science image
+    with fits.open(f"s3://{image_s3_fpath}", fsspec_kwargs={"anon": True}) as hdul:
+        # retrieve science data, which sits in the first extension
+        image_data = hdul[1].data
 
         # make wcs using header
-        coadd_wcs = wcs.WCS(hdul[0].header, naxis=2)
+        image_wcs = wcs.WCS(hdul[1].header)
 
-        return {'data': coadd_data, 'wcs': coadd_wcs, 'header': hdul[0].header}
+        # hand back the path too, so later cells can reuse it instead of
+        # asking the image search for it a second time
+        return {'data': image_data, 'wcs': image_wcs, 'header': hdul[1].header,
+                'fpath': image_s3_fpath}
 ```
 
-### Inspect a simulated Roman Coadd
+### Inspect a simulated Roman image
 
 +++
 
@@ -257,116 +187,89 @@ Choose a filter and a position within the survey footprint
 
 ```{code-cell} ipython3
 coord = SkyCoord(ra=9.6055383, dec=-44.1895542, unit="deg")
-filter_roman = 'H158' #F184, H158, J129, K213, and Y106 are available
+filter_roman = 'H158' #F184, H158, J129, K213, W146, and Y106 are available
 ```
 
-```{code-cell} ipython3
-print("This position falls in block (column, row):", find_block(coord))
-```
-
-Retrieve the data and header information from the simulated Roman coadd corresponding to the chosen position and filter.
+Retrieve the data and header information from the simulated Roman image corresponding to the chosen position and filter.
 
 ```{code-cell} ipython3
-coadd_roman = get_roman_coadd(coord, filter_roman)
+image_roman = get_roman_image(coord, filter_roman)
 ```
 
-With a real header in hand, we can confirm that the grid description above reproduces it exactly.
-
-```{code-cell} ipython3
-predicted = block_wcs(*find_block(coord))
-
-print("CRPIX from the file:", coadd_roman['header']['CRPIX1'], coadd_roman['header']['CRPIX2'])
-print("CRPIX from our grid:", *predicted.wcs.crpix)
-```
-
-### Understand the size of a simulated Roman coadd.
+### Understand the size of a simulated Roman image.
 
 ```{code-cell} ipython3
 # Number of pixels (Y, X)
-coadd_roman['data'].shape
+image_roman['data'].shape
 ```
 
 ```{code-cell} ipython3
 # Pixel size (scale Y, scale X) [degrees/pixel]
-coadd_roman['wcs'].proj_plane_pixel_scales()
+image_roman['wcs'].proj_plane_pixel_scales()
 ```
 
 ```{code-cell} ipython3
-# Coadd size (FOV Y, FOV X)
+# Image size (FOV Y, FOV X)
 [(num * size).to('arcsec') for num, size in zip(
-    coadd_roman['data'].shape, coadd_roman['wcs'].proj_plane_pixel_scales())]
+    image_roman['data'].shape, image_roman['wcs'].proj_plane_pixel_scales())]
 ```
 
-The field of view of Roman coadd is ~100 arcsec.
+The field of view of a Roman exposure is ~450 arcsec, at a pixel scale of ~0.11 arcsec.
 
 +++
 
-### Look at everything else the coadd file contains
+### Look at everything else the image file contains
 
-So far we have read a single plane out of the primary HDU, but each coadd file carries a good deal more than the science image. Listing the extensions shows the full picture. Finding each extension means walking past the one before it, so this cell has to work through the entire compressed file and takes a few seconds.
+So far we have read the science image, but each exposure carries more than that. Listing the extensions shows the full picture.
 
 ```{code-cell} ipython3
-with fits.open(f"s3://{get_roman_coadd_fpath(coord, filter_roman)}",
+with fits.open(f"s3://{image_roman['fpath']}",
                fsspec_kwargs={"anon": True}) as hdul:
     hdul.info()
 ```
 
-The Primary HDU for the coadded image is a cube with 15 layers, i.e., its shape is 1x15x2688x2688. The layers are as follows:
-
-0 = simulated "Science" image (Roman+Rubin simulation, units of e/(0.11 arcsec)^2/exposure)
-
-1 = lab noise: based on dark frames from the April 2023 test, masked at 3 e/p/s. Units: e/(0.11 arcsec)^2/s
-
-2 = GalSim stars, on HEALPix resolution 14 grid, normalized to total flux of 1
-
-3 = noisy stars, on HEALPix resolution 14 grid, normalized to total flux of 2.4e5 e with self-Poisson noise, including 86 e^2/input pixel background variance
-
-4 = stars, on HEALPix resolution 14 grid, total flux 1, but on in only one of the passes (to test transient response)
-
-5 = stars, on HEALPix resolution 14 grid, with total flux that varies by 5% from center to edge of the focal plane (to test what happens when the filter bandpass varies; 5% is highly exaggerated)
-
-6 = GalSim extended objects, on HEALPix resolution 14 grid, right now exponential profiles. The scale radius is log-distributed between 0.125 and 0.500 arcsec, and the ellipticity (g1,g2) is uniformly distributed in the disc of radius 0.5, i.e., g1^2+g2^2<0.5^2.
-
-7,8,9 = same objects as layer 6, but with applied shear of 0.02. The shear orientations are spaced by 60° in tangent vector space, so that in the (g1,g2)-space they are spaced by 120° and can be used for finite differences. Specifically, the directions are: layer 7 -> in East-West direction (shear PA = 270°). (g1,g2) = (0.02,0) layer 8 -> in NNW-SSE direction (shear PA = 330°). (g1,g2) = (-0.02/2,0.02√3/2) layer 9 -> in NNE-SSW direction (shear PA = 30°). (g1,g2) = (-0.02/2,-0.02√3/2)
-
-10 = coadded 1/f noise map, normalized to variance per ln f of 1
-
-11,12,13,14 = coadded white noise maps, different seeds, normalized to variance of 1 in each input pixel
-
-The remaining HDUs contain additional information:
-
-CONFIG = the configuration file
+Each file holds a primary header carrying the observation metadata but no pixels, followed by three 4088x4088 planes:
 
-INDATA = the input images used, as a binary table. The columns are: obsid (int32) -> observation ID sca (int16) -> SCA (1 through 18, inclusive) ra (float64) -> right ascension of pointing center in degrees dec (float64) -> declination of pointing center in degrees pa (float64) -> position angle of pointing in degrees valid (logical) -> input science data file is valid (should be True)
+SCI = the simulated science image, in electrons per second
 
-INWEIGHT = the mean input weights for how much each 1.25x1.25 arcsec postage stamp depends on each input exposure. The shape is 1 x Nin x 84 x 84, where Nin is the number of input images listed in INDATA. Note that each postage stamp is 32 output pixels, so 84x32=2688.
+ERR = the estimated uncertainty on each science pixel
 
-If summed on axis 1, this will normally be something close to 1. Deviations of ~10% are common, due to plate scale variations and the normalization issues introduced by diffraction spikes.
-
-INWTFLAT = a reshape of INWEIGHT suitable for display as a single image in a viewer such as DS9. The contributions from the Nin input exposures are rearranged into a 1 x 84 x (N_in*84) array.
-
-FIDELITY, SIGMA, INWTSUM, EFFCOVER = maps of U_alpha/C, S_alpha, sum_i T_{alpha i}, and the effective coverage as rescaled int16 maps. See Rowe et al. (2011) for details on the definitions of these quantities. The comment in the 'UNIT' keyword indicates how to rescale these.
+DQ = the data quality mask, flagging pixels that should not be trusted
 
 +++
 
 ### Use the WCS from the Roman simulated header to convert the specified coordinate into a pixel position.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def coord_to_xy(w, coord):
     return w.world_to_array_index(coord)[::-1] #reverse since 0th axis is y, 1st axis is x
+```
 
-coord_arr_idx = coord_to_xy(coadd_roman['wcs'], coord)
+```{code-cell} ipython3
+coord_arr_idx = coord_to_xy(image_roman['wcs'], coord)
 coord_arr_idx
 ```
 
-### Use matplotlib imshow to create a static visualization of the Roman simulated coadd and overplot the selected position.
+### Use matplotlib imshow to create a static visualization of the Roman simulated image and overplot the selected position.
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def stretch_color(data, clipPercent):
     return np.percentile(data, (0 + clipPercent, 100 - clipPercent))
+```
 
-plt.imshow(coadd_roman['data'], origin='lower', 
-           clim=stretch_color(coadd_roman['data'], 1)
+```{code-cell} ipython3
+plt.imshow(image_roman['data'], origin='lower', 
+           clim=stretch_color(image_roman['data'], 1)
            )
 
 plt.plot(*coord_arr_idx, 'r+', markersize=15)
@@ -374,7 +277,7 @@ plt.plot(*coord_arr_idx, 'r+', markersize=15)
 
 ## 3. Rubin Images
 
-OpenUniverse2024 includes simulated Rubin images in the following filters: u, g, r, i, z, y. These are individual visits rather than a fixed grid of mosaics, so instead of building a path we let the SIA service tell us which images cover our position. In this section, we define functions that retrieve a Rubin image for a chosen position and filter, returning data in the same structure as the functions we defined above for Roman.
+OpenUniverse2024 includes simulated Rubin images in the following filters: u, g, r, i, z, y. As with Roman, these are stored by visit and detector rather than by sky position, so we again let the SIA service tell us which images cover our position. In this section, we define functions that retrieve a Rubin image for a chosen position and filter, returning data in the same structure as the functions we defined above for Roman.
 
 +++
 
@@ -383,23 +286,43 @@ OpenUniverse2024 includes simulated Rubin images in the following filters: u, g,
 The image search returns one row per image, with the cloud location of each tucked inside a `cloud_access` JSON string. Many visits cover any given position, so we sort by observation time and take the earliest; that keeps the notebook reproducible from one run to the next.
 
 ```{code-cell} ipython3
-def get_rubin_image_fpath(coord, filter):
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
+def get_rubin_image_fpaths(coord, filters):
+    # One search covers every band, so ask once and sort the results out by band
+    # afterwards rather than searching again for each one.
     results = Irsa.query_sia(pos=(coord, SEARCH_RADIUS),
                              collection=OU_RUBIN_SIA_COLLECTION)
 
-    # Band names come back like "r_57", so compare only the part before the underscore.
-    in_band = results[[str(r['energy_bandpassname']).split('_')[0] == filter
-                       for r in results]]
-    if len(in_band) == 0:
-        raise ValueError(f"No Rubin {filter}-band images cover "
-                         f"{coord.to_string('hmsdms')}")
+    fpaths = {}
+    for filter in filters:
+        # Band names come back like "r_57", so compare only the part before the underscore.
+        in_band = results[[str(r['energy_bandpassname']).split('_')[0] == filter
+                           for r in results]]
+        if len(in_band) == 0:
+            raise ValueError(f"No Rubin {filter}-band images cover "
+                             f"{coord.to_string('hmsdms')}")
 
-    in_band.sort('t_min')
-    cloud_info = json.loads(in_band['cloud_access'][0])['aws']
-    return f"{cloud_info['bucket_name']}/{cloud_info['key']}"
+        in_band.sort('t_min')
+        cloud_info = json.loads(in_band['cloud_access'][0])['aws']
+        fpaths[filter] = f"{cloud_info['bucket_name']}/{cloud_info['key']}"
+
+    return fpaths
+
+
+def get_rubin_image_fpath(coord, filter):
+    return get_rubin_image_fpaths(coord, [filter])[filter]
 ```
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def get_rubin_image(coord, filter):
     image_s3_fpath = get_rubin_image_fpath(coord, filter)
 
@@ -410,7 +333,9 @@ def get_rubin_image(coord, filter):
         # make wcs using header
         image_wcs = wcs.WCS(hdul[1].header)
 
-        return {'data': image_data, 'wcs': image_wcs}
+        # hand back the path too, so later cells can reuse it instead of
+        # asking the image search for it a second time
+        return {'data': image_data, 'wcs': image_wcs, 'fpath': image_s3_fpath}
 ```
 
 ### Inspect a simulated Rubin image
@@ -442,7 +367,7 @@ image_rubin['wcs'].proj_plane_pixel_scales()
     image_rubin['data'].shape, image_rubin['wcs'].proj_plane_pixel_scales())]
 ```
 
-A single Rubin detector covers a far larger patch of sky than one Roman coadd block, at a coarser pixel scale.
+A single Rubin detector covers a somewhat larger patch of sky than one Roman exposure, at roughly twice the pixel scale.
 
 +++
 
@@ -460,6 +385,11 @@ plt.plot(*coord_to_xy(image_rubin['wcs'], coord), 'r+', markersize=15)
 Since the OpenUniverse2024 data is available through a public S3 bucket, we can access a given S3 file using HTTPS URL as follows:
 
 ```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
 def https_url(s3_fpath):
     s3_fpath_without_bucket = s3_fpath.split('/', 1)[1]
     return f"https://{BUCKET_NAME}.s3.amazonaws.com/{s3_fpath_without_bucket}"
@@ -468,7 +398,7 @@ def https_url(s3_fpath):
 Let's generate URL for the Rubin image we plotted above. Clicking on the returned URL will allow you to download this image locally.
 
 ```{code-cell} ipython3
-image_s3_fpath_rubin = get_rubin_image_fpath(coord, filter_rubin)
+image_s3_fpath_rubin = image_rubin['fpath']
 https_url(image_s3_fpath_rubin)
 ```
 
@@ -485,7 +415,7 @@ cutout_size = 50*u.arcsec
 ### Create the cutouts
 
 ```{code-cell} ipython3
-cutout_roman = Cutout2D(coadd_roman['data'], coord, size=cutout_size, wcs=coadd_roman['wcs'])
+cutout_roman = Cutout2D(image_roman['data'], coord, size=cutout_size, wcs=image_roman['wcs'])
 cutout_rubin = Cutout2D(image_rubin['data'], coord, size=cutout_size, wcs=image_rubin['wcs'])
 ```
 
@@ -513,7 +443,7 @@ plt.tight_layout(rect=[0, 0, 1, 0.97])
 
 ## 5. Use Firefly to interactively identify a blended source
 
-Clearly, the simulated Roman coadd has higher spatial resolution than the simulated Rubin image. Let's try to locate blended objects to compare in the simulated Rubin and Roman images. We will use Firefly's interactive visualization to make this task easier.
+The simulated Roman exposure has finer pixels and a sharper point spread function than the simulated Rubin image. Let's try to locate blended objects to compare in the simulated Rubin and Roman images. We will use Firefly's interactive visualization to make this task easier.
 
 +++
 
@@ -544,37 +474,28 @@ fc.show_fits(url=https_url(image_s3_fpath_rubin),
              )
 ```
 
-### Use ds9 region syntax to overplot the simulated Roman image blocks on the interactive display
+### Use ds9 region syntax to outline the Roman exposure on the interactive display
 
 The Firefly client includes several methods related to controlling ds9 region overlays. To 
 overlay a region layer on the loaded FITS images, we can use [`overlay_region_layer`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.overlay_region_layer).
 
 Region data is defined in ds9 region syntax that can be found [here](https://ds9.si.edu/doc/ref/region.html).
 
-Drawing all 1296 blocks would clutter the display, so we outline only the ones in the neighborhood of our target. Their centers come straight from the WCS of each block.
+The Rubin visit on display covers more sky than the Roman exposure does, so it is worth seeing where the Roman data actually falls. The corners of the Roman exposure come straight from its WCS, and we draw them as a polygon.
 
 ```{code-cell} ipython3
-def block_center(col, row):
-    """Sky position of the center of a given coadd block."""
-    return block_wcs(col, row).pixel_to_world((BLOCK_NPIX - 1) / 2,
-                                              (BLOCK_NPIX - 1) / 2)
-
-col, row = find_block(coord)
-nearby_blocks = [(c, r)
-                 for c in range(max(col - 3, 0), min(col + 4, N_BLOCKS))
-                 for r in range(max(row - 3, 0), min(row + 4, N_BLOCKS))]
+roman_corners = image_roman['wcs'].calc_footprint()
+roman_corners
 ```
 
 ```{code-cell} ipython3
-# mark the roman coadd blocks as boxes
-roman_regions = [
-    f'icrs;box {center.ra.deg}d {center.dec.deg}d {block_size.value}" {block_size.value}" 0d'
-    for center in (block_center(c, r) for c, r in nearby_blocks)
-]
+# outline the Roman exposure as a polygon
+corner_coords = ' '.join(f'{ra}d {dec}d' for ra, dec in roman_corners)
+roman_regions = [f'icrs;polygon {corner_coords} # text={{Roman {filter_roman}}}']
 
 roman_regions_id = 'roman_regions'
 fc.overlay_region_layer(region_data=roman_regions,
-                        title='Roman Mosiac', 
+                        title='Roman exposure', 
                         region_layer_id=roman_regions_id)
 ```
 
@@ -594,7 +515,7 @@ You can view the coordinates of your mouse pointer at the bottom left of the dis
 We have provided an example. You can change this based on your interests.
 
 ```{code-cell} ipython3
-coords_of_interest = SkyCoord('0h38m25.35s -44d00m10.1s', frame='icrs') # located and copied through UI
+coords_of_interest = SkyCoord('0h38m40.26s -44d04m22.26s', frame='icrs') # located and copied through UI
 coords_of_interest
 ```
 
@@ -613,17 +534,21 @@ fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id)
 
 ## 6. Plot cutouts of the identified blended source
 
+Two catalog sources sit 1.22 arcsec apart at this position. Roman separates them; Rubin merges them into one extended blob. Note that the Rubin image is a single visit, so faint detail washes out into the noise.
+
+The two cutouts are plotted in their own pixel frames, and because each telescope was pointed independently, north falls about 30 degrees apart between the panels. We leave them unrotated here, so expect the pair to appear tilted relative to each other rather than aligned.
+
 ```{code-cell} ipython3
-coadd_roman = get_roman_coadd(coords_of_interest, filter_roman)
+image_roman = get_roman_image(coords_of_interest, filter_roman)
 image_rubin = get_rubin_image(coords_of_interest, filter_rubin)
 ```
 
 ```{code-cell} ipython3
-cutout_size = 20*u.arcsec
+cutout_size = 8*u.arcsec
 ```
 
 ```{code-cell} ipython3
-cutout_roman = Cutout2D(coadd_roman['data'], coords_of_interest, size=cutout_size, wcs=coadd_roman['wcs'])
+cutout_roman = Cutout2D(image_roman['data'], coords_of_interest, size=cutout_size, wcs=image_roman['wcs'])
 cutout_rubin = Cutout2D(image_rubin['data'], coords_of_interest, size=cutout_size, wcs=image_rubin['wcs'])
 ```
 
@@ -635,13 +560,12 @@ axs[0].imshow(cutout_roman.data, origin='lower',
               )
 axs[0].set_title(f"ROMAN in filter {filter_roman}")
 
-# Let's also encircle the blended source we identified 
+# Let's also encircle the blended source we identified. The circle is 2 arcsec across
+# on both panels, so the two cutouts can be compared directly despite their different
+# pixel scales.
+roman_scale = cutout_roman.wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
 axs[0].add_patch(patches.Circle(coord_to_xy(cutout_roman.wcs, coords_of_interest), 
-                                radius=50, color='r', fill=False, linewidth=2))
-# and a bonus blended source that is close to it
-other_coords = SkyCoord(coords_of_interest.ra-9.2*u.arcsec, coords_of_interest.dec+5*u.arcsec)
-axs[0].add_patch(patches.Circle(coord_to_xy(cutout_roman.wcs, other_coords),
-                                radius=36, color='cyan', fill=False, linewidth=2))
+                                radius=2/roman_scale, color='r', fill=False, linewidth=2))
 
 
 axs[1].imshow(cutout_rubin.data, origin='lower', 
@@ -650,10 +574,9 @@ axs[1].imshow(cutout_rubin.data, origin='lower',
 axs[1].set_title(f"RUBIN in filter {filter_rubin}")
 
 # Let's also encircle the source we identified 
+rubin_scale = cutout_rubin.wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
 axs[1].add_patch(patches.Circle(coord_to_xy(cutout_rubin.wcs, coords_of_interest),
-                                radius=10, color='r', fill=False, linewidth=2))
-axs[1].add_patch(patches.Circle(coord_to_xy(cutout_rubin.wcs, other_coords),
-                                radius=8, color='cyan', fill=False, linewidth=2))
+                                radius=2/rubin_scale, color='r', fill=False, linewidth=2))
 
 
 fig.suptitle(f"Cutouts at ({coords_of_interest.ra:6f}, {coords_of_interest.dec:6f}) with {cutout_size} size", fontsize=14);
@@ -696,8 +619,6 @@ ra_bounds = (footprint[:, 0].min(), footprint[:, 0].max())
 dec_bounds = (footprint[:, 1].min(), footprint[:, 1].max())
 ```
 
-+++
-
 You can visualize catalogs interactively with Firefly using [`show_table`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_table). This capability can take many parameters. Here we will simply send our catalog to Firefly so that we can (a) see an interactive table; (b) see this table plotted over the image that we've already sent; and (c) use the GUI to quickly create exploratory plots. See if you can use the GUI to quickly determine approximately how many galaxies cover the Rubin image and what the redshift distribution of these galaxies is.
 
 ```{code-cell} ipython3
@@ -754,8 +675,9 @@ We already have a Rubin image with the catalog overlaid, let's make a 3 color im
 
 ```{code-cell} ipython3
 image_ff_id_rubin_3color = 'rubin-image-3color'
+rubin_rgb_fpaths = get_rubin_image_fpaths(coord, RUBIN_RGB_FILTERS)
 threeC = [
-    dict(url=https_url(get_rubin_image_fpath(coord, filter_name)),
+    dict(url=https_url(rubin_rgb_fpaths[filter_name]),
          Title="Rubin 3 color")
     for filter_name in RUBIN_RGB_FILTERS
 ]
@@ -780,39 +702,26 @@ point_region = f'icrs;point {high_z_gal_coords.ra.value}d {high_z_gal_coords.dec
 fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id)
 ```
 
-## 8. Plot 3-color Roman coadd containing your region of interest
-Let's inspect WCS of Roman coadd first
+## 8. Plot 3-color Roman image containing your region of interest
+Let's inspect WCS of the Roman image first
 
 ```{code-cell} ipython3
-coadd_roman['wcs']
+image_roman['wcs']
 ```
 
-### Prepare Roman coadds for displaying in Firefly
-
-Roman coadds have STG projection which cannot be read by Firefly yet. Firefly can display any FITS image but it needs to read the WCS for overlaying catalogs and other interactive features. So unlike the Rubin 3 color image where we directly passed URLs to Firefly, we will read Roman coadd files in Python, reproject them from STG to TAN, and write them back to FITS to pass them to Firefly.
+### Prepare Roman images for displaying in Firefly
 
-Let's first define functions to do so:
+The Roman exposures are gzipped on S3, so rather than hand Firefly a URL as we did for Rubin, we read each one in Python and pass it up as an in-memory FITS file.
 
 ```{code-cell} ipython3
-def reproject_to_TAN(coadd_roman):
-    # Define a new WCS with TAN projection (in CTYPE key)
-    output_wcs = coadd_roman['wcs'].deepcopy()
-    output_wcs.wcs.ctype = [ctype.replace('STG', 'TAN') for ctype in coadd_roman['wcs'].wcs.ctype]
-    
-    # Use reproject to convert a given data and wcs, to a desired wcs and shape
-    reprojected_data, _ = reproject_interp(
-        (coadd_roman['data'], coadd_roman['wcs']),
-        output_projection=output_wcs,
-        shape_out=coadd_roman['data'].shape
-    )
-
-    return {'data': reprojected_data, 'wcs': output_wcs}
-```
-
-```{code-cell} ipython3
-def get_fits_stream(coadd_roman):
-    # Create a FITS PrimaryHDU object with the coadd data
-    hdu = fits.PrimaryHDU(data=coadd_roman['data'], header=coadd_roman['wcs'].to_header())
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
+def get_fits_stream(image_roman):
+    # Create a FITS PrimaryHDU object with the image data
+    hdu = fits.PrimaryHDU(data=image_roman['data'], header=image_roman['wcs'].to_header())
 
     # Write the HDU to the in-memory stream (to save I/O time)
     fits_stream = BytesIO()
@@ -822,29 +731,23 @@ def get_fits_stream(coadd_roman):
     return fits_stream
 ```
 
-Then we perform all 3 operations we mentioned above for the RGB filters of Roman. This reads another coadd for each band, so the cell takes a couple of minutes:
+Then we read one exposure per RGB filter and turn each into a stream:
 
 ```{code-cell} ipython3
-coadds_rgb = []
-coadds_rgb_reprojected = []
-coadds_rgb_fits_stream = []
+images_rgb = []
+images_rgb_fits_stream = []
 
 for filter_name in ROMAN_RGB_FILTERS:
     print(f'\nFILTER: {filter_name}')
-    print('Retrieving Roman coadd...')
-    coadd_roman = get_roman_coadd(high_z_gal_coords, filter_name)
-    coadds_rgb.append(coadd_roman)
-
-    print('Reprojecting to TAN...')
-    coadd_roman_reprojected = reproject_to_TAN(coadd_roman)
-    coadds_rgb_reprojected.append(coadd_roman_reprojected)
+    print('Retrieving Roman image...')
+    image_roman = get_roman_image(high_z_gal_coords, filter_name)
+    images_rgb.append(image_roman)
 
-    print('Writing back to fits stream...')
-    coadd_roman_fits_stream = get_fits_stream(coadd_roman_reprojected)
-    coadds_rgb_fits_stream.append(coadd_roman_fits_stream)
+    print('Writing to fits stream...')
+    images_rgb_fits_stream.append(get_fits_stream(image_roman))
 ```
 
-### Use Firefly's show_fits_3color to create a 3 color image of Roman coadds
+### Use Firefly's show_fits_3color to create a 3 color image of Roman images
 
 Now we upload each fits stream (in-memory fits file) to firefly using [`upload_fits_data()`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.upload_fits_data) and prepare color params to pass to the [`show_fits_3color()`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits_3color).
 
@@ -852,17 +755,17 @@ Now we upload each fits stream (in-memory fits file) to firefly using [`upload_f
 three_color_params = [
     {
         'file': fc.upload_fits_data(fits_stream),
-        'Title': "Roman Coadd 3 color"
-    } for fits_stream in coadds_rgb_fits_stream]
+        'Title': "Roman 3 color"
+    } for fits_stream in images_rgb_fits_stream]
 ```
 
 ```{code-cell} ipython3
-coadd_ff_id_roman_3color = 'roman-coadd-3color-high_z_gal'
+image_ff_id_roman_3color = 'roman-3color-high_z_gal'
 fc.show_fits_3color(three_color_params=three_color_params,
-                    plot_id=coadd_ff_id_roman_3color)
+                    plot_id=image_ff_id_roman_3color)
 ```
 
-We can see 3 color image of Roman coadd containing the high-redshift galaxy sources. Try panning and zomming out, you can notice it spans over one block compared to the Rubin image which is much larger.
+We can see a 3 color image of the Roman exposures containing the high-redshift galaxy sources. Try panning and zooming out; you can notice it covers a smaller patch of sky than the Rubin image.
 
 ### Use Firefly's pan/zoom/align methods to locate high redshift sources
 Now, let's pan & zoom to the region where we located high-redshift galaxy sources. Also align & lock all images being displayed by WCS. For these operations we use these 3 methods (click on them to see their documentation):
@@ -871,8 +774,8 @@ Now, let's pan & zoom to the region where we located high-redshift galaxy source
 - [`align_images`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.align_images)
 
 ```{code-cell} ipython3
-fc.set_pan(plot_id=coadd_ff_id_roman_3color, x=high_z_gal_coords.ra.deg, y=high_z_gal_coords.dec.deg, coord='j2000')
-fc.set_zoom(plot_id=coadd_ff_id_roman_3color, factor=1)
+fc.set_pan(plot_id=image_ff_id_roman_3color, x=high_z_gal_coords.ra.deg, y=high_z_gal_coords.dec.deg, coord='j2000')
+fc.set_zoom(plot_id=image_ff_id_roman_3color, factor=1)
 fc.align_images(lock_match=True)
 ```
 
@@ -881,7 +784,7 @@ fc.align_images(lock_match=True)
 The image has a lot of noise that obscures our high redshift sources of interest. You can use the Firefly GUI to change the stretch of the image display. We identify that squared stretch from -2 to 10 sigma highlights the colors of our sources better. You can also use the Firefly client's [`set_stretch`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.set_stretch) to do this via Python. This is helpful for reproducibility and for scaling up to many images.
 
 ```{code-cell} ipython3
-fc.set_stretch(plot_id=coadd_ff_id_roman_3color, stype='sigma', algorithm='squared', 
+fc.set_stretch(plot_id=image_ff_id_roman_3color, stype='sigma', algorithm='squared', 
                band='ALL', lower_value=-2, upper_value=10)
 ```
 

From 6bedba5a20cc628f6c06e20f8213dfcf7443bd50 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 14 Aug 2026 20:59:48 -0400
Subject: [PATCH 6/9] Tighten the firefly notebook text and fix the
 high-redshift target

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../openuniverse2024_firefly_visualization.md | 400 +++++++++++++-----
 1 file changed, 295 insertions(+), 105 deletions(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
index a7429f74..122f804b 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
@@ -8,7 +8,7 @@ jupytext:
     extension: .md
     format_name: myst
     format_version: 0.13
-    jupytext_version: 1.16.3
+    jupytext_version: 1.19.5
 kernelspec:
   name: python3
   display_name: python3
@@ -71,11 +71,12 @@ If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart
 - json for reading the cloud location returned by the image search
 - firefly_client.FireflyClient for using the Firefly python client
 - astropy.nddata.Cutout2D for making image cutouts
+- reproject.reproject_interp for rotating cutouts so both telescopes' panels point the same way
 - io.BytesIO for writing a fits file to an in-memory stream
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-# !pip install numpy astropy matplotlib firefly_client astroquery hpgeom
+!pip install numpy astropy matplotlib firefly_client astroquery hpgeom reproject
 ```
 
 ```{code-cell} ipython3
@@ -91,6 +92,7 @@ import hpgeom
 import json
 from firefly_client import FireflyClient
 from astropy.nddata import Cutout2D
+from reproject import reproject_interp
 from io import BytesIO
 ```
 
@@ -98,7 +100,7 @@ from io import BytesIO
 
 The OpenUniverse2024 data are hosted in the cloud via Amazon Web Services (AWS). To access these data, you need to create a client to read data from Amazon's Simple Storage Service (s3) buckets, and you need to know some information about those buckets. OpenUniverse2024 contains simulations of the Roman Wide-Area Survey (WAS) and the Roman Time Domain Survey (TDS). In this tutorial, we will focus on the WAS.
 
-Both telescopes are reached the same way. Neither survey stores its images by sky position, so which file covers a given patch of sky is not something you can work out from a file path. Instead we ask IRSA's Simple Image Access (SIA) service, which answers exactly that question for Roman and Rubin alike.
+We ask IRSA's Simple Image Access (SIA) service to locate both the Roman and Rubin images by coordinate.
 
 ```{code-cell} ipython3
 BUCKET_NAME = "nasa-irsa-simulations"
@@ -123,13 +125,15 @@ SEARCH_RADIUS = 1 * u.arcsec
 
 The Nancy Grace Roman Space Telescope will carry out a wide-area survey (WAS) in the near infrared. OpenUniverse2024 simulates that survey as individual exposures, each covering one detector's worth of sky. Bands include F184, H158, J129, K213, W146, Y106.
 
-The survey stores those exposures by pointing and detector rather than by sky position, so in this section we define functions that ask the image search which exposure covers a position, and then read it.
+The survey stores those exposures by pointing and detector rather than by sky position, so in this section we define a function that asks the image search which exposure covers a position, and then reads it.
 
 +++
 
-### Define a function that returns the file path of a Roman image given a sky position and filter.
+### Define a function that retrieves a Roman simulated image given a sky position and filter.
+
+The Roman collection holds both the Wide-Area Survey and the Time Domain Survey, so we keep only the WAS exposures and then only the band we want. Many exposures cover any given position, so we sort by observation time and take the earliest.
 
-The Roman collection holds both the Wide-Area Survey and the Time Domain Survey, so we keep only the WAS exposures and then only the band we want. Many exposures cover any given position, so we sort by observation time and take the earliest; that keeps the notebook reproducible from one run to the next.
+Having found the file, we prefix its path with `s3://` and open it with astropy. The science image sits in the first extension, and we take the WCS from the same header as an `astropy.wcs.WCS` object. The function returns a dictionary of both, along with the header itself so we can look at it later, and the path so that later cells can reuse it instead of searching a second time.
 
 ```{code-cell} ipython3
 ---
@@ -137,7 +141,28 @@ jupyter:
   source_hidden: true
 tags: [hide-cell]
 ---
-def get_roman_image_fpath(coord, filter):
+def get_roman_image(coord, filter):
+    """Read the Roman WAS exposure covering a position in a given filter.
+
+    Parameters
+    ----------
+    coord : `~astropy.coordinates.SkyCoord`
+        Sky position to search around.
+    filter : str
+        Roman bandpass name, e.g. ``'H158'``.
+
+    Returns
+    -------
+    dict
+        Keys ``'data'`` (`~numpy.ndarray`), ``'wcs'`` (`~astropy.wcs.WCS`),
+        ``'header'`` (`~astropy.io.fits.Header`) and ``'fpath'`` (str, the S3
+        path of the exposure).
+
+    Raises
+    ------
+    ValueError
+        If no WAS exposure in that filter covers ``coord``.
+    """
     results = Irsa.query_sia(pos=(coord, SEARCH_RADIUS),
                              collection=OU_ROMAN_SIA_COLLECTION)
 
@@ -150,21 +175,7 @@ def get_roman_image_fpath(coord, filter):
 
     in_band.sort('t_min')
     cloud_info = json.loads(in_band['cloud_access'][0])['aws']
-    return f"{cloud_info['bucket_name']}/{cloud_info['key']}"
-```
-
-### Define a function that retrieves a Roman simulated image given a sky position and filter.
-
-Now we prefix that path with `s3://` and open it with astropy. The science image sits in the first extension, and we take the WCS from the same header as an `astropy.wcs.WCS` object. The function below returns a dictionary of both, along with the header itself so we can look at it later.
-
-```{code-cell} ipython3
----
-jupyter:
-  source_hidden: true
-tags: [hide-cell]
----
-def get_roman_image(coord, filter):
-    image_s3_fpath = get_roman_image_fpath(coord, filter)
+    image_s3_fpath = f"{cloud_info['bucket_name']}/{cloud_info['key']}"
 
     with fits.open(f"s3://{image_s3_fpath}", fsspec_kwargs={"anon": True}) as hdul:
         # retrieve science data, which sits in the first extension
@@ -173,8 +184,6 @@ def get_roman_image(coord, filter):
         # make wcs using header
         image_wcs = wcs.WCS(hdul[1].header)
 
-        # hand back the path too, so later cells can reuse it instead of
-        # asking the image search for it a second time
         return {'data': image_data, 'wcs': image_wcs, 'header': hdul[1].header,
                 'fpath': image_s3_fpath}
 ```
@@ -186,7 +195,7 @@ def get_roman_image(coord, filter):
 Choose a filter and a position within the survey footprint
 
 ```{code-cell} ipython3
-coord = SkyCoord(ra=9.6055383, dec=-44.1895542, unit="deg")
+coord = SkyCoord(ra=9.6205000, dec=-44.0641694, unit="deg")
 filter_roman = 'H158' #F184, H158, J129, K213, W146, and Y106 are available
 ```
 
@@ -247,6 +256,20 @@ jupyter:
 tags: [hide-cell]
 ---
 def coord_to_xy(w, coord):
+    """Convert a sky position into (x, y) pixel indices.
+
+    Parameters
+    ----------
+    w : `~astropy.wcs.WCS`
+        WCS of the image the indices refer to.
+    coord : `~astropy.coordinates.SkyCoord`
+        Sky position to convert.
+
+    Returns
+    -------
+    tuple of int
+        Pixel indices as ``(x, y)``.
+    """
     return w.world_to_array_index(coord)[::-1] #reverse since 0th axis is y, 1st axis is x
 ```
 
@@ -257,6 +280,8 @@ coord_arr_idx
 
 ### Use matplotlib imshow to create a static visualization of the Roman simulated image and overplot the selected position.
 
+Each telescope points independently, so neither survey's images arrive with north pointing up: this Roman exposure is turned about 170 degrees from north, and the Rubin visit in the next section about 160 degrees. Every figure in this notebook is therefore resampled onto a north-up grid before plotting, so that all of them, and both telescopes, can be read the same way. The rotation leaves the frame tilted inside a slightly larger box, with empty wedges in the corners.
+
 ```{code-cell} ipython3
 ---
 jupyter:
@@ -264,26 +289,106 @@ jupyter:
 tags: [hide-cell]
 ---
 def stretch_color(data, clipPercent):
-    return np.percentile(data, (0 + clipPercent, 100 - clipPercent))
-```
-
-```{code-cell} ipython3
-plt.imshow(image_roman['data'], origin='lower', 
-           clim=stretch_color(image_roman['data'], 1)
+    """Percentile clip limits for displaying an image.
+
+    Parameters
+    ----------
+    data : `~numpy.ndarray`
+        Image values; NaNs are ignored.
+    clipPercent : float
+        Percentage clipped off each end of the distribution.
+
+    Returns
+    -------
+    tuple of float
+        Lower and upper limits, for passing to ``clim``.
+    """
+    # nan-aware: rotating to north up leaves empty corners
+    return np.nanpercentile(data, (0 + clipPercent, 100 - clipPercent))
+
+
+def north_up(data, image_wcs, size, binning=1):
+    """Resample onto a north-up grid, keeping the original pixel scale.
+
+    Parameters
+    ----------
+    data : `~numpy.ndarray`
+        Image to resample.
+    image_wcs : `~astropy.wcs.WCS`
+        WCS of ``data``.
+    size : `~astropy.units.Quantity`
+        Angular width of the output grid.
+    binning : int, optional
+        Sample this many input pixels per output pixel. Values greater than 1
+        are enough for displaying a whole exposure and much quicker than
+        resampling every pixel.
+
+    Returns
+    -------
+    `~numpy.ndarray`
+        The resampled image.
+    `~astropy.wcs.WCS`
+        WCS of the resampled image.
+    """
+    scale = abs(image_wcs.proj_plane_pixel_scales()[0].to(u.deg).value) * binning
+    npix = int(round(size.to(u.deg).value / scale))
+
+    ny, nx = data.shape
+    center = image_wcs.pixel_to_world((nx - 1) / 2, (ny - 1) / 2)
+
+    target = wcs.WCS(naxis=2)
+    target.wcs.ctype = ['RA---TAN', 'DEC--TAN']
+    target.wcs.crpix = [(npix + 1) / 2, (npix + 1) / 2]
+    target.wcs.crval = [center.ra.deg, center.dec.deg]
+    target.wcs.cdelt = [-scale, scale]
+
+    return reproject_interp((data, image_wcs), target, shape_out=(npix, npix))[0], target
+
+
+def full_frame_north_up(image, binning=4):
+    """North-up view of a whole exposure, for display only.
+
+    Parameters
+    ----------
+    image : dict
+        Image dictionary with ``'data'`` and ``'wcs'`` keys, as returned by
+        `get_roman_image` or `get_rubin_image`.
+    binning : int, optional
+        Sample this many input pixels per output pixel.
+
+    Returns
+    -------
+    `~numpy.ndarray`
+        The resampled image.
+    `~astropy.wcs.WCS`
+        WCS of the resampled image.
+    """
+    ny, nx = image['data'].shape
+    fov = nx * image['wcs'].proj_plane_pixel_scales()[0].to(u.arcsec)
+    return north_up(image['data'], image['wcs'], 1.2 * fov, binning=binning)
+```
+
+```{code-cell} ipython3
+roman_display, roman_display_wcs = full_frame_north_up(image_roman)
+
+plt.imshow(roman_display, origin='lower', 
+           clim=stretch_color(roman_display, 1)
            )
 
-plt.plot(*coord_arr_idx, 'r+', markersize=15)
+plt.plot(*coord_to_xy(roman_display_wcs, coord), 'o', mfc='none', mec='r',
+         markersize=15, markeredgewidth=2)
+plt.axis('off')
 ```
 
 ## 3. Rubin Images
 
-OpenUniverse2024 includes simulated Rubin images in the following filters: u, g, r, i, z, y. As with Roman, these are stored by visit and detector rather than by sky position, so we again let the SIA service tell us which images cover our position. In this section, we define functions that retrieve a Rubin image for a chosen position and filter, returning data in the same structure as the functions we defined above for Roman.
+OpenUniverse2024 includes simulated Rubin images in the following filters: u, g, r, i, z, y. As with Roman, these are stored by visit and detector rather than by sky position, so we again let the SIA service tell us which images cover our position. In this section, we define functions that retrieve a Rubin image for a chosen position and filter, returning data in the same structure as the function we defined above for Roman.
 
 +++
 
 ### Retrieve Rubin images
 
-The image search returns one row per image, with the cloud location of each tucked inside a `cloud_access` JSON string. Many visits cover any given position, so we sort by observation time and take the earliest; that keeps the notebook reproducible from one run to the next.
+The image search returns one row per image, with the cloud location of each tucked inside a `cloud_access` JSON string. Many visits cover any given position, so we sort by observation time and take the earliest.
 
 ```{code-cell} ipython3
 ---
@@ -292,6 +397,25 @@ jupyter:
 tags: [hide-cell]
 ---
 def get_rubin_image_fpaths(coord, filters):
+    """Find the Rubin images covering a position in each of several filters.
+
+    Parameters
+    ----------
+    coord : `~astropy.coordinates.SkyCoord`
+        Sky position to search around.
+    filters : list of str
+        Rubin bandpass names, e.g. ``['g', 'r', 'i']``.
+
+    Returns
+    -------
+    dict
+        Filter name to S3 path of the image covering ``coord``.
+
+    Raises
+    ------
+    ValueError
+        If any requested filter has no image covering ``coord``.
+    """
     # One search covers every band, so ask once and sort the results out by band
     # afterwards rather than searching again for each one.
     results = Irsa.query_sia(pos=(coord, SEARCH_RADIUS),
@@ -311,10 +435,6 @@ def get_rubin_image_fpaths(coord, filters):
         fpaths[filter] = f"{cloud_info['bucket_name']}/{cloud_info['key']}"
 
     return fpaths
-
-
-def get_rubin_image_fpath(coord, filter):
-    return get_rubin_image_fpaths(coord, [filter])[filter]
 ```
 
 ```{code-cell} ipython3
@@ -324,7 +444,22 @@ jupyter:
 tags: [hide-cell]
 ---
 def get_rubin_image(coord, filter):
-    image_s3_fpath = get_rubin_image_fpath(coord, filter)
+    """Read the Rubin image covering a position in a given filter.
+
+    Parameters
+    ----------
+    coord : `~astropy.coordinates.SkyCoord`
+        Sky position to search around.
+    filter : str
+        Rubin bandpass name, e.g. ``'r'``.
+
+    Returns
+    -------
+    dict
+        Keys ``'data'`` (`~numpy.ndarray`), ``'wcs'`` (`~astropy.wcs.WCS`) and
+        ``'fpath'`` (str, the S3 path of the image).
+    """
+    image_s3_fpath = get_rubin_image_fpaths(coord, [filter])[filter]
 
     with fits.open(f"s3://{image_s3_fpath}", fsspec_kwargs={"anon": True}) as hdul:
         # retrieve science data, which sits in the first extension
@@ -373,33 +508,18 @@ A single Rubin detector covers a somewhat larger patch of sky than one Roman exp
 
 ### Use matplotlib imshow to create a static visualization of the Rubin simulated image and overplot the selected position.
 
-```{code-cell} ipython3
-plt.imshow(image_rubin['data'], origin='lower', 
-           clim=stretch_color(image_rubin['data'], 1)
-           )
-
-plt.plot(*coord_to_xy(image_rubin['wcs'], coord), 'r+', markersize=15)
-```
-
-### Define a function that returns the URL for a given S3 filepath
-Since the OpenUniverse2024 data is available through a public S3 bucket, we can access a given S3 file using HTTPS URL as follows:
+The two full-frame figures do not cover the same patch of sky. Each is centered on its own exposure and sized to its own field of view, so the Rubin panel is nearly twice as wide on the sky as the Roman one. What they share is the marked position, which falls inside both.
 
 ```{code-cell} ipython3
----
-jupyter:
-  source_hidden: true
-tags: [hide-cell]
----
-def https_url(s3_fpath):
-    s3_fpath_without_bucket = s3_fpath.split('/', 1)[1]
-    return f"https://{BUCKET_NAME}.s3.amazonaws.com/{s3_fpath_without_bucket}"
-```
+rubin_display, rubin_display_wcs = full_frame_north_up(image_rubin)
 
-Let's generate URL for the Rubin image we plotted above. Clicking on the returned URL will allow you to download this image locally.
+plt.imshow(rubin_display, origin='lower', 
+           clim=stretch_color(rubin_display, 1)
+           )
 
-```{code-cell} ipython3
-image_s3_fpath_rubin = image_rubin['fpath']
-https_url(image_s3_fpath_rubin)
+plt.plot(*coord_to_xy(rubin_display_wcs, coord), 'o', mfc='none', mec='r',
+         markersize=15, markeredgewidth=2)
+plt.axis('off')
 ```
 
 ## 4. Compare simulated Roman and Rubin cutouts for a selected position
@@ -415,8 +535,20 @@ cutout_size = 50*u.arcsec
 ### Create the cutouts
 
 ```{code-cell} ipython3
-cutout_roman = Cutout2D(image_roman['data'], coord, size=cutout_size, wcs=image_roman['wcs'])
-cutout_rubin = Cutout2D(image_rubin['data'], coord, size=cutout_size, wcs=image_rubin['wcs'])
+# Cut oversized: the rotation below needs data in the corners to draw from.
+cutout_roman = Cutout2D(image_roman['data'], coord, size=1.5*cutout_size, wcs=image_roman['wcs'])
+cutout_rubin = Cutout2D(image_rubin['data'], coord, size=1.5*cutout_size, wcs=image_rubin['wcs'])
+```
+
+### Rotate both cutouts so that north points up
+
+The full frames above were rotated for display; here we do the same to the cutouts, using the same helper. The two telescopes are turned about 10 degrees apart from each other, so without this the panels would not line up. Each keeps its own pixel scale, so Roman still shows the finer sampling.
+
+```{code-cell} ipython3
+# Rotating a same-size grid would leave empty corners, so each cutout above was cut
+# oversized and is resampled down to the size we actually want.
+roman_data, roman_wcs = north_up(cutout_roman.data, cutout_roman.wcs, cutout_size)
+rubin_data, rubin_wcs = north_up(cutout_rubin.data, cutout_rubin.wcs, cutout_size)
 ```
 
 ### Use matplotlib imshow to plot static side-by-side comparisons of the cutouts
@@ -425,25 +557,25 @@ cutout_rubin = Cutout2D(image_rubin['data'], coord, size=cutout_size, wcs=image_
 fig, axs = plt.subplots(1, 2, figsize=(12, 6))
 
 
-axs[0].imshow(cutout_roman.data, origin='lower', 
-              clim=stretch_color(cutout_roman.data, .5)
+axs[0].imshow(roman_data, origin='lower', 
+              clim=stretch_color(roman_data, .5)
               )
-axs[0].plot(*coord_to_xy(cutout_roman.wcs, coord), 'r+', markersize=15)
 axs[0].set_title(f"ROMAN in filter {filter_roman}")
 
-axs[1].imshow(cutout_rubin.data, origin='lower', 
-              clim=stretch_color(cutout_rubin.data, .5)
+axs[1].imshow(rubin_data, origin='lower', 
+              clim=stretch_color(rubin_data, .5)
               )
-axs[1].plot(*coord_to_xy(cutout_rubin.wcs, coord), 'r+', markersize=15)
 axs[1].set_title(f"RUBIN in filter {filter_rubin}")
 
 fig.suptitle(f"Cutouts at ({coord.ra}, {coord.dec}) with {cutout_size} size", fontsize=14)
 plt.tight_layout(rect=[0, 0, 1, 0.97])
 ```
 
-## 5. Use Firefly to interactively identify a blended source
+Both panels cover the same patch of sky, at the same size and with north up in each. The difference is stark: the simulated Roman exposure resolves a field full of galaxies, while the same patch in the simulated Rubin image is close to blank. A single Rubin visit is short, so only the brightest sources here rise above its noise, and from the ground its point spread function is broad enough to blend close neighbors together. Roman's finer pixels and sharper point spread function keep those sources separate. Rubin reaches its full depth by stacking many visits, which is not what we are looking at.
 
-The simulated Roman exposure has finer pixels and a sharper point spread function than the simulated Rubin image. Let's try to locate blended objects to compare in the simulated Rubin and Roman images. We will use Firefly's interactive visualization to make this task easier.
+## 5. Use Firefly to interactively explore the images
+
+Static cutouts only help once you know where to look. Choosing a position worth comparing means exploring the image interactively, and for that we hand it to Firefly. In this section we display the simulated Rubin visit, outline the area the Roman exposure covers, and pan and zoom around it to pick out a target.
 
 +++
 
@@ -454,14 +586,47 @@ You also need a Firefly server to communicate with your Firefly Python client. I
 
 ```{code-cell} ipython3
 # Uncomment when using within Jupyter Lab with jupyter_firefly_extensions installed
-# fc = FireflyClient.make_lab_client()
+fc = FireflyClient.make_lab_client()
 
 # Uncomment for contexts other than above 
-fc = FireflyClient.make_client(url="https://irsa.ipac.caltech.edu/irsaviewer")
+#fc = FireflyClient.make_client(url="https://irsa.ipac.caltech.edu/irsaviewer")
 
 fc.reinit_viewer() # to clean the state, if this cell ran earlier
 ```
 
+### Define a function that returns the URL for a given S3 filepath
+Since the OpenUniverse2024 data are in a public S3 bucket, any file can also be reached over HTTPS. Firefly loads images and catalogs by URL rather than by S3 path, so we use this helper throughout the rest of the notebook.
+
+```{code-cell} ipython3
+---
+jupyter:
+  source_hidden: true
+tags: [hide-cell]
+---
+def https_url(s3_fpath):
+    """Convert an S3 path into a public HTTPS URL.
+
+    Parameters
+    ----------
+    s3_fpath : str
+        Path of the form ``'bucket/key'``.
+
+    Returns
+    -------
+    str
+        HTTPS URL for the same file.
+    """
+    s3_fpath_without_bucket = s3_fpath.split('/', 1)[1]
+    return f"https://{BUCKET_NAME}.s3.amazonaws.com/{s3_fpath_without_bucket}"
+```
+
+Let's generate URL for the Rubin image we plotted above. Clicking on the returned URL will allow you to download this image locally.
+
+```{code-cell} ipython3
+image_s3_fpath_rubin = image_rubin['fpath']
+https_url(image_s3_fpath_rubin)
+```
+
 ### Send the simulated Rubin image to Firefly using show_fits.
 
 For displaying the FITS image of the Rubin visit in Firefly, we use [`show_fits`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.show_fits):
@@ -494,28 +659,31 @@ corner_coords = ' '.join(f'{ra}d {dec}d' for ra, dec in roman_corners)
 roman_regions = [f'icrs;polygon {corner_coords} # text={{Roman {filter_roman}}}']
 
 roman_regions_id = 'roman_regions'
+# Name the plot explicitly. Without plot_id the layer goes to whichever plot is
+# active, which is not necessarily this one if cells are run out of order.
 fc.overlay_region_layer(region_data=roman_regions,
                         title='Roman exposure', 
-                        region_layer_id=roman_regions_id)
+                        region_layer_id=roman_regions_id,
+                        plot_id=image_ff_id_rubin)
 ```
 
-### Use Firefly's pan and zoom capabilities to locate a region of interest (a blended source)
+### Use Firefly's pan and zoom capabilities to locate a region of interest
 
 You can view the coordinates of your mouse pointer at the bottom left of the display window. To copy the coordinates for a specific coordinate:
 
 - Toggle the "Click Lock" to "on" in the bottom right of the image display.
 - Click on the position of interest, and notice that the coordinate display is now frozen.
-- Click on "EQ-J2000", the coordinate label in the bottom left of the image display. In the dialog that opens, change copy options to "[Python] Astropy SkyCoord" so that we can directly work with them in python.
-- Close the dialog and click on the copy icon next to the coordinate values display.
+- Click on "EQ-J2000", the coordinate label in the bottom left of the image display. The pop-up that opens lets you choose the coordinate system, and the format used when copying — pick the Python format so the value can be pasted straight into a `SkyCoord`.
+- Close the pop-up and click the small clipboard icon next to the coordinate readout.
 
 +++
 
 ### Copy the coordinates from the coordinate display to the Python notebook
 
-We have provided an example. You can change this based on your interests.
+The position we have been using since Section 2 is the one marked below, so by default we carry it forward. Paste a position of your own here instead to compare a different target.
 
 ```{code-cell} ipython3
-coords_of_interest = SkyCoord('0h38m40.26s -44d04m22.26s', frame='icrs') # located and copied through UI
+coords_of_interest = coord  # or e.g. SkyCoord('0h38m28.92s -44d03m51.01s', frame='icrs')
 coords_of_interest
 ```
 
@@ -528,19 +696,21 @@ We can now use this [astropy `SkyCoord` object](https://docs.astropy.org/en/stab
 For this we use the id of the region layer we defined above, and add more region data using [`add_region_data`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.add_region_data).
 
 ```{code-cell} ipython3
-point_region = f'icrs;point {coords_of_interest.ra.value}d {coords_of_interest.dec.value}d # point=cross 15 text={{Blended source}}'
-fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id)
+point_region = f'icrs;point {coords_of_interest.ra.value}d {coords_of_interest.dec.value}d # point=cross 15 text={{Source of interest}}'
+fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id,
+                   plot_id=image_ff_id_rubin)
 ```
 
-## 6. Plot cutouts of the identified blended source
-
-Two catalog sources sit 1.22 arcsec apart at this position. Roman separates them; Rubin merges them into one extended blob. Note that the Rubin image is a single visit, so faint detail washes out into the noise.
+## 6. Plot cutouts of the identified source
 
-The two cutouts are plotted in their own pixel frames, and because each telescope was pointed independently, north falls about 30 degrees apart between the panels. We leave them unrotated here, so expect the pair to appear tilted relative to each other rather than aligned.
+This is a small group of galaxies, with a close neighbour 2.7 arcsec away. Roman shows them as separate objects; in Rubin they sit at the level of the noise. The Rubin image is a single visit rather than a stack, so faint sources like these are simply not detected in it, which is worth keeping in mind when comparing the two panels below.
 
 ```{code-cell} ipython3
-image_roman = get_roman_image(coords_of_interest, filter_roman)
-image_rubin = get_rubin_image(coords_of_interest, filter_rubin)
+# The images from Sections 2 and 3 already cover the default position, so only
+# re-read if you pasted a different one above.
+if coords_of_interest.separation(coord) > 0*u.arcsec:
+    image_roman = get_roman_image(coords_of_interest, filter_roman)
+    image_rubin = get_rubin_image(coords_of_interest, filter_rubin)
 ```
 
 ```{code-cell} ipython3
@@ -548,34 +718,38 @@ cutout_size = 8*u.arcsec
 ```
 
 ```{code-cell} ipython3
-cutout_roman = Cutout2D(image_roman['data'], coords_of_interest, size=cutout_size, wcs=image_roman['wcs'])
-cutout_rubin = Cutout2D(image_rubin['data'], coords_of_interest, size=cutout_size, wcs=image_rubin['wcs'])
+# As in section 4, cut oversized and then rotate both to north up.
+cutout_roman = Cutout2D(image_roman['data'], coords_of_interest, size=1.5*cutout_size, wcs=image_roman['wcs'])
+cutout_rubin = Cutout2D(image_rubin['data'], coords_of_interest, size=1.5*cutout_size, wcs=image_rubin['wcs'])
+
+roman_data, roman_wcs = north_up(cutout_roman.data, cutout_roman.wcs, cutout_size)
+rubin_data, rubin_wcs = north_up(cutout_rubin.data, cutout_rubin.wcs, cutout_size)
 ```
 
 ```{code-cell} ipython3
 fig, axs = plt.subplots(1, 2, figsize=(12, 6))
 
-axs[0].imshow(cutout_roman.data, origin='lower', 
-              clim=stretch_color(cutout_roman.data, .5)
+axs[0].imshow(roman_data, origin='lower', 
+              clim=stretch_color(roman_data, .5)
               )
 axs[0].set_title(f"ROMAN in filter {filter_roman}")
 
-# Let's also encircle the blended source we identified. The circle is 2 arcsec across
+# Let's also encircle the source we identified. The circle is 2 arcsec across
 # on both panels, so the two cutouts can be compared directly despite their different
 # pixel scales.
-roman_scale = cutout_roman.wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
-axs[0].add_patch(patches.Circle(coord_to_xy(cutout_roman.wcs, coords_of_interest), 
+roman_scale = roman_wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
+axs[0].add_patch(patches.Circle(coord_to_xy(roman_wcs, coords_of_interest), 
                                 radius=2/roman_scale, color='r', fill=False, linewidth=2))
 
 
-axs[1].imshow(cutout_rubin.data, origin='lower', 
-              clim=stretch_color(cutout_rubin.data, .5)
+axs[1].imshow(rubin_data, origin='lower', 
+              clim=stretch_color(rubin_data, .5)
               )
 axs[1].set_title(f"RUBIN in filter {filter_rubin}")
 
 # Let's also encircle the source we identified 
-rubin_scale = cutout_rubin.wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
-axs[1].add_patch(patches.Circle(coord_to_xy(cutout_rubin.wcs, coords_of_interest),
+rubin_scale = rubin_wcs.proj_plane_pixel_scales()[0].to(u.arcsec).value
+axs[1].add_patch(patches.Circle(coord_to_xy(rubin_wcs, coords_of_interest),
                                 radius=2/rubin_scale, color='r', fill=False, linewidth=2))
 
 
@@ -591,6 +765,8 @@ OpenUniverse2024 includes the input truth files that were used to create the sim
 
 The catalogs are split by HEALPix sky region (nside=32, RING ordering) and the region index appears in each filename, so we can convert our sky position into that index with [`hpgeom`](https://hpgeom.readthedocs.io/en/latest/) rather than hunting through the directory listing. This is the same approach the [Quickstart](openuniverse2024_quickstart) tutorial takes.
 
+Each region comes as two separate catalogs: a point source catalog holding the stars, and a galaxy catalog. We load both below and overlay each on the image.
+
 ```{code-cell} ipython3
 region = hpgeom.angle_to_pixel(32, coord.ra.deg, coord.dec.deg, lonlat=True, nest=False)
 print(f"HEALPix region for our position: {region}")
@@ -650,13 +826,13 @@ For each row in the table you can notice a marker in the image. Selecting a row
 
 ### Use Firefly's apply_table_filters to show only high-redshift galaxies
 
-High redshift galaxies are the most interesting, so let's filter the table we sent to Firefly to only include z>3 galaxies. Notice how the table display and image overlay change. Notice how the chart becomes a scatterplot from a heatmap because the sources reduce. You can remove this filter or add new ones through the GUI.
+High redshift galaxies are the most interesting, so let's filter the table we sent to Firefly to only include z>2.5 galaxies. Notice how the table display and image overlay change. Notice how the chart becomes a scatterplot from a heatmap because the sources reduce. You can remove this filter or add new ones through the GUI.
 
 For filtering, we will use [`apply_table_filters`](https://caltech-ipac.github.io/firefly_client/api/firefly_client.FireflyClient.html#firefly_client.FireflyClient.apply_table_filters) method on the galaxy table we loaded above.
 
 ```{code-cell} ipython3
 fc.apply_table_filters(tbl_id=gal_cat_tbl_id,
-                       filters=" AND ".join(cat_filters+['"redshift" > 3']))
+                       filters=" AND ".join(cat_filters+['"redshift" > 2.5']))
 ```
 
 You can play with the filters directly from the UI as well. Try removing adding more filters in the tables and see how markers change.
@@ -692,14 +868,15 @@ For example, we found a region of the sky that seems to have a high number of hi
 
 ```{code-cell} ipython3
 # located and copied through UI
-high_z_gal_coords = SkyCoord('0h38m00.77s -44d12m10.2s', frame='icrs')
+high_z_gal_coords = SkyCoord('0h38m28.62s -44d00m50.52s', frame='icrs')
 high_z_gal_coords
 ```
 
 ```{code-cell} ipython3
 # let's also mark it in our region layer, so that it's easy to pinpoint later
-point_region = f'icrs;point {high_z_gal_coords.ra.value}d {high_z_gal_coords.dec.value}d # point=cross 15 text={{z>3 mock galaxies}}'
-fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id)
+point_region = f'icrs;point {high_z_gal_coords.ra.value}d {high_z_gal_coords.dec.value}d # point=cross 15 text={{z>2.5 mock galaxies}}'
+fc.add_region_data(region_data=point_region, region_layer_id=roman_regions_id,
+                   plot_id=image_ff_id_rubin)
 ```
 
 ## 8. Plot 3-color Roman image containing your region of interest
@@ -720,6 +897,19 @@ jupyter:
 tags: [hide-cell]
 ---
 def get_fits_stream(image_roman):
+    """Write an image to an in-memory FITS file for handing to Firefly.
+
+    Parameters
+    ----------
+    image_roman : dict
+        Image dictionary with ``'data'`` and ``'wcs'`` keys, as returned by
+        `get_roman_image`.
+
+    Returns
+    -------
+    `~io.BytesIO`
+        FITS file positioned at the start, ready to be read.
+    """
     # Create a FITS PrimaryHDU object with the image data
     hdu = fits.PrimaryHDU(data=image_roman['data'], header=image_roman['wcs'].to_header())
 
@@ -792,6 +982,6 @@ fc.set_stretch(plot_id=image_ff_id_roman_3color, stype='sigma', algorithm='squar
 
 ## About This Notebook
 
-**Updated:** 2026-08-05
+**Updated:** 2026-08-11
 
 **Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems.

From 195665cd2a6f55488964565f284f51e040c7b1b3 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 14 Aug 2026 21:09:58 -0400
Subject: [PATCH 7/9] Restore the testing overrides and add reproject to the
 install line

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../openuniverse2024_firefly_visualization.md          | 10 +++++-----
 1 file changed, 5 insertions(+), 5 deletions(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
index 122f804b..2447f7dc 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
@@ -8,7 +8,7 @@ jupytext:
     extension: .md
     format_name: myst
     format_version: 0.13
-    jupytext_version: 1.19.5
+    jupytext_version: 1.16.3
 kernelspec:
   name: python3
   display_name: python3
@@ -76,7 +76,7 @@ If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-!pip install numpy astropy matplotlib firefly_client astroquery hpgeom reproject
+# !pip install numpy astropy matplotlib firefly_client astroquery hpgeom reproject
 ```
 
 ```{code-cell} ipython3
@@ -586,10 +586,10 @@ You also need a Firefly server to communicate with your Firefly Python client. I
 
 ```{code-cell} ipython3
 # Uncomment when using within Jupyter Lab with jupyter_firefly_extensions installed
-fc = FireflyClient.make_lab_client()
+# fc = FireflyClient.make_lab_client()
 
 # Uncomment for contexts other than above 
-#fc = FireflyClient.make_client(url="https://irsa.ipac.caltech.edu/irsaviewer")
+fc = FireflyClient.make_client(url="https://irsa.ipac.caltech.edu/irsaviewer")
 
 fc.reinit_viewer() # to clean the state, if this cell ran earlier
 ```
@@ -982,6 +982,6 @@ fc.set_stretch(plot_id=image_ff_id_roman_3color, stype='sigma', algorithm='squar
 
 ## About This Notebook
 
-**Updated:** 2026-08-11
+**Updated:** 2026-08-14
 
 **Contact:** [the IRSA Helpdesk](https://irsa.ipac.caltech.edu/docs/help_desk.html) with questions or reporting problems.

From b402b7f4aaf80b59e5c9e9defb061bcbbf8dadd6 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 14 Aug 2026 21:09:58 -0400
Subject: [PATCH 8/9] Ignore the local PR description note

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .gitignore | 1 +
 1 file changed, 1 insertion(+)

diff --git a/.gitignore b/.gitignore
index d9bb91cf..3f7c7834 100644
--- a/.gitignore
+++ b/.gitignore
@@ -58,3 +58,4 @@ __pycache__
 # Local working notes, not part of the tutorials
 CLAUDE.md
 OU_PR_REVIEW.md
+PR_DESCRIPTION.md

From 99ba4fa6833883ec3db0e1450a495388c587f465 Mon Sep 17 00:00:00 2001
From: Jessica Krick <jkrick@caltech.edu>
Date: Fri, 14 Aug 2026 21:11:03 -0400
Subject: [PATCH 9/9] Add s3fs back to the install line

Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
---
 .../OpenUniverse2024/openuniverse2024_firefly_visualization.md  | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
index 2447f7dc..71400cf2 100644
--- a/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
+++ b/tutorials/simulated-data/OpenUniverse2024/openuniverse2024_firefly_visualization.md
@@ -76,7 +76,7 @@ If you are new to OpenUniverse2024, the [Quickstart](openuniverse2024_quickstart
 
 ```{code-cell} ipython3
 # Uncomment the next line to install dependencies if needed.
-# !pip install numpy astropy matplotlib firefly_client astroquery hpgeom reproject
+# !pip install numpy astropy s3fs matplotlib firefly_client astroquery hpgeom reproject
 ```
 
 ```{code-cell} ipython3