DM-54397: Add task for computing a rough shapelet expansion of isolated stars. by TallJimbo · Pull Request #467 · lsst/meas_algorithms

TallJimbo · 2026-03-17T16:30:38Z

No description provided.

laurenam

Just some flyby comments while perusing…

laurenam · 2026-03-17T17:58:09Z

include/lsst/meas/algorithms/SpanSetMoments.h

+ * `lsst.afw.geom.SpanSet`.
+ *
+ * This class provides low-level pixel-processing code for
+ * `ComputeRoughtPsfShapeletsTask`.


Drop the “t” in Rought

laurenam · 2026-03-17T18:00:47Z

include/lsst/meas/algorithms/SpanSetMoments.h

+    /// The pixels actually used to compute the moments.
+    std::shared_ptr<afw::geom::SpanSet> spans;
+
+    /// Flag if there were too many bad pixels to compute the moments.


“Flag if” —> “Flag set if” to match those below?

laurenam · 2026-03-17T18:04:17Z

include/lsst/meas/algorithms/SpanSetMoments.h

+    /// Flag set if there was a bad pixel too close to the center.
+    bool bad_pixel_in_center = false;
+
+    /// Flag set if the second moments do not resolve to an ellipse.


Tense is changed from above, i.e. “do not” vs. “did not” & “was”.

laurenam · 2026-03-17T18:07:29Z

include/lsst/meas/algorithms/SpanSetMoments.h

+    /// Flag set if the second moments do not resolve to an ellipse.
+    bool singular_second_moments = false;
+
+    /// Test whether any failure falg is set.


“falg” —> “flag”

laurenam · 2026-03-17T18:08:27Z

include/lsst/meas/algorithms/SpanSetMoments.h

+     *
+     * @param[in] spans          Pixel region to use.
+     * @param[in] masked_image   Image to measure.
+     * @param[in] bad_bitmask    Mask of bad pixels to remove from `spans` to


Extra “to”?

laurenam · 2026-03-17T18:09:23Z

include/lsst/meas/algorithms/SpanSetMoments.h

+     *                                        before giving up.
+     * @param[in] bad_pixel_exlusion_radius   Radius around the estimated
+     *                                        center where bad pixels will
+     *                                        cause the algorithm to give up.


I don’t quite follow this description…

It's a radius that tries to act analogously to a cut on base_PixelFlags_[bad_bitmask]Center, but about the center the algorithm is computing itself. I've reworded slightly.

laurenam

A few more fly-by comments...

laurenam · 2026-03-17T18:23:34Z

python/lsst/meas/algorithms/computePsfChiSq.py

+    bad_mask_planes = ListField(
+        "Mask planes to identify pixels to drop from the calculation.",
+        dtype=str,
+        default=["SAT", "SUSPECT"],


Any reason not to include BAD and/or NO_DATA in the list? (I always forget what SUSPECT includes.)

laurenam · 2026-03-17T18:30:56Z

python/lsst/meas/algorithms/computePsfChiSq.py

+    )
+
+
+class ComputePsfChiSqTask(Task):


Oops, the commit that adds this task was not supposed to be in the PR, as I didn't think it was useful in the end. I deleted it, and then it accidentally came back in a rebase when I was transferring the code between USDF and my laptop. I'll go delete it again.

laurenam · 2026-03-17T18:41:26Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+        default=2.0,
+    )
+    max_footprint_area = Field(
+        "Footprints larger than this threshold are dropped before computing moments.",


“Footprints whose area has a pixel count larger…”?

laurenam · 2026-03-17T21:58:01Z

python/lsst/meas/algorithms/spanSetMoments.cc

+                cls.def_readonly("flux", &SpanSetMoments::flux);
+                cls.def_readonly("variance", &SpanSetMoments::variance);
+                // def_readonly doesn't ever seem to do the right thing w.r.t.
+                // memory management, even when explicit return value policies.


Sentence (after the comma) needs some rewording.

laurenam · 2026-03-17T22:00:45Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+    """Exception raised when ComputeRawPsfMomentsTask fails to find any usable
+    stars.
+    """
+


Unnecessary blank line?

No, I think there's always supposed to be a blank line after the docstring. It might be one of those things that I'm used to having ruff complain about or fix while flake8 doesn't care.

laurenam · 2026-03-26T19:13:07Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+                min_count,
+                name,
+                kind,
+                values[sorter[min_count - 1]],


In adding ComputeRoughPsfShapeletsTask as subtask in calibrateImage over on DM-54482, I got this in running the tests:

File "/sdf/data/rubin/user/laurenma/LSST/meas_algorithms/python/lsst/meas/algorithms/computeRoughPsfShapelets.py", line 709, in _threshold_with_bounds values[sorter[min_count - 1]], ~~~~~~^^^^^^^^^^^^^^^ IndexError: index 9 is out of bounds for axis 0 with size 4

I think values[sorter[min_count - 1]] needs to be something like values[sorter[min(min_count, len(values)) - 1]]. Do we want to add a hard minimum for this task to pass?

The min_count values were supposed to be that hard minimum, but as you noticed the logic was bad. Fixed in a new commit.

This means we'll raise an AlgorithmError in these cases, and that might mean CalibrateImageTask should catch that in a lower-level try block (rather than the big one in runQuantum) if it wants to proceed anyway when this new task can't do its job. Whether that's a good idea or not depends on whether this failure mode happens on good (but e.g. crowded or super-shallow-u) data.

This reverts commit 0346200.

TallJimbo

I've just deleted a file that wasn't supposed to be here. Sorry!

TallJimbo · 2026-03-27T15:19:23Z

python/lsst/meas/algorithms/computePsfChiSq.py

+    )
+
+
+class ComputePsfChiSqTask(Task):


