Synthetic Source Injection with LSST Code

notebooks/ssi_with_slslim.ipynb

@@ -0,0 +1,302 @@
+{
+ "cells": [
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "4e80a982-6590-457a-ab5d-c5727c17847f",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import os\n",
+    "import sys\n",
+    "import numpy as np\n",
+    "import subprocess\n",
+    "from astropy.table import Table\n",
+    "from astropy.time import Time\n",
+    "import lsst.daf.butler as dafButler\n",
+    "from lsst.pipe.base import Pipeline\n",
+    "\n",
+    "%matplotlib ipympl"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "aaeba882-202a-48bc-82e6-6b28e348ebfb",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import lsst.afw.display as afwDisplay\n",
+    "\n",
+    "afwDisplay.setDefaultBackend(\"firefly\")\n",
+    "display1 = afwDisplay.Display(frame=1)\n",
+    "display2 = afwDisplay.Display(frame=2)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "b9af51d8-42fd-4125-b210-e335d0e41ec3",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Catalog generation parameters\n",
+    "ra_min, ra_max = 52, 54\n",
+    "dec_min, dec_max = -29, -27\n",
+    "n_galaxies = 1000\n",
+    "sky_area_value = 5.0\n",
+    "\n",
+    "# Butler configuration\n",
+    "repo_path = \"/repo/main\"\n",
+    "instrument = \"LSSTComCam\"\n",
+    "\n",
+    "# Specific exposure and detector to process\n",
+    "visit = 2024110800246\n",
+    "detector = 5\n",
+    "\n",
+    "# Collection with pvi images\n",
+    "collection = \"LSSTComCam/runs/DRP/DP1/w_2025_10/DM-49359\"\n",
+    "\n",
+    "# Collections for processing results\n",
+    "inject_collection = \"u/rea3/slsim/injection_test10_calexp\"\n",
+    "output_collection = \"u/rea3/slsim/test10_calexp\"\n",
+    "\n",
+    "# Filter visits to use for ingestion\n",
+    "filters = [\"r\"]\n",
+    "\n",
+    "# Time based selection of visits\n",
+    "start_time = \"2024-11-01T00:00:00\"\n",
+    "end_time = \"2024-12-29T00:00:00\"\n",
+    "\n",
+    "ngal = 10000\n",
+    "sky_area = 0.5\n",
+    "output_file = \"catalog.fits\"\n",
+    "\n",
+    "# path to slsim\n",
+    "slsim_path = \"../slsim/slsim/LsstSciencePipeline/\""
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "c1395cc9-cf4c-4698-a255-3d0498dde63d",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "gen_cmd = (\n",
+    "    f\"python {slsim_path}/generate_catalog.py --ra-min {ra_min} --ra-max {ra_max} \"\n",
+    "    f\" --dec-min {dec_min} --dec-max {dec_max} --start-time {start_time} --end-time {end_time}\"\n",
+    "    f\" --n-galaxies {ngal} --sky-area {sky_area} --repo-path {repo_path}\"\n",
+    "    f\" --collection {collection} --instrument {instrument}\"\n",
+    "    f\" --filters {\" \".join(filters)} --output-file {output_file}\"\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "39d287e6-4cb5-4fb8-9509-1d8e157fb7a4",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# os.system(gen_cmd)\n",
+    "print(gen_cmd)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "9ecece51-459b-4726-ae2d-b89f5ba4dd6a",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Use the ingest_injection_catalog command line tool to ingest the previously created catalog\n",
+    "injest_cmd = (\n",
+    "    f\"ingest_injection_catalog -b {repo_path} -o {inject_collection}\"\n",
+    "    f\" --injection-catalog {output_file} {\" \".join(filters)} -t injection_slsim\"\n",
+    "    f\" --format fits\"\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "4645d58b-be42-42aa-ad83-08a6af5eb6a0",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# os.system(injest_cmd)\n",
+    "print(injest_cmd)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "00db276e-e6a2-4b9b-a56e-c2b108d52c84",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# Create a pipeline configuration file for SLSim injection with DRP processing\n",
+    "pipeline_yaml = \"config_wfake.yaml\"\n",
+    "\n",
+    "with open(pipeline_yaml, \"w\") as f:\n",
+    "    f.write(\n",
+    "        \"\"\"\n",
+    "description: Strong lensing source injection pipeline with DRP processing\n",
+    "instrument: lsst.obs.lsst.LsstComCam\n",
+    "\n",
+    "tasks:\n",
+    "  inject_visit_slsim:\n",
+    "    class: slsim.LsstSciencePipeline.inject_slsim_visit.VisitInjectSLSimTask\n",
+    "    config:\n",
+    "      external_psf: false\n",
+    "      external_photo_calib: false\n",
+    "      external_wcs: false\n",
+    "      connections.input_exposure: pvi\n",
+    "      connections.output_exposure: injected_slsim_pvi\n",
+    "      selection: \"np.isin(injection_catalog['visit'], {visit})\"\n",
+    "  subtractImages:\n",
+    "    class: lsst.ip.diffim.subtractImages.AlardLuptonSubtractTask\n",
+    "    config:\n",
+    "    - connections.coaddName: goodSeeing\n",
+    "      connections.template: goodSeeingDiff_templateExp\n",
+    "      connections.science: injected_slsim_pvi\n",
+    "      connections.difference: injected_slsim_goodSeeingDiff_differenceTempExp\n",
+    "      connections.matchedTemplate: injected_slsim_goodSeeingDiff_matchedExp\n",
+    "  detectAndMeasureDiaSources:\n",
+    "    class: lsst.ip.diffim.detectAndMeasure.DetectAndMeasureTask\n",
+    "    config:\n",
+    "    - connections.coaddName: goodSeeing\n",
+    "      connections.diaSources: injected_slsim_goodSeeingDiff_diaSrc\n",
+    "      connections.subtractedMeasuredExposure: injected_slsim_goodSeeingDiff_differenceExp\n",
+    "      connections.science: injected_slsim_pvi\n",
+    "      connections.difference: injected_slsim_goodSeeingDiff_differenceTempExp\n",
+    "      connections.matchedTemplate: injected_slsim_goodSeeingDiff_matchedExp\n",
+    "\n",
+    "subsets:\n",
+    "  injected_DRP:\n",
+    "    subset:\n",
+    "      - inject_visit_slsim\n",
+    "      - subtractImages\n",
+    "      - detectAndMeasureDiaSources\n",
+    "    description: >\n",
+    "      Pipeline for injecting strong lensing sources and running difference imaging\n",
+    "\"\"\"\n",
+    "    )\n",
+    "\n",
+    "print(f\"Created pipeline configuration in {pipeline_yaml}\")\n",
+    "print(\n",
+    "    \"The pipeline includes strong lensing injection followed by subtraction and detection tasks, exactly matching test.yaml\"\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "26f8ff66-f605-4e15-a8f8-681ff01475b1",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "run_cmd = (\n",
+    "    f\"pipetask run -b {repo_path}\"\n",
+    "    f' -d \"exposure={visit} and detector={detector}\" -p {pipeline_yaml}#injected_DRP'\n",
+    "    f\" -i {collection},{inject_collection}\"\n",
+    "    f\" -o {output_collection} --register-dataset-types\"\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "182106f5-287a-4767-9f6b-78aced648ccf",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "# os.system(run_cmd)\n",
+    "print(run_cmd)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "46c94c87-6352-4b98-bfcd-4429af277c35",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "butler = dafButler.Butler(repo_path, collections=output_collection)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "6ec30488-5a32-4d11-a72f-743d4bb87bed",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "injected_pvi = butler.get(\"injected_slsim_pvi\", visit=visit, detector=detector)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "2fad76e0-b704-4318-8f46-8e1df3b9ffb7",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "diff_pvi = butler.get(\n",
+    "    \"injected_slsim_goodSeeingDiff_differenceTempExp\", visit=visit, detector=detector\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "91fee920-382b-4951-98dc-3abedad7301b",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "display1.mtv(injected_pvi)\n",
+    "display2.mtv(diff_pvi)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "d29c7a4f-e101-42e1-ad57-d5ae07292f97",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "output_collection"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "fbfb23d6-7c95-45d9-b008-f91723c605b0",
+   "metadata": {},
+   "outputs": [],
+   "source": []
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "LSST",
+   "language": "python",
+   "name": "lsst"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.12.9"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}