DM-37727: Refactor and develop further the ref cat matching code

pipelines/coaddQualityCore.yaml

@@ -33,9 +33,33 @@ tasks:
         from lsst.analysis.tools.contexts import *
         from lsst.analysis.tools.actions.plot import *
   catalogMatchTract:
-    class: lsst.analysis.tools.tasks.catalogMatch.CatalogMatchTask
+    class: lsst.analysis.tools.tasks.astrometricCatalogMatch.AstrometricCatalogMatchTask
+    config:
+      connections.refCatalog: gaia_dr2_20200414
   refCatObjectTract:
     class: lsst.analysis.tools.tasks.refCatObjectAnalysis.RefCatObjectAnalysisTask
+    config:
+      atools.astromDiffRAScatterPlot: TargetRefCatDeltaRAScatterPlot
+      atools.astromDiffDecScatterPlot: TargetRefCatDeltaDecScatterPlot
+      atools.astromDiffRASkyPlot: TargetRefCatDeltaRASkyPlot
+      atools.astromDiffDecSkyPlot: TargetRefCatDeltaDecSkyPlot
+      python: |
+        from lsst.analysis.tools.atools import *
+
+    ## Perform photometric comparison to the reference catalog
+  photometricCatalogMatch:
+    class: lsst.analysis.tools.tasks.photometricCatalogMatch.PhotometricCatalogMatchTask
+    config:
+      connections.refCatalog: ps1_pv3_3pi_20170110
+  photometricRefCatObjectTract:
+    class: lsst.analysis.tools.tasks.refCatObjectPhotometricAnalysis.RefCatObjectPhotometricAnalysisTask
+    config:
+      atools.targetRefCatDeltaPsfScatterPlot: TargetRefCatDeltaPsfScatterPlot
+      atools.targetRefCatDeltaCModelScatterPlot: TargetRefCatDeltaCModelScatterPlot
+      atools.targetRefCatDeltaPsfSkyPlot: TargetRefCatDeltaPsfSkyPlot
+      atools.targetRefCatDeltaCModelSkyPlot: TargetRefCatDeltaCModelSkyPlot
+      python: |
+        from lsst.analysis.tools.atools import *
   plotPropertyMapTract:
     class: lsst.analysis.tools.tasks.PropertyMapTractAnalysisTask
     config:

pipelines/visitQualityCore.yaml

@@ -12,3 +12,28 @@ tasks:
       python: |
         from lsst.analysis.tools.atools import *
         from lsst.analysis.tools.contexts import *
+  catalogMatchVisit:
+    class: lsst.analysis.tools.tasks.astrometricCatalogMatch.AstrometricCatalogMatchVisitTask
+    config:
+      connections.refCatalog: gaia_dr2_20200414
+  astrometricRefCatSourceVisit:
+    class: lsst.analysis.tools.tasks.refCatSourceAnalysis.RefCatSourceAnalysisTask
+    config:
+      atools.astromDiffRASkyVisitPlot: TargetRefCatDeltaRASkyVisitPlot
+      atools.astromDiffDecSkyVisitPlot: TargetRefCatDeltaDecSkyVisitPlot
+      atools.astromDiffRAScatterVisitPlot: TargetRefCatDeltaRAScatterVisitPlot
+      atools.astromDiffDecScatterVisitPlot: TargetRefCatDeltaDecScatterVisitPlot
+      python: |
+        from lsst.analysis.tools.atools import *
+  photometricMatchVisit:
+    class: lsst.analysis.tools.tasks.photometricCatalogMatch.PhotometricCatalogMatchVisitTask
+    config:
+      extraColumns: ["x", "y", "ap09Flux", "ap09FluxErr"]
+  photometricRefCatSourceVisit:
+    class: lsst.analysis.tools.tasks.refCatSourcePhotometricAnalysis.RefCatSourcePhotometricAnalysisTask
+    config:
+      atools.photomDiffPsfSkyVisitPlot: TargetRefCatDeltaPsfSkyVisitPlot
+      atools.photomDiffAp09SkyVisitPlot: TargetRefCatDeltaAp09SkyVisitPlot
+      atools.photoDiffPsfScatterVisitPlot: TargetRefCatDeltaPsfScatterVisitPlot
+      atools.photoDiffCModelScatterVisitPlot: TargetRefCatDeltaAp09ScatterVisitPlot
+      python: from lsst.analysis.tools.atools import *

pipelines/visitQualityExtended.yaml

@@ -2,6 +2,6 @@ description: |
   Tier2 plots and metrics to assess visit quality
 tasks:
   catalogMatchVisit:
