Update not enough data behavior in stellarLocus fit

We set a minimum number of objects surviving all cuts to be deemed
suitable for including in the fit (defaults to 3).

python/lsst/analysis/tools/actions/keyedData/stellarLocusFit.py

@@ -60,6 +60,9 @@ def stellarLocusFit(xs, ys, mags, paramDict):
             The hardwired gradient for the fit (`float`).
         ``"bHw"``
             The hardwired intercept of the fit (`float`).
+        ``"minObjectForFit"``
+            Minimum number of objects surviving cuts to attempt fit.  If not
+            met, return NANs for values in ``fitParams`` (`int`).
 
     Returns
     -------
@@ -106,6 +109,10 @@ def stellarLocusFit(xs, ys, mags, paramDict):
         & (ys > paramDict["yMin"])
         & (ys < paramDict["yMax"])
     )
+    if sum(fitPoints) < paramDict["minObjectForFit"]:
+        fitParams = _setFitParamsNans(fitParams, fitPoints, paramDict)
+        return fitParams
+
     linear = scipyODR.polynomial(1)
 
     fitData = scipyODR.Data(xs[fitPoints], ys[fitPoints])
@@ -138,6 +145,9 @@ def stellarLocusFit(xs, ys, mags, paramDict):
 
         # Use these perpendicular lines to choose the data and refit
         fitPoints = (ys > mPerp0 * xs + bPerpMin) & (ys < mPerp0 * xs + bPerpMax)
+        if sum(fitPoints) < paramDict["minObjectForFit"]:
+            fitParams = _setFitParamsNans(fitParams, fitPoints, paramDict)
+            return fitParams
         p1 = np.array([0, bODR0])
         p2 = np.array([-bODR0 / mODR0, 0])
 
@@ -147,6 +157,9 @@ def stellarLocusFit(xs, ys, mags, paramDict):
         allDists = np.array(perpDistance(p1, p2, zip(xs, ys))) * 1000
         keep = np.logical_not(np.abs(allDists) > clippedStats.clipValue)
         fitPoints &= keep
+        if sum(fitPoints) < paramDict["minObjectForFit"]:
+            fitParams = _setFitParamsNans(fitParams, fitPoints, paramDict)
+            return fitParams
         fitData = scipyODR.Data(xs[fitPoints], ys[fitPoints])
         odr = scipyODR.ODR(fitData, linear, beta0=[bODR0, mODR0])
         params = odr.run()
@@ -167,6 +180,19 @@ def stellarLocusFit(xs, ys, mags, paramDict):
     return fitParams
 
 
