Skip to content

Commit 6933a21

Browse files
christinawilliamschristinawilliams
committed
updating image stamps nb
1 parent 396ff06 commit 6933a21

1 file changed

Lines changed: 86 additions & 76 deletions

File tree

‎DP1/100_How_to_Use_RSP_Tools/103_Image_access_and_display/103_5_Image_stamps.ipynb‎

Lines changed: 86 additions & 76 deletions
Original file line numberDiff line numberDiff line change
@@ -22,7 +22,7 @@
2222
"Data Release: <a href=\"https://dp1.lsst.io/\">Data Preview 1</a> <br>\n",
2323
"Container Size: large <br>\n",
2424
"LSST Science Pipelines version: r29.2.0 <br>\n",
25-
"Last verified to run: 2026-01-31 <br>\n",
25+
"Last verified to run: 2026-02-10 <br>\n",
2626
"Repository: <a href=\"https://github.com/lsst/tutorial-notebooks\">github.com/lsst/tutorial-notebooks</a> <br>"
2727
]
2828
},
@@ -33,9 +33,9 @@
3333
"source": [
3434
"**Learning objective:** An overview of candidate transient (Ia Supernovae) identification in DP1.\n",
3535
"\n",
36-
"**LSST data products:** `DiaObject`, `ForcedSourceOnDiaObject`, `Visit`, `visit_image`, `template_coadd`, `difference_image`\n",
36+
"**LSST data products:** `DiaObject`, `Visit`, `visit_image`, `template_coadd`, `difference_image`\n",
3737
"\n",
38-
"**Packages:** `matplotlib`, `numpy`, `lsst.rsp`,`lsst.afw`, `lsst.geom`, `lsst.resources`, `lsst.utils.plotting`, `pyvo.dal.adhoc`\n",
38+
"**Packages:** `matplotlib`, `numpy`, `lsst.rsp`,`lsst.afw`, `lsst.geom`, `lsst.resources`, `lsst.utils.plotting`, `pyvo.dal.adhoc`, `reproject`\n",
3939
"\n",
4040
"**Credit:**\n",
4141
"Originally developed by the Rubin Community Science team with feedback from Eric Bellm. This notebook also utilizes a transient detection from LSSTComCam first identified by Dan Taranu.\n",
@@ -56,9 +56,9 @@
5656
"source": [
5757
"## 1. Introduction\n",
5858
"\n",
59-
"This notebook demonstrates how to use the light curve statistics in the DiaObject table to identify candidate transients. The specific case in this notebook is for Ia Supernova candidates.\n",
59+
"This notebook demonstrates how to produce bulk image cutouts using the Rubin cutout service.\n",
6060
"\n",
61-
"**Related tutorials**: There are 200-level tutorials on `difference_images` as well as the `DiaSource`, `DiaObject` and `ForcedSourceOnDiaObject` catalogs. Use of the image cutout service is covered in the 100-level tutorials. There is another notebook in this 306 series (306.1) that presents how to make a light curve of an extragalactic transient."
61+
"**Related tutorials**: This is the second of two 100-level tutorials demonstrating use of the image cutout service. Tutorial notebook 103.4 demonstrates the return of cutouts of Rubin LSST image type ExposureF and MaskedImageF."
6262
]
6363
},
6464
{
@@ -82,7 +82,8 @@
8282
"\n",
8383
"From the `lsst` package, import modules for accessing the Table Access Protocol (TAP) service, image display functions, and the Simple Image Access (SIA) service from the LSST Science Pipelines (<a href=\"https://pipelines.lsst.io/\">pipelines.lsst.io</a>).\n",
8484
"\n",
85-
"Functions from `pyvo` are imported to perform and parse results from Datalink and SODA queries. This is used to obtain image products cutouts from the SIA service."
85+
"Functions from `pyvo` are imported to perform and parse results from Datalink and SODA queries. This is used to obtain image products cutouts from the SIA service.\n",
86+
"\n"
8687
]
8788
},
8889
{
@@ -93,6 +94,8 @@
9394
"outputs": [],
9495
"source": [
9596
"import matplotlib.pyplot as plt\n",
97+
"import math\n",
98+
"\n",
9699
"import numpy as np\n",
97100
"from astropy import units as u\n",
98101
"\n",
@@ -104,8 +107,7 @@
104107
" get_multiband_plot_symbols)\n",
105108
"\n",
106109
"import lsst.afw.display as afwDisplay\n",
107-
"from lsst.afw.image import ExposureF, ImageF\n",
108-
"from lsst.afw.math import Warper, WarperConfig\n",
110+
"from lsst.afw.image import ImageF\n",
109111
"from lsst.afw.fits import MemFileManager\n",
110112
"import lsst.geom as geom\n",
111113
"import lsst.resources\n",
@@ -115,9 +117,9 @@
115117
"from astropy.wcs import WCS\n",
116118
"from astropy.io import fits\n",
117119
"import io\n",
118-
"import tempfile\n",
119120
"import getpass\n",
120121
"import glob\n",
122+
"import os\n",
121123
"\n",
122124
"from reproject import reproject_interp\n",
123125
"\n",
@@ -351,6 +353,34 @@
351353
" return ImageF(mem), hdul"
352354
]
353355
},
356+
{
357+
"cell_type": "code",
358+
"execution_count": null,
359+
"id": "413fc23e-6c5e-47b6-997e-45d704d709a8",
360+
"metadata": {},
361+
"outputs": [],
362+
"source": [
363+
"def make_subplot_grid(n_subplots, figsize_per_plot=(4, 3)):\n",
364+
" \"\"\"\n",
365+
" Create an optimal grid of subplots for n_subplots.\n",
366+
"\n",
367+
" Returns\n",
368+
" -------\n",
369+
" fig : matplotlib.figure.Figure\n",
370+
" axes : list of matplotlib.axes.Axes\n",
371+
" \"\"\"\n",
372+
" n_cols = math.ceil(math.sqrt(n_subplots))\n",
373+
" n_rows = math.ceil(n_subplots / n_cols)\n",
374+
"\n",
375+
" figsize = (figsize_per_plot[0] * n_cols,\n",
376+
" figsize_per_plot[1] * n_rows)\n",
377+
"\n",
378+
" fig, axes = plt.subplots(n_rows, n_cols, figsize=figsize)\n",
379+
" axes = axes.flatten() if n_subplots > 1 else [axes]\n",
380+
"\n",
381+
" return fig, axes[:n_subplots]"
382+
]
383+
},
354384
{
355385
"cell_type": "markdown",
356386
"id": "b1091fba-dbf4-4bdf-9970-f1a69951a4f2",
@@ -373,10 +403,8 @@
373403
"metadata": {},
374404
"outputs": [],
375405
"source": [
376-
"SNcandID = 611255759837069401\n",
377-
"#ra = DiaObjs[DiaObjs['diaObjectId'] == SNcandID]['ra'][0]\n",
378-
"#dec = DiaObjs[DiaObjs['diaObjectId'] == SNcandID]['dec'][0]\n",
379-
"ra = 53.124767650110215 \n",
406+
"#SNcandID = 611255759837069401\n",
407+
"ra = 53.124767650110215\n",
380408
"dec = -27.739814591611168\n",
381409
"spherePoint = lsst.geom.SpherePoint(ra*geom.degrees, dec*geom.degrees)"
382410
]
@@ -505,8 +533,8 @@
505533
"outputs": [],
506534
"source": [
507535
"print('begin date:', scitab['t_max'][0], difftab['t_max'][0])\n",
508-
"print('end date:', scitab['t_max'][len(scitab['t_max'])-1])#, difftab['t_max'][0])\n",
509-
"print('number of visits:', len(scitab['t_max']))\n"
536+
"print('end date:', scitab['t_max'][len(scitab['t_max'])-1])\n",
537+
"print('number of visits:', len(scitab['t_max']))"
510538
]
511539
},
512540
{
@@ -550,9 +578,7 @@
550578
"metadata": {},
551579
"outputs": [],
552580
"source": [
553-
"#stamps could be smaller? \n",
554-
"\n",
555-
"fov = 0.003\n",
581+
"fov = 0.002\n",
556582
"sci, scihdul = get_cutout(dl_result_sci, spherePoint, get_pyvo_auth(), fov)\n",
557583
"ref, refhdul = get_cutout(dl_result_ref, spherePoint, get_pyvo_auth(), fov)\n",
558584
"diff, diffhdul = get_cutout(dl_result_diff, spherePoint, get_pyvo_auth(), fov)"
@@ -634,30 +660,54 @@
634660
"\n",
635661
"ref_primary_header = refhdul[0].header\n",
636662
"ref_hdu = refhdul[1]\n",
637-
"ref_header =refhdul[1].header\n",
663+
"ref_header = refhdul[1].header\n",
638664
"refdata = refhdul[1].data\n",
639665
"refwcs = WCS(refhdul[1].header)\n"
640666
]
641667
},
668+
{
669+
"cell_type": "markdown",
670+
"id": "02974977-001b-467e-ae0d-dc44f105872c",
671+
"metadata": {},
672+
"source": [
673+
"Use the `reproject_interp` function to put the science visit image onto the WCS of the template image."
674+
]
675+
},
676+
{
677+
"cell_type": "code",
678+
"execution_count": null,
679+
"id": "635013a0-e962-4958-bb09-007ac1cb169e",
680+
"metadata": {},
681+
"outputs": [],
682+
"source": [
683+
"reproj_scidata, footprint = reproject_interp((scidata, sciwcs), refwcs)"
684+
]
685+
},
686+
{
687+
"cell_type": "markdown",
688+
"id": "6b3360a9-ccc6-47bb-bc45-b87f208dacda",
689+
"metadata": {},
690+
"source": [
691+
"Plot the reprojected science visit image together with the template image to compare."
692+
]
693+
},
642694
{
643695
"cell_type": "code",
644696
"execution_count": null,
645697
"id": "d9734ade-5848-4e44-8e38-ea70aebd6a47",
646698
"metadata": {},
647699
"outputs": [],
648700
"source": [
649-
"\n",
650-
"reproj_scidata, footprint = reproject_interp((scidata, sciwcs), refwcs)\n",
651-
"\n",
652-
"ax1 = plt.subplot(1,2,1, projection=WCS(ref_hdu.header))\n",
701+
"ax1 = plt.subplot(1, 2, 1, projection=WCS(ref_hdu.header))\n",
653702
"ax1.imshow(reproj_scidata, origin='lower', vmin=np.nanpercentile(reproj_scidata, 1),\n",
654703
" vmax=np.nanpercentile(reproj_scidata, 99))\n",
655704
"ax1.coords['ra'].set_axislabel('Right Ascension')\n",
656705
"ax1.coords['dec'].set_axislabel('Declination')\n",
657706
"ax1.set_title('Reprojected Visit Image')\n",
658707
"\n",
659-
"ax2 = plt.subplot(1,2,2, projection=WCS(ref_hdu.header))\n",
660-
"ax2.imshow(refdata, origin='lower',vmin=np.nanpercentile(refdata, 1), vmax=np.nanpercentile(refdata, 99))\n",
708+
"ax2 = plt.subplot(1, 2, 2, projection=WCS(ref_hdu.header))\n",
709+
"ax2.imshow(refdata, origin='lower', vmin=np.nanpercentile(refdata, 1),\n",
710+
" vmax=np.nanpercentile(refdata, 99))\n",
661711
"ax2.coords['ra'].set_axislabel('Right Ascension')\n",
662712
"ax2.coords['dec'].set_axislabel('Declination')\n",
663713
"ax2.coords['dec'].set_axislabel_position('r')\n",
@@ -683,41 +733,6 @@
683733
"Bulk image cutouts will be a key visualization tool for DIA objects. Below, demonstrate how to retrieve (and save), and display many image cutouts."
684734
]
685735
},
686-
{
687-
"cell_type": "code",
688-
"execution_count": null,
689-
"id": "413fc23e-6c5e-47b6-997e-45d704d709a8",
690-
"metadata": {},
691-
"outputs": [],
692-
"source": [
693-
"import math\n",
694-
"import matplotlib.pyplot as plt\n",
695-
"\n",
696-
"def make_subplot_grid(n_subplots, figsize_per_plot=(4, 3)):\n",
697-
" \"\"\"\n",
698-
" Create an optimal grid of subplots for n_subplots.\n",
699-
"\n",
700-
" Returns\n",
701-
" -------\n",
702-
" fig : matplotlib.figure.Figure\n",
703-
" axes : list of matplotlib.axes.Axes\n",
704-
" \"\"\"\n",
705-
" n_cols = math.ceil(math.sqrt(n_subplots))\n",
706-
" n_rows = math.ceil(n_subplots / n_cols)\n",
707-
"\n",
708-
" figsize = (figsize_per_plot[0] * n_cols,\n",
709-
" figsize_per_plot[1] * n_rows)\n",
710-
"\n",
711-
" fig, axes = plt.subplots(n_rows, n_cols, figsize=figsize)\n",
712-
" axes = axes.flatten() if n_subplots > 1 else [axes]\n",
713-
"\n",
714-
" # Hide unused axes\n",
715-
" for ax in axes[n_subplots:]:\n",
716-
" ax.set_visible(False)\n",
717-
"\n",
718-
" return fig, axes[:n_subplots]"
719-
]
720-
},
721736
{
722737
"cell_type": "code",
723738
"execution_count": null,
@@ -734,30 +749,26 @@
734749
" j = i + 5\n",
735750
" ax.set_title(f\"MJD {np.round(scitab['t_max'][j], 4)}\")\n",
736751
"\n",
737-
" dl_result_sci = DatalinkResults.from_result_url(scitab['access_url'][j], session=get_pyvo_auth())\n",
738-
" #print(dl_result_sci)\n",
752+
" dl_result_sci = DatalinkResults.from_result_url(scitab['access_url'][j],\n",
753+
" session=get_pyvo_auth())\n",
754+
"\n",
739755
" sci, scihdul = get_cutout(dl_result_sci, spherePoint, get_pyvo_auth(), fov)\n",
740756
" sci_hdu = scihdul[1]\n",
741757
" sci_header = scihdul[1].header\n",
742758
" scidata = scihdul[1].data\n",
743759
" sciwcs = WCS(scihdul[1].header)\n",
744760
" reproj_scidata, footprint = reproject_interp((scidata, sciwcs), refwcs)\n",
745761
"\n",
746-
" ax.imshow(reproj_scidata, origin='lower', vmin=np.nanpercentile(reproj_scidata, 1), vmax=np.nanpercentile(reproj_scidata, 99))\n",
747-
"\n",
762+
" ax.imshow(reproj_scidata, origin='lower', vmin=np.nanpercentile(reproj_scidata, 1),\n",
763+
" vmax=np.nanpercentile(reproj_scidata, 99))\n",
748764
"\n",
749765
" #figname = os.path.join(tempdir, 'cutout_' + str(i) + '.png')\n",
750766
" #if os.path.isfile(figname):\n",
751767
" # os.remove(figname)\n",
752768
" #plt.savefig(figname)\n",
753769
"\n",
754-
"\n",
755770
"plt.tight_layout()\n",
756-
"plt.show()\n",
757-
"\n",
758-
"#for i in range(len(num_visits)):\n",
759-
"\n",
760-
"# "
771+
"plt.show()"
761772
]
762773
},
763774
{
@@ -787,28 +798,27 @@
787798
"outputs": [],
788799
"source": [
789800
"for i in range(num_visits):\n",
790-
" #ax.plot([0, 1, 2], [i, i + 1, i + 2])\n",
791801
" j = i #+ 5\n",
792802
" plt.title(f\"MJD {np.round(scitab['t_max'][j], 4)}\")\n",
793803
"\n",
794-
" dl_result_sci = DatalinkResults.from_result_url(scitab['access_url'][j], session=get_pyvo_auth())\n",
795-
" #print(dl_result_sci)\n",
804+
" dl_result_sci = DatalinkResults.from_result_url(scitab['access_url'][j],\n",
805+
" session=get_pyvo_auth())\n",
806+
"\n",
796807
" sci, scihdul = get_cutout(dl_result_sci, spherePoint, get_pyvo_auth(), fov)\n",
797808
" sci_hdu = scihdul[1]\n",
798809
" sci_header = scihdul[1].header\n",
799810
" scidata = scihdul[1].data\n",
800811
" sciwcs = WCS(scihdul[1].header)\n",
801812
" reproj_scidata, footprint = reproject_interp((scidata, sciwcs), refwcs)\n",
802813
"\n",
803-
" plt.imshow(reproj_scidata, origin='lower', vmin=np.nanpercentile(reproj_scidata, 1), vmax=np.nanpercentile(reproj_scidata, 99))\n",
804-
"\n",
814+
" plt.imshow(reproj_scidata, origin='lower', vmin=np.nanpercentile(reproj_scidata, 1),\n",
815+
" vmax=np.nanpercentile(reproj_scidata, 99))\n",
805816
"\n",
806817
" figname = os.path.join(tempdir, 'cutout_' + str(i) + '.png')\n",
807818
" if os.path.isfile(figname):\n",
808819
" os.remove(figname)\n",
809820
" plt.savefig(figname)\n",
810821
"\n",
811-
"\n",
812822
"#plt.tight_layout()\n",
813823
"#plt.show()"
814824
]

0 commit comments

Comments
 (0)