|
96 | 96 | "import matplotlib.pyplot as plt\n", |
97 | 97 | "\n", |
98 | 98 | "from lsst.daf.butler import Butler\n", |
99 | | - "\n", |
100 | 99 | "from lsst.rsp import get_tap_service\n", |
101 | 100 | "import lsst.geom as geom\n", |
102 | | - "\n", |
103 | 101 | "import lsst.afw.display as afwDisplay\n", |
104 | 102 | "from lsst.afw.image import MultibandExposure\n", |
105 | 103 | "\n", |
106 | 104 | "import lsst.scarlet.lite as sl\n", |
107 | | - "import lsst.meas.extensions.scarlet as mes\n" |
| 105 | + "import lsst.meas.extensions.scarlet as mes" |
108 | 106 | ] |
109 | 107 | }, |
110 | 108 | { |
|
176 | 174 | "id": "16eeb034-6728-406e-bfc8-7b875c3a9416", |
177 | 175 | "metadata": {}, |
178 | 176 | "source": [ |
179 | | - "## 2. Identify a blend¶\n", |
| 177 | + "## 2. Identify a blend\n", |
180 | 178 | "\n", |
181 | 179 | "This first section will identify a blend in DP1, query the `Object` table for the deblended child objects, and visualize the blend.\n" |
182 | 180 | ] |
|
223 | 221 | " \"FROM dp1.Object \" + \\\n", |
224 | 222 | " \"WHERE CONTAINS(POINT('ICRS', coord_ra, coord_dec ), \" + \\\n", |
225 | 223 | " \"CIRCLE('ICRS', \" + str(ra) + \", \" + str(dec) + \", .1)) = 1 \" + \\\n", |
226 | | - " \"AND parentObjectId = \" + str(parentId)\n" |
| 224 | + " \"AND parentObjectId = \" + str(parentId)" |
227 | 225 | ] |
228 | 226 | }, |
229 | 227 | { |
|
280 | 278 | "id": "4288252b-6639-42de-82aa-70db8252f284", |
281 | 279 | "metadata": {}, |
282 | 280 | "source": [ |
283 | | - "To visualize the blend, first define a function to generate an image cutout. " |
| 281 | + "To visualize the blend, plot a `deep_coadd` image subset that contains the blend with children identified. First, find the patch and tract containing the blend, and set the band to be i-band." |
284 | 282 | ] |
285 | 283 | }, |
286 | 284 | { |
287 | 285 | "cell_type": "code", |
288 | 286 | "execution_count": null, |
289 | | - "id": "a1bf0635-8a88-42ae-b48d-79e69213f710", |
| 287 | + "id": "245bec12-edcf-4564-811e-9de0d06fd3e5", |
290 | 288 | "metadata": {}, |
291 | 289 | "outputs": [], |
292 | 290 | "source": [ |
293 | | - "def cutout_coadd(butler, ra, dec, band='i', datasetType='deepCoadd',\n", |
294 | | - " skymap=None, cutoutSideLength=51, **kwargs):\n", |
295 | | - " \"\"\"\n", |
296 | | - " Produce a cutout from a coadd at the given ra, dec position.\n", |
297 | | - "\n", |
298 | | - " Adapted from DC2 tutorial notebook by Michael Wood-Vasey.\n", |
299 | | - "\n", |
300 | | - " Parameters\n", |
301 | | - " ----------\n", |
302 | | - " butler: lsst.daf.persistence.Butler\n", |
303 | | - " Helper object providing access to a data repository\n", |
304 | | - " ra: float\n", |
305 | | - " Right ascension of the center of the cutout, in degrees\n", |
306 | | - " dec: float\n", |
307 | | - " Declination of the center of the cutout, in degrees\n", |
308 | | - " band: string\n", |
309 | | - " Filter of the image to load\n", |
310 | | - " datasetType: string ['deepCoadd']\n", |
311 | | - " Which type of coadd to load. Doesn't support 'calexp'\n", |
312 | | - " skymap: lsst.afw.skyMap.SkyMap [optional]\n", |
313 | | - " Pass in to avoid the Butler read. Useful if you have lots of them.\n", |
314 | | - " cutoutSideLength: float [optional]\n", |
315 | | - " Size of the cutout region in pixels.\n", |
316 | | - "\n", |
317 | | - " Returns\n", |
318 | | - " -------\n", |
319 | | - " MaskedImage\n", |
320 | | - " \"\"\"\n", |
321 | | - " radec = geom.SpherePoint(ra, dec, geom.degrees)\n", |
322 | | - " cutoutSize = geom.ExtentI(cutoutSideLength, cutoutSideLength)\n", |
323 | | - "\n", |
324 | | - " if skymap is None:\n", |
325 | | - " skymap = butler.get(\"skyMap\")\n", |
| 291 | + "radec = geom.SpherePoint(children['coord_ra'][0],children['coord_dec'][0],geom.degrees)\n", |
| 292 | + "skymap = butler.get(\"skyMap\")\n", |
326 | 293 | "\n", |
327 | | - " tractInfo = skymap.findTract(radec)\n", |
328 | | - " patchInfo = tractInfo.findPatch(radec)\n", |
329 | | - " xy = geom.PointI(tractInfo.getWcs().skyToPixel(radec))\n", |
330 | | - " bbox = geom.BoxI(xy - cutoutSize // 2, cutoutSize)\n", |
331 | | - " print(\"bbox = \", bbox, \"xy = \", xy, \"cutoutSize = \", cutoutSize)\n", |
332 | | - " patch = tractInfo.getSequentialPatchIndex(patchInfo)\n", |
333 | | - " tract = tractInfo.getId()\n", |
| 294 | + "tractInfo = skymap.findTract(radec)\n", |
334 | 295 | "\n", |
335 | | - " parameters = {'bbox': bbox}\n", |
336 | | - "\n", |
337 | | - " query = \"\"\"band.name = '{}' AND patch = {} AND tract = {}\n", |
338 | | - " \"\"\".format(band, patch, tract)\n", |
339 | | - " print(query)\n", |
340 | | - "\n", |
341 | | - " dataset_refs = butler.query_datasets(\"deep_coadd\", where=query)\n", |
342 | | - "\n", |
343 | | - " cutout_image = butler.get(dataset_refs[0], parameters=parameters)\n", |
344 | | - "\n", |
345 | | - " return cutout_image" |
| 296 | + "tract = tractInfo.getId()\n", |
| 297 | + "patch = tractInfo.findPatch(radec).getSequentialIndex()\n", |
| 298 | + "band='i'" |
346 | 299 | ] |
347 | 300 | }, |
348 | 301 | { |
349 | 302 | "cell_type": "markdown", |
350 | | - "id": "5546408c-6009-4230-9ef4-2faf1c302940", |
| 303 | + "id": "80f1ea90-55d5-43b4-b6da-0256baf8702c", |
351 | 304 | "metadata": {}, |
352 | 305 | "source": [ |
353 | | - "Make a cutout of the `deep_coadd` imaging, using the `cutout_coadd` function to visualize blend and its child objects. Set the size of the cutout to be 5x the square root of the total `footprintArea` of the blend.\n", |
354 | | - " " |
| 306 | + "Next, define the dimensions of the image subset to be 5x the square root of the total `footprintArea` of the blend." |
355 | 307 | ] |
356 | 308 | }, |
357 | 309 | { |
358 | 310 | "cell_type": "code", |
359 | 311 | "execution_count": null, |
360 | | - "id": "b2f6c287-7502-4343-89e3-e76cb3e98e87", |
| 312 | + "id": "f845cc95-c731-40d2-a202-eea9e2681afe", |
361 | 313 | "metadata": {}, |
362 | 314 | "outputs": [], |
363 | 315 | "source": [ |
364 | 316 | "total_footprint_area = np.sum(children['footprintArea'])\n", |
365 | | - "\n", |
366 | | - "cutout = cutout_coadd(butler, children['coord_ra'][0],\n", |
367 | | - " children['coord_dec'][0],\n", |
368 | | - " band='i', datasetType='deepCoadd',\n", |
369 | | - " cutoutSideLength=5.0\n", |
370 | | - " * np.sqrt(total_footprint_area))" |
| 317 | + "cutoutSideLength=5.0 * np.sqrt(total_footprint_area)" |
| 318 | + ] |
| 319 | + }, |
| 320 | + { |
| 321 | + "cell_type": "markdown", |
| 322 | + "id": "a1922530-e638-4514-bce0-7087f6696c32", |
| 323 | + "metadata": {}, |
| 324 | + "source": [ |
| 325 | + "Define the bounding box for the image subset. " |
| 326 | + ] |
| 327 | + }, |
| 328 | + { |
| 329 | + "cell_type": "code", |
| 330 | + "execution_count": null, |
| 331 | + "id": "ab876741-2e6a-491e-bdec-de5ca72646bf", |
| 332 | + "metadata": {}, |
| 333 | + "outputs": [], |
| 334 | + "source": [ |
| 335 | + "cutoutSize = geom.ExtentI(cutoutSideLength, cutoutSideLength)\n", |
| 336 | + "xy = geom.PointI(tractInfo.getWcs().skyToPixel(radec))\n", |
| 337 | + "bbox = geom.BoxI(xy - cutoutSize // 2, cutoutSize)\n", |
| 338 | + "parameters = {'bbox': bbox}" |
| 339 | + ] |
| 340 | + }, |
| 341 | + { |
| 342 | + "cell_type": "markdown", |
| 343 | + "id": "0ffb1405-2389-41a4-b91c-3a65729f21fb", |
| 344 | + "metadata": {}, |
| 345 | + "source": [ |
| 346 | + "Finally, send the bounding box to the butler and request the image subset from the `deep_coadd` to visualize blend and its child objects. " |
| 347 | + ] |
| 348 | + }, |
| 349 | + { |
| 350 | + "cell_type": "code", |
| 351 | + "execution_count": null, |
| 352 | + "id": "8fdd4188-33ae-4b30-8513-52d82e0dce9d", |
| 353 | + "metadata": {}, |
| 354 | + "outputs": [], |
| 355 | + "source": [ |
| 356 | + "query = \"\"\"band.name = '{}' AND patch = {} AND tract = {}\n", |
| 357 | + " \"\"\".format(band, patch, tract)\n", |
| 358 | + "dataset_refs = butler.query_datasets(\"deep_coadd\", where=query)\n", |
| 359 | + "cutout = butler.get(dataset_refs[0], parameters=parameters)" |
371 | 360 | ] |
372 | 361 | }, |
373 | 362 | { |
|
428 | 417 | "registry = butler.registry\n", |
429 | 418 | "\n", |
430 | 419 | "for dt in sorted(registry.queryDatasetTypes('*scarlet*')):\n", |
431 | | - " print(dt)\n" |
| 420 | + " print(dt)" |
432 | 421 | ] |
433 | 422 | }, |
434 | 423 | { |
|
448 | 437 | "source": [ |
449 | 438 | "query = \"\"\"band.name = '{}' AND patch.region OVERLAPS POINT({}, {})\n", |
450 | 439 | " \"\"\".format('i', children['coord_ra'][0], children['coord_dec'][0])\n", |
451 | | - "print(query)\n" |
| 440 | + "print(query)" |
452 | 441 | ] |
453 | 442 | }, |
454 | 443 | { |
|
467 | 456 | "outputs": [], |
468 | 457 | "source": [ |
469 | 458 | "dataset_refs = butler.query_datasets(\"deep_coadd\", where=query)\n", |
470 | | - "dataId = dataset_refs[0].dataId\n" |
| 459 | + "dataId = dataset_refs[0].dataId" |
471 | 460 | ] |
472 | 461 | }, |
473 | 462 | { |
|
533 | 522 | "from numpy.typing import DTypeLike\n", |
534 | 523 | "from lsst.scarlet.lite import Box, Blend, Observation\n", |
535 | 524 | "\n", |
536 | | - "\n", |
537 | 525 | "def minimal_data_to_blend(cls, model_psf: np.ndarray, dtype: DTypeLike) -> Blend:\n", |
538 | 526 | "\n", |
539 | 527 | " \"\"\"Convert the storage data model into a scarlet lite blend\n", |
|
0 commit comments