-
Notifications
You must be signed in to change notification settings - Fork 5
Commit
This commit does not belong to any branch on this repository, and may belong to a fork outside of the repository.
Merge pull request #166 from LSSTDESC/tqz/add_blending_notebook
add blending notebook
- Loading branch information
Showing
1 changed file
with
283 additions
and
0 deletions.
There are no files selected for viewing
283 changes: 283 additions & 0 deletions
283
examples/creation_examples/blending_degrader_demo.ipynb
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Original file line number | Diff line number | Diff line change |
---|---|---|
@@ -0,0 +1,283 @@ | ||
{ | ||
"cells": [ | ||
{ | ||
"cell_type": "markdown", | ||
"id": "2610a0f0-0c71-4401-896f-734442bcd66d", | ||
"metadata": {}, | ||
"source": [ | ||
"## Blending Degrader demo\n", | ||
"\n", | ||
"author: Shuang Liang\n", | ||
"\n", | ||
"This notebook demonstrate the use of `rail.creation.degradation.unrec_bl_model`, which uses Friends of Friends to finds sources close to each other and merge them into unrecognized blends" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "f7a6adc3-68e8-4a1d-842f-bfb0960a1c4a", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"from rail.creation.degraders.unrec_bl_model import UnrecBlModel\n", | ||
"\n", | ||
"from rail.core.data import PqHandle\n", | ||
"from rail.core.stage import RailStage\n", | ||
"\n", | ||
"import matplotlib.pyplot as plt\n", | ||
"import pandas as pd, numpy as np" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "6912a740-31ea-4987-b06d-81ff17cd895a", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"DS = RailStage.data_store\n", | ||
"DS.__class__.allow_overwrite = True\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "markdown", | ||
"id": "a282c2ed-141b-4507-8254-dc8fbc9864dc", | ||
"metadata": {}, | ||
"source": [ | ||
"### Create a random catalog with ugrizy+YJHF bands as the the true input" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "1078bc2a-fc54-41c3-bd30-6c447bb971d4", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"data = np.random.normal(23, 3, size = (1000,12))\n", | ||
"data[:, 0] = np.random.uniform(low=0, high=0.03, size=1000)\n", | ||
"data[:, 1] = np.random.uniform(low=0, high=0.03, size=1000)\n", | ||
"\n", | ||
"data_df = pd.DataFrame(data=data, # values\n", | ||
" columns=['ra', 'dec', 'u', 'g', 'r', 'i', 'z', 'y', 'Y', 'J', 'H', 'F'])\n", | ||
"\n", | ||
"data_truth_handle = DS.add_data('input', data_df, PqHandle)\n", | ||
"data_truth = data_truth_handle.data" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "33c99a4d-8375-4003-9a9a-70fa85a3eb82", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"#data_df.to_parquet('bl_test.pq')" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "a5636721-a734-4746-bd93-8101bc558b6e", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"plt.scatter(data_truth['ra'], data_truth['dec'], s=5)\n", | ||
"plt.xlabel(\"Ra [Deg]\", fontsize=14)\n", | ||
"plt.ylabel(\"Dec [Deg]\", fontsize=14)\n", | ||
"plt.show()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "markdown", | ||
"id": "1da27deb-d167-4f38-8c59-f270184d6ab3", | ||
"metadata": {}, | ||
"source": [ | ||
"### The blending model" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "a07f72a0-e24c-4844-90f0-d5a49ac4362b", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"## model configuration; linking length is in arcsecs\n", | ||
"\n", | ||
"blModel = UnrecBlModel.make_stage(name='unrec_bl_model', ra_label='ra', dec_label='dec', linking_lengths=1.0, \\\n", | ||
" bands='ugrizy')\n", | ||
"blModel.get_config_dict()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "e5f4862a-0621-46d4-8901-7e84b461c424", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"# run the model\n", | ||
"\n", | ||
"outputs = blModel(data_truth)\n", | ||
"\n", | ||
"samples_w_bl = outputs['output'].data\n", | ||
"component_ind = outputs['compInd'].data\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "bc5158dd-f474-4731-b847-b4a7358656b9", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"fig, ax = plt.subplots(figsize=(6, 5), dpi=100)\n", | ||
"\n", | ||
"ax.scatter(data_truth['ra'], data_truth['dec'], s=10, facecolors='none', edgecolors='b', label='Original')\n", | ||
"ax.scatter(samples_w_bl['ra'], samples_w_bl['dec'], s=5, c='r', label='w. Unrec-BL')\n", | ||
"\n", | ||
"ax.legend(loc=2, fontsize=12)\n", | ||
"ax.set_xlabel(\"Ra [Deg]\", fontsize=14)\n", | ||
"ax.set_ylabel(\"Dec [Deg]\", fontsize=14)\n", | ||
"\n", | ||
"plt.show()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "268b3d37-b7fd-4ac1-8457-2104a87c9e6d", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"b = 'r'\n", | ||
"plt.hist(data_truth[b], bins=np.linspace(10, 30, 20), label='Original')\n", | ||
"plt.hist(samples_w_bl[b], bins=np.linspace(10, 30, 20), fill=False, label='w. Unrec-BL')\n", | ||
"\n", | ||
"plt.xlabel(fr'Magnitude ${b}$', fontsize=14)\n", | ||
"plt.legend(fontsize=12)\n", | ||
"plt.show()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "a1d51c15-1e04-4b22-9abb-9b267965dbeb", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"flux = 10**(-(data_truth[b]-28.10)/2.5) # r band zp for lsst is 28.10\n", | ||
"flux_bl = 10**(-(samples_w_bl[b]-28.10)/2.5)\n", | ||
"\n", | ||
"plt.hist(flux, bins=np.linspace(0, 10000, 40), label='Original')\n", | ||
"plt.hist(flux_bl, bins=np.linspace(0, 10000, 40), fill=False, label='w. Unrec-BL')\n", | ||
"\n", | ||
"plt.xlabel(fr'Flux ${b}$', fontsize=14)\n", | ||
"plt.yscale('log')\n", | ||
"plt.legend(fontsize=12)\n", | ||
"plt.show()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "markdown", | ||
"id": "f3ba003e-da62-4bfc-b70e-c07c1112efc0", | ||
"metadata": {}, | ||
"source": [ | ||
"### Study one BL case" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "4d3fbd87-b227-43bf-b712-e8a069b51a54", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"## find a source with more than 1 truth component\n", | ||
"\n", | ||
"group_size = 1\n", | ||
"while group_size==1:\n", | ||
"\n", | ||
" rand_ind = np.random.randint(len(samples_w_bl))\n", | ||
" this_bl = samples_w_bl.iloc[rand_ind]\n", | ||
" group_id = this_bl['group_id']\n", | ||
" \n", | ||
" FoF_group = component_ind.query(f\"group_id == {group_id}\")\n", | ||
" group_size = len(FoF_group)\n", | ||
"\n", | ||
"truth_comp = data_truth.iloc[FoF_group.index]\n", | ||
"\n", | ||
"print('Truth RA / Blended RA:')\n", | ||
"print(truth_comp['ra'].to_numpy(), '/', this_bl['ra'])\n", | ||
"print(\"\")\n", | ||
"\n", | ||
"print('Truth DEC / Blended DEC:')\n", | ||
"print(truth_comp['dec'].to_numpy(), '/', this_bl['dec'])\n", | ||
"print(\"\")\n", | ||
"\n", | ||
"for b in 'ugrizy':\n", | ||
" print(f'Truth mag {b} / Blended mag {b}:')\n", | ||
" print(truth_comp[b].to_numpy(), '/', this_bl[b])\n", | ||
" print(\"\")\n", | ||
" " | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "8dacb910-dd26-404f-ba61-4278094b6355", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [ | ||
"\n", | ||
"fig, ax = plt.subplots(figsize=(6, 5), dpi=100)\n", | ||
"\n", | ||
"ax.scatter(this_bl['ra']*3600, this_bl['dec']*3600, s=1e4, c='r')\n", | ||
"ax.scatter(truth_comp['ra']*3600, truth_comp['dec']*3600, s=1e4, facecolors='none', edgecolors='b')\n", | ||
"\n", | ||
"ax.scatter([], [], s=1e2, facecolors='none', edgecolors='b', label='Truth Components')\n", | ||
"ax.scatter([], [], s=1e2, c='r', label='Merged Source')\n", | ||
"\n", | ||
"fig_size = 1 ## in arcsecs\n", | ||
"ax.set_xlim(this_bl['ra']*3600-fig_size, this_bl['ra']*3600+fig_size)\n", | ||
"ax.set_ylim(this_bl['dec']*3600-fig_size, this_bl['dec']*3600+fig_size)\n", | ||
"\n", | ||
"ax.legend(fontsize=12)\n", | ||
"ax.set_xlabel(\"Ra [arcsecs]\", fontsize=14)\n", | ||
"ax.set_ylabel(\"Dec [arcsecs]\", fontsize=14)\n", | ||
"\n", | ||
"plt.show()\n" | ||
] | ||
}, | ||
{ | ||
"cell_type": "code", | ||
"execution_count": null, | ||
"id": "5fc4b38b-55d1-43ff-9039-ee9c49c54f4d", | ||
"metadata": {}, | ||
"outputs": [], | ||
"source": [] | ||
} | ||
], | ||
"metadata": { | ||
"kernelspec": { | ||
"display_name": "Python 3 (ipykernel)", | ||
"language": "python", | ||
"name": "python3" | ||
}, | ||
"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.3" | ||
} | ||
}, | ||
"nbformat": 4, | ||
"nbformat_minor": 5 | ||
} |