Rebuilt how profiles of objects are detected and fit. by cmccully · Pull Request #101 · LCOGT/banzai-floyds

cmccully · 2026-09-01T23:41:32Z

This PR closes #88

This is a major rewrite of fitting the center shape of the profile of objects. The procedure is now, we remove some of the background with a rolling median filter, and the run a Gaussian match filter in an overlapping part of both the blue and red orders. This gives us the location of point sources in the slit. We choose one of the point sources to extract. We normally pick the brightest object in the slit. If there are multiple points sources with a few tens of percent of each other, we choose the one closer to the center of the slit because we assume that's where we tried to do the acquisition. This will not solve every extraction problem for every user, but does a good job for an automated pipeline. We will still need to provide the user UI tools to re-extract. After we have a chosen point source, we step both directions along the slit, tracing the object using similar match filter centroiding like the wavelength solution. We then fit the trace with a Legendre polynomial. With the centroids in hand, we estimate the FWHM in an empirical way: we subctract a linear backgroun by taking a median of regions +-3-5 sigma on each side of the object, find the peak flux (and half that) and then interpolate to the width where the profile reaches the half-max. We iterate this procedure so that the background region is taken with the current fwhm in mind. For objects that are bright enough and are isolated enough, we fit a voigt profile to characterize the shape of the PSF. Otherwise we adopt the most recent.

I had to make some edits to extraction and background subtraction to get good enough results to compare different profile fitting procedures so they are included here but are not the main point of this PR so don't spend much time on it. I'm going to clean that code up in an upcoming pr.

- Fix the trace point count guard to count points per order instead of cumulatively across orders. Previously a faint order with only 2-3 trace points could be fit with a degree-7 polynomial (rank deficient, wildly oscillating) because the bright order's points satisfied the guard. - Bound the per-bin Gaussian fits (center on the grid, 0.5 < sigma < order_height/2, normalization >= 0), guard curve_fit against RuntimeError/ValueError, and reject fits with non-finite covariance or center errors > 2 pixels. - Add iterative MAD-based sigma clipping of trace points before the polynomial fits. Clipping iterations use unweighted fits so a bad point with a spuriously small error can't drag the polynomial through itself and escape the clip. - Fix the weight convention in Legendre.fit calls (here and in the background fit): numpy applies weights to the unsquared residual, so inverse-variance weighting is w=1/sigma, not 1/sigma**2. - Fix UnboundLocalError when an order has no usable bins (sigma was only assigned inside the loop); the width fallback now always uses the initial FWHM guess. - Only interpolate each wavelength bin onto the part of the common grid it covers; masked pixels could shrink a bin's y range below the chunk's and raise ValueError. Sort and de-duplicate y values per bin first. - Use an integrated-normalization initial guess for the Gaussian fit instead of the peak flux. - Fix the trace table schema (sigma/sigma_error instead of vestigial width/width_error columns) and add a 'used' flag for clipped points. - Log warnings on profile fallbacks and record QC headers: L1PRNP<order> (trace points used) and L1PRFB<order> (fallback flag).

The degree-7 Legendre fit of the trace center could swing wildly off the slit, which silently ruined the background windows and the extraction. Each 25 column chunk picked its center from a global argmax over the whole slit with no continuity constraint, and the polynomial was a bare weighted fit with no outlier rejection, no degree adaptation, and no check that the result stayed anywhere near the data. - Detect the object once per order at high s/n. The order is stacked in blocks of a few hundred columns and a coarse low order trace is fit through the brightest peak of each block. Each chunk then searches only a window around that coarse trace, updated with a running offset from the last few accepted chunks and swept outward from the middle of the order. - Project the smooth slit illumination out of the matched filter template so the detection measures the object above the sky rather than the sky. Stacked over hundreds of columns the sky is millions of counts, and a small error in the shape of the sky model was otherwise a very significant detection of nothing. - Fit the trace polynomials in two robust passes: a low order trend rejects points several pixels off a smooth curve, then the survivors are fit at a degree the data supports, reduced until the polynomial stays on the slit and within a profile width of the measured range over the whole domain. - Fall back progressively instead of silently jumping to the order center: reduced degree, constant at the measured position, the position and width from the other order, then the order center. Log every step and record L1PRFB (fallback level), L1PRNP, L1PRDG, and L1PRSN per order. - Count trace points per order rather than cumulatively, so a faint order can't be fit at degree 7 through the bright order's points. - Clamp profile widths when building the profile image, and warn when extraction drops bins that missed the profile. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>