+def _setFitParamsNans(fitParams, fitPoints, paramDict):
+    nan = float("nan")
+    fitParams["bPerpMin"] = nan
+    fitParams["bPerpMax"] = nan
+    fitParams["mODR"] = nan
+    fitParams["bODR"] = nan
+    fitParams["mPerp"] = nan
+    fitParams["goodPoints"] = nan
+    fitParams["fitPoints"] = fitPoints
+    fitParams["paramDict"] = paramDict
+    return fitParams
+
+
 def perpDistance(p1, p2, points):
     """Calculate the perpendicular distance to a line from a point.
 
@@ -257,11 +283,22 @@ class StellarLocusFitAction(KeyedDataAction):
     r"""Determine Stellar Locus fit parameters from given input `Vector`\ s."""
 
     stellarLocusFitDict = DictField[str, float](
-        doc="The parameters to use for the stellar locus fit. The default parameters are examples and are "
-        "not useful for any of the fits. The dict needs to contain xMin/xMax/yMin/yMax which are the "
-        "limits of the initial box for fitting the stellar locus, mHW and bHW are the initial "
-        "intercept and gradient for the fitting.",
-        default={"xMin": 0.1, "xMax": 0.2, "yMin": 0.1, "yMax": 0.2, "mHW": 0.5, "bHW": 0.0},
+        doc="The parameters to use for the stellar locus fit. Except for minObjectForFit, the "
+        "default parameters are examples and are not useful for any of the fits, so should be "
+        "updated in the PlotAction definition in the atools directory. The dict needs to "
+        "contain xMin/xMax/yMin/yMax which are the limits of the initial box for fitting the "
+        "stellar locus, mHW and bHW are the initial intercept and gradient for the fitting. "
+        "minObjectForFit sets a minimum number of points deemed suitable for inclusion in the "
+        "fit in order to bother attempting the fit.",
+        default={
+            "xMin": 0.1,
+            "xMax": 0.2,
+            "yMin": 0.1,
+            "yMax": 0.2,
+            "mHW": 0.5,
+            "bHW": 0.0,
+            "minObjectForFit": 3,
+        },
     )
 
     def getInputSchema(self) -> KeyedDataSchema:

python/lsst/analysis/tools/actions/plot/colorColorFitPlot.py

@@ -62,6 +62,9 @@ class ColorColorFitPlot(PlotAction):
         doc="Minimum number of valid objects to bother attempting a fit.",
         default=5,
         min=1,
+        deprecated="This field is no longer used. The value should go as an "
+        "entry to the paramsDict keyed as minObjectForFit.  Will be removed "
+        "after v27.",
     )
 
     xLims = ListField[float](
@@ -228,6 +231,20 @@ def makePlot(
         fitPoints = data.pop("fitPoints")
         # Points with finite values for x, y, and mag
         goodPoints = data.pop("goodPoints")
+
+        # TODO: Make a no data fig function and use here
+        if sum(fitPoints) < paramDict["minObjectForFit"]:
+            fig = plt.figure(dpi=120)
+            ax = fig.add_axes([0.12, 0.25, 0.43, 0.62])
+            ax.tick_params(labelsize=7)
+            noDataText = (
+                "Number of objects after cuts ({})\nis less than the minimum required\nby "
+                "paramDict[minObjectForFit] ({})".format(sum(fitPoints), int(paramDict["minObjectForFit"]))
+            )
+            plt.text(0.5, 0.5, noDataText, ha="center", va="center", fontsize=8)
+            fig = addPlotInfo(plt.gcf(), plotInfo)
+            return fig
+
         # Define new colormaps
         newBlues = mkColormap(["darkblue", "paleturquoise"])
         newGrays = mkColormap(["lightslategray", "white"])
@@ -245,17 +262,6 @@ def makePlot(
         ys = cast(Vector, data["y"])
         mags = cast(Vector, data["mag"])
 
-        # TODO: Make a no data fig function and use here
-        if len(fitPoints) < self.minPointsForFit:
-            fig = plt.figure(dpi=120)
-            noDataText = (
-                "Number of objects after cuts ({}) is less than the\nminimum required by "
-                "minPointsForFit ({})".format(len(fitPoints), self.minPointsForFit)
-            )
-            fig.text(0.5, 0.5, noDataText, ha="center", va="center")
-            fig = addPlotInfo(plt.gcf(), plotInfo)
-            return fig
-
         # Plot the initial fit box
         (initialBox,) = ax.plot(
             [paramDict["xMin"], paramDict["xMax"], paramDict["xMax"], paramDict["xMin"], paramDict["xMin"]],

python/lsst/analysis/tools/atools/stellarLocus.py

@@ -94,14 +94,16 @@ def setDefaults(self):
         self.process.buildActions.mag = ConvertFluxToMag(vectorKey=colorBands[1] + "_" + fluxType)
 
         self.process.calculateActions.wPerp_psfFlux = StellarLocusFitAction()
-        self.process.calculateActions.wPerp_psfFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 0.28,
             "xMax": 1.0,
             "yMin": 0.02,
             "yMax": 0.48,
             "mHW": 0.483,
             "bHW": -0.042,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.wPerp_psfFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "wPerp_psfFlux_sigmaMAD": "mmag",
@@ -138,14 +140,16 @@ def setDefaults(self):
         self.process.buildActions.mag = ConvertFluxToMag(vectorKey=colorBands[1] + "_" + fluxType)
 
         self.process.calculateActions.wPerp_cModelFlux = StellarLocusFitAction()
-        self.process.calculateActions.wPerp_cModelFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 0.28,
             "xMax": 1.0,
             "yMin": 0.02,
             "yMax": 0.48,
             "mHW": 0.483,
             "bHW": -0.042,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.wPerp_cModelFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "wPerp_cModelFlux_sigmaMAD": "mmag",
@@ -170,14 +174,16 @@ def setDefaults(self):
         fluxType = "psfFlux"
 
         self.process.calculateActions.xPerp_psfFlux = StellarLocusFitAction()
-        self.process.calculateActions.xPerp_psfFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 1.05,
             "xMax": 1.55,
             "yMin": 0.78,
             "yMax": 1.62,
             "mHW": 60.0,
             "bHW": -75.0,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.xPerp_psfFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "xPerp_psfFlux_sigmaMAD": "mmag",
@@ -192,14 +198,16 @@ def setDefaults(self):
         fluxType = "cModelFlux"
 
         self.process.calculateActions.xPerp_cModelFlux = StellarLocusFitAction()
-        self.process.calculateActions.xPerp_cModelFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 1.05,
             "xMax": 1.55,
             "yMin": 0.78,
             "yMax": 1.62,
             "mHW": 60.0,
             "bHW": -75.0,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.xPerp_cModelFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "xPerp_cModelFlux_sigmaMAD": "mmag",
@@ -231,14 +239,16 @@ def setDefaults(self):
         self.process.buildActions.mag = ConvertFluxToMag(vectorKey=colorBands[1] + "_" + fluxType)
 
         self.process.calculateActions.yPerp_psfFlux = StellarLocusFitAction()
-        self.process.calculateActions.yPerp_psfFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 0.82,
             "xMax": 2.01,
             "yMin": 0.37,
             "yMax": 0.90,
             "mHW": 0.385,
             "bHW": 0.064,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.yPerp_psfFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "yPerp_psfFlux_sigmaMAD": "mmag",
@@ -275,14 +285,16 @@ def setDefaults(self):
         self.process.buildActions.mag = ConvertFluxToMag(vectorKey=colorBands[1] + "_" + fluxType)
 
         self.process.calculateActions.yPerp_cModelFlux = StellarLocusFitAction()
-        self.process.calculateActions.yPerp_cModelFlux.stellarLocusFitDict = {
+        stellarLocusFitDict = {
             "xMin": 0.82,
             "xMax": 2.01,
             "yMin": 0.37,
             "yMax": 0.90,
             "mHW": 0.385,
             "bHW": 0.064,
         }
+        for item in stellarLocusFitDict.items():
+            self.process.calculateActions.yPerp_cModelFlux.stellarLocusFitDict.update([item])
 
         self.produce.metric.units = {
             "yPerp_cModelFlux_sigmaMAD": "mmag",

tests/test_plotUtils.py

@@ -38,7 +38,7 @@ def testFitLine(self):
         mags = np.arange(17.5, 22, 0.5)
 
         # Define an initial fit box that encompasses the points
-        testParams = {"xMin": 0, "xMax": 11, "yMin": 0, "yMax": 11, "mHW": 1, "bHW": 0}
+        testParams = {"xMin": 0, "xMax": 11, "yMin": 0, "yMax": 11, "mHW": 1, "bHW": 0, "minObjectForFit": 3}
         paramsOut = stellarLocusFit(xs, ys, mags, testParams)
 
         # stellarLocusFit performs two iterations of fitting and also