DM-26545: Add spline linearizer #56

czwa · 2020-10-08T17:49:52Z

This adds the code to fit a spline linearizer to the data contained within the PTC curve.

mfisherlevine · 2020-10-08T22:07:36Z

python/lsst/cp/pipe/linearity.py

+    fitOrder = pexConfig.Field(
        dtype=int,
-        doc="Degree of polynomial to fit.",
+        doc="Degree of polynomial to fit or number of spline knots to use",


I see the reasoning here on the rename (have both the number of spline knots and the poly order controlled by a single variable), but it does have the downside that 3 is a reasonable one for the poly order, but 3 spline knots is likely not a good default. I think there's certainly precedent elsewhere for having two different params here, so that one could change the type and still get a reasonable behaviour. I think it's your call, but just wanted to point that out as a choice so that it's conscious.

You can, of course, recombine them around L213 to use them in the same way, so the code there doesn't change at all, but just have 2 variables here so that you can set different defaults.

mfisherlevine · 2020-10-08T22:08:06Z

python/lsst/cp/pipe/linearity.py

            "LookupTable": "Create a lookup table solution.",
            "Polynomial": "Create an arbitrary polynomial solution.",
            "Squared": "Create a single order squared solution.",
+            "Spline": "Create a spline based solution.",
            "None": "Create a dummy solution.",


Is the departure from ALLCAPS for options like this a Gen3 thing?

This type needs to match the value listed in the LinearizeBase subclasses like at
https://github.com/lsst/ip_isr/blob/master/python/lsst/ip/isr/linearize.py#L539

mfisherlevine · 2020-10-08T22:09:15Z

python/lsst/cp/pipe/linearity.py

+    maxLinearAdu = pexConfig.Field(
+        dtype=float,
+        doc="Maximum DN value to use to estimate linear term.",
+        default=20000.0,
+    )
+    minLinearAdu = pexConfig.Field(
+        dtype=float,
+        doc="Minimum DN value to use to estimate linear term.",
+        default=10000.0,


I'm a little surprised that these values would be so high for bias-subtracted images, but you and @plazas are likely more up to date on that than I am, so again, whatever you say, just making sure that's a decent default.

I've dropped the value of minLinearAdu to 2000, as that's more representative of the data I've been testing on.

mfisherlevine · 2020-10-08T22:26:25Z

python/lsst/cp/pipe/linearity.py

@@ -31,12 +33,11 @@
 from .utils import (fitLeastSq, funcPolynomial)


Not part of this ticket necessarily, but the docstring on funcPolynomial is pretty weird, claims to return a covariance. Fancy cleaning it up while you're in here? (Also random comment on the return line too saying # C_00)

OK, I'll file a ticket on this, because I do think it's weird and a quick fix.

mfisherlevine · 2020-10-08T22:35:03Z

python/lsst/cp/pipe/linearity.py

+
+            # Exclude low end outliers
+            threshold = self.config.nSigmaClipLinear * np.sqrt(linearOrdinate)
+            fluxMask = np.abs(inputAbscissa - linearOrdinate) < threshold


I was totally following up to here, but wasn't quite sure about inputAbscissa - linearOrdinate but I'll trust that it's right, just writing this to save getting stuck here for any longer.

linearOrdinate is now a "flux-like" array, so it compares directly with inputAbscissa, the input flux array. The additional confusing point is that this fluxMask is masking values to retain, not exclude.

This was incorrect, and now uses np.abs(inputOrdinate - linearOrdinate). I confused myself (in many ways, but also) by converting to using the linear form of the flux (linearOrdinate) as the new X value in the fits, but without being clear about that.

mfisherlevine · 2020-10-15T23:50:17Z

python/lsst/cp/pipe/linearity.py

+
+                self.debugFit('splineFit', binCenters, np.abs(values), values, None, ampName)
+                interp = afwMath.makeInterpolate(binCenters.tolist(), values.tolist(),
+                                                 afwMath.stringToInterpStyle("AKIMA_SPLINE"))


Do we want AKIMA_SPLINE hard-coded here? This isn't the same as the not-a-knot, right? Don't be also want/need to be able to do that? I'd have thought that smooth handling of the boundaries was more important than dealing with fast-changing 2nd derivatives, but I could well be wrong.

The boundaries should work, I think. The low end will switch to a first order extrapolation (in the difference of input and ideal) based on the padding I put in below. The high end of the model should work up until nearly saturation, at which the linearity is never going to be accurate.

mfisherlevine · 2020-10-15T23:52:41Z

python/lsst/cp/pipe/linearity.py

+                # If we exclude a lot of points, we may end up with
+                # less than fitOrder points.  Pad out the low-flux end
+                # to ensure equal lengths.
+                if len(binCenters) != self.config.fitOrder:


A < would be more defensive, in that, if it's more, you've got binning/other problems, and they won't be handled by this block. Then again, having it fail on the pad line is probably no bad thing, so maybe ignore this.

mfisherlevine · 2020-10-15T23:54:35Z

python/lsst/cp/pipe/linearity.py

+                    binCenters = np.pad(binCenters, (padN, 0), 'linear_ramp',
+                                        end_values=(binCenters.min() - 1.0, ))
+                    # This stores the correction, which is zero at low values.
+                    values = np.pad(values, (padN, 0))


I don't have quite enough of this in my head right now to be confident that this padding is legit, but I certainly trust you, but just noting that I didn't think this through hard.

mfisherlevine · 2020-10-15T23:56:46Z

python/lsst/cp/pipe/linearity.py

+                    values = np.pad(values, (padN, 0))
+
+                # Pack the spline into a single array.
+                linearityFit = np.concatenate((binCenters.tolist(), values.tolist())).tolist()


This is an interesting data format, but it certainly at least explains the split I saw elsewhere!

It's the problem of packing things into the same single vector coefficients instead of adding a bunch of new attributes that are only optionally filled. I went with this because we pack and unpack crosstalk coefficients in an analogous way.

mfisherlevine · 2020-10-16T00:03:34Z

python/lsst/cp/pipe/linearity.py


+            image = afwImage.ImageF(len(inputOrdinate), 1)
+            image.getArray()[:, :] = inputAbscissa


I get that they're the same length, I just found it odd to declare it with the length of inputOrdinate but then assign it to inputAbscissa instead (unless this is a very-late-in-the-day form of cross checking)

- Add step to transform the input abscissa (exp time or MONDIODE value) into a linear flux measurement, based on the config options. - Filter the data based on the linear fit to exclude significant outliers (needed to exclude BOT neutral density data). - Add iteratively reweighted least squares solver to reduce the impact of remaining outlier points. - Update debugFit method to be more informative, and to provide information during the fitting stage. - Add gen2 implementation Task, and associated bin.src script. - Ensure spline data is padded to fitOrder length. - Properly write out output. - Set hasLinearity.

Be pedantic on updateMetadata to avoid issues with ingest.

Add spline algorithm comments to make the code less confusing.

czwa force-pushed the tickets/DM-26545 branch from 044087b to b618cc5 Compare October 8, 2020 19:26

czwa requested a review from mfisherlevine October 8, 2020 20:28

czwa force-pushed the tickets/DM-26545 branch from b618cc5 to 67f154a Compare October 12, 2020 19:51

mfisherlevine reviewed Oct 16, 2020

View reviewed changes

czwa force-pushed the tickets/DM-26545 branch 3 times, most recently from 0e6221a to 06b948f Compare October 20, 2020 23:23

czwa added 2 commits October 22, 2020 14:38

Use common metadata handler.

0c34e06

Be pedantic on updateMetadata to avoid issues with ingest.

Clarify variable names. Split fitOrder config options.

e996c76

Add spline algorithm comments to make the code less confusing.

czwa force-pushed the tickets/DM-26545 branch from 06b948f to e996c76 Compare October 22, 2020 19:39

czwa merged commit c512666 into master Oct 23, 2020

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

DM-26545: Add spline linearizer #56

DM-26545: Add spline linearizer #56

czwa commented Oct 8, 2020

mfisherlevine Oct 8, 2020

mfisherlevine Oct 8, 2020

mfisherlevine Oct 8, 2020

czwa Oct 19, 2020

mfisherlevine Oct 8, 2020

czwa Oct 19, 2020

mfisherlevine Oct 8, 2020

mfisherlevine Oct 21, 2020 •

edited

mfisherlevine Oct 21, 2020

mfisherlevine Oct 8, 2020

czwa Oct 19, 2020

czwa Oct 20, 2020

mfisherlevine Oct 15, 2020

czwa Oct 20, 2020

mfisherlevine Oct 15, 2020

mfisherlevine Oct 15, 2020

mfisherlevine Oct 15, 2020

czwa Oct 19, 2020

mfisherlevine Oct 16, 2020

		@@ -31,12 +33,11 @@
		from .utils import (fitLeastSq, funcPolynomial)


		image = afwImage.ImageF(len(inputOrdinate), 1)
		image.getArray()[:, :] = inputAbscissa

DM-26545: Add spline linearizer #56

DM-26545: Add spline linearizer #56

Conversation

czwa commented Oct 8, 2020

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mfisherlevine Oct 21, 2020 • edited

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mfisherlevine Oct 21, 2020 •

edited