Copilot

🟡 Changes recommended

The no-object early-return paths leave required frame products (BACKGROUND/EXTRACTED) unset while downstream stages still assume they exist, which can crash the pipeline and prevent writing the intended 2D-only output.

Once you've addressed the issues Copilot identified, you can request another Copilot review.

Pull request overview

Major rewrite of the object profile detection and fitting pipeline to make source detection, tracing, and PSF characterization more robust (per Issue #88), with follow-on updates to background modeling and extraction behavior.

Changes:

Replaces profile fitting with a detection → source selection → trace → empirical FWHM estimation → Voigt-wing shape fitting workflow (plus persistence of shape metadata).
Reworks background fitting to jointly fit object + sky per wavelength bin and adapt polynomial degree based on object width.
Updates extraction/output behavior and expands test coverage substantially for the new profile/background logic.

File summaries

File	Description
uv.lock	Adds locked dependencies for scikit-learn (and transitive deps).
pyproject.toml	Adds scikit-learn runtime dependency.
characterization_testing/README.md	Documents new characterization scripts to run.
characterization_testing/fringe_correction_results.py	Makes stats panel fontsize configurable.
characterization_testing/extraction_results.py	Refactors reduction workflow to use `reduce_to_stage`.
CHANGES.md	Updates release notes to mention profile fitting hardening.
banzai_floyds/utils/wavelength_utils.py	Removes `lsf_sigma` convenience property usage.
banzai_floyds/utils/telluric_utils.py	Removes unused telluric scaling helpers (now relying on estimate path).
banzai_floyds/utils/profile_utils.py	Switches profile model to Voigt + seeing scaling; normalizes per-column; updates header round-trip.
banzai_floyds/utils/fitting_utils.py	Adds Voigt model, robust nonlinear fitting, clamped Legendre, and robust linear-fit refactor.
banzai_floyds/tests/utils.py	Extends fake-frame generator to support multi-source/wavelength-range/profile width scenarios.
banzai_floyds/tests/test_wavelengths.py	Updates wavelength tests for removed `lsf_sigma`.
banzai_floyds/tests/test_telluric.py	Makes telluric test robust to stochastic noise amplification in bands.
banzai_floyds/tests/test_profile.py	Adds extensive tests for detection, tracing, width/shape estimation, and header round-trips.
banzai_floyds/tests/test_orders.py	Asserts order solution has no proprietary period.
banzai_floyds/tests/test_fringing.py	Updates expectations around super-fringe hole filling and overlap masking.
banzai_floyds/tests/test_frames.py	Adds tests for suppressing 1D output when no object is detected.
banzai_floyds/tests/test_extract.py	Updates tests to new profile tuple structure (centers + FWHM + shape).
banzai_floyds/tests/test_background.py	Updates tests for new joint-fit background model and adaptive degree behavior.
banzai_floyds/profile.py	Replaces profile stage implementation with detection/trace/FWHM/shape fitting pipeline.
banzai_floyds/orders.py	Ensures order-solution products are immediately public.
banzai_floyds/fringe.py	Simplifies fringe mask semantics (no longer flags interpolated-as-usable pixels).
banzai_floyds/frames.py	Updates profile serialization, extracted/2D product handling, and binned-data profile columns.
banzai_floyds/extract.py	Adjusts extraction window default and adds no-object handling hooks.
banzai_floyds/dbs.py	Adds ProfileShape DB model + accessors for recent shape fallback.
banzai_floyds/background.py	Replaces background-window approach with joint object+sky fit per wavelength bin.

Review details

Files reviewed: 25/26 changed files
Comments generated: 6
Review effort level: Lite

💡 Add a code-review agent skill or configure MCP servers for context-aware, tailored reviews. Learn more in the docs.

+        # Without a profile the binned data has no profile_sigma or profile_gamma_ratio column, and reading
+        # one raises. banzai catches that by dropping the frame from the reduction entirely, so a
+        # frame with no object in the slit used to produce no product at all rather than an
+        # unextracted one. There is no object to fit a sky around here, so hand the frame back.
+        if image.profile_fits is None:
+            logger.warning('No object was detected, so there is no profile to fit the sky around.',
+                           image=image)
+            return image


    def do_stage(self, image):
+        # Nothing was found in the slit. Hand the frame back untouched rather than raising on the
+        # missing profile columns, which banzai would turn into the frame being dropped from the
+        # reduction with no product written at all.
+        if image.profile_fits is None:
+            logger.warning('No object was detected, so there is nothing to extract.', image=image)
+            return image


+        # An order the profile stage could not trace, and could not borrow a trace for, carries a
+        # placeholder profile at the center of the order. Leaving its window empty is what keeps that
+        # placeholder from being quietly extracted as if it were a measurement: the order drops out
+        # of the extracted spectrum, which is visible, instead of appearing at the wrong position,
+        # which is not.
+        if not image.meta.get(f'L1PRTR{order_id}', True):
+            logger.warning(f'No trace was measured for order {order_id}, so it is not extracted.',
+                           image=image)
+            continue


Copilot

🟡 Changes recommended

Missing database and characterization files plus several profile and no-object control-flow defects prevent reliable deployment.

Once you've addressed the issues Copilot identified, you can request another Copilot review.

Review details

Suppressed comments (6)

Previously missed (5) — in code that hasn't changed since the last review.

banzai_floyds/profile.py:390

This threshold counts every measured center, including points rejected by max_center_error. If at least seven centers were measured but none (or too few) are fittable, the following min/max operates on an empty array or fits an underdetermined degree-5 polynomial. Apply the minimum to the retained points and ensure there are at least polynomial_order + 1.
banzai_floyds/profile.py:394
Restricting the Legendre domain to measured wavelengths does not restrict evaluation: NumPy extrapolates the degree-5 polynomial outside that interval. The profile setter later evaluates this center over the entire order, so faint traces that stop early can swing out of the slit and lose extraction bins. Wrap the fit with ClampedLegendre and reconstruct that behavior in load_profile_fits so persisted profiles remain clamped.
characterization_testing/extraction_results.py:44
reduction_utils.py is not present under characterization_testing, and this is the only reference to it, so running this script now fails immediately with ModuleNotFoundError. Include the shared helper in this PR or retain the local reduction implementation.
banzai_floyds/frames.py:148
The persisted profile schema changed here, but docs/banzai_floyds/data-products.rst:72-79 still documents in_background plus per-point sigma/sigma_error, which this implementation no longer emits, and omits PROFFWHM, PROFGAM, and profile_gamma_ratio. Update the data-product documentation so consumers can interpret new files correctly.
characterization_testing/README.md:13
Neither referenced script exists in characterization_testing, so these newly documented commands cannot run. Add both scripts to the PR or remove the commands until they are available.

banzai_floyds/background.py:178

When no source is found, this return leaves the frame without a BACKGROUND extension. The next configured stage is CosmicRayDetector, which unconditionally evaluates image.background (settings.py:21-25, cosmics.py:87-92), so the intended 2D-only path fails before output generation. Provide a sky-only fallback background (or coordinate a no-profile bypass in the downstream stages) so no-object frames can complete reduction.

        if image.profile_fits is None:
            logger.warning('No object was detected, so there is no profile to fit the sky around.',
                           image=image)
            return image

Files reviewed: 25/26 changed files
Comments generated: 2
Review effort level: Balanced

+class ProfileShape(Base):
+    __tablename__ = 'profileshape'
+    id = Column(Integer, primary_key=True, autoincrement=True)
+    instrument_id = Column(Integer, ForeignKey("instruments.id"), index=True)
+    filename = Column(String(100), unique=True)


+        if any(trace is None for trace in profile_center):
+            logger.warning('The object was not traced in every order, so no profile was fit.', image=image)
+            image.meta['L1OBJDET'] = (False, 'Was an object detected in the slit?')
+            return image


sfoale · 2026-09-02T10:16:27Z

+    return peaks
+
+
+def choose_source_to_extract(point_sources: list[dict], snr_ratio: float = 0.6) -> dict | None:


the default snr ratio here is a magic number - could it be a class constant?

sfoale · 2026-09-02T10:19:13Z

+    return model


 def interp_with_errors(x, y, yerr, x_new):


I think this is dead code now with no callers

sfoale · 2026-09-02T10:29:14Z

+        for chunks in chunks_from_detection(order_data, point_source['detection_wavelength'], chunk_size):
+            center_guess = point_source['center']
+            for chunk_low, chunk_high in chunks:
+                stacked_y, stacked_flux, stacked_flux_error = stack_slit_profile(


as far as I can tell these stack_slit_profile calls are repeated with identical arguments later on in measure_chunk_fwhms so consider storing the results from here to use there

sfoale · 2026-09-02T10:30:16Z

+            for chunk_low, chunk_high in chunks:
+                wavelength = 0.5 * (chunk_low + chunk_high)
+                center = float(trace(wavelength))
+                stacked_y, stacked_flux, stacked_flux_error = stack_slit_profile(


repeated stack_slit_profile call - see comment at line 368

sfoale

Copilot did a better job than me and those comments should be addressed. I made a couple of style/efficiency comments.

jchate6

A few comments/questions and a typo or two.

jchate6 · 2026-09-14T20:41:09Z



+def set_up_profile(frame):
+    """Give a fake frame the profile the background stage needs, from the values it was built with."""


The phrasing of this is a little awkward. What exactly does "from the values it was built with" mean?
Even just clarifying the antecedent for "it" would help.

jchate6 · 2026-09-14T20:49:30Z

+    # prepare_fringe_data cuts whole columns at the cutoff while this mask is per pixel, and the
+    # order is tilted, so a few pixels clear 6000 A in columns that were never padded. There is
+    # nothing to interpolate from there.
+    overlap = np.logical_and(fake_frame.wavelengths.data[order_region][2:-2] >= 6000.0,


It feels like 6000 should be a variable set somewhere, rather than just a number included repeatedly.

jchate6 · 2026-09-14T21:25:40Z

+
+def fit_background(data, background_order=3, minimum_fit_pixels=MINIMUM_FIT_PIXELS):
+    """
+    Fit the sky in each wavelength bin, with the object in the model.


I'm not sure what "with the object in the model" means

jchate6 · 2026-09-14T21:30:22Z

-        background = fit_background(image.binned_data)
+        # Without a profile the binned data has no profile_sigma or profile_gamma_ratio column, and reading
+        # one raises. banzai catches that by dropping the frame from the reduction entirely, so a
+        # frame with no object in the slit used to produce no product at all rather than an


"used to produce"? Is this no longer the case?
If we miss the target, does no reduction happen?

jchate6 · 2026-09-14T21:46:16Z

+    Our detection algorithm is to combine about a hundred pixels of a wavelength region,
+    do a median filter along the y-axis to remove any smooth background component, and
+    then run a match filter to a Gaussian with provided fwhm to detect objects. This was found to me more
+    stable than trying to simultaneously fit a background with a polynomial do a match filter. The median


typo? do -> to?

jchate6 · 2026-09-14T21:50:03Z

+def choose_source_to_extract(sources_by_order: dict[int, list[dict]],
+                             snr_ratio: float = 0.8) -> dict[int, dict]:
+    """Pick which object to extract in each order. We choose the brightest object unless the top two objects
+    are within a few tens of percent of each other, then we choose the closest to center of the slit.


"...closer of the two to the center..."

jchate6 · 2026-09-14T21:53:50Z

+    -----
+    Acquisition puts the requested coordinates at the center of the slit,
+    so we choose that one if the sources are close to the same brightness (Set by the `snr_ratio` parameter).
+    Both orders are detected over the same wavelengths, so choosing in each of them independently lands


Is it possible for the distance from the center to be shifted between orders?
I'm concerned for the scenario where there are 2 traces of similar brightness straddling the center, where the top one is closer in one order, and the bottom is closer in the other.

jchate6 · 2026-09-14T22:09:11Z

+            clip_sigma=self.N_SIGMA_CLIP, max_chunk_shift=self.MAX_CHUNK_SHIFT
+        )
+        if any(trace is None for trace in profile_center):
+            logger.warning('The object was not traced in every order, so no profile was fit.', image=image)


We will not trace objects that are too faint, correct? How do we handle a very red object that is very obvious in one order but indistinguishable from background in the other?

Copilot

Copilot review overview

🟡 Changes recommended

Critical missing modules, migration and lockfile issues, plus unresolved profile, compatibility, and no-object processing defects, block approval.

Get a fresh assessment by requesting another Copilot review.

Review effort: Balanced
Findings: 12 · 2 · 2

Open (16)

Add migration for the profileshape table and instrument index · New
Add missing Gaia utilities module or remove profile imports · New
Match detections across orders to select one shared source · New
Support legacy PROFILEFITS sigma headers · New
Add the missing characterization reduction utility · New
Add missing characterization reduction and reporting helpers · New
Add missing characterization reduction and reporting helpers · New
Regenerate uv.lock to include astroquery · New
Returning when either trace is missing discards a valid trace in the other order, labels the whole…
The new ORM table has no corresponding Alembic revision (the only revision creates orderheights).…
Extractor.do_stage returns early when image.profile_fits is None, but later stages (e.g.…
BackgroundFitter returns early when image.profile_fits is None, but CosmicRayDetector runs…
Clip fitted profile widths before extraction and background modeling · New
Clamp gamma to prevent invalid Voigt widths · New
Add or correct the documented characterization entry point · New
L1PRTR{order} is checked to decide whether to skip extracting an order, but this header keyword…

Resolved since last review (3)

Docstring typo: “emperical” → “empirical”.
Docstring typo: “requried” → “required”.
Docstring grammar: “This was found to me more stable…” reads as a typo/grammatical error and is…

+            wavelengths = self.binned_data['wavelength'][in_order]
+            center = centers[order - 1](wavelengths)
+            self.binned_data['y_profile'][in_order] = self.binned_data['y_order'][in_order] - center
+            self.binned_data['profile_sigma'][in_order] = profile_sigmas(wavelengths, fwhms[order - 1])


+            center = centers[order - 1](wavelengths)
+            self.binned_data['y_profile'][in_order] = self.binned_data['y_order'][in_order] - center
+            self.binned_data['profile_sigma'][in_order] = profile_sigmas(wavelengths, fwhms[order - 1])
+            self.binned_data['profile_gamma_ratio'][in_order] = gamma_ratios[order - 1](wavelengths)


cmccully · 2026-09-23T23:02:56Z

profile_cross_sections.pdf
trace_overlay_results.pdf

cmccully and others added 14 commits August 4, 2026 20:59

WIP. Working to disentangle code spaghetti.

a5f7f9c

Fix to comment.

cfd6637

WIP: Refactoring profile fitting to be more readable.

44ac180

WIP. I think our point source detection works now.

cdeacdc

WIP: Choose source to extract done.

1a55805

WIP: Tracing should now work.

b739cc1

FWHM measurement code should be better now.

124f611

WIP: Hopefully profile.py is about done now.

3d7f219

WIP: One step closer.

d4941b0

I think we have a fully working profile fitting scheme?

abe04a6

Fixed some comments.

eda41e3

Merge with upstream.

151cefd

cmccully requested review from drhaz, jchate6 and sfoale and a lite review from Copilot September 1, 2026 23:41

Copilot started reviewing on behalf of cmccully September 1, 2026 23:42 View session

Update to orders comment.

2d3d89f

Copilot AI reviewed Sep 1, 2026

View reviewed changes

cmccully requested a balanced review from Copilot September 2, 2026 00:49

Copilot started reviewing on behalf of cmccully September 2, 2026 00:50 View session

Copilot AI reviewed Sep 2, 2026

View reviewed changes

sfoale reviewed Sep 2, 2026

View reviewed changes

sfoale approved these changes Sep 2, 2026

View reviewed changes

Added my characterization scripts.

e4ee678

jchate6 approved these changes Sep 14, 2026

View reviewed changes

Stablized profile fitting. Addressed issues from pr comments.

c13e43c

cmccully requested a balanced review from Copilot September 23, 2026 20:48

Copilot started reviewing on behalf of cmccully September 23, 2026 20:48 View session

Copilot AI reviewed Sep 23, 2026

View reviewed changes

Fixes based on copilot notes.

058d170

cmccully self-assigned this Sep 23, 2026

cmccully merged commit f818555 into main Sep 24, 2026
7 checks passed

cmccully deleted the profile-rebuild branch September 24, 2026 19:37

		return peaks


		def choose_source_to_extract(point_sources: list[dict], snr_ratio: float = 0.6) -> dict \| None:



		def set_up_profile(frame):
		"""Give a fake frame the profile the background stage needs, from the values it was built with."""

Conversation

cmccully commented Sep 1, 2026

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

🟡 Changes recommended

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

🟡 Changes recommended

Suppressed comments (6)

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

sfoale left a comment

Choose a reason for hiding this comment

Uh oh!

jchate6 left a comment

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Copilot AI left a comment

Choose a reason for hiding this comment

Copilot review overview

🟡 Changes recommended

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

Uh oh!

cmccully commented Sep 23, 2026

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

4 participants