Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

PIPE2D-1592: fitPfsFluxReference: Don't use low-SNR flux standards #369

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
61 changes: 61 additions & 0 deletions python/pfs/drp/stella/fitPfsFluxReference.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,7 @@
import copy
import dataclasses
import math
import warnings

from typing import Literal, overload
from collections.abc import Generator, Mapping, Sequence
Expand Down Expand Up @@ -146,6 +147,22 @@ class FitPfsFluxReferenceConfig(PipelineTaskConfig, pipelineConnections=FitPfsFl
doc="Right ends of wavelength ranges ignored (because e.g. of strong atmospheric absorption)"
" when comparing middle-resolution observation spectra to models.",
)
cutoffSNR = Field(
dtype=float,
default=10,
doc="Minimally required pixel-wise S/N of observed spectra."
" Spectra with S/N worse than this value will be discarded.",
)
cutoffSNRRangeLeft = Field(
dtype=float,
default=840,
doc="Left edge of the wavelength range in which S/N is averaged to be compared with ``cutoffSNR``",
)
cutoffSNRRangeRight = Field(
dtype=float,
default=880,
doc="Right edge of the wavelength range in which S/N is averaged to be compared with ``cutoffSNR``",
)
badMask = ListField(
dtype=str, default=["BAD", "SAT", "CR", "NO_DATA", "SUSPECT"], doc="Mask planes for bad pixels"
)
Expand Down Expand Up @@ -402,6 +419,14 @@ def selectPfsConfig(pfsConfig: PfsConfig, flagName: str, isGood: Sequence[bool])
return goodConfig

pfsConfig = selectPfsConfig(pfsConfig, "ABSENT_FIBER", np.isin(pfsConfig.fiberId, pfsMerged.fiberId))

pfsMerged = removeBadSpectra(
pfsMerged,
self.config.cutoffSNR,
(self.config.cutoffSNRRangeLeft, self.config.cutoffSNRRangeRight),
)
pfsConfig = selectPfsConfig(pfsConfig, "LOW_SNR_FIBER", np.isin(pfsConfig.fiberId, pfsMerged.fiberId))

pfsConfig = selectPfsConfig(pfsConfig, "BAD_FIBER", (pfsConfig.fiberStatus == FiberStatus.GOOD))
pfsConfig = selectPfsConfig(
pfsConfig,
Expand Down Expand Up @@ -1996,6 +2021,42 @@ def removeBadFluxes(
pfsConfig.totalFluxErr = totalFluxErrArrays


def removeBadSpectra(
pfsMerged: PfsMerged,
cutoffSNR: float,
cutoffSNRRange: tuple[float, float],
) -> PfsMerged:
"""Remove bad spectra from ``pfsMerged``.

Spectra are bad if their average pixel-wise S/N is worse than ``cutoffSNR``.

Parameters
----------
pfsMerged : `PfsMerged`
Merged spectra from exposure.
cutoffSNR : `float`
Cut-off S/N.
cutoffSNRRange : `tuple` [`float`, `float`]
Wavelength range (nm) in which pixel-wise S/N is averaged to be compared
with ``cutoffSNR``.

Returns
-------
pfsMerged : `PfsMerged`
PfsMerged from which bad spectra have been removed.
"""
left, right = cutoffSNRRange
good = np.zeros(shape=(len(pfsMerged),), dtype=bool)
with warnings.catch_warnings():
warnings.simplefilter("ignore")
for i in range(len(pfsMerged)):
sampleIndex = (left <= pfsMerged.wavelength[i]) & (pfsMerged.wavelength[i] <= right)
snr = pfsMerged.flux[i, sampleIndex] / np.sqrt(pfsMerged.covar[i, 0, sampleIndex])
good[i] = np.nanmedian(snr) >= cutoffSNR

return pfsMerged[good]


def promoteSimpleSpectrumToFiberArray(spectrum: PfsSimpleSpectrum, snr: float) -> PfsFiberArray:
"""Promote an instance of PfsSimpleSpectrum to PfsFiberArray.

Expand Down