{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Identifying sources\n",
    "\n",
    "`phoptic` uses `photutils` to find sources in images. This notebook will demonstrate how to define source finders for use with `phoptic`, as well as explain `phoptic`'s default behaviour when no source finder is specified.\n",
    "\n",
    "## Test Image\n",
    "\n",
    "First thing's first, let's open an image that contains some sources. For this example, I'll use one of the images from the [Basic Usage tutorial](basic_usage.ipynb):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:12.294136Z",
     "iopub.status.busy": "2026-01-23T11:20:12.294058Z",
     "iopub.status.idle": "2026-01-23T11:20:15.246842Z",
     "shell.execute_reply": "2026-01-23T11:20:15.246576Z"
    }
   },
   "outputs": [],
   "source": [
    "from pathlib import Path\n",
    "\n",
    "import phoptic\n",
    "\n",
    "out_dir = Path('out')\n",
    "\n",
    "phoptic.generate_observations(\n",
    "    out_directory=out_dir / 'data',\n",
    "    n_images=3,\n",
    "    circular_aperture=False,\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.248056Z",
     "iopub.status.busy": "2026-01-23T11:20:15.247801Z",
     "iopub.status.idle": "2026-01-23T11:20:15.271707Z",
     "shell.execute_reply": "2026-01-23T11:20:15.271357Z"
    }
   },
   "outputs": [],
   "source": [
    "from astropy.io import fits\n",
    "import numpy as np\n",
    "import os\n",
    "\n",
    "files = os.listdir(out_dir / 'data')\n",
    "file = files[0]\n",
    "\n",
    "with fits.open(out_dir / 'data' / file) as hdul:\n",
    "    print(repr(hdul[0].header))\n",
    "    image = np.array(hdul[0].data)\n",
    "    binning_factor = int(hdul[0].header['BINNING'][0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.272792Z",
     "iopub.status.busy": "2026-01-23T11:20:15.272577Z",
     "iopub.status.idle": "2026-01-23T11:20:15.411831Z",
     "shell.execute_reply": "2026-01-23T11:20:15.411494Z"
    }
   },
   "outputs": [],
   "source": [
    "from astropy.visualization import simple_norm\n",
    "from matplotlib import pyplot as plt\n",
    "\n",
    "fig, ax = plt.subplots(tight_layout=True)\n",
    "\n",
    "im = ax.imshow(\n",
    "    image,\n",
    "    norm=simple_norm(image, stretch=\"sqrt\"),\n",
    "    origin=\"lower\",\n",
    "    cmap=\"Greys\",\n",
    "    )\n",
    "\n",
    "ax.set_xlabel(\"X\")\n",
    "ax.set_ylabel(\"Y\")\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Default Finder\n",
    "\n",
    "`phoptic` implements two default source finders: `Finder` (default) and `CrowdedFinder` (better for crowded fields, but more expensive). `Finder` does not implement any source deblending, while `CrowdedFinder` does. Both of these source finders are wrappers for `photutils.segmentation.SourceFinder` with some added convenience tailored to OPTICAM.\n",
    "\n",
    "Let's use `phoptic.Finder` to identify the sources in the above image:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.413020Z",
     "iopub.status.busy": "2026-01-23T11:20:15.412879Z",
     "iopub.status.idle": "2026-01-23T11:20:15.464127Z",
     "shell.execute_reply": "2026-01-23T11:20:15.463845Z"
    }
   },
   "outputs": [],
   "source": [
    "from phoptic import DefaultFinder\n",
    "\n",
    "from photutils.segmentation import detect_threshold\n",
    "\n",
    "# default value for n_pixels is 128 // binning_factor**2\n",
    "# default value for border_width is 1/16th of the image width\n",
    "default_finder = DefaultFinder(n_pixels=128 // 4**2, border_width=image.shape[0] // 16)\n",
    "\n",
    "default_tbl = default_finder(\n",
    "    image,\n",
    "    threshold=detect_threshold(image, n_sigma=5),  # 5 sigma detection threshold\n",
    "    )\n",
    "print(type(default_tbl))\n",
    "default_tbl.pprint_all()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "When calling a `DefaultFinder()` instance, an `astropy.table.QTable` instance is returned that is sorted in descending order of source flux. This is how catalogs are defined in `phoptic`. We can use the table to visualise which sources have been identified and how their fluxes are ranked:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.465184Z",
     "iopub.status.busy": "2026-01-23T11:20:15.465047Z",
     "iopub.status.idle": "2026-01-23T11:20:15.599475Z",
     "shell.execute_reply": "2026-01-23T11:20:15.599188Z"
    }
   },
   "outputs": [],
   "source": [
    "from matplotlib.patches import Circle\n",
    "\n",
    "fig, ax = plt.subplots(\n",
    "    tight_layout=True,\n",
    "    )\n",
    "\n",
    "im = ax.imshow(\n",
    "    image,\n",
    "    norm=simple_norm(image, stretch=\"log\"),\n",
    "    origin=\"lower\",\n",
    "    cmap=\"Greys\",\n",
    "    )\n",
    "\n",
    "for i, row in enumerate(default_tbl):\n",
    "    \n",
    "    xc = row['x_centroid']\n",
    "    yc = row['y_centroid']\n",
    "    \n",
    "    ax.add_patch(\n",
    "        Circle(\n",
    "            xy=(xc, yc),\n",
    "            radius=5 * row['semimajor_axis'].value,\n",
    "            facecolor='none',\n",
    "            edgecolor='cyan',\n",
    "            lw=1,\n",
    "        )\n",
    "    )\n",
    "    \n",
    "    ax.text(\n",
    "        xc * 1.05,\n",
    "        yc * 1.05,\n",
    "        f'{i + 1}',\n",
    "        color='cyan',\n",
    "    )\n",
    "\n",
    "ax.set_xlabel(\"X\")\n",
    "ax.set_ylabel(\"Y\")\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Custom Source Finders\n",
    "\n",
    "### Defining the Source Finder\n",
    "\n",
    "Let's now define a custom source finder. A custom source finder must be a `callable` that takes two parameters: `image` and `threshold`, and returns a `QTable` instance. `image` should be an `NDArray` containing the image data. `threshold` defines the threshold for source detection in ADU. The simplest way to define a custom source finder is by using the `photutils.segmentation.SourceFinder` class, which combines source detection and deblending:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.600468Z",
     "iopub.status.busy": "2026-01-23T11:20:15.600283Z",
     "iopub.status.idle": "2026-01-23T11:20:15.602425Z",
     "shell.execute_reply": "2026-01-23T11:20:15.602223Z"
    }
   },
   "outputs": [],
   "source": [
    "%%writefile out/custom_routines.py\n",
    "from photutils.segmentation import SourceCatalog, SourceFinder\n",
    "\n",
    "def custom_finder(image, threshold):\n",
    "    finder = SourceFinder(\n",
    "        n_pixels=128 // 4**2,  # same as DefaultFinder\n",
    "        deblend=True,  # same as DefaultFinder\n",
    "        n_levels=256,  # higher than DefaultFinder (more deblending)\n",
    "        contrast=0,  # lower than DefaultFinder (more deblending)\n",
    "        progress_bar=False,  # disable progress bar (same as DefaultFinder)\n",
    "        )\n",
    "    segment_img = finder(image, threshold)\n",
    "    \n",
    "    tbl = SourceCatalog(image, segment_img).to_table()\n",
    "    tbl.sort('segment_flux', reverse=True)  # sort in descending order\n",
    "    \n",
    "    return tbl"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Due to compatibility issues, routines defined inside of IPython notebooks cannot be used by `multiprocessing` (used by `phoptic`). To overcome this issue, we have used the `%%writefile` magic command to write our `custom_finder()` function to a temporary module called `custom_routines.py` in our `out` directory. This is required in Python $\\geq$ 3.14.\n",
    "\n",
    "To understand how `SourceFinder` works, and what parameters can be weaked, I refer to the excellent `photutils` documentation: https://photutils.readthedocs.io/en/stable/user_guide/segmentation.html.\n",
    "\n",
    "Under-the-hood, `phoptic.DefaultFinder` also uses `photutils.segmentation.SourceFinder`, but here we have defined our custom source finder to use different parameter values for better source deblending (though deblending isn't really needed for our example images). In addition to the different parameter values, our custom finder differs from `DefaultFinder` in that `DefaultFinder` automatically omits sources that are close to the edge of an image (where \"close\" is defined as 1/16th of the image width, the same size as the [default background pixels](backgrounds.ipynb)). For this example, this will not make a difference but it can be important in practise. \n",
    "\n",
    "`SourceFinder` returns a `SegmentationImage`, but `phoptic` assumes an `astropy.table.QTable` is returned. Here we have converted the returned `SegmentationImage` to a `QTable` by first converting it to a `photutils.segmentation.SourceCatalog`, and then calling the `to_table()` method of `SourceCatalog`.\n",
    "\n",
    "Let's now initialise this custom source finder and use it to identify sources in the above image."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.603362Z",
     "iopub.status.busy": "2026-01-23T11:20:15.603286Z",
     "iopub.status.idle": "2026-01-23T11:20:15.651415Z",
     "shell.execute_reply": "2026-01-23T11:20:15.651197Z"
    }
   },
   "outputs": [],
   "source": [
    "from out.custom_routines import custom_finder\n",
    "\n",
    "\n",
    "custom_tbl = custom_finder(\n",
    "    image,\n",
    "    threshold=detect_threshold(image, n_sigma=5),  # 5 sigma detection threshold\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.652562Z",
     "iopub.status.busy": "2026-01-23T11:20:15.652482Z",
     "iopub.status.idle": "2026-01-23T11:20:15.742438Z",
     "shell.execute_reply": "2026-01-23T11:20:15.742178Z"
    }
   },
   "outputs": [],
   "source": [
    "fig, ax = plt.subplots(\n",
    "    tight_layout=True,\n",
    "    )\n",
    "\n",
    "im = ax.imshow(\n",
    "    image,\n",
    "    norm=simple_norm(image, stretch=\"log\"),\n",
    "    origin=\"lower\",\n",
    "    cmap=\"Greys\",\n",
    "    )\n",
    "\n",
    "for i, row in enumerate(custom_tbl):\n",
    "    \n",
    "    xc = row['x_centroid']\n",
    "    yc = row['y_centroid']\n",
    "    \n",
    "    ax.add_patch(\n",
    "        Circle(\n",
    "            xy=(xc, yc),\n",
    "            radius=5 * row['semimajor_axis'].value,\n",
    "            facecolor='none',\n",
    "            edgecolor='cyan',\n",
    "            lw=1,\n",
    "        )\n",
    "    )\n",
    "    \n",
    "    ax.text(\n",
    "        xc * 1.05,\n",
    "        yc * 1.05,\n",
    "        f'{i + 1}',\n",
    "        color='cyan',\n",
    "    )\n",
    "\n",
    "ax.set_xlabel(\"X\")\n",
    "ax.set_ylabel(\"Y\")\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we can see, we have once again recovered all six sources. Let's now see how to use this custom source finder with `phoptic.Reducer`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:15.744181Z",
     "iopub.status.busy": "2026-01-23T11:20:15.744093Z",
     "iopub.status.idle": "2026-01-23T11:20:20.245694Z",
     "shell.execute_reply": "2026-01-23T11:20:20.245467Z"
    }
   },
   "outputs": [],
   "source": [
    "custom_reducer = phoptic.Reducer(\n",
    "    out_directory=out_dir / 'reduced' / 'custom_finder',\n",
    "    data_directory=out_dir / 'data',\n",
    "    finder=custom_finder,  # use custom source finder\n",
    ")\n",
    "\n",
    "custom_reducer.create_catalogs()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Our custom source finder has been used to successfully identify all six sources. Let's compare this to using `DefaultFinder`:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "execution": {
     "iopub.execute_input": "2026-01-23T11:20:20.246849Z",
     "iopub.status.busy": "2026-01-23T11:20:20.246708Z",
     "iopub.status.idle": "2026-01-23T11:20:24.523589Z",
     "shell.execute_reply": "2026-01-23T11:20:24.523355Z"
    }
   },
   "outputs": [],
   "source": [
    "default_reducer = phoptic.Reducer(\n",
    "    out_directory=out_dir / 'reduced' / 'default_finder',\n",
    "    data_directory=out_dir / 'data',\n",
    ")\n",
    "\n",
    "default_reducer.create_catalogs()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Unsurprisingly, the default source finder has also identified all six sources. Admittedly, this rather simple example does a poor job of demonstrating why defining a custom source finder may be useful, but hopefully it is clear how custom source finders can be implemented. For more information on defining custom source finders, I again refer to the excellent `photutils` documentation: https://photutils.readthedocs.io/en/stable/user_guide/segmentation.html. Note that `phoptic` requires that source finder routines return an `astropy.table.QTable` instance when called, while `photutils` returns a `SourceCatalog` instance. To convert the `SourceCatalog` to a `QTable`, you can use the [`to_table()`](https://photutils.readthedocs.io/en/stable/api/photutils.segmentation.SourceCatalog.html#photutils.segmentation.SourceCatalog.to_table) method.\n",
    "\n",
    "That concludes the source finder tutorial for `phoptic`!"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "opticam7",
   "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.14.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}