-    class: lsst.analysis.tools.tasks.catalogMatch.CatalogMatchVisitTask
+    class: lsst.analysis.tools.tasks.astrometricCatalogMatch.AstrometricCatalogMatchVisitTask
   refCatSourceVisit:
     class: lsst.analysis.tools.tasks.refCatSourceAnalysis.RefCatSourceAnalysisTask

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

@@ -508,7 +508,7 @@ def _getPanelRange(self, data, panel, mads=None, meds=None):
             minMed = np.nanmin(meds)
             panel_range = [minMed - lowerRange * maxMad, maxMed + upperRange * maxMad]
             if panel_range[1] - panel_range[0] == 0:
-                self.log.info(
+                log.info(
                     "NOTE: panel_range for {} based on med/sigMad was 0. Computing using "
                     "percentile range instead.".format(panel)
                 )

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

@@ -123,7 +123,10 @@ def generateSummaryStats(data, skymap, plotInfo):
         # Once the objectTable_tract catalogues are using gen 3 patches
         # this will go away
         onPatch = data["patch"] == patch
-        stat = np.nanmedian(data[yCol][onPatch])
+        if sum(onPatch) == 0:
+            stat = np.nan
+        else:
+            stat = np.nanmedian(data[yCol][onPatch])
         try:
             patchTuple = (int(patch.split(",")[0]), int(patch.split(",")[-1]))
             patchInfo = tractInfo.getPatchInfo(patchTuple)

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

@@ -178,7 +178,7 @@ class ScatterPlotWithTwoHists(PlotAction):
 
     addSummaryPlot = Field[bool](
         doc="Add a summary plot to the figure?",
-        default=False,
+        default=True,
     )
 
     _stats = ("median", "sigmaMad", "count", "approxMag")
@@ -673,9 +673,9 @@ def _makeTopHistogram(
         # Top histogram
         topHist = figure.add_subplot(gs[0, :-1], sharex=ax)
 
-        if "all" in self.plotTypes:
-            x_all = f"x{self._datatypes['all'].suffix_xy}"
-            keys_notall = [x for x in self.plotTypes if x != "all"]
+        if "any" in self.plotTypes:
+            x_all = data[f"x{self._datatypes['any'].suffix_xy}"]
+            keys_notall = [x for x in self.plotTypes if x != "any"]
         else:
             x_all = np.concatenate([data[f"x{self._datatypes[key].suffix_xy}"] for key in self.plotTypes])
             keys_notall = self.plotTypes
@@ -712,16 +712,16 @@ def _makeSideHistogram(
         # Side histogram
         sideHist = figure.add_subplot(gs[1:, -1], sharey=ax)
 
-        if "all" in self.plotTypes:
-            y_all = f"y{self._datatypes['all'].suffix_xy}"
-            keys_notall = [x for x in self.plotTypes if x != "all"]
+        if "any" in self.plotTypes:
+            y_all = data[f"y{self._datatypes['any'].suffix_xy}"]
+            keys_notall = [x for x in self.plotTypes if x != "any"]
         else:
             y_all = np.concatenate([data[f"y{self._datatypes[key].suffix_xy}"] for key in self.plotTypes])
             keys_notall = self.plotTypes
 
         y_min, y_max = ax.get_ylim()
         bins = np.linspace(y_min, y_max, 100)
-        sideHist.hist(y_all, bins=bins, color="grey", alpha=0.3, orientation="horizontal", log=True)
+        sideHist.hist(np.array(y_all), bins=bins, color="grey", alpha=0.3, orientation="horizontal", log=True)
         kwargs_hist = dict(
             bins=bins,
             histtype="step",

python/lsst/analysis/tools/actions/scalar/scalarActions.py

@@ -36,7 +36,12 @@ def getInputSchema(self) -> KeyedDataSchema:
 
     def __call__(self, data: KeyedData, **kwargs) -> Scalar:
         mask = self.getMask(**kwargs)
-        return cast(Scalar, float(np.nanmedian(cast(Vector, data[self.vectorKey.format(**kwargs)])[mask])))
+        if len(data[self.vectorKey.format(**kwargs)][mask]) != 0:
+            med = cast(Scalar, float(np.nanmedian(cast(Vector, data[self.vectorKey.format(**kwargs)])[mask])))
+        else:
+            med = np.NaN
+
+        return med
 
 
 class MeanAction(ScalarAction):
@@ -49,7 +54,12 @@ def getInputSchema(self) -> KeyedDataSchema:
 
     def __call__(self, data: KeyedData, **kwargs) -> Scalar:
         mask = self.getMask(**kwargs)
-        return cast(Scalar, float(np.nanmean(cast(Vector, data[self.vectorKey.format(**kwargs)])[mask])))
+        if len(data[self.vectorKey.format(**kwargs)][mask]) != 0:
+            mean = cast(Scalar, float(np.nanmean(cast(Vector, data[self.vectorKey.format(**kwargs)])[mask])))
+        else:
+            mean = np.NaN
+
+        return mean
 
 
 class StdevAction(ScalarAction):

python/lsst/analysis/tools/actions/vector/selectors.py

@@ -127,6 +127,15 @@ class CoaddPlotFlagSelector(FlagSelector):
     def getInputSchema(self) -> KeyedDataSchema:
         yield from super().getInputSchema()
 
+    def refMatchContext(self):
+        self.selectWhenFalse = [
+            "{band}_psfFlux_flag_target",
+            "{band}_pixelFlags_saturatedCenter_target",
+            "{band}_extendedness_flag_target",
+            "xy_flag_target",
+        ]
+        self.selectWhenTrue = ["detect_isPatchInner_target", "detect_isDeblendedSource_target"]
+
     def __call__(self, data: KeyedData, **kwargs) -> Vector:
         result: Optional[Vector] = None
         bands: tuple[str, ...]
@@ -162,9 +171,19 @@ class VisitPlotFlagSelector(FlagSelector):
     (i.e., using sourceTable_visit catalogs).
     """
 
+    catalogSuffix = Field[str](doc="The suffix to apply to all the keys.", default="")
+
     def getInputSchema(self) -> KeyedDataSchema:
         yield from super().getInputSchema()
 
+    def refMatchContext(self):
+        self.selectWhenFalse = [
+            "psfFlux_flag_target",
+            "pixelFlags_saturatedCenter_target",
+            "extendedness_flag_target",
+            "centroid_flag_target",
+        ]
+
     def __call__(self, data: KeyedData, **kwargs) -> Vector:
         result: Optional[Vector] = None
         temp = super().__call__(data, **kwargs)
@@ -227,8 +246,11 @@ class SnSelector(SelectorBase):
     )
 
     def getInputSchema(self) -> KeyedDataSchema:
-        yield (fluxCol := self.fluxType), Vector
-        yield f"{fluxCol}{self.uncertaintySuffix}", Vector
+        fluxCol = self.fluxType
+        fluxInd = fluxCol.find("lux") + len("lux")
+        errCol = f"{fluxCol}"[:fluxInd] + f"{self.uncertaintySuffix}" + f"{fluxCol}"[fluxInd:]
+        yield fluxCol, Vector
+        yield errCol, Vector
 
     def __call__(self, data: KeyedData, **kwargs) -> Vector:
         """Makes a mask of objects that have S/N greater than
@@ -260,8 +282,11 @@ def __call__(self, data: KeyedData, **kwargs) -> Vector:
                 bands = ("",)
         for band in bands:
             fluxCol = self.fluxType.format(**(kwargs | dict(band=band)))
-            errCol = f"{fluxCol}{self.uncertaintySuffix.format(**kwargs)}"
-            vec = cast(Vector, data[fluxCol]) / cast(Vector, data[errCol])
+            fluxInd = fluxCol.find("lux") + len("lux")
+            errCol = (
+                f"{fluxCol}"[:fluxInd] + f"{self.uncertaintySuffix.format(**kwargs)}" + f"{fluxCol}"[fluxInd:]
+            )
+            vec = cast(Vector, data[fluxCol]) / data[errCol]
             temp = (vec > self.threshold) & (vec < self.maxSN)
             if mask is not None:
                 mask &= temp  # type: ignore

python/lsst/analysis/tools/actions/vector/vectorActions.py

@@ -289,9 +289,9 @@ def getInputSchema(self) -> KeyedDataSchema:
         return ((self.decKey, Vector), (self.raKey, Vector))
 
     def __call__(self, data: KeyedData, **kwargs) -> Vector:
-        ra = data[self.raKey]
-        dec = data[self.decKey]
-        return ra.to_numpy() * np.cos((dec.to_numpy() * u.degree).to(u.radian).value)
+        ra = np.array(data[self.raKey])
+        dec = np.array(data[self.decKey])
+        return ra * np.cos((dec * u.degree).to(u.radian).value)
 
 
 # Statistical vectorActions

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

@@ -1,5 +1,4 @@
 from .astrometricRepeatability import *
-from .astrometryWithReference import *
 from .coveragePlots import *
 from .diaSolarSystemObjectMetrics import *
 from .diaSourceMetrics import *
@@ -15,6 +14,7 @@
 from .photometricRepeatability import *
 from .photometry import *
 from .propertyMap import *
+from .refCatMatchPlots import *
 from .seeingMetric import *
 from .shapes import *
 from .simpleDiaPlot import *