Oops, the commit that adds this task was not supposed to be in the PR, as I didn't think it was useful in the end. I deleted it, and then it accidentally came back in a rebase when I was transferring the code between USDF and my laptop. I'll go delete it again.

TallJimbo · 2026-03-27T15:21:39Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+    """Exception raised when ComputeRawPsfMomentsTask fails to find any usable
+    stars.
+    """
+


No, I think there's always supposed to be a blank line after the docstring. It might be one of those things that I'm used to having ruff complain about or fix while flake8 doesn't care.

TallJimbo · 2026-03-27T15:30:25Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+                min_count,
+                name,
+                kind,
+                values[sorter[min_count - 1]],


The min_count values were supposed to be that hard minimum, but as you noticed the logic was bad. Fixed in a new commit.

This means we'll raise an AlgorithmError in these cases, and that might mean CalibrateImageTask should catch that in a lower-level try block (rather than the big one in runQuantum) if it wants to proceed anyway when this new task can't do its job. Whether that's a good idea or not depends on whether this failure mode happens on good (but e.g. crowded or super-shallow-u) data.

ctslater

Looks good, just a few questions on things I didn't understand.

ctslater · 2026-03-27T18:13:51Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+        peaks, _ = scipy.signal.find_peaks(scores)
+        if not peaks.size:
+            raise NoStarsForShapeletsError("Radius distribute has no mode.")
+        return sorted_radii[peaks.min()], kde


I don't understand how peaks can be used to index into sorted_radii, since the scores were created from the original radii?

Yup, that's a bug; good catch! scores should have been created from sorted_radii.

I think that says that this logic isn't super important, and that sort of makes sense, because this is about keeping galaxies and blends out of the star sample. I think at S/N=50 we usually mostly get stars without trying very hard, and the rest of the robust statistics can kind of take it from there.

include/lsst/meas/algorithms/SpanSetMoments.h

ctslater · 2026-03-27T22:19:59Z

src/SpanSetMoments.cc

+            ndarray::Array<float, 2, 2>(ndarray::allocate(n_coeff, n_pixels)).transpose();
+    design_matrix.deep() = 0.0;
+    std::size_t start = 0;
+    for (auto const& source_moments : moments) {


It took me a long time to decode that this is a loop over sources. The naming is tricky here because the relevant values per-source are moments, and so a vector of those can fairly be called "moments", but that doesn't leave any room in the naming scheme for saying what distinguishes these elements from each other. I'm tempted to suggest the input parameter be named sources? I'm not sure it's much better but it's the best alternative I can think of.

I've renamed the argument to sources and the loop variable to source; adding "moments" to either doesn't add much information for the reader, even if it's correct.

ctslater · 2026-03-27T23:13:54Z

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

+        default=False,
+    )
+    radius_mode_min = Field(
+        "Minimum radius at which to start searching for the first mode. "


Units are pixels?

👍, will update.

python/lsst/meas/algorithms/computeRoughPsfShapelets.py

Co-authored-by: Colin Slater <ctslater@uw.edu>

TallJimbo added 2 commits March 13, 2026 12:34

Add simple task for computing PSF model Chi^2 over an image.

0346200

Add task for computing a rough shapelet expansion of isolated stars.

aa9b3dc

TallJimbo force-pushed the tickets/DM-54397 branch from 54502cd to aa9b3dc Compare March 17, 2026 16:33

laurenam reviewed Mar 17, 2026

View reviewed changes

Doc fixes from review.

1410a75

laurenam reviewed Mar 26, 2026

View reviewed changes

Revert "Add simple task for computing PSF model Chi^2 over an image."

b8de2c3

This reverts commit 0346200.

TallJimbo commented Mar 27, 2026

View reviewed changes

TallJimbo added 2 commits March 27, 2026 11:33

Fail in computeRoughPsfShapelets when there aren't enough stars.

64bb1cc

More doc fixes from review.

598d141

TallJimbo force-pushed the tickets/DM-54397 branch from 7bbd9dc to 598d141 Compare March 27, 2026 15:33

ctslater approved these changes Mar 30, 2026

View reviewed changes

TallJimbo and others added 2 commits March 31, 2026 09:46

WILL SQUASH: fix typos spotted in review

ce2428b

Co-authored-by: Colin Slater <ctslater@uw.edu>

WILL SQUASH: bug and naming fixes from review

71bd3ca

TallJimbo force-pushed the tickets/DM-54397 branch from 697b79b to 71bd3ca Compare March 31, 2026 13:50

TallJimbo merged commit 6bc04b0 into main Mar 31, 2026
4 checks passed

TallJimbo deleted the tickets/DM-54397 branch March 31, 2026 15:41

Conversation

TallJimbo commented Mar 17, 2026

Uh oh!

laurenam 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!

laurenam 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!

Choose a reason for hiding this comment

Uh oh!

TallJimbo 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!

ctslater 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!

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!

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants