From 8c027567e5b38496ddaf5ab48023dbc2bea98086 Mon Sep 17 00:00:00 2001 From: Trys McCann Date: Thu, 21 Mar 2024 08:21:10 -0700 Subject: [PATCH] Add QA plots for source injection to analysis_tools --- pipelines/injectedCoaddQualityCore.yaml | 16 + .../tools/actions/keyedData/__init__.py | 1 + .../tools/actions/keyedData/magPercentiles.py | 96 ++++++ .../analysis/tools/actions/plot/__init__.py | 1 + .../tools/actions/plot/completenessPlot.py | 183 ++++++++++++ python/lsst/analysis/tools/atools/__init__.py | 1 + .../tools/atools/sourceInjectionPlots.py | 281 ++++++++++++++++++ python/lsst/analysis/tools/math.py | 7 + python/lsst/analysis/tools/tasks/__init__.py | 1 + .../tools/tasks/injectedObjectAnalysis.py | 64 ++++ 10 files changed, 651 insertions(+) create mode 100644 pipelines/injectedCoaddQualityCore.yaml create mode 100644 python/lsst/analysis/tools/actions/keyedData/magPercentiles.py create mode 100644 python/lsst/analysis/tools/actions/plot/completenessPlot.py create mode 100644 python/lsst/analysis/tools/atools/sourceInjectionPlots.py create mode 100644 python/lsst/analysis/tools/tasks/injectedObjectAnalysis.py diff --git a/pipelines/injectedCoaddQualityCore.yaml b/pipelines/injectedCoaddQualityCore.yaml new file mode 100644 index 000000000..6e4d30ea4 --- /dev/null +++ b/pipelines/injectedCoaddQualityCore.yaml @@ -0,0 +1,16 @@ +description: | + Tier1 plots and metrics to assess injected coadd quality +tasks: + injectedObjectAnalysis: + class: lsst.analysis.tools.tasks.injectedObjectAnalysis.InjectedObjectAnalysisTask + config: + atools.completenessHist: CompletenessPurityTool + atools.targetInjectedCatDeltaRAScatterPlot: TargetInjectedCatDeltaRAScatterPlot + atools.targetInjectedCatDeltaDecScatterPlot: TargetInjectedCatDeltaDecScatterPlot + atools.targetInjectedCatDeltaPsfScatterPlot: TargetInjectedCatDeltaPsfScatterPlot + atools.injectedMatchDiffMetrics: TargetInjectedCatDeltaMetrics + atools.injectedMatchDiffMetrics.applyContext: CoaddContext + bands: ["g", "r", "i", "z", "y"] + python: | + from lsst.analysis.tools.atools import * + from lsst.analysis.tools.contexts import * diff --git a/python/lsst/analysis/tools/actions/keyedData/__init__.py b/python/lsst/analysis/tools/actions/keyedData/__init__.py index 50450b7e6..2ba32241d 100644 --- a/python/lsst/analysis/tools/actions/keyedData/__init__.py +++ b/python/lsst/analysis/tools/actions/keyedData/__init__.py @@ -1,3 +1,4 @@ from .calcDistances import * from .keyedDataActions import * +from .magPercentiles import * from .stellarLocusFit import * diff --git a/python/lsst/analysis/tools/actions/keyedData/magPercentiles.py b/python/lsst/analysis/tools/actions/keyedData/magPercentiles.py new file mode 100644 index 000000000..a842de8d4 --- /dev/null +++ b/python/lsst/analysis/tools/actions/keyedData/magPercentiles.py @@ -0,0 +1,96 @@ +# This file is part of analysis_tools. +# +# Developed for the LSST Data Management System. +# This product includes software developed by the LSST Project +# (https://www.lsst.org). +# See the COPYRIGHT file at the top-level directory of this distribution +# for details of code ownership. +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +from __future__ import annotations + +__all__ = ("MagPercentileAction",) + +import logging + +import numpy as np + +# import pandas as pd +from astropy import units as u + +# from astropy.coordinates import SkyCoord +from lsst.pex.config import Field, ListField + +from ...interfaces import KeyedData, KeyedDataSchema, Scalar, Vector, VectorAction +from ...math import fluxToMag, isPercent # divide, fluxToMag, isPercent, log10 + +# from typing import Optional, cast + +# from lsst.pex.config.configurableActions import ConfigurableActionField, ConfigurableActionStructField + +# from .selectors import VectorSelector + +_LOG = logging.getLogger(__name__) + + +class MagPercentileAction(VectorAction): + """Calculates the magnitude at the given percentile for completeness""" + + matchDistanceKey = Field[str]("Match distance Vector") + vectorKey = Field[str](doc="Key of vector which should be loaded") + fluxUnits = Field[str](doc="Units for the column.", default="nanojansky") + percentiles = ListField[float]( + doc="The percentiles to find the magnitude at.", default=[16.0, 50.0, 84.0], itemCheck=isPercent + ) + + def getInputSchema(self) -> KeyedDataSchema: + return ( + (self.matchDistanceKey, Vector), + (self.vectorKey, Vector), + ) + + def getOutputSchema(self) -> KeyedDataSchema: + result = [] + for pct in self.percentiles: + name = self.getPercentileName(pct) + result.append((name, Scalar)) + return result + + def getPercentileName(self, percentile: float) -> str: + return f"mag_{percentile:.2f}" + + def __call__(self, data: KeyedData, **kwargs) -> KeyedData: + matched = np.isfinite(data[self.matchDistanceKey]) + fluxValues = data[self.vectorKey.format(**kwargs)] + values = fluxToMag(fluxValues, flux_unit=u.Unit(self.fluxUnits)) + nInput, bins = np.histogram( + values, + range=(np.nanmin(values), np.nanmax(values)), + bins=100, + ) + nOutput, _ = np.histogram( + values[matched], + range=(np.nanmin(values[matched]), np.nanmax(values[matched])), + bins=bins, + ) + # Find bin where the fraction recovered first falls below a percentile. + mags: KeyedData = {} + for pct in self.percentiles: + name = self.getPercentileName(pct) + belowPercentile = np.where((nOutput / nInput < pct / 100))[0] + if len(belowPercentile) == 0: + mags[name] = np.nan + else: + mags[name] = np.min(bins[belowPercentile]) + return mags diff --git a/python/lsst/analysis/tools/actions/plot/__init__.py b/python/lsst/analysis/tools/actions/plot/__init__.py index 4eee8fd60..af77f1712 100644 --- a/python/lsst/analysis/tools/actions/plot/__init__.py +++ b/python/lsst/analysis/tools/actions/plot/__init__.py @@ -1,6 +1,7 @@ from .barPlots import * from .calculateRange import * from .colorColorFitPlot import * +from .completenessPlot import * from .diaSkyPlot import * from .focalPlanePlot import * from .gridPlot import * diff --git a/python/lsst/analysis/tools/actions/plot/completenessPlot.py b/python/lsst/analysis/tools/actions/plot/completenessPlot.py new file mode 100644 index 000000000..7455aee94 --- /dev/null +++ b/python/lsst/analysis/tools/actions/plot/completenessPlot.py @@ -0,0 +1,183 @@ +# This file is part of analysis_tools. +# +# Developed for the LSST Data Management System. +# This product includes software developed by the LSST Project +# (https://www.lsst.org). +# See the COPYRIGHT file at the top-level directory of this distribution +# for details of code ownership. +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . + + +from typing import Mapping + +import matplotlib.pyplot as plt +import numpy as np +from lsst.pex.config import Field, ListField +from matplotlib.figure import Figure + +from ...interfaces import KeyedData, KeyedDataSchema, PlotAction, Scalar, ScalarType, Vector +from .plotUtils import addPlotInfo + +__all__ = ("CompletenessHist",) + + +class CompletenessHist(PlotAction): + """Makes a scatter plot of the data with a marginal + histogram for each axis. + """ + + magKey = Field[str](doc="Name of the magnitude column.", default="mag") + matchDistanceKey = Field[str](doc="Name of the match distance column.", default="matchDistance") + xAxisLabel = Field[str](doc="Label for the x axis.", default="Input Magnitude (mag)") + inputLabel = Field[str](doc="Label for the input source histogram.", default="Synthetic Inputs") + outputLabel = Field[str](doc="Label for the recovered source histogram.", default="Synthetic Recovered") + numBins = Field[int](doc="Number of bins to use for the histograms.", default=100) + completenessPercentiles = ListField[float]( + doc="Record the magnitudes at these percentiles", default=[84.0, 50.0, 16.0] + ) + + def getInputSchema(self) -> KeyedDataSchema: + base: list[tuple[str, type[Vector] | ScalarType]] = [] + base.append((self.magKey, Vector)) + base.append((self.matchDistanceKey, Vector)) + return base + + def __call__(self, data: KeyedData, **kwargs) -> Mapping[str, Figure] | Figure: + self._validateInput(data, **kwargs) + return self.makePlot(data, **kwargs) + + def _validateInput(self, data: KeyedData, **kwargs) -> None: + """NOTE currently can only check that something is not a Scalar, not + check that the data is consistent with Vector + """ + needed = self.getFormattedInputSchema(**kwargs) + if remainder := {key.format(**kwargs) for key, _ in needed} - { + key.format(**kwargs) for key in data.keys() + }: + raise ValueError(f"Task needs keys {remainder} but they were not found in input") + for name, typ in needed: + isScalar = issubclass((colType := type(data[name.format(**kwargs)])), Scalar) + if isScalar and typ != Scalar: + raise ValueError(f"Data keyed by {name} has type {colType} but action requires type {typ}") + + def makePlot(self, data, plotInfo, **kwargs): + """Makes a plot showing the fraction of injected sources recovered by + input magnitude. + + Parameters + ---------- + data : `KeyedData` + All the data + plotInfo : `dict` + A dictionary of information about the data being plotted with keys: + ``camera`` + The camera used to take the data (`lsst.afw.cameraGeom.Camera`) + ``"cameraName"`` + The name of camera used to take the data (`str`). + ``"filter"`` + The filter used for this data (`str`). + ``"ccdKey"`` + The ccd/dectector key associated with this camera (`str`). + ``"visit"`` + The visit of the data; only included if the data is from a + single epoch dataset (`str`). + ``"patch"`` + The patch that the data is from; only included if the data is + from a coadd dataset (`str`). + ``"tract"`` + The tract that the data comes from (`str`). + ``"photoCalibDataset"`` + The dataset used for the calibration, e.g. "jointcal" or "fgcm" + (`str`). + ``"skyWcsDataset"`` + The sky Wcs dataset used (`str`). + ``"rerun"`` + The rerun the data is stored in (`str`). + + Returns + ------ + ``fig`` + The figure to be saved (`matplotlib.figure.Figure`). + + Notes + ----- + Makes a histogram showing the fraction recovered in each magnitude + bin with the number input and recovered overplotted. + """ + + # Make plot showing the fraction recovered in magnitude bins + fig, axLeft = plt.subplots(dpi=300) + axLeft.tick_params(axis="y", labelcolor="C0") + axLeft.set_xlabel(self.xAxisLabel) + axLeft.set_ylabel("Fraction Recovered", color="C0") + axRight = axLeft.twinx() + axRight.set_ylabel("Number of Sources") + matched = np.isfinite(data[self.matchDistanceKey]) + nInput, bins, _ = axRight.hist( + data[self.magKey], + range=(np.nanmin(data[self.magKey]), np.nanmax(data[self.magKey])), + bins=self.numBins, + log=True, + histtype="step", + label=self.inputLabel, + color="black", + ) + nOutput, _, _ = axRight.hist( + data[self.magKey][matched], + range=(np.nanmin(data[self.magKey][matched]), np.nanmax(data[self.magKey][matched])), + bins=bins, + log=True, + histtype="step", + label=self.outputLabel, + color="grey", + ) + + # Find bin where the fraction recovered falls below a given percentile. + percentileInfo = [] + xlims = plt.gca().get_xlim() + for pct in self.completenessPercentiles: + pct /= 100 + magArray = np.where((nOutput / nInput < pct))[0] + if len(magArray) == 0: + mag = np.nan + else: + mag = np.min(bins[magArray]) + axLeft.plot([xlims[0], mag], [pct, pct], ls=":", color="grey") + axLeft.plot([mag, mag], [0, pct], ls=":", color="grey") + percentileInfo.append("Magnitude at {}% recovered: {:0.2f}".format(pct * 100, mag)) + plt.xlim(xlims) + axLeft.set_ylim(0, 1.05) + axRight.legend(loc="lower left", ncol=2) + axLeft.axhline(1, color="grey", ls="--") + axLeft.bar( + bins[:-1], + nOutput / nInput, + width=np.diff(bins), + align="edge", + color="C0", + alpha=0.5, + zorder=10, + ) + + # Add useful information to the plot + fig = plt.gcf() + addPlotInfo(fig, plotInfo) + statsText = "" + for info in percentileInfo: + statsText += f"{info}\n" + bbox = dict(edgecolor="grey", linestyle=":", facecolor="none") + fig.text(0.7, 0.075, statsText[:-1], bbox=bbox, transform=fig.transFigure, fontsize=6) + fig.subplots_adjust(bottom=0.2) + return fig diff --git a/python/lsst/analysis/tools/atools/__init__.py b/python/lsst/analysis/tools/atools/__init__.py index 0ec8bc027..a1aeab570 100644 --- a/python/lsst/analysis/tools/atools/__init__.py +++ b/python/lsst/analysis/tools/atools/__init__.py @@ -33,6 +33,7 @@ from .skyFluxStatisticMetrics import * from .skyObject import * from .skySource import * +from .sourceInjectionPlots import * from .sources import * from .stellarLocus import * from .wholeSkyPlotTool import * diff --git a/python/lsst/analysis/tools/atools/sourceInjectionPlots.py b/python/lsst/analysis/tools/atools/sourceInjectionPlots.py new file mode 100644 index 000000000..092512c4d --- /dev/null +++ b/python/lsst/analysis/tools/atools/sourceInjectionPlots.py @@ -0,0 +1,281 @@ +# This file is part of analysis_tools. +# +# Developed for the LSST Data Management System. +# This product includes software developed by the LSST Project +# (https://www.lsst.org). +# See the COPYRIGHT file at the top-level directory of this distribution +# for details of code ownership. +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +from __future__ import annotations + +__all__ = ( + "CompletenessPurityTool", + "TargetInjectedCatDelta", + "TargetInjectedCatDeltaScatterAstrom", + "TargetInjectedCatDeltaRAScatterPlot", + "TargetInjectedCatDeltaDecScatterPlot", + "TargetInjectedCatDeltaScatterPhotom", + "TargetInjectedCatDeltaPsfScatterPlot", + "TargetInjectedCatDeltaCModelScatterPlot", + "TargetInjectedCatDeltaMetrics", +) + +from lsst.pex.config import Field, ListField + +from ..actions.keyedData import MagPercentileAction +from ..actions.plot.completenessPlot import CompletenessHist +from ..actions.plot.scatterplotWithTwoHists import ScatterPlotStatsAction, ScatterPlotWithTwoHists +from ..actions.scalar.scalarActions import MedianAction, SigmaMadAction +from ..actions.vector import ( + ConvertFluxToMag, + ConvertUnits, + DownselectVector, + LoadVector, + MagDiff, + RAcosDec, + RangeSelector, + SnSelector, + SubtractVector, +) +from ..contexts import CoaddContext +from ..interfaces import AnalysisTool + + +class CompletenessPurityTool(AnalysisTool): + """Plot the completeness or purity of injected sources by magnitude.""" + + completenessPercentiles = ListField[float](doc="tmp", default=[84.0, 50.0, 16.0]) + + def setDefaults(self): + super().setDefaults() + + self.process.buildActions.mag = ConvertFluxToMag(vectorKey="ref_{band}_flux") + self.process.buildActions.matchDistance = LoadVector(vectorKey="match_distance") + self.process.calculateActions.magPercentiles = MagPercentileAction( + matchDistanceKey="match_distance", + vectorKey="ref_{band}_flux", + percentiles=self.completenessPercentiles, + ) + + self.produce.plot = CompletenessHist() + self.produce.metric.units = {} + self.produce.metric.newNames = {} + + def finalize(self): + super().finalize() + units = {} + newNames = {} + for percentile in self.completenessPercentiles: + name = self.process.calculateActions.magPercentiles.getPercentileName(percentile) + units[name] = "mag" + newNames[name] = "{band}_" + name + self.produce.metric.units = units + self.produce.metric.newNames = newNames + + +class TargetInjectedCatDelta(AnalysisTool): + """Plot the difference between a target catalog and an + injected catalog for the quantity set in `setDefaults`. + """ + + parameterizedBand = Field[bool]( + doc="Does this AnalysisTool support band as a name parameter", default=True + ) + + def coaddContext(self) -> None: + """Apply coadd options for the ref cat plots. + Applies the coadd plot flag selector and sets + flux types. + """ + self.prep.selectors.snSelector = SnSelector(fluxType="{band}_psfFlux", threshold=200) + self.process.calculateActions.stars = ScatterPlotStatsAction( + vectorKey="yStars", + fluxType="{band}_psfFlux", + lowSNSelector=SnSelector(fluxType="{band}_psfFlux", threshold=300), + highSNSelector=SnSelector(fluxType="{band}_psfFlux", threshold=2700), + ) + self.process.buildActions.starStatMask = SnSelector(fluxType="{band}_psfFlux") + self.process.buildActions.patch = LoadVector(vectorKey="patch") + + +class TargetInjectedCatDeltaScatterAstrom(TargetInjectedCatDelta): + """Plot the difference in milliseconds between a target catalog and an + injected catalog for the coordinate set in `setDefaults`. Plot it on + a scatter plot. + """ + + def setDefaults(self): + super().setDefaults() + + self.process.buildActions.yStars = ConvertUnits( + buildAction=SubtractVector, inUnit="degree", outUnit="milliarcsecond" + ) + self.process.buildActions.xStars = ConvertFluxToMag(vectorKey="{band}_psfFlux") + self.process.calculateActions.stars = ScatterPlotStatsAction( + vectorKey="yStars", fluxType="{band}_psfFlux" + ) + self.produce = ScatterPlotWithTwoHists( + plotTypes=["stars"], magLabel="PSF Magnitude (mag)", xAxisLabel="PSF Magnitude (mag)" + ) + self.applyContext(CoaddContext) + + +class TargetInjectedCatDeltaRAScatterPlot(TargetInjectedCatDeltaScatterAstrom): + """Plot the difference in milliseconds between the RA of a target catalog + and an injected catalog + """ + + def setDefaults(self): + super().setDefaults() + self.process.buildActions.yStars.buildAction.actionA = RAcosDec(raKey="coord_ra", decKey="coord_dec") + self.process.buildActions.yStars.buildAction.actionB = RAcosDec(raKey="ref_ra", decKey="ref_dec") + + self.produce.yAxisLabel = "RA$_{{output}}$ - RA$_{{input}}$ (marcsec)" + + +class TargetInjectedCatDeltaDecScatterPlot(TargetInjectedCatDeltaScatterAstrom): + """Plot the difference in milliseconds between the Dec of a target catalog + and an injected catalog + """ + + def setDefaults(self): + super().setDefaults() + self.process.buildActions.yStars.buildAction.actionA = LoadVector(vectorKey="coord_dec") + self.process.buildActions.yStars.buildAction.actionB = LoadVector(vectorKey="ref_dec") + + self.produce.yAxisLabel = "Dec$_{{output}}$ - Dec$_{{input}}$ (marcsec)" + + +class TargetInjectedCatDeltaScatterPhotom(TargetInjectedCatDelta): + """Plot the difference in millimags between a target catalog and an + injected catalog for the flux type set in `setDefaults`. + """ + + def setDefaults(self): + super().setDefaults() + + self.process.buildActions.yStars = MagDiff(col2="ref_{band}_flux") + self.process.buildActions.xStars = ConvertFluxToMag(vectorKey="{band}_psfFlux") + self.produce = ScatterPlotWithTwoHists( + plotTypes=["stars"], + magLabel="PSF Magnitude (mag)", + xAxisLabel="PSF Magnitude (mag)", + yAxisLabel="Output Mag - Input Mag (mmag)", + ) + self.applyContext(CoaddContext) + + +class TargetInjectedCatDeltaPsfScatterPlot(TargetInjectedCatDeltaScatterPhotom): + """Plot the difference in millimags between the PSF flux + of a target catalog and an injected catalog + """ + + def setDefaults(self): + super().setDefaults() + + self.process.buildActions.yStars.col1 = "{band}_psfFlux" + + +class TargetInjectedCatDeltaCModelScatterPlot(TargetInjectedCatDeltaScatterPhotom): + """Plot the difference in millimags between the CModel flux + of a target catalog and an injected catalog. + """ + + def setDefaults(self): + super().setDefaults() + + self.process.buildActions.yStars.col1 = "{band}_cModelFlux" + + +class TargetInjectedCatDeltaMetrics(AnalysisTool): + """Calculate the diff metric and the sigma MAD from the difference between + the target and injected catalog coordinates and photometry. + """ + + def coaddContext(self) -> None: + """Apply coadd options for the metrics. Applies the coadd plot flag + selector and sets flux types. + """ + + self.process.buildActions.mags.vectorKey = "{band}_psfFlux" + + self.produce.metric.newNames = { + "injected_RA_diff_median": "injected_{band}_RA_diff_median_coadd", + "injected_RA_diff_sigmaMad": "injected_{band}_RA_diff_sigmaMad_coadd", + "injected_Dec_diff_median": "injected_{band}_Dec_diff_median_coadd", + "injected_Dec_diff_sigmaMad": "injected_{band}_Dec_diff_sigmaMad_coadd", + "injected_phot_diff_median": "injected_{band}_mag_diff_median_coadd", + "injected_phot_diff_sigmaMad": "injected_{band}_mag_diff_sigmaMad_coadd", + } + + def setDefaults(self): + super().setDefaults() + + # Calculate difference in RA + self.process.buildActions.astromDiffRA = ConvertUnits( + buildAction=SubtractVector, inUnit="degree", outUnit="milliarcsecond" + ) + self.process.buildActions.astromDiffRA.buildAction.actionA = RAcosDec( + raKey="coord_ra", decKey="coord_dec" + ) + self.process.buildActions.astromDiffRA.buildAction.actionB = RAcosDec( + raKey="ref_ra", decKey="ref_dec" + ) + + # Calculate difference in Dec + self.process.buildActions.astromDiffDec = ConvertUnits( + buildAction=SubtractVector, inUnit="degree", outUnit="milliarcsecond" + ) + self.process.buildActions.astromDiffDec.buildAction.actionA = LoadVector(vectorKey="coord_dec") + self.process.buildActions.astromDiffDec.buildAction.actionB = LoadVector(vectorKey="ref_dec") + + # Calculate difference in photometry + self.process.buildActions.photDiff = MagDiff(col1="{band}_psfFlux", col2="ref_{band}_flux") + + # Filter down to only objects with mag 17-21.5 + self.process.buildActions.mags = ConvertFluxToMag() + self.process.filterActions.brightStarsRA = DownselectVector(vectorKey="astromDiffRA") + self.process.filterActions.brightStarsRA.selector = RangeSelector( + vectorKey="mags", minimum=17, maximum=21.5 + ) + self.process.filterActions.brightStarsDec = DownselectVector(vectorKey="astromDiffDec") + self.process.filterActions.brightStarsDec.selector = RangeSelector( + vectorKey="mags", minimum=17, maximum=21.5 + ) + self.process.filterActions.brightStarsPhot = DownselectVector(vectorKey="photDiff") + self.process.filterActions.brightStarsPhot.selector = RangeSelector( + vectorKey="mags", minimum=17, maximum=21.5 + ) + + # Calculate median and sigmaMad + self.process.calculateActions.injected_RA_diff_median = MedianAction(vectorKey="brightStarsRA") + self.process.calculateActions.injected_RA_diff_sigmaMad = SigmaMadAction(vectorKey="brightStarsRA") + + self.process.calculateActions.injected_Dec_diff_median = MedianAction(vectorKey="brightStarsDec") + self.process.calculateActions.injected_Dec_diff_sigmaMad = SigmaMadAction(vectorKey="brightStarsDec") + + self.process.calculateActions.injected_phot_diff_median = MedianAction(vectorKey="brightStarsPhot") + self.process.calculateActions.injected_phot_diff_sigmaMad = SigmaMadAction( + vectorKey="brightStarsPhot" + ) + + self.produce.metric.units = { + "injected_RA_diff_median": "mas", + "injected_RA_diff_sigmaMad": "mas", + "injected_Dec_diff_median": "mas", + "injected_Dec_diff_sigmaMad": "mas", + "injected_phot_diff_median": "mmag", + "injected_phot_diff_sigmaMad": "mmag", + } diff --git a/python/lsst/analysis/tools/math.py b/python/lsst/analysis/tools/math.py index f7aadedc0..c6bf748d4 100644 --- a/python/lsst/analysis/tools/math.py +++ b/python/lsst/analysis/tools/math.py @@ -22,6 +22,7 @@ __all__ = ( "cos", "divide", + "isPercent", "fluxToMag", "nanMax", "nanMean", @@ -110,6 +111,12 @@ def fluxToMag( return mag +def isPercent(value: Scalar) -> bool: + """Return true if the value is between 0-100""" + result = 0.0 <= value <= 100.0 + return result + + def log(values: Scalar | Vector) -> Scalar | Vector: """Return the natural logarithm of values.""" with warnings.catch_warnings(): diff --git a/python/lsst/analysis/tools/tasks/__init__.py b/python/lsst/analysis/tools/tasks/__init__.py index 5da8e2395..c87d11d02 100644 --- a/python/lsst/analysis/tools/tasks/__init__.py +++ b/python/lsst/analysis/tools/tasks/__init__.py @@ -15,6 +15,7 @@ from .gatherResourceUsage import * from .makeMetricTable import * from .metricAnalysis import * +from .injectedObjectAnalysis import * from .objectTableSurveyAnalysis import * from .objectTableTractAnalysis import * from .photometricCatalogMatch import * diff --git a/python/lsst/analysis/tools/tasks/injectedObjectAnalysis.py b/python/lsst/analysis/tools/tasks/injectedObjectAnalysis.py new file mode 100644 index 000000000..e0c299ccc --- /dev/null +++ b/python/lsst/analysis/tools/tasks/injectedObjectAnalysis.py @@ -0,0 +1,64 @@ +# This file is part of analysis_tools. +# +# Developed for the LSST Data Management System. +# This product includes software developed by the LSST Project +# (https://www.lsst.org). +# See the COPYRIGHT file at the top-level directory of this distribution +# for details of code ownership. +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +from __future__ import annotations + +__all__ = ("InjectedObjectAnalysisConfig", "InjectedObjectAnalysisTask") + +from lsst.pipe.base import connectionTypes as ct +from lsst.skymap import BaseSkyMap + +from ..interfaces import AnalysisBaseConfig, AnalysisBaseConnections, AnalysisPipelineTask + + +class InjectedObjectAnalysisConnections( + AnalysisBaseConnections, + dimensions=("skymap", "tract"), + defaultTemplates={ + "outputName": "matched_injected_deepCoadd_catalog_tract_injected_objectTable_tract", + }, +): + data = ct.Input( + doc="Tract based object table to load from the butler", + name="{outputName}", + storageClass="ArrowAstropy", + deferLoad=True, + dimensions=("skymap", "tract"), + ) + + skymap = ct.Input( + doc="The skymap that covers the tract that the data is from.", + name=BaseSkyMap.SKYMAP_DATASET_TYPE_NAME, + storageClass="SkyMap", + dimensions=("skymap",), + ) + + +class InjectedObjectAnalysisConfig(AnalysisBaseConfig, pipelineConnections=InjectedObjectAnalysisConnections): + pass + + +class InjectedObjectAnalysisTask(AnalysisPipelineTask): + """Make plots and metrics using a table of objects matched to reference + catalog sources. + """ + + ConfigClass = InjectedObjectAnalysisConfig + _DefaultName = "injectedObjectAnalysisTask"