Skip to content

Eac/mst matching #27

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
eacharles merged 5 commits into from
Nov 4, 2026
Merged

Eac/mst matching #27

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
8196a8a
Add MST-based alternative initial grouper
eacharles Nov 4, 2026
48950c7
Add RunMatchingMST.ipynb — MST-matcher version of RunMatching
eacharles Nov 4, 2026
f77c6a1
Replace pixel-zoom splitting with MST-native subdivision
eacharles Nov 4, 2026
4f61e07
Add showMSTCluster and showMSTObjects visualization functions
eacharles Nov 4, 2026
60f752f
Add VizMatchesMST.ipynb using MST visualization functions
eacharles Nov 4, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
388 changes: 388 additions & 0 deletions examples/rubin_dp2/RunMatchingMST.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,388 @@
{
"cells": [
{
"cell_type": "markdown",
"id": "60f22d4f-5b19-4316-9636-2eb93d7b8bbb",
"metadata": {},
"source": [
"# (re-)Run the matchers — MST version\n",
"\n",
"This is a drop-in replacement for `RunMatching.ipynb` that uses the\n",
"MST-based matchers (`MSTShearMatch`, `MSTWcsMatch`) instead of the\n",
"footprint-based ones.\n",
"\n",
"The key difference is that initial source grouping is done by\n",
"single-linkage clustering at a configurable angular separation\n",
"(`mst_cut`, in arcseconds) rather than by detecting connected pixels\n",
"in a pixelised counts map. This removes the dependence on pixel scale\n",
"choices for the initial grouping step; the subsequent hierarchical\n",
"confusion-resolution is unchanged.\n",
"\n",
"**Tuning parameters added here vs `RunMatching.ipynb`:**\n",
"- `mst_cut_shear` — grouping radius for shear matchers (arcsec)\n",
"- `mst_cut_wcs` — grouping radius for WCS meta-matcher (arcsec)"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "5a074eb6-d9aa-40e1-a442-52262d0fd08b",
"metadata": {},
"outputs": [],
"source": [
"# Standard imports\n",
"import tables_io\n",
"import hpmcm\n",
"import matplotlib.pyplot as plt\n",
"import numpy as np"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "e03fee38-b519-4708-a051-82b9b6c82db1",
"metadata": {},
"outputs": [],
"source": [
"# Set the tract\n",
"tract = 8006\n",
"# This will need to change with tract\n",
"ref_dir=(340.5000000000, -9.6694214876)\n",
"our_geom = hpmcm.shear_utils.ShearCellGeometry(\n",
" ref_dir=ref_dir,\n",
")"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "eddb05fe-f39d-4de4-b040-81a6116c3356",
"metadata": {},
"outputs": [],
"source": [
"# If this is True only run on a single cell\n",
"test_only = False"
]
},
{
"cell_type": "code",
"execution_count": 4,
"id": "a1b2c3d4-e5f6-7890-abcd-ef1234567890",
"metadata": {},
"outputs": [],
"source": [
"# MST grouping radii (arcseconds)\n",
"# mst_cut_shear: for shear catalogs (0.2 arcsec/pixel native scale)\n",
"# comparable to pixel_match_scale=2 footprint grouping (~0.4-0.8 arcsec)\n",
"mst_cut_shear = 0.8\n",
"\n",
"# mst_cut_wcs: for the WCS meta-matcher (0.2 arcsec/pixel, pixel_match_scale=4)\n",
"# comparable to pixel_match_scale=4 footprint grouping (~0.8-1.6 arcsec)\n",
"mst_cut_wcs = 1.5"
]
},
{
"cell_type": "code",
"execution_count": 5,
"id": "3edac63a-e803-48a0-8ad1-cb09c5e76b27",
"metadata": {},
"outputs": [],
"source": [
"# Files and catalogs ids for rubin metadetect\n",
"rubin_shear_match_files = [\n",
" f'{tract}/shear_object_reduced_uncleaned_{tract}_ns.parq',\n",
" f'{tract}/shear_object_reduced_uncleaned_{tract}_1m.parq',\n",
" f'{tract}/shear_object_reduced_uncleaned_{tract}_1p.parq',\n",
" f'{tract}/shear_object_reduced_uncleaned_{tract}_2m.parq',\n",
" f'{tract}/shear_object_reduced_uncleaned_{tract}_2p.parq', \n",
"]\n",
"rubin_shear_match_ids = list(range(len(rubin_shear_match_files)))\n",
"\n",
"# Files and catalogs ids for desc metadetect\n",
"desc_shear_match_files = [\n",
" f'{tract}/mdet-run-reduced_uncleaned_{tract}_ns.parq',\n",
" f'{tract}/mdet-run-reduced_uncleaned_{tract}_1m.parq',\n",
" f'{tract}/mdet-run-reduced_uncleaned_{tract}_1p.parq',\n",
"]\n",
"desc_shear_match_ids = list(range(len(desc_shear_match_files)))\n",
"\n",
"# Files and catalogs ids for meta-matching\n",
"all_files = [\n",
" f'{tract}/object_reduced_{tract}.parq',\n",
" f'{tract}/anacal_reduced_{tract}.parq',\n",
" f'{tract}/desc_mst_md_match_object_joined.pq',\n",
" f'{tract}/rubin_mst_md_match_object_joined.pq',\n",
"]\n",
"all_ids = list(range(len(all_files)))"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "ffd22736-c485-4bf7-acc9-6c0736513b7d",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"......... 10!\n",
"......... 20!\n",
"......... 30!\n",
"......... 40!\n",
".."
]
}
],
"source": [
"# Allow for re-running various matchers\n",
"rematch_rubin_shear = True\n",
"rematch_desc_shear = True\n",
"rematch_objects = True\n",
"\n",
"if rematch_rubin_shear:\n",
" rubin_shear_match = hpmcm.MSTShearMatch.createShearMatch(\n",
" geometry=our_geom,\n",
" mst_cut=mst_cut_shear,\n",
" pixel_match_scale=2,\n",
" pixel_r2_cut=16,\n",
" deshear=-0.01,\n",
" shear_names=['ns', '1m', '1p', '2m', '2p'], \n",
" )\n",
" rubin_shear_match.reduceData(rubin_shear_match_files, rubin_shear_match_ids)\n",
" if test_only:\n",
" rubin_shear_match.analysisLoop(range(100, 101), range(100, 101))\n",
" else:\n",
" rubin_shear_match.analysisLoop()\n",
" rubin_stats = rubin_shear_match.extractStats()\n",
" rubin_stats.update(**rubin_shear_match.extractShearStats())\n",
" tables_io.write(rubin_stats, f'{tract}/rubin_mst_md_match_')\n",
" rubin_shear_match.save(f'{tract}/rubin_mst_md_match')\n",
" \n",
"\n",
"if rematch_desc_shear:\n",
" desc_shear_match = hpmcm.MSTShearMatch.createShearMatch(\n",
" geometry=our_geom,\n",
" mst_cut=mst_cut_shear,\n",
" deshear=-0.01,\n",
" pixel_match_scale=2,\n",
" pixel_r2_cut=16, \n",
" shear_names=['ns', '1m', '1p'],\n",
" )\n",
" desc_shear_match.reduceData(desc_shear_match_files, desc_shear_match_ids)\n",
" if test_only: \n",
" desc_shear_match.analysisLoop(range(100, 101), range(100, 101))\n",
" else:\n",
" desc_shear_match.analysisLoop()\n",
" desc_stats = desc_shear_match.extractStats()\n",
" desc_stats.update(**desc_shear_match.extractShearStats())\n",
" tables_io.write(desc_stats, f'{tract}/desc_mst_md_match_') \n",
" desc_shear_match.save(f'{tract}/desc_mst_md_match')\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "608b0387-feb2-42e1-aa47-bf93ceb39030",
"metadata": {},
"outputs": [],
"source": [
"shear_rubin_obj_stats = tables_io.read(f'{tract}/rubin_mst_md_match_object_stats.pq')\n",
"shear_rubin_obj_assoc = tables_io.read(f'{tract}/rubin_mst_md_match_object_assoc.pq')\n",
"shear_rubin_obj_shear = tables_io.read(f'{tract}/rubin_mst_md_match_object_shear.pq')\n",
"\n",
"joined_rubin_shear_obj = hpmcm.output_tables.buildJoinedObjectTable(\n",
" shear_rubin_obj_stats,\n",
" shear_rubin_obj_assoc,\n",
" rubin_shear_match_files,\n",
" rubin_shear_match_ids,\n",
" object_shear=f'{tract}/rubin_mst_md_match_object_shear.pq',\n",
")\n",
"stats_cols = [f\"{band}_gaussFlux\" for band in 'griz']\n",
"stats_cols += [f\"{band}_gaussFluxErr\" for band in 'griz']\n",
"drop_cols = []\n",
"for col_prefix in [\n",
" 'tract', 'ra', 'dec', 'gauss_snr', 'patch_x', 'patch_y', \n",
" 'cell_idx_x', 'cell_idx_y', 'shear', 'meta_step',\n",
"]:\n",
" drop_cols += [f\"{col_prefix}_{id}\" for id in rubin_shear_match_ids]\n",
"\n",
"rd_joined_rubin_shear_obj = hpmcm.reduceJoinedTable(\n",
" joined_rubin_shear_obj,\n",
" catalog_ids=rubin_shear_match_ids,\n",
" drop_cols=drop_cols,\n",
" stats_cols=stats_cols,\n",
")\n",
"tables_io.write(rd_joined_rubin_shear_obj, f'{tract}/rubin_mst_md_match_object_joined.pq')"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "0ccc0274-d134-434b-9865-31e2a9416003",
"metadata": {},
"outputs": [],
"source": [
"shear_desc_obj_stats = tables_io.read(f'{tract}/desc_mst_md_match_object_stats.pq')\n",
"shear_desc_obj_assoc = tables_io.read(f'{tract}/desc_mst_md_match_object_assoc.pq')\n",
"shear_desc_obj_shear = tables_io.read(f'{tract}/desc_mst_md_match_object_shear.pq')\n",
"\n",
"joined_desc_shear_obj = hpmcm.output_tables.buildJoinedObjectTable(\n",
" shear_desc_obj_stats,\n",
" shear_desc_obj_assoc,\n",
" desc_shear_match_files,\n",
" desc_shear_match_ids,\n",
" object_shear=f'{tract}/desc_mst_md_match_object_shear.pq',\n",
")\n",
"stats_cols = [f\"flux_{band}\" for band in 'riz']\n",
"stats_cols += [f\"flux_err_{band}\" for band in 'riz']\n",
"drop_cols = []\n",
"for col_prefix in [\n",
" 'ra', 'dec', 'gauss_snr', 'patch_x', 'patch_y', \n",
" 'cell_idx_x', 'cell_idx_y', 'shear', 'meta_step',\n",
"]:\n",
" drop_cols += [f\"{col_prefix}_{id}\" for id in rubin_shear_match_ids]\n",
"\n",
"rd_joined_desc_shear_obj = hpmcm.reduceJoinedTable(\n",
" joined_desc_shear_obj,\n",
" catalog_ids=desc_shear_match_ids,\n",
" drop_cols=drop_cols,\n",
" stats_cols=stats_cols,\n",
")\n",
"tables_io.write(rd_joined_desc_shear_obj, f'{tract}/desc_mst_md_match_object_joined.pq')"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "408c4f7b-3b36-46ac-b744-b5340c2e60c1",
"metadata": {},
"outputs": [],
"source": [
"# Draw a cluster from the rubin shear matching\n",
"# countsMap() re-computes the pixel map from source positions for visualisation\n",
"if rematch_rubin_shear:\n",
" cell_data = rubin_shear_match.cell_dict[20200]\n",
" cluster_id = list(cell_data.cluster_dict.keys())[1]\n",
" cluster_data = cell_data.cluster_dict[cluster_id]\n",
" img = cell_data.countsMap()\n",
" _ = hpmcm.viz_utils.showCluster(\n",
" img,\n",
" cluster_data,\n",
" cell_data,\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "3328791e-876d-4a4d-8982-4f04afaaab29",
"metadata": {},
"outputs": [],
"source": [
"# Draw a cluster from the desc shear matching\n",
"if rematch_desc_shear:\n",
" cell_data = desc_shear_match.cell_dict[20200]\n",
" cluster_id = list(cell_data.cluster_dict.keys())[1]\n",
" cluster_data = cell_data.cluster_dict[cluster_id]\n",
" img = cell_data.countsMap()\n",
" _ = hpmcm.viz_utils.showCluster(\n",
" img,\n",
" cluster_data,\n",
" cell_data,\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "5bd02341-0983-4dc5-85ab-368230a0ac41",
"metadata": {},
"outputs": [],
"source": [
"if rematch_objects:\n",
" pixel_scale = 0.2\n",
" n_pix_region = 30000\n",
" cell_size = 150\n",
" region_size_deg = n_pix_region*pixel_scale/3600\n",
" mst_wcs_match = hpmcm.MSTWcsMatch.create(\n",
" ref_dir=ref_dir,\n",
" region_size=(region_size_deg, region_size_deg),\n",
" pixel_size=pixel_scale/3600,\n",
" mst_cut=mst_cut_wcs,\n",
" cell_size=cell_size,\n",
" cell_buffer=5,\n",
" cell_max_object=10000,\n",
" pixel_match_scale=4,\n",
" pixel_r2_cut=16, \n",
" )\n",
" \n",
" mst_wcs_match.reduceData(all_files, all_ids)\n",
" if test_only:\n",
" mst_wcs_match.analysisLoop(range(100, 101), range(100, 101))\n",
" else:\n",
" mst_wcs_match.analysisLoop()\n",
" stats = mst_wcs_match.extractStats()\n",
" tables_io.write(stats, f'{tract}/all_mst_match_') \n",
" mst_wcs_match.save(f'{tract}/all_mst_match')"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "c900bb1b-8605-4cf2-92bf-5d1991a925d8",
"metadata": {},
"outputs": [],
"source": [
"# Draw a cluster from the meta matching\n",
"if rematch_objects:\n",
" cell_data = mst_wcs_match.cell_dict[20200]\n",
" cluster_id = list(cell_data.cluster_dict.keys())[0]\n",
" cluster_data = cell_data.cluster_dict[cluster_id]\n",
" img = cell_data.countsMap()\n",
" _ = hpmcm.viz_utils.showCluster(\n",
" img,\n",
" cluster_data,\n",
" cell_data,\n",
" )"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "ee7ec80d-7681-44b5-92bf-5f5a4a82e6fa",
"metadata": {},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": null,
"id": "c984950d-375a-4de5-a1c7-0affd74cc6f7",
"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.12"
}
},
"nbformat": 4,
"nbformat_minor": 5
}
Loading
Loading