From 2771ecab8c52dca9749a06759f49d7aae9d53332 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sun, 22 Jun 2025 16:32:18 -0400 Subject: [PATCH 01/24] Implemented partial loading for spectra stored as `TargetSpectra`. --- bin/testPfsTargetSpectraPartialIO.py | 76 +++++++++++++++++ python/pfs/datamodel/pfsTargetSpectra.py | 102 +++++++++++++++++------ tests/test_PfsTargetSpectra.py | 11 +++ 3 files changed, 163 insertions(+), 26 deletions(-) create mode 100755 bin/testPfsTargetSpectraPartialIO.py diff --git a/bin/testPfsTargetSpectraPartialIO.py b/bin/testPfsTargetSpectraPartialIO.py new file mode 100755 index 00000000..4f7777cf --- /dev/null +++ b/bin/testPfsTargetSpectraPartialIO.py @@ -0,0 +1,76 @@ +#!/usr/bin/env python + +import os +from time import time +import resource +from argparse import ArgumentParser +from pfs.datamodel import PfsCalibrated + + +def get_bytes_read(pid): + """Return the number of bytes read so far for the given process ID.""" + + with open(f"/proc/{pid}/io", "r") as f: + for line in f: + if line.startswith("rchar:"): + return int(line.split()[1]) + + return 0 + + +def main(filename, **kwargs): + # Remove None values from kwargs + kwargs = {k: v for k, v in kwargs.items() if v is not None} + + # Get the PID of the current process + pid = os.getpid() + + usage = resource.getrusage(resource.RUSAGE_SELF) + bytes_start = get_bytes_read(pid) + time_start = time() + print('Bytes read so far:', bytes_start) + print('Max RSS:', usage.ru_maxrss, 'kB') + print('Minor page faults so far:', usage.ru_minflt) + print('Major page faults so far:', usage.ru_majflt) + + # Load a subset of the spectra + spectra = PfsCalibrated.readFits(filename, **kwargs) + print(f"Loaded {len(spectra)} spectra with kwargs: {kwargs}") + + usage = resource.getrusage(resource.RUSAGE_SELF) + bytes_partial = get_bytes_read(pid) + time_partial = time() + print('Bytes read during loading:', bytes_partial - bytes_start) + print('Time taken to load a subset of spectra:', time_partial - time_start) + print('Max RSS:', usage.ru_maxrss, 'kB') + print('Minor page faults so far:', usage.ru_minflt) + print('Major page faults so far:', usage.ru_majflt) + + # Load the entire file + spectra = PfsCalibrated.readFits(filename) + print(f"Loaded {len(spectra)} spectra.") + + usage = resource.getrusage(resource.RUSAGE_SELF) + bytes_all = get_bytes_read(pid) + time_all = time() + print('Bytes read during loading:', bytes_all - bytes_partial) + print('Time taken to load all spectra:', time_all - time_partial) + print('Max RSS:', usage.ru_maxrss, 'kB') + print('Minor page faults so far:', usage.ru_minflt) + print('Major page faults so far:', usage.ru_majflt) + + +if __name__ == "__main__": + parser = ArgumentParser(description=""" +Transitioning to gen3 Butler resulted in combining target spectra into large FITS files which +contain many PfsSingle or PfsObject spectra. This script can be used to benchmark the IO cost +saved by partially loading PfsCalibrated or PfsCoadd spectra, rather than loading the +entire file. +""") + parser.add_argument("filename", help="Name of the file to load.") + parser.add_argument("--targetId", type=int, help="Target IDs to load (default: all).") + parser.add_argument("--targetType", type=int, help="TargetType to load (default: all).") + parser.add_argument("--catId", type=int, help="Objects with catId to load (default: all).") + args = parser.parse_args() + main(args.filename, + targetId=args.targetId) diff --git a/python/pfs/datamodel/pfsTargetSpectra.py b/python/pfs/datamodel/pfsTargetSpectra.py index a62dcc3c..c746163c 100644 --- a/python/pfs/datamodel/pfsTargetSpectra.py +++ b/python/pfs/datamodel/pfsTargetSpectra.py @@ -119,13 +119,18 @@ def __contains__(self, target) -> bool: return target in self.spectra @classmethod - def readFits(cls, filename: str) -> "PfsTargetSpectra": - """Read from FITS file + def readFits(cls, filename: str, **kwargs) -> "PfsTargetSpectra": + """Read a collection of spectra from a FITS file, optionally + filtered by keyword arguments. Parameters ---------- filename : `str` Filename of FITS file. + **kwargs : `dict` + Keyword arguments to filter the spectra. The keys supported depend on the + subclass of `PfsTargetSpectra`. For example, the `PfsFiberArray` class supports + filtering by `targetType`, `catId`, and `objId`. Returns ------- @@ -133,7 +138,7 @@ def readFits(cls, filename: str) -> "PfsTargetSpectra": Constructed instance, from FITS file. """ - def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: + def readData(hduName: str, dtype: type, mask=None) -> Union[np.ndarray, List[np.ndarray]]: """Read data from FITS file Parameters @@ -142,6 +147,8 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: Name of the HDU. dtype : `type` Data type. + mask : `np.ndarray`, optional + Mask to apply to the data. If provided, the data will be read as a masked array. Returns ------- @@ -150,35 +157,78 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: """ hdu = fits[hduName] if isinstance(hdu, (ImageHDU, CompImageHDU)): - return hdu.data.astype(dtype) - numRows = len(hdu.data.dtype) - if "row_0" in hdu.data.dtype.names: - data: List[np.ndarray] = [] - for ii in range(len(hdu.data)): - data.append(np.array([hdu.data[f"row_{jj}"][ii] for jj in range(numRows)], dtype=dtype)) - return data - if numRows == 1: - return [row["value"].astype(dtype) for row in hdu.data] - raise RuntimeError(f"Cannot read data from HDU {hduName}") + if mask is None: + return hdu.data.astype(dtype) + else: + if hdu.data.ndim == 1: + # Same for all spectra, such as wavelength, do not filter + return hdu.data.astype(dtype) + else: + # Different for each spectrum, such as flux, do filter + return hdu.data[mask].astype(dtype) + elif isinstance(hdu, BinTableHDU): + # This is a special case of storing arrays in a table. + # Rows of the arrays are stored as columns of a table. + numRows = len(hdu.data.dtype) + if "row_0" in hdu.data.dtype.names: + data: List[np.ndarray] = [] + for ii in range(len(hdu.data)): + data.append(np.array([hdu.data[f"row_{jj}"][ii] for jj in range(numRows)], dtype=dtype)) + return data + if numRows == 1: + # Single row, such as wavelength fixed wavelength + return [row["value"].astype(dtype) for row in hdu.data] + else: + # One row per spectrum, such as flux + data: List[np.ndarray] = [] + if mask is not None: + index = np.where(mask)[0] + else: + index = range(len(hdu.data)) + + for ii in index: + data.append(np.array( + [hdu.data[f"row_{jj}"][ii] for jj in range(numRows)], + dtype=dtype)) + + return data + # raise RuntimeError(f"Cannot read data from HDU {hduName}") + else: + raise TypeError(f"HDU '{hduName}' is not an ImageHDU, CompImageHDU or BinTableHDU.") spectra = [] with astropy.io.fits.open(filename) as fits: - header = astropyHeaderToDict(fits[0].header) + # Read the main headers and the list of targets. + header = fits[0].header targetHdu = fits["TARGET"].data + + # If any keyword arguments are provided, filter the spectra + mask = None + for key, value in kwargs.items(): + if key not in targetHdu.columns.names: + raise KeyError(f"Keyword argument '{key}' not found in TARGET HDU.") + m = targetHdu[key] == value + mask = m if mask is None else mask & m + targetFluxHdu = fits["TARGETFLUX"].data observationsHdu = fits["OBSERVATIONS"].data - wavelengthData = readData("WAVELENGTH", np.float64) - fluxData = readData("FLUX", np.float32) - maskData = readData("MASK", np.int32) - skyData = readData("SKY", np.float32) - covarData = readData("COVAR", np.float32) - covar2Hdu = fits["COVAR2"].data if "COVAR2" in fits else None + wavelengthData = readData("WAVELENGTH", np.float64, mask=mask) + fluxData = readData("FLUX", np.float32, mask=mask) + maskData = readData("MASK", np.int32, mask=mask) + skyData = readData("SKY", np.float32, mask=mask) + covarData = readData("COVAR", np.float32, mask=mask) + covar2Data = readData("COVAR2", np.float32, mask=mask) if "COVAR2" in fits else None metadataHdu = fits["METADATA"].data fluxTableHdu = fits["FLUXTABLE"].data notesTable = cls.NotesClass.readHdu(fits) - for ii, row in enumerate(targetHdu): - targetId = row["targetId"] + if mask is None: + targets = enumerate(targetHdu.targetId) + else: + targets = zip(np.where(mask)[0], targetHdu.targetId[mask]) + + for ii, (index, targetId) in enumerate(targets): + row = targetHdu[index] select = targetFluxHdu.targetId == targetId fiberFlux = dict( zip( @@ -216,7 +266,7 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: expTime, ) - metadataRow = metadataHdu[ii] + metadataRow = metadataHdu[index] assert metadataRow["targetId"] == targetId metadata = yaml.load( @@ -226,7 +276,7 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: ) flags = MaskHelper.fromFitsHeader(metadata, strip=True) - fluxTableRow = fluxTableHdu[ii] + fluxTableRow = fluxTableHdu[index] assert fluxTableRow["targetId"] == targetId fluxTable = FluxTable( fluxTableRow["wavelength"], @@ -237,7 +287,7 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: ) notes = cls.PfsFiberArrayClass.NotesClass( - **{col.name: getattr(notesTable, col.name)[ii] for col in notesTable.schema} + **{col.name: getattr(notesTable, col.name)[index] for col in notesTable.schema} ) # Wavelength: we might write a single array if everything has the same length and value @@ -255,7 +305,7 @@ def readData(hduName: str, dtype: type) -> Union[np.ndarray, List[np.ndarray]]: maskData[ii], skyData[ii], covarData[ii], - covar2Hdu[ii] if covar2Hdu is not None else [], + covar2Data[ii] if covar2Data is not None else [], flags, metadata, fluxTable, diff --git a/tests/test_PfsTargetSpectra.py b/tests/test_PfsTargetSpectra.py index 710c65ed..458ed1fa 100644 --- a/tests/test_PfsTargetSpectra.py +++ b/tests/test_PfsTargetSpectra.py @@ -349,6 +349,7 @@ def testLookup(self): def testIO(self): """Test I/O functionality""" + spectra = self.makePfsTargetSpectra(list(range(123))) with lsst.utils.tests.getTempFilePath(".fits") as filename: spectra.writeFits(filename) @@ -359,6 +360,16 @@ def testIO(self): with astropy.io.fits.open(filename) as hdus: self.assertIn("DAMD_VER", hdus[0].header) + def testIOWithFilters(self): + """Test I/O functionality with filters""" + + spectra = self.makePfsTargetSpectra(list(range(123))) + with lsst.utils.tests.getTempFilePath(".fits") as filename: + spectra.writeFits(filename) + + copy = PfsCoadd.readFits(filename, objId=22) + self.assertPfsTargetSpectraEqual(copy, PfsCoadd([spectra[22]])) + def testDifferentLengths(self): """Test I/O when the spectra have different lengths""" spectra = self.makePfsTargetSpectra(list(range(123)), length=itertools.cycle([999, 1000, 1001])) From d71beab55a63e1a4fff61a24f08054859c4b0f91 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Mon, 23 Jun 2025 08:53:27 -0400 Subject: [PATCH 02/24] Fixed partial loading of spectra with variable array length. --- python/pfs/datamodel/pfsTargetSpectra.py | 13 ++++++---- tests/test_PfsTargetSpectra.py | 33 +++++++++++++++++++++++- 2 files changed, 40 insertions(+), 6 deletions(-) diff --git a/python/pfs/datamodel/pfsTargetSpectra.py b/python/pfs/datamodel/pfsTargetSpectra.py index c746163c..00bcdc23 100644 --- a/python/pfs/datamodel/pfsTargetSpectra.py +++ b/python/pfs/datamodel/pfsTargetSpectra.py @@ -167,8 +167,8 @@ def readData(hduName: str, dtype: type, mask=None) -> Union[np.ndarray, List[np. # Different for each spectrum, such as flux, do filter return hdu.data[mask].astype(dtype) elif isinstance(hdu, BinTableHDU): - # This is a special case of storing arrays in a table. - # Rows of the arrays are stored as columns of a table. + # This is a special case of storing variable length arrays in a table. + # Rows of 2D arrays are stored as columns of a table (e.g. covariance) numRows = len(hdu.data.dtype) if "row_0" in hdu.data.dtype.names: data: List[np.ndarray] = [] @@ -176,10 +176,13 @@ def readData(hduName: str, dtype: type, mask=None) -> Union[np.ndarray, List[np. data.append(np.array([hdu.data[f"row_{jj}"][ii] for jj in range(numRows)], dtype=dtype)) return data if numRows == 1: - # Single row, such as wavelength fixed wavelength - return [row["value"].astype(dtype) for row in hdu.data] + # 1D array stored as a single column + if mask is None: + return [row["value"].astype(dtype) for row in hdu.data] + else: + return [row.astype(dtype) for row in hdu.data["value"][mask]] else: - # One row per spectrum, such as flux + # 2D array stored as multiple columns data: List[np.ndarray] = [] if mask is not None: index = np.where(mask)[0] diff --git a/tests/test_PfsTargetSpectra.py b/tests/test_PfsTargetSpectra.py index 458ed1fa..9933d6bc 100644 --- a/tests/test_PfsTargetSpectra.py +++ b/tests/test_PfsTargetSpectra.py @@ -363,13 +363,20 @@ def testIO(self): def testIOWithFilters(self): """Test I/O functionality with filters""" - spectra = self.makePfsTargetSpectra(list(range(123))) + spectra = self.makePfsTargetSpectra(list(range(123)), catId=itertools.cycle([12345, 12346, 12347])) with lsst.utils.tests.getTempFilePath(".fits") as filename: spectra.writeFits(filename) + # Filter single spectrum by objId copy = PfsCoadd.readFits(filename, objId=22) self.assertPfsTargetSpectraEqual(copy, PfsCoadd([spectra[22]])) + # Filter multiple spectra by catId + copy = PfsCoadd.readFits(filename, catId=12345) + for (_, left), right in zip(copy.items(), [s for t, s in spectra.items() if t.catId == 12345]): + self.assertTargetEqual(left.target, right.target) + self.assertPfsFiberArrayEqual(left, right) + def testDifferentLengths(self): """Test I/O when the spectra have different lengths""" spectra = self.makePfsTargetSpectra(list(range(123)), length=itertools.cycle([999, 1000, 1001])) @@ -378,6 +385,30 @@ def testDifferentLengths(self): copy = PfsCoadd.readFits(filename) self.assertPfsTargetSpectraEqual(copy, spectra) + def testDifferentLengthsWithFilters(self): + """Test I/O when the spectra have different lengths and filters are used""" + spectra = self.makePfsTargetSpectra( + list(range(123)), + catId=itertools.cycle([12345, 12346, 12347]), + length=itertools.cycle([999, 1000, 1001, 1002, 1003]) + ) + with lsst.utils.tests.getTempFilePath(".fits") as filename: + spectra.writeFits(filename) + + # Filter single spectrum by objId + for ii in range(22, 27): + copy = PfsCoadd.readFits(filename, objId=ii) + self.assertPfsTargetSpectraEqual(copy, PfsCoadd([spectra[ii]])) + + # Filter multiple spectra by catId + copy = PfsCoadd.readFits(filename, catId=12345) + for (_, left), right in zip( + copy.items(), + [s for t, s in spectra.items() if t.catId == 12345]): + + self.assertTargetEqual(left.target, right.target) + self.assertPfsFiberArrayEqual(left, right) + class TestMemory(lsst.utils.tests.MemoryTestCase): pass From ffe5bb8c374325ff99c632b3a9de0e1b13e1a118 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 24 May 2024 15:37:44 -0400 Subject: [PATCH 03/24] Minor updates to work with gapipe --- python/pfs/datamodel/pfsFiberArray.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/python/pfs/datamodel/pfsFiberArray.py b/python/pfs/datamodel/pfsFiberArray.py index e8929f02..256d2fd7 100644 --- a/python/pfs/datamodel/pfsFiberArray.py +++ b/python/pfs/datamodel/pfsFiberArray.py @@ -44,6 +44,7 @@ class PfsFiberArray(PfsSimpleSpectrum): """ filenameFormat = None # Subclasses should override NotesClass: Type[Notes] # Subclasses should override + FluxTableClass = FluxTable # Subclasses may override def __init__( self, @@ -115,7 +116,7 @@ def _readImpl(cls, fits): data["covar"] = fits["COVAR"].data.astype(np.float32) data["covar2"] = fits["COVAR2"].data.astype(np.float32) try: - fluxTable = FluxTable.fromFits(fits) + fluxTable = cls.FluxTableClass.fromFits(fits) except KeyError as exc: # Only want to catch "Extension XXX not found." if not exc.args[0].startswith("Extension"): @@ -152,3 +153,4 @@ def _writeImpl(self, fits): self.notes.writeFits(fits) if self.fluxTable is not None: self.fluxTable.toFits(fits) + return header \ No newline at end of file From 90c14810b9dc179f3ca2a85b9f9ba02a1c67b570 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 31 May 2024 18:42:15 -0400 Subject: [PATCH 04/24] Added first version of GA datamodel with unit tests. --- python/pfs/datamodel/__init__.py | 1 + python/pfs/datamodel/fluxTable.py | 21 +- python/pfs/datamodel/ga.py | 308 ++++++++++++++++++++++++++ python/pfs/datamodel/pfsFiberArray.py | 2 +- tests/test_pfsGAObject.py | 183 +++++++++++++++ 5 files changed, 509 insertions(+), 6 deletions(-) create mode 100644 python/pfs/datamodel/ga.py create mode 100644 tests/test_pfsGAObject.py diff --git a/python/pfs/datamodel/__init__.py b/python/pfs/datamodel/__init__.py index 04c97a92..da81f310 100644 --- a/python/pfs/datamodel/__init__.py +++ b/python/pfs/datamodel/__init__.py @@ -22,3 +22,4 @@ from .pfsFiberNorms import * from .pfsZCandidates import * from .pfsCoZCandidates import * +from .ga import * diff --git a/python/pfs/datamodel/fluxTable.py b/python/pfs/datamodel/fluxTable.py index c60c1cd4..48fde0d6 100644 --- a/python/pfs/datamodel/fluxTable.py +++ b/python/pfs/datamodel/fluxTable.py @@ -38,11 +38,11 @@ class FluxTable: _hduName = "FLUX_TABLE" # HDU name to use def __init__(self, wavelength, flux, error, mask, flags): - dims = np.array([len(wavelength.shape), len(flux.shape), len(error.shape), len(mask.shape)]) - lengths = set([wavelength.shape, flux.shape, error.shape, mask.shape]) - if np.any(dims != 1) or len(lengths) > 1: - raise RuntimeError("Bad shapes for wavelength,flux,error,mask: %s,%s,%s,%s" % - (wavelength.shape, flux.shape, error.shape, mask.shape)) + self.checkShapes(wavelength=wavelength, + flux=flux, + error=error, + mask=mask) + self.wavelength = wavelength self.flux = flux self.error = error @@ -52,6 +52,17 @@ def __init__(self, wavelength, flux, error, mask, flags): def __len__(self): """Return number of elements""" return len(self.wavelength) + + def checkShapes(self, **kwargs): + keys = list(sorted(kwargs.keys())) + dims = np.array([ len(kwargs[k].shape) for k in keys ]) + lengths = set([ kwargs[k].shape for k in keys ]) + + if np.any(dims != 1) or len(lengths) > 1: + names = ','.join(keys) + shapes = ','.join([str(kwargs[k].shape) for k in keys]) + raise RuntimeError("Bad shapes for %s: %s" % + (names, shapes)) def toFits(self, fits): """Write to a FITS file diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py new file mode 100644 index 00000000..c65d94f1 --- /dev/null +++ b/python/pfs/datamodel/ga.py @@ -0,0 +1,308 @@ +from typing import Type +import numpy as np + +from .notes import makeNotesClass, Notes +from .pfsFiberArray import PfsFiberArray +from .fluxTable import FluxTable +from .pfsTable import PfsTable, Column +from .utils import inheritDocstrings +from .utils import astropyHeaderToDict, astropyHeaderFromDict +from .masks import MaskHelper + +__all__ = [ + "VelocityCorrections", + "StellarParams", + "Abundances", + "GAFluxTable", + "PfsGAObjectNotes", + "PfsGAObject", +] + +class VelocityCorrections(PfsTable): + """A table of velocity corrections applied to the individual visits.""" + + damdVer = 2 + schema = [ + Column("visit", np.int32, "ID of the visit these corrections apply for", -1), + Column("JD", np.float32, "Julian date of the visit", -1), + Column("helio", np.float32, "Heliocentric correction", np.nan), + Column("bary", np.float32, "Barycentric correction", np.nan), + ] + fitsExtName = 'VELCORR' + +class StellarParams(PfsTable): + """List of measured stellar parameters for a target.""" + + damdVer = 2 + schema = [ + Column("method", str, "Line-of-sight velocity measurement method", ""), + Column("frame", str, "Reference frame of velocity: helio, bary", ""), + Column("param", str, "Stellar parameter: v_los, M_H, T_eff, log_g, a_M", ""), + Column("covarId", np.uint8, "Param position within covariance matrix", -1), + Column("unit", str, "Physical unit of parameter", ""), + Column("value", np.float32, "Stellar parameter value", np.nan), + Column("valueErr", np.float32, "Stellar parameter error", np.nan), + # TODO: add quantiles or similar for MCMC results + Column("flag", bool, "Measurement flag (true means bad)", False), + Column("status", str, "Measurement flags", ""), + ] + fitsExtName = 'STELLARPARAM' + +class Abundances(PfsTable): + """List of measured abundance parameters for stellar targets.""" + + damdVer = 2 + schema = [ + Column("method", str, "Abundance measurement method", ""), + Column("element", str, "Chemical element the abundance is measured for", ""), + Column("covarId", np.uint8, "Param position within covariance matrix", -1), + Column("value", np.float32, "Abundance value", np.nan), + Column("valueErr", np.float32, "Abundance error", np.nan), + ] + fitsExtName = 'ABUND' + +class GAFluxTable(FluxTable): + """Table of coadded fluxes at near-original sampling and model fits + + Merged and coadded spectra have been resampled to a standard wavelength + sampling. This representation provides coadded fluxes at approximately the + native wavelength sampling, for those that want the data with a minimum of + resampling. This is mostly of use for single exposures and coadds made from + back-to-back exposures with the same top-end configuration. For coadds made + from exposures with different top-end configurations, the different + wavelength samplings obtained from the different fibers means there's no + single native wavelength sampling, and so this is less useful. + + This is like a `pfs.datamodel.PfsSimpleSpectrum`, except that it includes a + variance array, and is written to a FITS HDU rather than a file (so it can + be incorporated within a `pfs.datamodel.PfsSpectrum`). + + Parameters + ---------- + wavelength : `numpy.ndarray` of `float` + Array of wavelengths. + flux : `numpy.ndarray` of `float` + Array of fluxes. + error : `numpy.ndarray` of `float` + Array of flux errors. + model : `numpy.ndarray` of `float` + Array of best-fit model flux. + cont : `numpy.ndarray` of `float` + Array of continuum model. + norm_flux : `numpy.ndarray` of `float` + Array of continuum-normalized flux. + norm_error : `numpy.ndarray` of `float` + Array of continuum-normalized flux error. + norm_model : `numpy.ndarray` of `float` + Array of continuum-normalized model. + mask : `numpy.ndarray` of `int` + Array of mask pixels. + flags : `pfs.datamodel.MaskHelper` + Helper for dealing with symbolic names for mask values. + """ + _hduName = "FLUX_TABLE" # HDU name to use + + def __init__(self, wavelength, flux, error, model, cont, norm_flux, norm_error, norm_model, mask, flags): + self.checkShapes(wavelength=wavelength, + flux=flux, + error=error, + model=model, + cont=cont, + norm_flux=norm_flux, + norm_error=norm_error, + norm_model=norm_model, + mask=mask) + + self.wavelength = wavelength + self.flux = flux + self.error = error + self.model = model + self.cont = cont + self.norm_flux = norm_flux + self.norm_error = norm_error + self.norm_model = norm_model + self.mask = mask + self.flags = flags + + def toFits(self, fits): + """Write to a FITS file + + Parameters + ---------- + fits : `astropy.io.fits.HDUList` + Opened FITS file. + """ + # NOTE: When making any changes to this method that modify the output + # format, increment the DAMD_VER header value and record the change in + # the versions.txt file. + from astropy.io.fits import BinTableHDU, Column + header = astropyHeaderFromDict(self.flags.toFitsHeader()) + header['DAMD_VER'] = (1, "GAFluxTable datamodel version") + hdu = BinTableHDU.from_columns([ + Column("wavelength", "D", array=self.wavelength), + Column("flux", "E", array=self.flux), + Column("error", "E", array=self.error), + Column("model", "E", array=self.model), + Column("cont", "E", array=self.cont), + Column("norm_flux", "E", array=self.norm_flux), + Column("norm_error", "E", array=self.norm_error), + Column("norm_model", "E", array=self.norm_model), + Column("mask", "K", array=self.mask), + ], header=header, name=self._hduName) + fits.append(hdu) + + @classmethod + def fromFits(cls, fits): + """Construct from a FITS file + + Parameters + ---------- + fits : `astropy.io.fits.HDUList` + Opened FITS file. + + Returns + ------- + self : `FluxTable` + Constructed `FluxTable`. + """ + hdu = fits[cls._hduName] + header = astropyHeaderToDict(hdu.header) + flags = MaskHelper.fromFitsHeader(header) + return cls(hdu.data["wavelength"].astype(float), + hdu.data["flux"].astype(np.float32), + hdu.data["error"].astype(np.float32), + hdu.data["model"].astype(np.float32), + hdu.data["cont"].astype(np.float32), + hdu.data["norm_flux"].astype(np.float32), + hdu.data["norm_error"].astype(np.float32), + hdu.data["norm_model"].astype(np.float32), + hdu.data["mask"].astype(np.int32), + flags) + +PfsGAObjectNotes = makeNotesClass( + "PfsGAObjectNotes", + [] +) + +@inheritDocstrings +class PfsGAObject(PfsFiberArray): + """Coadded spectrum of a GA target with derived quantities. + + Produced by ˙˙gapipe`` + + Parameters + ---------- + target : `pfs.datamodel.Target` + Target information. + observations : `pfs.datamodel.Observations` + Observations of the target. + wavelength : `numpy.ndarray` of `float` + Array of wavelengths. + flux : `numpy.ndarray` of `float` + Array of fluxes. + mask : `numpy.ndarray` of `int` + Array of mask pixels. + sky : `numpy.ndarray` of `float` + Array of sky values. + covar : `numpy.ndarray` of `float` + Near-diagonal (diagonal and either side) part of the covariance matrix. + covar2 : `numpy.ndarray` of `float` + Low-resolution non-sparse covariance estimate. + flags : `MaskHelper` + Helper for dealing with symbolic names for mask values. + metadata : `dict` (`str`: POD), optional + Keyword-value pairs for the header. + fluxTable : `pfs.datamodel.GAFluxTable`, optional + Table of coadded fluxes and continuum-normalized flux from contributing observations. + stellarParams: `pfs.datamodel.StellarParams`, optional + Table of measured stellar parameters. + velocityCorrections: `pfs.datamodel.VelocityCorrections`, optional + Table of velocity corrections applied to the individual visits. + abundances: `pfs.datamodel.Abundances`, optional + Table of measured abundance parameters. + paramsCovar: `numpy.ndarray` of `float`, optional + Covariance matrix for stellar parameters. + abundCovar: `numpy.ndarray` of `float`, optional + Covariance matrix for abundance parameters. + notes : `Notes`, optional + Reduction notes. + """ + + filenameFormat = ("pfsGAObject-%(catId)05d-%(tract)05d-%(patch)s-%(objId)016x" + "-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") + filenameRegex = r"^pfsGAObject-(\d{5})-(\d{5})-(.*)-([0-9a-f]{16})-(\d{3})-0x([0-9a-f]{16})\.fits.*$" + filenameKeys = [("catId", int), ("tract", int), ("patch", str), ("objId", int), + ("nVisit", int), ("pfsVisitHash", int)] + NotesClass = PfsGAObjectNotes + FluxTableClass = GAFluxTable + + StellarParamsFitsExtName = "STELLARCOVAR" + AbundancesFitsExtName = "ABUNDCOVAR" + + def __init__( + self, + target, + observations, + wavelength, + flux, + mask, + sky, + covar, + covar2, + flags, + metadata=None, + fluxTable=None, + stellarParams=None, + velocityCorrections=None, + abundances=None, + paramsCovar=None, + abundCovar=None, + notes: Notes = None, + ): + super().__init__(target, observations, wavelength, flux, mask, sky, covar, covar2, flags, metadata=metadata, fluxTable=fluxTable, notes=notes) + + self.stellarParams = stellarParams + self.velocityCorrections = velocityCorrections + self.abundances = abundances + self.paramsCovar = paramsCovar + self.abundCovar = abundCovar + + def validate(self): + """Validate that all the arrays are of the expected shape""" + super().validate() + + # TODO: write any validation code + + @classmethod + def _readImpl(cls, fits): + data = super()._readImpl(fits) + + # TODO: handle missing extensions + + data["stellarParams"] = StellarParams.readHdu(fits) + data["velocityCorrections"] = VelocityCorrections.readHdu(fits) + data["abundances"] = Abundances.readHdu(fits) + if cls.StellarParamsFitsExtName in fits: + data["paramsCovar"] = fits[cls.StellarParamsFitsExtName].data.astype(np.float32) + if cls.AbundancesFitsExtName in fits: + data["abundCovar"] = fits[cls.AbundancesFitsExtName].data.astype(np.float32) + + return data + + def _writeImpl(self, fits): + from astropy.io.fits import ImageHDU + + header = super()._writeImpl(fits) + + if self.stellarParams is not None: + self.stellarParams.writeHdu(fits) + if self.velocityCorrections is not None: + self.velocityCorrections.writeHdu(fits) + if self.abundances is not None: + self.abundances.writeHdu(fits) + if self.paramsCovar is not None: + fits.append(ImageHDU(self.paramsCovar.astype(np.float32), header=header, name=self.StellarParamsFitsExtName)) + if self.abundCovar is not None: + fits.append(ImageHDU(self.abundCovar.astype(np.float32), header=header, name=self.AbundancesFitsExtName)) + + return header \ No newline at end of file diff --git a/python/pfs/datamodel/pfsFiberArray.py b/python/pfs/datamodel/pfsFiberArray.py index 256d2fd7..7c2cba90 100644 --- a/python/pfs/datamodel/pfsFiberArray.py +++ b/python/pfs/datamodel/pfsFiberArray.py @@ -44,7 +44,7 @@ class PfsFiberArray(PfsSimpleSpectrum): """ filenameFormat = None # Subclasses should override NotesClass: Type[Notes] # Subclasses should override - FluxTableClass = FluxTable # Subclasses may override + FluxTableClass: Type[FluxTable] = FluxTable # Subclasses may override def __init__( self, diff --git a/tests/test_pfsGAObject.py b/tests/test_pfsGAObject.py new file mode 100644 index 00000000..f720dc8b --- /dev/null +++ b/tests/test_pfsGAObject.py @@ -0,0 +1,183 @@ +import os +import re +import numpy as np +from unittest import TestCase + +from pfs.datamodel import Target, TargetType +from pfs.datamodel import Observations +from pfs.datamodel import MaskHelper +from pfs.datamodel import GAFluxTable +from pfs.datamodel import PfsGAObject, StellarParams, Abundances, VelocityCorrections + +class PfsGAObjectTestCase(TestCase): + """ Check the format of example datamodel files are + consistent with that specified in the corresponding + datamodel classes. + """ + + def makePfsGAObject(self): + """Construct a PfsGAObject with dummy values for testing.""" + + catId = 12345 + tract = 1 + patch = '1,1' + objId = 123456789 + ra = -100.63654 + dec = -68.591576 + targetType = TargetType.SCIENCE + + target = Target(catId, tract, patch, objId, ra, dec, targetType) + + visit = np.array([ 83219, 83220 ]) + arm = np.array([ 'b', 'm', ]) + spectrograph = np.array([1, 1]) + pfsDesignId = np.array([8854764194165386399, 8854764194165386400]) + fiberId = np.array([476, 476]) + pfiNominal = np.array([[ra, dec], [ra, dec]]) + pfiCenter = np.array([[ra, dec], [ra, dec]]) + + observations = Observations(visit, arm, spectrograph, pfsDesignId, fiberId, pfiNominal, pfiCenter) + + npix = 4096 + wavelength = np.concatenate([ + np.linspace(380, 650, npix, dtype=np.float32), + np.linspace(710, 885, npix, dtype=np.float32) + ]) + flux = np.zeros_like(wavelength) + error = np.zeros_like(wavelength) + model = np.zeros_like(wavelength) + cont = np.zeros_like(wavelength) + norm_flux = np.zeros_like(wavelength) + norm_error = np.zeros_like(wavelength) + norm_model = np.zeros_like(wavelength) + mask = np.zeros_like(wavelength, dtype=np.int32) + sky = np.zeros_like(wavelength) + covar = np.zeros((3, wavelength.size), dtype=np.float32) # Tridiagonal covariance matrix of flux + covar2 = np.zeros((1, 1), dtype=np.float32) # ? + + flags = MaskHelper() # {'BAD': 0, 'BAD_FIBERTRACE': 11, 'BAD_FLAT': 9, 'BAD_FLUXCAL': 13, 'BAD_SKY': 12, 'CR': 3, 'DETECTED': 5, 'DETECTED_NEGATIVE': 6, 'EDGE': 4, 'FIBERTRACE': 10, 'INTRP': 2, 'IPC': 14, 'NO_DATA': 8, 'REFLINE': 15, 'SAT': 1, 'SUSPECT': 7, 'UNMASKEDNAN': 16}) + metadata = {} # Key-value pairs to put in the header + fluxTable = GAFluxTable(wavelength, flux, error, model, cont, + norm_flux, norm_error, norm_model, + mask, flags) + + stellarParams = StellarParams( + method=np.array(['rvfit', 'rvfit', 'rvfit', 'rvfit', 'rvfit']), + frame=np.array(['helio', '', '', '', '']), + param=np.array(['v_los', 'Fe_H', 'T_eff', 'log_g', 'a_Fe']), + covarId=np.array([0, 1, 2, 3, 4]), + unit=np.array(['km s-1', 'dex', 'K', '', 'dex']), + value=np.array([0.0, 0.0, 0.0, 0.0, 0.0]), + valueErr=np.array([0.0, 0.0, 0.0, 0.0, 0.0]), + flag=np.array([False, False, False, False, False]), + status=np.array(['', '', '', '', '']), + ) + + velocityCorrections = VelocityCorrections( + visit=visit, + JD=np.zeros_like(visit, dtype=np.float32), + helio=np.zeros_like(visit, dtype=np.float32), + bary=np.zeros_like(visit, dtype=np.float32), + ) + + abundances = Abundances( + method=np.array(['chemfit', 'chemfit', 'chemfit']), + element=np.array(['Mg', 'Ti', 'Si']), + covarId=np.array([0, 1, 2]), + value=np.array([0.0, 0.0, 0.0]), + valueErr=np.array([0.0, 0.0, 0.0]), + ) + + paramsCovar = np.eye(3, dtype=np.float32) + abundCovar = np.eye(4, dtype=np.float32) + notes = None + + return PfsGAObject(target, observations, + wavelength, flux, mask, sky, covar, covar2, + flags, metadata, + fluxTable, + stellarParams, + velocityCorrections, + abundances, + paramsCovar, + abundCovar, + notes) + + def assertPfsGAObject(self, lhs, rhs): + np.testing.assert_array_equal(lhs.observations.visit, rhs.observations.visit) + + # TODO: add more tests here + + def extractAttributes(self, cls, fileName): + matches = re.search(cls.filenameRegex, fileName) + if not matches: + self.fail( + "Unable to parse filename: {} using regex {}" + .format(fileName, cls.filenameRegex)) + + # Cannot use algorithm in PfsSpectra._parseFilename(), + # specifically cls.filenameKeys, due to ambiguity in parsing + # integers in hex format (eg objId). Need to parse cls.filenameFormat + ff = re.search(r'^[a-zA-Z]+(.*)\.fits', cls.filenameFormat)[1] + cmps = re.findall(r'-{0,1}(0x){0,1}\%\((\w+)\)\d*(\w)', ff) + fmts = [(kk, tt) for ox, kk, tt in cmps] + + d = {} + for (kk, tt), vv in zip(fmts, matches.groups()): + if tt == 'd': + ii = int(vv) + elif tt == 'x': + ii = int(vv, 16) + elif tt == 's': + ii = vv + d[kk] = ii + return d + + def test_filenameRegex(self): + d = self.extractAttributes( + PfsGAObject, + 'pfsGAObject-07621-01234-2,2-02468ace1234abcd-003-0x0123456789abcdef.fits') + self.assertEqual(d['catId'], 7621) + self.assertEqual(d['tract'], 1234) + self.assertEqual(d['patch'], '2,2') + self.assertEqual(d['objId'], 163971054118939597) + self.assertEqual(d['nVisit'], 3) + self.assertEqual(d['pfsVisitHash'], 81985529216486895) + + def test_getIdentity(self): + """Construct a PfsGAObject and get its identity.""" + + pfsGAObject = self.makePfsGAObject() + identity = pfsGAObject.getIdentity() + filename = pfsGAObject.filenameFormat % identity + + self.assertEqual('pfsGAObject-12345-00001-1,1-00000000075bcd15-002-0x05a95bc24d8ce16f.fits', filename) + + def test_validate(self): + """Construct a PfsGAObject and run validation.""" + + pfsGAObject = self.makePfsGAObject() + pfsGAObject.validate() + + def test_writeFits_fromFits(self): + """Construct a PfsGAObject and save it to a FITS file.""" + + pfsGAObject = self.makePfsGAObject() + + dirName = os.path.splitext(__file__)[0] + if not os.path.exists(dirName): + os.makedirs(dirName) + + id = pfsGAObject.getIdentity() + filename = os.path.join(dirName, pfsGAObject.filenameFormat % id) + if os.path.exists(filename): + os.unlink(filename) + + try: + pfsGAObject.writeFits(filename) + other = PfsGAObject.readFits(filename) + self.assertPfsGAObject(pfsGAObject, other) + except Exception as e: + raise + finally: + os.unlink(filename) From 50f308911ef7edf8c268d6e9904c35d15838f29a Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 31 May 2024 18:43:18 -0400 Subject: [PATCH 05/24] Do not fail on members defining class type of other members. --- tests/test_docstrings.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/tests/test_docstrings.py b/tests/test_docstrings.py index 4e590fd9..88eb14f3 100644 --- a/tests/test_docstrings.py +++ b/tests/test_docstrings.py @@ -5,18 +5,20 @@ import lsst.utils.tests from pfs.datamodel.drp import PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject +from pfs.datamodel.ga import PfsGAObject class DocstringsTestCase(unittest.TestCase): def testDocstrings(self): - for cls in (PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject): + for cls in (PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject, PfsGAObject): for name, attr in inspect.getmembers(cls): if not hasattr(attr, "__doc__") or not attr.__doc__: continue docstring = attr.__doc__ for base in cls.__mro__[1:-1]: - self.assertNotIn(base.__name__, docstring, - f"{cls.__name__}.{name}.__doc__ contains {base.__name__}: {docstring}") + if not name.endswith("Class"): + self.assertNotIn(base.__name__, docstring, + f"{cls.__name__}.{name}.__doc__ contains {base.__name__}: {docstring}") class TestMemory(lsst.utils.tests.MemoryTestCase): From ff25b1954f7613f9040f19435ad8a0f2c8c7a4c0 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 14 Jun 2024 15:50:42 -0400 Subject: [PATCH 06/24] Flags and status added to abundances, reorganized order of HDUs. --- python/pfs/datamodel/ga.py | 15 +++++++++------ tests/test_pfsGAObject.py | 2 ++ 2 files changed, 11 insertions(+), 6 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index c65d94f1..e4bf10f7 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -58,6 +58,9 @@ class Abundances(PfsTable): Column("covarId", np.uint8, "Param position within covariance matrix", -1), Column("value", np.float32, "Abundance value", np.nan), Column("valueErr", np.float32, "Abundance error", np.nan), + # TODO: will we have systematic errors? + Column("flag", bool, "Measurement flag (true means bad)", False), + Column("status", str, "Measurement flags", ""), ] fitsExtName = 'ABUND' @@ -279,11 +282,11 @@ def _readImpl(cls, fits): # TODO: handle missing extensions - data["stellarParams"] = StellarParams.readHdu(fits) data["velocityCorrections"] = VelocityCorrections.readHdu(fits) - data["abundances"] = Abundances.readHdu(fits) + data["stellarParams"] = StellarParams.readHdu(fits) if cls.StellarParamsFitsExtName in fits: data["paramsCovar"] = fits[cls.StellarParamsFitsExtName].data.astype(np.float32) + data["abundances"] = Abundances.readHdu(fits) if cls.AbundancesFitsExtName in fits: data["abundCovar"] = fits[cls.AbundancesFitsExtName].data.astype(np.float32) @@ -294,14 +297,14 @@ def _writeImpl(self, fits): header = super()._writeImpl(fits) - if self.stellarParams is not None: - self.stellarParams.writeHdu(fits) if self.velocityCorrections is not None: self.velocityCorrections.writeHdu(fits) - if self.abundances is not None: - self.abundances.writeHdu(fits) + if self.stellarParams is not None: + self.stellarParams.writeHdu(fits) if self.paramsCovar is not None: fits.append(ImageHDU(self.paramsCovar.astype(np.float32), header=header, name=self.StellarParamsFitsExtName)) + if self.abundances is not None: + self.abundances.writeHdu(fits) if self.abundCovar is not None: fits.append(ImageHDU(self.abundCovar.astype(np.float32), header=header, name=self.AbundancesFitsExtName)) diff --git a/tests/test_pfsGAObject.py b/tests/test_pfsGAObject.py index f720dc8b..4ec10261 100644 --- a/tests/test_pfsGAObject.py +++ b/tests/test_pfsGAObject.py @@ -86,6 +86,8 @@ def makePfsGAObject(self): covarId=np.array([0, 1, 2]), value=np.array([0.0, 0.0, 0.0]), valueErr=np.array([0.0, 0.0, 0.0]), + flag=np.array([False, False, False]), + status=np.array(['', '', '']), ) paramsCovar = np.eye(3, dtype=np.float32) From a4126e63c433d3100cab1cc70c7fe21bc857b77a Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 14 Jun 2024 15:51:33 -0400 Subject: [PATCH 07/24] Added docs on GA data model. --- datamodel.txt | 58 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 58 insertions(+) diff --git a/datamodel.txt b/datamodel.txt index 9cc4ce98..0ec5ddc3 100644 --- a/datamodel.txt +++ b/datamodel.txt @@ -699,6 +699,64 @@ The format will be identical to that of the pfsObject files. -=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=- +Physical parameter measurements for Galactic Archeology targets including flux-calibrated +combined spectra and synthetic stellar template fits. + + "pfsGAObject-%05d-%05d-%s-%016x-%03d-0x%016x.fits" + % (catId, tract, patch, objId, nVisit % 1000, pfsVisitHash) + +The format is similar to pfsObject files, but with additional HDUs and a different fluxTable format. + +HDU #j VELCORR Velocity corrections used at each visit [BINARY FITS TABLE] +HDU #i STELLARPARAMS Fundamental stellar parameter measurements [BINARY FITS TABLE] +HDU #l STELLARCOVAR Covariance matrix of stellar parameters [BINARY FITS TABLE] n*n +HDU #k ABUND Single element abundances [BINARY FITS TABLE] +HDU #m ABUNDCOVAR Covarinace matrix of single element abundances [32-bit FLOAT] n*n + +In the data tables outlined below, we allow for the possibility of multiple measurements of the same +physical parameters indicated by the `method` field but we provide the full covariance matrix for +the primary method only. The column `covarId` indicates the position of the parameter within the +corresponding covariance matrix. + +The VELCORR table lists the velocity corrections used at each visit: + + visit visit identifier 32-bit int + JD Julian Date of the visit 32-bit float + helio Heliocentric correction at the visit [km s-1] 32-bit float + bary Barycentric correction at the visit [km s-2] 32-bit float + +The STELLARPARAMS table lists the measured fundamental parameters + + method Method of measuring the parameters string + frame Reference frame for measuring the velocity string + param Stellar parameter string + covarId Index of parameter within the covariance matrix 8-bit uint + unit Physical unit of the parameter string + value Measured value of the parameter 32-bit float + valueErr Uncertainty of the measured value 32-bit float + flag Measurement flag (true means bad) bool + status Flags describing the quality of the measurement string + +The ABUND table lists the measured single element abundances + + method Method of measuring the parameters string + element Element name string + covarId Index of parameter within the covariance matrix 8-bit uint + value Measured abundance 32-bit float + valueErr Uncertainty of the measured abundance 32-bit float + flag Measurement flag (true means bad) bool + status Flags describing the quality of the measurement string + +The FLUXTABLE of pfsGAObject files is extended and includes the following additional columns: + + model Best fit fluxed template model 32-bit float + cont Continuum model 32-bit float + norm_flux Continuum-normalized flux 32-bit float + norm_error Error of continuum-normalized flux 32-bit float + norm_model Best fit continuum-normalized model 32-bit float + +-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=- + 2D spectral line measurements The 2D pipeline measures the position and flux of sky and arc lines as part of the reduction process. From 7882a2f37e728eea917aefb65b71a966bd3f1af8 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 14 Jun 2024 17:45:24 -0400 Subject: [PATCH 08/24] Fixed all flake8 errors. --- python/pfs/datamodel/fluxTable.py | 8 +++--- python/pfs/datamodel/ga.py | 26 ++++++++++++++------ python/pfs/datamodel/pfsFiberArray.py | 2 +- tests/test_docstrings.py | 5 ++-- tests/test_pfsGAObject.py | 35 ++++++++++++++------------- 5 files changed, 44 insertions(+), 32 deletions(-) diff --git a/python/pfs/datamodel/fluxTable.py b/python/pfs/datamodel/fluxTable.py index 48fde0d6..192827f7 100644 --- a/python/pfs/datamodel/fluxTable.py +++ b/python/pfs/datamodel/fluxTable.py @@ -42,7 +42,7 @@ def __init__(self, wavelength, flux, error, mask, flags): flux=flux, error=error, mask=mask) - + self.wavelength = wavelength self.flux = flux self.error = error @@ -52,11 +52,11 @@ def __init__(self, wavelength, flux, error, mask, flags): def __len__(self): """Return number of elements""" return len(self.wavelength) - + def checkShapes(self, **kwargs): keys = list(sorted(kwargs.keys())) - dims = np.array([ len(kwargs[k].shape) for k in keys ]) - lengths = set([ kwargs[k].shape for k in keys ]) + dims = np.array([len(kwargs[k].shape) for k in keys]) + lengths = set([kwargs[k].shape for k in keys]) if np.any(dims != 1) or len(lengths) > 1: names = ','.join(keys) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index e4bf10f7..69ad3a70 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -1,4 +1,3 @@ -from typing import Type import numpy as np from .notes import makeNotesClass, Notes @@ -18,6 +17,7 @@ "PfsGAObject", ] + class VelocityCorrections(PfsTable): """A table of velocity corrections applied to the individual visits.""" @@ -30,6 +30,7 @@ class VelocityCorrections(PfsTable): ] fitsExtName = 'VELCORR' + class StellarParams(PfsTable): """List of measured stellar parameters for a target.""" @@ -48,6 +49,7 @@ class StellarParams(PfsTable): ] fitsExtName = 'STELLARPARAM' + class Abundances(PfsTable): """List of measured abundance parameters for stellar targets.""" @@ -64,6 +66,7 @@ class Abundances(PfsTable): ] fitsExtName = 'ABUND' + class GAFluxTable(FluxTable): """Table of coadded fluxes at near-original sampling and model fits @@ -115,7 +118,7 @@ def __init__(self, wavelength, flux, error, model, cont, norm_flux, norm_error, norm_error=norm_error, norm_model=norm_model, mask=mask) - + self.wavelength = wavelength self.flux = flux self.error = error @@ -182,17 +185,19 @@ def fromFits(cls, fits): hdu.data["mask"].astype(np.int32), flags) + PfsGAObjectNotes = makeNotesClass( "PfsGAObjectNotes", [] ) + @inheritDocstrings class PfsGAObject(PfsFiberArray): """Coadded spectrum of a GA target with derived quantities. - + Produced by ˙˙gapipe`` - + Parameters ---------- target : `pfs.datamodel.Target` @@ -262,7 +267,8 @@ def __init__( abundCovar=None, notes: Notes = None, ): - super().__init__(target, observations, wavelength, flux, mask, sky, covar, covar2, flags, metadata=metadata, fluxTable=fluxTable, notes=notes) + super().__init__(target, observations, wavelength, flux, mask, sky, + covar, covar2, flags, metadata=metadata, fluxTable=fluxTable, notes=notes) self.stellarParams = stellarParams self.velocityCorrections = velocityCorrections @@ -302,10 +308,14 @@ def _writeImpl(self, fits): if self.stellarParams is not None: self.stellarParams.writeHdu(fits) if self.paramsCovar is not None: - fits.append(ImageHDU(self.paramsCovar.astype(np.float32), header=header, name=self.StellarParamsFitsExtName)) + fits.append(ImageHDU(self.paramsCovar.astype(np.float32), + header=header, + name=self.StellarParamsFitsExtName)) if self.abundances is not None: self.abundances.writeHdu(fits) if self.abundCovar is not None: - fits.append(ImageHDU(self.abundCovar.astype(np.float32), header=header, name=self.AbundancesFitsExtName)) + fits.append(ImageHDU(self.abundCovar.astype(np.float32), + header=header, + name=self.AbundancesFitsExtName)) - return header \ No newline at end of file + return header diff --git a/python/pfs/datamodel/pfsFiberArray.py b/python/pfs/datamodel/pfsFiberArray.py index 7c2cba90..9c08c9e9 100644 --- a/python/pfs/datamodel/pfsFiberArray.py +++ b/python/pfs/datamodel/pfsFiberArray.py @@ -153,4 +153,4 @@ def _writeImpl(self, fits): self.notes.writeFits(fits) if self.fluxTable is not None: self.fluxTable.toFits(fits) - return header \ No newline at end of file + return header diff --git a/tests/test_docstrings.py b/tests/test_docstrings.py index 88eb14f3..8b1553eb 100644 --- a/tests/test_docstrings.py +++ b/tests/test_docstrings.py @@ -17,8 +17,9 @@ def testDocstrings(self): docstring = attr.__doc__ for base in cls.__mro__[1:-1]: if not name.endswith("Class"): - self.assertNotIn(base.__name__, docstring, - f"{cls.__name__}.{name}.__doc__ contains {base.__name__}: {docstring}") + self.assertNotIn( + base.__name__, docstring, + f"{cls.__name__}.{name}.__doc__ contains {base.__name__}: {docstring}") class TestMemory(lsst.utils.tests.MemoryTestCase): diff --git a/tests/test_pfsGAObject.py b/tests/test_pfsGAObject.py index 4ec10261..2353c64c 100644 --- a/tests/test_pfsGAObject.py +++ b/tests/test_pfsGAObject.py @@ -9,6 +9,7 @@ from pfs.datamodel import GAFluxTable from pfs.datamodel import PfsGAObject, StellarParams, Abundances, VelocityCorrections + class PfsGAObjectTestCase(TestCase): """ Check the format of example datamodel files are consistent with that specified in the corresponding @@ -28,8 +29,8 @@ def makePfsGAObject(self): target = Target(catId, tract, patch, objId, ra, dec, targetType) - visit = np.array([ 83219, 83220 ]) - arm = np.array([ 'b', 'm', ]) + visit = np.array([83219, 83220]) + arm = np.array(['b', 'm']) spectrograph = np.array([1, 1]) pfsDesignId = np.array([8854764194165386399, 8854764194165386400]) fiberId = np.array([476, 476]) @@ -55,7 +56,7 @@ def makePfsGAObject(self): covar = np.zeros((3, wavelength.size), dtype=np.float32) # Tridiagonal covariance matrix of flux covar2 = np.zeros((1, 1), dtype=np.float32) # ? - flags = MaskHelper() # {'BAD': 0, 'BAD_FIBERTRACE': 11, 'BAD_FLAT': 9, 'BAD_FLUXCAL': 13, 'BAD_SKY': 12, 'CR': 3, 'DETECTED': 5, 'DETECTED_NEGATIVE': 6, 'EDGE': 4, 'FIBERTRACE': 10, 'INTRP': 2, 'IPC': 14, 'NO_DATA': 8, 'REFLINE': 15, 'SAT': 1, 'SUSPECT': 7, 'UNMASKEDNAN': 16}) + flags = MaskHelper() # metadata = {} # Key-value pairs to put in the header fluxTable = GAFluxTable(wavelength, flux, error, model, cont, norm_flux, norm_error, norm_model, @@ -95,16 +96,16 @@ def makePfsGAObject(self): notes = None return PfsGAObject(target, observations, - wavelength, flux, mask, sky, covar, covar2, - flags, metadata, - fluxTable, - stellarParams, - velocityCorrections, - abundances, - paramsCovar, - abundCovar, - notes) - + wavelength, flux, mask, sky, covar, covar2, + flags, metadata, + fluxTable, + stellarParams, + velocityCorrections, + abundances, + paramsCovar, + abundCovar, + notes) + def assertPfsGAObject(self, lhs, rhs): np.testing.assert_array_equal(lhs.observations.visit, rhs.observations.visit) @@ -134,11 +135,11 @@ def extractAttributes(self, cls, fileName): ii = vv d[kk] = ii return d - + def test_filenameRegex(self): d = self.extractAttributes( - PfsGAObject, - 'pfsGAObject-07621-01234-2,2-02468ace1234abcd-003-0x0123456789abcdef.fits') + PfsGAObject, + 'pfsGAObject-07621-01234-2,2-02468ace1234abcd-003-0x0123456789abcdef.fits') self.assertEqual(d['catId'], 7621) self.assertEqual(d['tract'], 1234) self.assertEqual(d['patch'], '2,2') @@ -180,6 +181,6 @@ def test_writeFits_fromFits(self): other = PfsGAObject.readFits(filename) self.assertPfsGAObject(pfsGAObject, other) except Exception as e: - raise + raise e finally: os.unlink(filename) From f925cc2c93551c39abed52356f7f8e503934c39d Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Tue, 22 Apr 2025 15:07:23 -0400 Subject: [PATCH 09/24] Added GA catalog data model first version. --- python/pfs/datamodel/ga.py | 258 ++++++++++++++++++++++++++++++++++++- tests/test_PfsGACatalog.py | 143 ++++++++++++++++++++ 2 files changed, 397 insertions(+), 4 deletions(-) create mode 100644 tests/test_PfsGACatalog.py diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 69ad3a70..c6e51492 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -1,3 +1,4 @@ +import os import numpy as np from .notes import makeNotesClass, Notes @@ -6,7 +7,11 @@ from .pfsTable import PfsTable, Column from .utils import inheritDocstrings from .utils import astropyHeaderToDict, astropyHeaderFromDict +from .utils import calculatePfsVisitHash, wraparoundNVisit from .masks import MaskHelper +from .observations import Observations + +GA_DAMD_VER = 2 __all__ = [ "VelocityCorrections", @@ -21,7 +26,7 @@ class VelocityCorrections(PfsTable): """A table of velocity corrections applied to the individual visits.""" - damdVer = 2 + damdVer = GA_DAMD_VER schema = [ Column("visit", np.int32, "ID of the visit these corrections apply for", -1), Column("JD", np.float32, "Julian date of the visit", -1), @@ -34,7 +39,7 @@ class VelocityCorrections(PfsTable): class StellarParams(PfsTable): """List of measured stellar parameters for a target.""" - damdVer = 2 + damdVer = GA_DAMD_VER schema = [ Column("method", str, "Line-of-sight velocity measurement method", ""), Column("frame", str, "Reference frame of velocity: helio, bary", ""), @@ -53,7 +58,7 @@ class StellarParams(PfsTable): class Abundances(PfsTable): """List of measured abundance parameters for stellar targets.""" - damdVer = 2 + damdVer = GA_DAMD_VER schema = [ Column("method", str, "Abundance measurement method", ""), Column("element", str, "Chemical element the abundance is measured for", ""), @@ -143,7 +148,7 @@ def toFits(self, fits): # the versions.txt file. from astropy.io.fits import BinTableHDU, Column header = astropyHeaderFromDict(self.flags.toFitsHeader()) - header['DAMD_VER'] = (1, "GAFluxTable datamodel version") + header['DAMD_VER'] = (GA_DAMD_VER, "GAFluxTable datamodel version") hdu = BinTableHDU.from_columns([ Column("wavelength", "D", array=self.wavelength), Column("flux", "E", array=self.flux), @@ -319,3 +324,248 @@ def _writeImpl(self, fits): name=self.AbundancesFitsExtName)) return header + +class GACatalogTable(PfsTable): + """Catalog of GA objects with associated parameters.""" + + damdVer = GA_DAMD_VER + schema = [ + Column("catId", np.int32, "PFS catalog identifier", -1), + Column("objId", np.int32, "Object identifier.", -1), + Column("gaiaId", np.int64, "GAIA identifier", -1), + Column("ps1Id", np.int64, "PS1 identifier", -1), + Column("hscId", np.int64, "HSC identifier", -1), + Column("miscId", np.int64, "Miscellaneous identifier", -1), + Column("ra", np.float32, "Right ascension ICRS [deg]", np.nan), + Column("dec", np.float32, "Declination ICRS [deg]", np.nan), + Column("epoch", np.float32, "Coordinate epoch [Jyr]", 2016.0), + Column("pmRa", np.float32, "Proper motion pmracosdec [mas/yr]", np.nan), + Column("pmDec", np.float32, "Proper motion dec [mas/yr]", np.nan), + Column("parallax", np.float32, "Parallax [mas]", np.nan), + Column("targetType", np.int16, "Target type.", -1), + Column("proposalId", str, "Proposal ID", ""), + Column("obCode", str, "Observing Block", ""), + + Column("nVisit_b", np.int16, "Number of visits in B", -1), + Column("nVisit_m", np.int16, "Number of visits in M", -1), + Column("nVisit_r", np.int16, "Number of visits in R", -1), + Column("nVisit_n", np.int16, "Number of visits in N", -1), + + Column("expTimeEff_b", np.float32, "Effective exposure time in B [s]", np.nan), + Column("expTimeEff_m", np.float32, "Effective exposure time in M [s]", np.nan), + Column("expTimeEff_r", np.float32, "Effective exposure time in R [s]", np.nan), + Column("expTimeEff_n", np.float32, "Effective exposure time in N [s]", np.nan), + + Column("rv", np.float32, "Radial velocity [km/s]", np.nan), + Column("rvErr", np.float32, "Radial velocity error [km/s]", np.nan), + Column("tEff", np.float32, "Effective temperature [K]", np.nan), + Column("tEffErr", np.float32, "Effective temperature error [K]", np.nan), + Column("m_h", np.float32, "Metallicity [dex]", np.nan), + Column("m_hErr", np.float32, "Metallicity error [dex]", np.nan), + Column("logg", np.float32, "log g", np.nan), + Column("loggErr", np.float32, "log g error", np.nan), + + Column("flag", bool, "Measurement flag (true means bad)", False), + Column("status", str, "Measurement flags", ""), + ] + fitsExtName = 'GACATALOG' + +PfsGACatalogNotes = makeNotesClass( + "PfsGACatalogNotes", + [] +) + +class PfsGACatalog(): + filenameFormat = ("pfsGACatalog-%(catId)05d-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") + filenameRegex = r"^pfsGACatalog-(\d{5})-([0-9]{3})-0x([0-9a-f]{16})\.fits.*$" + filenameKeys = [("catId", int), ("nVisit", int), ("pfsVisitHash", int)] + + NotesClass = PfsGACatalogNotes + + def __init__( + self, + catId, + observations: Observations, + catalog: GACatalogTable, + metadata=None, + notes: Notes = None): + + self.catId = catId + self.observations = observations + self.nVisit = wraparoundNVisit(len(observations.visit)) + self.pfsVisitHash = calculatePfsVisitHash(observations.visit) + self.catalog = catalog + self.metadata = metadata if metadata is not None else {} + self.notes = notes if notes is not None else self.NotesClass() + self.length = len(catalog) + + self.validate() + + def getIdentity(self): + """Return the identity of the catalog + + Returns + ------- + identity : `dict` + Key-value pairs that identify this catalog + """ + identity = dict( + catId=self.catId, + visit=self.observations.visit, + nVisit=self.nVisit, + pfsVisitHash=self.pfsVisitHash + ) + return identity + + def validate(self): + + # TODO: Make sure catId is the same for each object + + pass + + def __len__(self): + """Return the length of the arrays""" + return self.length + + @classmethod + def _readImpl(cls, fits): + data = {} + + version = fits[0].header["DAMD_VER"] + + if version >= GA_DAMD_VER: + data["notes"] = cls.NotesClass.readFits(fits) + + try: + observations = Observations.fromFits(fits) + data['observations'] = observations + except KeyError as exc: + # Only want to catch "Extension XXX not found." + if not exc.args[0].startswith("Extension"): + raise + data['observations'] = None + + try: + catalog = GACatalogTable.readHdu(fits) + data['catalog'] = catalog + data['catId'] = catalog.catId[0] + except KeyError as exc: + # Only want to catch "Extension XXX not found." + if not exc.args[0].startswith("Extension"): + raise + data['catalog'] = None + + return data + + @classmethod + def readFits(cls, filename): + """ + Read from FITS file + + This API is intended for use by the LSST data butler, which handles + translating the desired identity into a filename. + + Parameters + ---------- + filename : `str` + Filename of FITS file. + + Returns + ------- + self : ``cls`` + Constructed instance, from FITS file. + """ + + import astropy.io.fits + with astropy.io.fits.open(filename) as fd: + data = cls._readImpl(fd) + return cls(**data) + + @classmethod + def read(cls, catId, pfsVisitHash, dirName="."): + """Read from file + + This API is intended for use by science users, as it allows selection + of the correct file from parameters that make sense, such as which + catId, objId, etc. + + Parameters + ---------- + identity : `dict` + Keyword-value pairs identifying the data of interest. Common keywords + include ``catId``, ``tract``, ``patch``, ``objId``. + dirName : `str`, optional + Directory from which to read. + + Returns + ------- + self : ``cls`` + Spectrum read from file. + """ + filename = os.path.join(dirName, cls.filenameFormat % (catId, pfsVisitHash)) + return cls.readFits(filename) + + def _writeImpl(self, fits): + """ + Implementation for writing to FITS file + + Parameters + ---------- + fits : `astropy.io.fits.HDUList` + List of FITS HDUs. This has a Primary HDU already, the header of + which may be supplemented with additional keywords. + """ + + from astropy.io.fits import ImageHDU, Header + + header = Header() + + header['DAMD_VER'] = GA_DAMD_VER + + if self.observations is not None: + # Override catId from class + self.observations.catId = np.array(self.observations.num * [self.catId]) + self.observations.toFits(fits) + if self.catalog is not None: + self.catalog.writeHdu(fits) + if self.metadata is not None: + header.extend(astropyHeaderFromDict(self.metadata)) + if self.notes is not None: + self.notes.writeFits(fits) + + return header + + def writeFits(self, filename): + """Write to FITS file + + This API is intended for use by the LSST data butler, which handles + translating the desired identity into a filename. + + Parameters + ---------- + filename : `str` + Filename of FITS file. + """ + from astropy.io.fits import HDUList, PrimaryHDU + fits = HDUList() + header = self._writeImpl(fits) + fits[0].header.update(header) + + with open(filename, "wb") as fd: + fits.writeto(fd) + + def write(self, dirName="."): + """Write to file + + This API is intended for use by science users, as it allows setting the + correct filename from parameters that make sense, such as which + catId, objId, etc. + + Parameters + ---------- + dirName : `str`, optional + Directory to which to write. + """ + identity = self.getIdentity() + filename = os.path.join(dirName, self.filenameFormat % identity) + return self.writeFits(filename) \ No newline at end of file diff --git a/tests/test_PfsGACatalog.py b/tests/test_PfsGACatalog.py new file mode 100644 index 00000000..1146ca98 --- /dev/null +++ b/tests/test_PfsGACatalog.py @@ -0,0 +1,143 @@ +import os +import re +import numpy as np +from unittest import TestCase + +from pfs.datamodel import TargetType, Observations +from pfs.datamodel.ga import PfsGACatalog, PfsGACatalogNotes, GACatalogTable + + +class PfsGACatalogTestCase(TestCase): + """ Check the format of example datamodel files are + consistent with that specified in the corresponding + datamodel classes. + """ + + def makePfsGACatalog(self): + """Construct a PfsGACatalog with dummy values for testing.""" + + catId = 12345 + + observations = Observations( + visit=np.array([10000, 10001, 10002, 10005]), + arm=np.array(['bmn', 'bmn', 'bmn', 'bmn']), + spectrograph=np.array(4 * [-1]), + pfsDesignId=np.array([100001, 100002, 100003, 100004]), + fiberId=np.array(4 * [-1]), + pfiNominal=np.array(4 * [2 * [np.nan]]), + pfiCenter=np.array(4 * [2 * [np.nan]]), + obsTime=np.array(['2025-01-10 12:00:00', '2025-01-10 12:15:00', '2025-01-10 12:30:00', '2025-01-10 13:00:00']), + expTime=np.array([1800.0, 1800.0, 1800.0, 1800.0]), + ) + + catalog = GACatalogTable( + catId=np.array([10001, 10001, 10001]), + objId=np.array([1, 2, 3]), + gaiaId=np.array([11, 12, 13]), + ps1Id=np.array([21, 22, 23]), + hscId=np.array([21, 22, 23]), + miscId=np.array([-1, -1, 33]), + ra=np.array([210.01, 210.02, 210.03]), + dec=np.array([67.1, 67.2, 67.3]), + epoch=np.array([2016.0, 2016.0, 2016.0]), + pmRa=np.array([0.0, 0.0, 0.0]), + pmDec=np.array([0.0, 0.0, 0.0]), + parallax=np.array([0.0, 0.0, 0.0]), + targetType=np.array([TargetType.SCIENCE, TargetType.SCIENCE, TargetType.SCIENCE]), + proposalId=np.array(['A', 'A', 'A']), + obCode=np.array(['A', 'B', 'C']), + + nVisit_b=np.array([3, 3, 3]), + nVisit_m=np.array([3, 3, 3]), + nVisit_r=np.array([3, 3, 3]), + nVisit_n=np.array([3, 3, 3]), + + expTimeEff_b=np.array([5400, 5400, 5400]), + expTimeEff_m=np.array([5400, 5400, 5400]), + expTimeEff_r=np.array([5400, 5400, 5400]), + expTimeEff_n=np.array([5400, 5400, 5400]), + + rv=np.array([0.0, 0.0, 0.0]), + rvErr=np.array([0.0, 0.0, 0.0]), + tEff=np.array([0.0, 0.0, 0.0]), + tEffErr=np.array([0.0, 0.0, 0.0]), + m_h=np.array([0.0, 0.0, 0.0]), + m_hErr=np.array([0.0, 0.0, 0.0]), + logg=np.array([0.0, 0.0, 0.0]), + loggErr=np.array([0.0, 0.0, 0.0]), + + flag=np.array([False, False, False]), + status=np.array(['', '', '']), + ) + + return PfsGACatalog( + catId, observations, catalog, metadata=None, notes=None + ) + + def assertPfsGACatalog(self, lhs, rhs): + np.testing.assert_array_equal(lhs.observations.visit, rhs.observations.visit) + + # TODO: add more tests here + + def extractAttributes(self, cls, fileName): + matches = re.search(cls.filenameRegex, fileName) + if not matches: + self.fail( + "Unable to parse filename: {} using regex {}" + .format(fileName, cls.filenameRegex)) + + # Cannot use algorithm in PfsSpectra._parseFilename(), + # specifically cls.filenameKeys, due to ambiguity in parsing + # integers in hex format (eg objId). Need to parse cls.filenameFormat + ff = re.search(r'^[a-zA-Z]+(.*)\.fits', cls.filenameFormat)[1] + cmps = re.findall(r'-{0,1}(0x){0,1}\%\((\w+)\)\d*(\w)', ff) + fmts = [(kk, tt) for ox, kk, tt in cmps] + + d = {} + for (kk, tt), vv in zip(fmts, matches.groups()): + if tt == 'd': + ii = int(vv) + elif tt == 'x': + ii = int(vv, 16) + elif tt == 's': + ii = vv + d[kk] = ii + return d + + def test_filenameRegex(self): + d = self.extractAttributes( + PfsGACatalog, + 'pfsGACatalog-07621-003-0x0123456789abcdef.fits') + self.assertEqual(d['catId'], 7621) + self.assertEqual(d['nVisit'], 3) + self.assertEqual(d['pfsVisitHash'], 81985529216486895) + + def test_validate(self): + """Construct a PfsGACatalog and run validation.""" + + pfsGACatalog = self.makePfsGACatalog() + pfsGACatalog.validate() + + def test_writeFits_fromFits(self): + """Construct a PfsGACatalog and save it to a FITS file.""" + + pfsGACatalog = self.makePfsGACatalog() + + dirName = os.path.splitext(__file__)[0] + if not os.path.exists(dirName): + os.makedirs(dirName) + + id = pfsGACatalog.getIdentity() + filename = os.path.join(dirName, pfsGACatalog.filenameFormat % id) + if os.path.exists(filename): + os.unlink(filename) + + try: + pfsGACatalog.writeFits(filename) + other = PfsGACatalog.readFits(filename) + self.assertPfsGACatalog(pfsGACatalog, other) + except Exception as e: + raise e + finally: + if os.path.exists(filename): + os.unlink(filename) From f7b0d45bec920649fb044a923d9f2bff5967500e Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Tue, 22 Apr 2025 15:21:25 -0400 Subject: [PATCH 10/24] Fixed flake8 issues. --- python/pfs/datamodel/ga.py | 19 +++++++++++-------- tests/test_PfsGACatalog.py | 7 ++++--- 2 files changed, 15 insertions(+), 11 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index c6e51492..abfacb2f 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -325,6 +325,7 @@ def _writeImpl(self, fits): return header + class GACatalogTable(PfsTable): """Catalog of GA objects with associated parameters.""" @@ -370,11 +371,13 @@ class GACatalogTable(PfsTable): ] fitsExtName = 'GACATALOG' + PfsGACatalogNotes = makeNotesClass( "PfsGACatalogNotes", [] ) + class PfsGACatalog(): filenameFormat = ("pfsGACatalog-%(catId)05d-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") filenameRegex = r"^pfsGACatalog-(\d{5})-([0-9]{3})-0x([0-9a-f]{16})\.fits.*$" @@ -389,7 +392,7 @@ def __init__( catalog: GACatalogTable, metadata=None, notes: Notes = None): - + self.catId = catId self.observations = observations self.nVisit = wraparoundNVisit(len(observations.visit)) @@ -418,7 +421,7 @@ def getIdentity(self): return identity def validate(self): - + # TODO: Make sure catId is the same for each object pass @@ -516,12 +519,12 @@ def _writeImpl(self, fits): which may be supplemented with additional keywords. """ - from astropy.io.fits import ImageHDU, Header - + from astropy.io.fits import Header + header = Header() header['DAMD_VER'] = GA_DAMD_VER - + if self.observations is not None: # Override catId from class self.observations.catId = np.array(self.observations.num * [self.catId]) @@ -546,11 +549,11 @@ def writeFits(self, filename): filename : `str` Filename of FITS file. """ - from astropy.io.fits import HDUList, PrimaryHDU + from astropy.io.fits import HDUList fits = HDUList() header = self._writeImpl(fits) fits[0].header.update(header) - + with open(filename, "wb") as fd: fits.writeto(fd) @@ -568,4 +571,4 @@ def write(self, dirName="."): """ identity = self.getIdentity() filename = os.path.join(dirName, self.filenameFormat % identity) - return self.writeFits(filename) \ No newline at end of file + return self.writeFits(filename) diff --git a/tests/test_PfsGACatalog.py b/tests/test_PfsGACatalog.py index 1146ca98..f0844ee3 100644 --- a/tests/test_PfsGACatalog.py +++ b/tests/test_PfsGACatalog.py @@ -4,7 +4,7 @@ from unittest import TestCase from pfs.datamodel import TargetType, Observations -from pfs.datamodel.ga import PfsGACatalog, PfsGACatalogNotes, GACatalogTable +from pfs.datamodel.ga import PfsGACatalog, GACatalogTable class PfsGACatalogTestCase(TestCase): @@ -26,7 +26,8 @@ def makePfsGACatalog(self): fiberId=np.array(4 * [-1]), pfiNominal=np.array(4 * [2 * [np.nan]]), pfiCenter=np.array(4 * [2 * [np.nan]]), - obsTime=np.array(['2025-01-10 12:00:00', '2025-01-10 12:15:00', '2025-01-10 12:30:00', '2025-01-10 13:00:00']), + obsTime=np.array(['2025-01-10 12:00:00', '2025-01-10 12:15:00', + '2025-01-10 12:30:00', '2025-01-10 13:00:00']), expTime=np.array([1800.0, 1800.0, 1800.0, 1800.0]), ) @@ -46,7 +47,7 @@ def makePfsGACatalog(self): targetType=np.array([TargetType.SCIENCE, TargetType.SCIENCE, TargetType.SCIENCE]), proposalId=np.array(['A', 'A', 'A']), obCode=np.array(['A', 'B', 'C']), - + nVisit_b=np.array([3, 3, 3]), nVisit_m=np.array([3, 3, 3]), nVisit_r=np.array([3, 3, 3]), From 5afda9d6ef67957100fa7db1353510473ea58ede Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Tue, 22 Apr 2025 15:22:14 -0400 Subject: [PATCH 11/24] Fixed file name capitalization --- tests/{test_pfsGAObject.py => test_PfsGAObject.py} | 0 1 file changed, 0 insertions(+), 0 deletions(-) rename tests/{test_pfsGAObject.py => test_PfsGAObject.py} (100%) diff --git a/tests/test_pfsGAObject.py b/tests/test_PfsGAObject.py similarity index 100% rename from tests/test_pfsGAObject.py rename to tests/test_PfsGAObject.py From a6bdfbcbf5d1e782236300dd6b5db56e7347667a Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 25 Apr 2025 06:29:35 -0400 Subject: [PATCH 12/24] Working on the catalog table data model. --- python/pfs/datamodel/ga.py | 23 ++++++++++++++--------- 1 file changed, 14 insertions(+), 9 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index abfacb2f..07db88ff 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -20,6 +20,9 @@ "GAFluxTable", "PfsGAObjectNotes", "PfsGAObject", + "GACatalogTable", + "PfsGACatalogNotes", + "PfsGACatalog", ] @@ -339,7 +342,7 @@ class GACatalogTable(PfsTable): Column("miscId", np.int64, "Miscellaneous identifier", -1), Column("ra", np.float32, "Right ascension ICRS [deg]", np.nan), Column("dec", np.float32, "Declination ICRS [deg]", np.nan), - Column("epoch", np.float32, "Coordinate epoch [Jyr]", 2016.0), + Column("epoch", str, "Coordinate epoch [Jyr]", ''), Column("pmRa", np.float32, "Proper motion pmracosdec [mas/yr]", np.nan), Column("pmDec", np.float32, "Proper motion dec [mas/yr]", np.nan), Column("parallax", np.float32, "Parallax [mas]", np.nan), @@ -357,14 +360,16 @@ class GACatalogTable(PfsTable): Column("expTimeEff_r", np.float32, "Effective exposure time in R [s]", np.nan), Column("expTimeEff_n", np.float32, "Effective exposure time in N [s]", np.nan), - Column("rv", np.float32, "Radial velocity [km/s]", np.nan), - Column("rvErr", np.float32, "Radial velocity error [km/s]", np.nan), - Column("tEff", np.float32, "Effective temperature [K]", np.nan), - Column("tEffErr", np.float32, "Effective temperature error [K]", np.nan), - Column("m_h", np.float32, "Metallicity [dex]", np.nan), - Column("m_hErr", np.float32, "Metallicity error [dex]", np.nan), - Column("logg", np.float32, "log g", np.nan), - Column("loggErr", np.float32, "log g error", np.nan), + # TODO: add SNR + + Column("v_los", np.float32, "Radial velocity [km/s]", np.nan), + Column("v_losErr", np.float32, "Radial velocity error [km/s]", np.nan), + Column("T_eff", np.float32, "Effective temperature [K]", np.nan), + Column("T_effErr", np.float32, "Effective temperature error [K]", np.nan), + Column("M_H", np.float32, "Metallicity [dex]", np.nan), + Column("M_HErr", np.float32, "Metallicity error [dex]", np.nan), + Column("log_g", np.float32, "log g", np.nan), + Column("log_gErr", np.float32, "log g error", np.nan), Column("flag", bool, "Measurement flag (true means bad)", False), Column("status", str, "Measurement flags", ""), From 64a6fab996656d268f08751a3ecfa77bcfa178e6 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Wed, 11 Jun 2025 11:26:23 -0400 Subject: [PATCH 13/24] Added tempfit error flags. --- python/pfs/datamodel/ga.py | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 07db88ff..18debab8 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -1,6 +1,8 @@ import os import numpy as np +from enum import IntFlag +from .pfsConfig import DocEnum from .notes import makeNotesClass, Notes from .pfsFiberArray import PfsFiberArray from .fluxTable import FluxTable @@ -14,6 +16,7 @@ GA_DAMD_VER = 2 __all__ = [ + "TempFitFlag", "VelocityCorrections", "StellarParams", "Abundances", @@ -25,6 +28,21 @@ "PfsGACatalog", ] +class TempFitFlag(IntFlag): + """Flags for template and RV fitting.""" + + # NOTE: Keep this consistent with the flags in pfs.ga.pfsspec.tempfit.TempFitFlag + + OK = 0 # No flags + BADINIT = 1 << 0 # "Bad initial values" + NOCONVERGE = 1 << 1 # "Optimization did not convergence" + BADCONVERGE = 1 << 2 # "Bad convergence" + MAXITER = 1 << 3 # "Maximum iterations reached" + PARAMEDGE = 1 << 4 # "Template parameters at the bounds" + BADERROR = 1 << 5 # "Errors could not be calculated" + BADCOV = 1 << 6 # "Covariance matrix could not be calculated" + UNLIKELYPRIOR = 1 << 7 # "Unlikely prior" + class VelocityCorrections(PfsTable): """A table of velocity corrections applied to the individual visits.""" From bfe271ce49be8288e630ca1b1b350ebb0974e79e Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sat, 30 Aug 2025 09:09:27 -0400 Subject: [PATCH 14/24] Handle variable length string columns --- python/pfs/datamodel/pfsTargetSpectra.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/python/pfs/datamodel/pfsTargetSpectra.py b/python/pfs/datamodel/pfsTargetSpectra.py index 00bcdc23..e03257fb 100644 --- a/python/pfs/datamodel/pfsTargetSpectra.py +++ b/python/pfs/datamodel/pfsTargetSpectra.py @@ -210,7 +210,13 @@ def readData(hduName: str, dtype: type, mask=None) -> Union[np.ndarray, List[np. for key, value in kwargs.items(): if key not in targetHdu.columns.names: raise KeyError(f"Keyword argument '{key}' not found in TARGET HDU.") - m = targetHdu[key] == value + + if targetHdu.dtype[key].kind == "V": + # Handle variable length string columns + m = np.array([''.join(v) for v in targetHdu[key]]) == value + else: + m = np.array(targetHdu[key]) == value + mask = m if mask is None else mask & m targetFluxHdu = fits["TARGETFLUX"].data From ca1be79e2a32bb97731afc379efb42c4f78cfe58 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sat, 18 Oct 2025 09:37:07 -0400 Subject: [PATCH 15/24] Update stellarParams method to use 'tempfit' for consistency in PfsGAObject tests --- tests/test_PfsGAObject.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_PfsGAObject.py b/tests/test_PfsGAObject.py index 2353c64c..63e42464 100644 --- a/tests/test_PfsGAObject.py +++ b/tests/test_PfsGAObject.py @@ -63,7 +63,7 @@ def makePfsGAObject(self): mask, flags) stellarParams = StellarParams( - method=np.array(['rvfit', 'rvfit', 'rvfit', 'rvfit', 'rvfit']), + method=np.array(['tempfit', 'tempfit', 'tempfit', 'tempfit', 'tempfit']), frame=np.array(['helio', '', '', '', '']), param=np.array(['v_los', 'Fe_H', 'T_eff', 'log_g', 'a_Fe']), covarId=np.array([0, 1, 2, 3, 4]), From 08a979e8957132f15444d0c9144f5c6c867fb957 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sat, 18 Oct 2025 14:27:57 -0400 Subject: [PATCH 16/24] Refactor GA data model to Star terminology and add corresponding test cases --- python/pfs/datamodel/ga.py | 65 +++++++++------- .../{test_PfsGAObject.py => test_PfsStar.py} | 67 ++++++++-------- ...PfsGACatalog.py => test_PfsStarCatalog.py} | 78 +++++++++++-------- tests/test_docstrings.py | 4 +- 4 files changed, 118 insertions(+), 96 deletions(-) rename tests/{test_PfsGAObject.py => test_PfsStar.py} (76%) rename tests/{test_PfsGACatalog.py => test_PfsStarCatalog.py} (65%) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 18debab8..defc6680 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -20,12 +20,12 @@ "VelocityCorrections", "StellarParams", "Abundances", - "GAFluxTable", - "PfsGAObjectNotes", - "PfsGAObject", - "GACatalogTable", - "PfsGACatalogNotes", - "PfsGACatalog", + "StarFluxTable", + "PfsStarNotes", + "PfsStar", + "StarCatalogTable", + "PfsStarCatalogNotes", + "PfsStarCatalog", ] class TempFitFlag(IntFlag): @@ -93,7 +93,7 @@ class Abundances(PfsTable): fitsExtName = 'ABUND' -class GAFluxTable(FluxTable): +class StarFluxTable(FluxTable): """Table of coadded fluxes at near-original sampling and model fits Merged and coadded spectra have been resampled to a standard wavelength @@ -169,7 +169,7 @@ def toFits(self, fits): # the versions.txt file. from astropy.io.fits import BinTableHDU, Column header = astropyHeaderFromDict(self.flags.toFitsHeader()) - header['DAMD_VER'] = (GA_DAMD_VER, "GAFluxTable datamodel version") + header['DAMD_VER'] = (GA_DAMD_VER, "StarFluxTable datamodel version") hdu = BinTableHDU.from_columns([ Column("wavelength", "D", array=self.wavelength), Column("flux", "E", array=self.flux), @@ -212,14 +212,14 @@ def fromFits(cls, fits): flags) -PfsGAObjectNotes = makeNotesClass( - "PfsGAObjectNotes", +PfsStarNotes = makeNotesClass( + "PfsStarNotes", [] ) @inheritDocstrings -class PfsGAObject(PfsFiberArray): +class PfsStar(PfsFiberArray): """Coadded spectrum of a GA target with derived quantities. Produced by ˙˙gapipe`` @@ -246,7 +246,7 @@ class PfsGAObject(PfsFiberArray): Helper for dealing with symbolic names for mask values. metadata : `dict` (`str`: POD), optional Keyword-value pairs for the header. - fluxTable : `pfs.datamodel.GAFluxTable`, optional + fluxTable : `pfs.datamodel.StarFluxTable`, optional Table of coadded fluxes and continuum-normalized flux from contributing observations. stellarParams: `pfs.datamodel.StellarParams`, optional Table of measured stellar parameters. @@ -262,13 +262,13 @@ class PfsGAObject(PfsFiberArray): Reduction notes. """ - filenameFormat = ("pfsGAObject-%(catId)05d-%(tract)05d-%(patch)s-%(objId)016x" + filenameFormat = ("pfsStar-%(catId)05d-%(tract)05d-%(patch)s-%(objId)016x" "-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") - filenameRegex = r"^pfsGAObject-(\d{5})-(\d{5})-(.*)-([0-9a-f]{16})-(\d{3})-0x([0-9a-f]{16})\.fits.*$" + filenameRegex = r"^pfsStar-(\d{5})-(\d{5})-(.*)-([0-9a-f]{16})-(\d{3})-0x([0-9a-f]{16})\.fits.*$" filenameKeys = [("catId", int), ("tract", int), ("patch", str), ("objId", int), ("nVisit", int), ("pfsVisitHash", int)] - NotesClass = PfsGAObjectNotes - FluxTableClass = GAFluxTable + NotesClass = PfsStarNotes + FluxTableClass = StarFluxTable StellarParamsFitsExtName = "STELLARCOVAR" AbundancesFitsExtName = "ABUNDCOVAR" @@ -347,8 +347,8 @@ def _writeImpl(self, fits): return header -class GACatalogTable(PfsTable): - """Catalog of GA objects with associated parameters.""" +class StarCatalogTable(PfsTable): + """Catalog of stellar objects with associated parameters.""" damdVer = GA_DAMD_VER schema = [ @@ -378,41 +378,50 @@ class GACatalogTable(PfsTable): Column("expTimeEff_r", np.float32, "Effective exposure time in R [s]", np.nan), Column("expTimeEff_n", np.float32, "Effective exposure time in N [s]", np.nan), - # TODO: add SNR + Column("snr_b", np.float32, "Signal-to-noise ratio in B", np.nan), + Column("snr_m", np.float32, "Signal-to-noise ratio in M", np.nan), + Column("snr_r", np.float32, "Signal-to-noise ratio in R", np.nan), + Column("snr_n", np.float32, "Signal-to-noise ratio in N", np.nan), Column("v_los", np.float32, "Radial velocity [km/s]", np.nan), Column("v_losErr", np.float32, "Radial velocity error [km/s]", np.nan), + Column("EBV", np.float32, "E(B-V) extinction", np.nan), + Column("EBVErr", np.float32, "E(B-V) extinction error", np.nan), Column("T_eff", np.float32, "Effective temperature [K]", np.nan), Column("T_effErr", np.float32, "Effective temperature error [K]", np.nan), Column("M_H", np.float32, "Metallicity [dex]", np.nan), Column("M_HErr", np.float32, "Metallicity error [dex]", np.nan), + Column("a_M", np.float32, "Alpha abundance [dex]", np.nan), + Column("a_MErr", np.float32, "Alpha abundance error [dex]", np.nan), + Column("C", np.float32, "Relative carbon", np.nan), + Column("CErr", np.float32, "Relative carbon error", np.nan), Column("log_g", np.float32, "log g", np.nan), Column("log_gErr", np.float32, "log g error", np.nan), Column("flag", bool, "Measurement flag (true means bad)", False), Column("status", str, "Measurement flags", ""), ] - fitsExtName = 'GACATALOG' + fitsExtName = 'STARCATALOG' -PfsGACatalogNotes = makeNotesClass( - "PfsGACatalogNotes", +PfsStarCatalogNotes = makeNotesClass( + "PfsStarCatalogNotes", [] ) -class PfsGACatalog(): - filenameFormat = ("pfsGACatalog-%(catId)05d-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") - filenameRegex = r"^pfsGACatalog-(\d{5})-([0-9]{3})-0x([0-9a-f]{16})\.fits.*$" +class PfsStarCatalog(): + filenameFormat = ("pfsStarCatalog-%(catId)05d-%(nVisit)03d-0x%(pfsVisitHash)016x.fits") + filenameRegex = r"^pfsStarCatalog-(\d{5})-([0-9]{3})-0x([0-9a-f]{16})\.fits.*$" filenameKeys = [("catId", int), ("nVisit", int), ("pfsVisitHash", int)] - NotesClass = PfsGACatalogNotes + NotesClass = PfsStarCatalogNotes def __init__( self, catId, observations: Observations, - catalog: GACatalogTable, + catalog: StarCatalogTable, metadata=None, notes: Notes = None): @@ -472,7 +481,7 @@ def _readImpl(cls, fits): data['observations'] = None try: - catalog = GACatalogTable.readHdu(fits) + catalog = StarCatalogTable.readHdu(fits) data['catalog'] = catalog data['catId'] = catalog.catId[0] except KeyError as exc: diff --git a/tests/test_PfsGAObject.py b/tests/test_PfsStar.py similarity index 76% rename from tests/test_PfsGAObject.py rename to tests/test_PfsStar.py index 63e42464..78da8cbe 100644 --- a/tests/test_PfsGAObject.py +++ b/tests/test_PfsStar.py @@ -1,23 +1,24 @@ import os import re + import numpy as np -from unittest import TestCase + +import lsst.utils.tests from pfs.datamodel import Target, TargetType from pfs.datamodel import Observations from pfs.datamodel import MaskHelper -from pfs.datamodel import GAFluxTable -from pfs.datamodel import PfsGAObject, StellarParams, Abundances, VelocityCorrections - +from pfs.datamodel import StarFluxTable +from pfs.datamodel import PfsStar, StellarParams, Abundances, VelocityCorrections -class PfsGAObjectTestCase(TestCase): +class PfsStarTestCase(lsst.utils.tests.TestCase): """ Check the format of example datamodel files are consistent with that specified in the corresponding datamodel classes. """ - def makePfsGAObject(self): - """Construct a PfsGAObject with dummy values for testing.""" + def makePfsStar(self): + """Construct a PfsStar with dummy values for testing.""" catId = 12345 tract = 1 @@ -58,9 +59,9 @@ def makePfsGAObject(self): flags = MaskHelper() # metadata = {} # Key-value pairs to put in the header - fluxTable = GAFluxTable(wavelength, flux, error, model, cont, - norm_flux, norm_error, norm_model, - mask, flags) + fluxTable = StarFluxTable(wavelength, flux, error, model, cont, + norm_flux, norm_error, norm_model, + mask, flags) stellarParams = StellarParams( method=np.array(['tempfit', 'tempfit', 'tempfit', 'tempfit', 'tempfit']), @@ -95,7 +96,7 @@ def makePfsGAObject(self): abundCovar = np.eye(4, dtype=np.float32) notes = None - return PfsGAObject(target, observations, + return PfsStar(target, observations, wavelength, flux, mask, sky, covar, covar2, flags, metadata, fluxTable, @@ -106,7 +107,7 @@ def makePfsGAObject(self): abundCovar, notes) - def assertPfsGAObject(self, lhs, rhs): + def assertPfsStar(self, lhs, rhs): np.testing.assert_array_equal(lhs.observations.visit, rhs.observations.visit) # TODO: add more tests here @@ -136,10 +137,10 @@ def extractAttributes(self, cls, fileName): d[kk] = ii return d - def test_filenameRegex(self): + def testFilenameRegex(self): d = self.extractAttributes( - PfsGAObject, - 'pfsGAObject-07621-01234-2,2-02468ace1234abcd-003-0x0123456789abcdef.fits') + PfsStar, + 'pfsStar-07621-01234-2,2-02468ace1234abcd-003-0x0123456789abcdef.fits') self.assertEqual(d['catId'], 7621) self.assertEqual(d['tract'], 1234) self.assertEqual(d['patch'], '2,2') @@ -147,39 +148,39 @@ def test_filenameRegex(self): self.assertEqual(d['nVisit'], 3) self.assertEqual(d['pfsVisitHash'], 81985529216486895) - def test_getIdentity(self): - """Construct a PfsGAObject and get its identity.""" + def testGetIdentity(self): + """Construct a PfsStar and get its identity.""" - pfsGAObject = self.makePfsGAObject() - identity = pfsGAObject.getIdentity() - filename = pfsGAObject.filenameFormat % identity + pfsStar = self.makePfsStar() + identity = pfsStar.getIdentity() + filename = pfsStar.filenameFormat % identity - self.assertEqual('pfsGAObject-12345-00001-1,1-00000000075bcd15-002-0x05a95bc24d8ce16f.fits', filename) + self.assertEqual('pfsStar-12345-00001-1,1-00000000075bcd15-002-0x05a95bc24d8ce16f.fits', filename) - def test_validate(self): - """Construct a PfsGAObject and run validation.""" + def testValidate(self): + """Construct a PfsStar and run validation.""" - pfsGAObject = self.makePfsGAObject() - pfsGAObject.validate() + pfsStar = self.makePfsStar() + pfsStar.validate() - def test_writeFits_fromFits(self): - """Construct a PfsGAObject and save it to a FITS file.""" + def testWriteFitsFromFits(self): + """Construct a PfsStar and save it to a FITS file.""" - pfsGAObject = self.makePfsGAObject() + pfsStar = self.makePfsStar() dirName = os.path.splitext(__file__)[0] if not os.path.exists(dirName): os.makedirs(dirName) - id = pfsGAObject.getIdentity() - filename = os.path.join(dirName, pfsGAObject.filenameFormat % id) + id = pfsStar.getIdentity() + filename = os.path.join(dirName, pfsStar.filenameFormat % id) if os.path.exists(filename): os.unlink(filename) try: - pfsGAObject.writeFits(filename) - other = PfsGAObject.readFits(filename) - self.assertPfsGAObject(pfsGAObject, other) + pfsStar.writeFits(filename) + other = PfsStar.readFits(filename) + self.assertPfsStar(pfsStar, other) except Exception as e: raise e finally: diff --git a/tests/test_PfsGACatalog.py b/tests/test_PfsStarCatalog.py similarity index 65% rename from tests/test_PfsGACatalog.py rename to tests/test_PfsStarCatalog.py index f0844ee3..5dcff43d 100644 --- a/tests/test_PfsGACatalog.py +++ b/tests/test_PfsStarCatalog.py @@ -1,20 +1,21 @@ import os import re + import numpy as np -from unittest import TestCase -from pfs.datamodel import TargetType, Observations -from pfs.datamodel.ga import PfsGACatalog, GACatalogTable +import lsst.utils.tests +from pfs.datamodel import TargetType, Observations +from pfs.datamodel import PfsStarCatalog, StarCatalogTable -class PfsGACatalogTestCase(TestCase): +class PfsStarCatalogTestCase(lsst.utils.tests.TestCase): """ Check the format of example datamodel files are consistent with that specified in the corresponding datamodel classes. """ - def makePfsGACatalog(self): - """Construct a PfsGACatalog with dummy values for testing.""" + def makePfsStarCatalog(self): + """Construct a PfsStarCatalog with dummy values for testing.""" catId = 12345 @@ -31,7 +32,7 @@ def makePfsGACatalog(self): expTime=np.array([1800.0, 1800.0, 1800.0, 1800.0]), ) - catalog = GACatalogTable( + catalog = StarCatalogTable( catId=np.array([10001, 10001, 10001]), objId=np.array([1, 2, 3]), gaiaId=np.array([11, 12, 13]), @@ -40,7 +41,7 @@ def makePfsGACatalog(self): miscId=np.array([-1, -1, 33]), ra=np.array([210.01, 210.02, 210.03]), dec=np.array([67.1, 67.2, 67.3]), - epoch=np.array([2016.0, 2016.0, 2016.0]), + epoch=np.array(['J2016.0', 'J2016.0', 'J2016.0']), pmRa=np.array([0.0, 0.0, 0.0]), pmDec=np.array([0.0, 0.0, 0.0]), parallax=np.array([0.0, 0.0, 0.0]), @@ -58,24 +59,35 @@ def makePfsGACatalog(self): expTimeEff_r=np.array([5400, 5400, 5400]), expTimeEff_n=np.array([5400, 5400, 5400]), - rv=np.array([0.0, 0.0, 0.0]), - rvErr=np.array([0.0, 0.0, 0.0]), - tEff=np.array([0.0, 0.0, 0.0]), - tEffErr=np.array([0.0, 0.0, 0.0]), - m_h=np.array([0.0, 0.0, 0.0]), - m_hErr=np.array([0.0, 0.0, 0.0]), - logg=np.array([0.0, 0.0, 0.0]), - loggErr=np.array([0.0, 0.0, 0.0]), + snr_b=np.array([0.0, 0.0, 0.0]), + snr_m=np.array([0.0, 0.0, 0.0]), + snr_r=np.array([0.0, 0.0, 0.0]), + snr_n=np.array([0.0, 0.0, 0.0]), + + v_los=np.array([0.0, 0.0, 0.0]), + v_losErr=np.array([0.0, 0.0, 0.0]), + EBV=np.array([0.0, 0.0, 0.0]), + EBVErr=np.array([0.0, 0.0, 0.0]), + T_eff=np.array([0.0, 0.0, 0.0]), + T_effErr=np.array([0.0, 0.0, 0.0]), + M_H=np.array([0.0, 0.0, 0.0]), + M_HErr=np.array([0.0, 0.0, 0.0]), + a_M=np.array([0.0, 0.0, 0.0]), + a_MErr=np.array([0.0, 0.0, 0.0]), + C=np.array([0.0, 0.0, 0.0]), + CErr=np.array([0.0, 0.0, 0.0]), + log_g=np.array([0.0, 0.0, 0.0]), + log_gErr=np.array([0.0, 0.0, 0.0]), flag=np.array([False, False, False]), status=np.array(['', '', '']), ) - return PfsGACatalog( + return PfsStarCatalog( catId, observations, catalog, metadata=None, notes=None ) - def assertPfsGACatalog(self, lhs, rhs): + def assertPfsStarCatalog(self, lhs, rhs): np.testing.assert_array_equal(lhs.observations.visit, rhs.observations.visit) # TODO: add more tests here @@ -105,38 +117,38 @@ def extractAttributes(self, cls, fileName): d[kk] = ii return d - def test_filenameRegex(self): + def testFilenameRegex(self): d = self.extractAttributes( - PfsGACatalog, - 'pfsGACatalog-07621-003-0x0123456789abcdef.fits') + PfsStarCatalog, + 'pfsStarCatalog-07621-003-0x0123456789abcdef.fits') self.assertEqual(d['catId'], 7621) self.assertEqual(d['nVisit'], 3) self.assertEqual(d['pfsVisitHash'], 81985529216486895) - def test_validate(self): - """Construct a PfsGACatalog and run validation.""" + def testValidate(self): + """Construct a PfsStarCatalog and run validation.""" - pfsGACatalog = self.makePfsGACatalog() - pfsGACatalog.validate() + pfsStarCatalog = self.makePfsStarCatalog() + pfsStarCatalog.validate() - def test_writeFits_fromFits(self): - """Construct a PfsGACatalog and save it to a FITS file.""" + def testWriteFitsFromFits(self): + """Construct a PfsStarCatalog and save it to a FITS file.""" - pfsGACatalog = self.makePfsGACatalog() + pfsStarCatalog = self.makePfsStarCatalog() dirName = os.path.splitext(__file__)[0] if not os.path.exists(dirName): os.makedirs(dirName) - id = pfsGACatalog.getIdentity() - filename = os.path.join(dirName, pfsGACatalog.filenameFormat % id) + id = pfsStarCatalog.getIdentity() + filename = os.path.join(dirName, pfsStarCatalog.filenameFormat % id) if os.path.exists(filename): os.unlink(filename) try: - pfsGACatalog.writeFits(filename) - other = PfsGACatalog.readFits(filename) - self.assertPfsGACatalog(pfsGACatalog, other) + pfsStarCatalog.writeFits(filename) + other = PfsStarCatalog.readFits(filename) + self.assertPfsStarCatalog(pfsStarCatalog, other) except Exception as e: raise e finally: diff --git a/tests/test_docstrings.py b/tests/test_docstrings.py index 8b1553eb..564372e0 100644 --- a/tests/test_docstrings.py +++ b/tests/test_docstrings.py @@ -5,12 +5,12 @@ import lsst.utils.tests from pfs.datamodel.drp import PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject -from pfs.datamodel.ga import PfsGAObject +from pfs.datamodel import PfsStar, PfsStarCatalog class DocstringsTestCase(unittest.TestCase): def testDocstrings(self): - for cls in (PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject, PfsGAObject): + for cls in (PfsArm, PfsMerged, PfsReference, PfsSingle, PfsObject, PfsStar, PfsStarCatalog): for name, attr in inspect.getmembers(cls): if not hasattr(attr, "__doc__") or not attr.__doc__: continue From 83f14cc8feb8a6884d224c8f4f79c3aa02259268 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Wed, 22 Oct 2025 12:11:03 -0400 Subject: [PATCH 17/24] Added measurement flags for each algorithm we run. Added some more stellar params to the catalog. --- python/pfs/datamodel/ga.py | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index defc6680..455242f9 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -20,6 +20,7 @@ "VelocityCorrections", "StellarParams", "Abundances", + "MeasurementFlags", "StarFluxTable", "PfsStarNotes", "PfsStar", @@ -92,6 +93,16 @@ class Abundances(PfsTable): ] fitsExtName = 'ABUND' +class MeasurementFlags(PfsTable): + """List of measurement flags for a target, for each algorithm.""" + + damdVer = GA_DAMD_VER + schema = [ + Column("method", str, "Measurement method", ""), + Column("flag", bool, "Measurement flag (true means bad)", False), + Column("status", str, "Measurement flags", ""), + ] + fitsExtName = 'MEASFLAG' class StarFluxTable(FluxTable): """Table of coadded fluxes at near-original sampling and model fits @@ -258,6 +269,8 @@ class PfsStar(PfsFiberArray): Covariance matrix for stellar parameters. abundCovar: `numpy.ndarray` of `float`, optional Covariance matrix for abundance parameters. + measurementFlags: `pfs.datamodel.MeasurementFlags`, optional + Table of measurement flags for a target, for each algorithm. notes : `Notes`, optional Reduction notes. """ @@ -291,6 +304,7 @@ def __init__( abundances=None, paramsCovar=None, abundCovar=None, + measurementFlags=None, notes: Notes = None, ): super().__init__(target, observations, wavelength, flux, mask, sky, @@ -301,6 +315,7 @@ def __init__( self.abundances = abundances self.paramsCovar = paramsCovar self.abundCovar = abundCovar + self.measurementFlags = measurementFlags def validate(self): """Validate that all the arrays are of the expected shape""" @@ -321,6 +336,7 @@ def _readImpl(cls, fits): data["abundances"] = Abundances.readHdu(fits) if cls.AbundancesFitsExtName in fits: data["abundCovar"] = fits[cls.AbundancesFitsExtName].data.astype(np.float32) + data["measurementFlags"] = MeasurementFlags.readHdu(fits) return data @@ -343,6 +359,8 @@ def _writeImpl(self, fits): fits.append(ImageHDU(self.abundCovar.astype(np.float32), header=header, name=self.AbundancesFitsExtName)) + if self.measurementFlags is not None: + self.measurementFlags.writeHdu(fits) return header @@ -357,6 +375,7 @@ class StarCatalogTable(PfsTable): Column("gaiaId", np.int64, "GAIA identifier", -1), Column("ps1Id", np.int64, "PS1 identifier", -1), Column("hscId", np.int64, "HSC identifier", -1), + Column("sdssId", np.int64, "SDSS identifier", -1), Column("miscId", np.int64, "Miscellaneous identifier", -1), Column("ra", np.float32, "Right ascension ICRS [deg]", np.nan), Column("dec", np.float32, "Declination ICRS [deg]", np.nan), @@ -385,18 +404,25 @@ class StarCatalogTable(PfsTable): Column("v_los", np.float32, "Radial velocity [km/s]", np.nan), Column("v_losErr", np.float32, "Radial velocity error [km/s]", np.nan), + Column("v_losStatus", str, "Radial velocity flag", ""), Column("EBV", np.float32, "E(B-V) extinction", np.nan), Column("EBVErr", np.float32, "E(B-V) extinction error", np.nan), + Column("EBVStatus", str, "E(B-V) extinction flag", ""), Column("T_eff", np.float32, "Effective temperature [K]", np.nan), Column("T_effErr", np.float32, "Effective temperature error [K]", np.nan), + Column("T_effStatus", str, "Effective temperature flag", ""), Column("M_H", np.float32, "Metallicity [dex]", np.nan), Column("M_HErr", np.float32, "Metallicity error [dex]", np.nan), + Column("M_HStatus", str, "Metallicity flag", ""), Column("a_M", np.float32, "Alpha abundance [dex]", np.nan), Column("a_MErr", np.float32, "Alpha abundance error [dex]", np.nan), + Column("a_MStatus", str, "Alpha abundance flag", ""), Column("C", np.float32, "Relative carbon", np.nan), Column("CErr", np.float32, "Relative carbon error", np.nan), + Column("CStatus", str, "Relative carbon flag", ""), Column("log_g", np.float32, "log g", np.nan), Column("log_gErr", np.float32, "log g error", np.nan), + Column("log_gStatus", str, "log g flag", ""), Column("flag", bool, "Measurement flag (true means bad)", False), Column("status", str, "Measurement flags", ""), From 6bf062c48070cec50bc87658bded026d76dd3b23 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Mon, 3 Nov 2025 15:57:51 -0500 Subject: [PATCH 18/24] Added new flags to mark multimodal log L --- python/pfs/datamodel/ga.py | 20 +++++++++++--------- 1 file changed, 11 insertions(+), 9 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 455242f9..15f22105 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -34,15 +34,17 @@ class TempFitFlag(IntFlag): # NOTE: Keep this consistent with the flags in pfs.ga.pfsspec.tempfit.TempFitFlag - OK = 0 # No flags - BADINIT = 1 << 0 # "Bad initial values" - NOCONVERGE = 1 << 1 # "Optimization did not convergence" - BADCONVERGE = 1 << 2 # "Bad convergence" - MAXITER = 1 << 3 # "Maximum iterations reached" - PARAMEDGE = 1 << 4 # "Template parameters at the bounds" - BADERROR = 1 << 5 # "Errors could not be calculated" - BADCOV = 1 << 6 # "Covariance matrix could not be calculated" - UNLIKELYPRIOR = 1 << 7 # "Unlikely prior" + OK = 0 # No flags + BADINIT = 1 << 0 # "Bad initial values" + NOCONVERGE = 1 << 1 # "Optimization did not convergence" + BADCONVERGE = 1 << 2 # "Bad convergence" + MAXITER = 1 << 3 # "Maximum iterations reached" + PARAMEDGE = 1 << 4 # "Template parameters at the bounds" + BADERROR = 1 << 5 # "Errors could not be calculated" + BADCOV = 1 << 6 # "Covariance matrix could not be calculated" + UNLIKELYPRIOR = 1 << 7 # "Unlikely prior" + RVGUESSMULTIMODAL = 1 << 8 # "log L(RV) during guess is multimodal" + RVFITMULTIMODAL = 1 << 9 # "log L(RV) around best fit is multimodal" class VelocityCorrections(PfsTable): From d09145a31f9651740ea4e75e08e7188b76d582cd Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Mon, 3 Nov 2025 15:58:56 -0500 Subject: [PATCH 19/24] Renamed pfsGAObject to pfsStar --- datamodel.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/datamodel.txt b/datamodel.txt index 0ec5ddc3..de6092c9 100644 --- a/datamodel.txt +++ b/datamodel.txt @@ -702,7 +702,7 @@ The format will be identical to that of the pfsObject files. Physical parameter measurements for Galactic Archeology targets including flux-calibrated combined spectra and synthetic stellar template fits. - "pfsGAObject-%05d-%05d-%s-%016x-%03d-0x%016x.fits" + "pfsStar-%05d-%05d-%s-%016x-%03d-0x%016x.fits" % (catId, tract, patch, objId, nVisit % 1000, pfsVisitHash) The format is similar to pfsObject files, but with additional HDUs and a different fluxTable format. @@ -747,7 +747,7 @@ The ABUND table lists the measured single element abundances flag Measurement flag (true means bad) bool status Flags describing the quality of the measurement string -The FLUXTABLE of pfsGAObject files is extended and includes the following additional columns: +The FLUXTABLE of pfsStar files is extended and includes the following additional columns: model Best fit fluxed template model 32-bit float cont Continuum model 32-bit float From 1d613a01317cd35317bae4839d7c775e6e41c66d Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sat, 21 Feb 2026 13:46:23 -0500 Subject: [PATCH 20/24] Added fiberId and cobraId to catalog. --- python/pfs/datamodel/ga.py | 3 +++ tests/test_PfsStarCatalog.py | 4 ++++ 2 files changed, 7 insertions(+) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 15f22105..e4ccd653 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -389,6 +389,9 @@ class StarCatalogTable(PfsTable): Column("proposalId", str, "Proposal ID", ""), Column("obCode", str, "Observing Block", ""), + Column("fiberId", np.int16, "Fiber identifier", -1), + Column("cobraId", np.int16, "Cobra identifier", -1), + Column("nVisit_b", np.int16, "Number of visits in B", -1), Column("nVisit_m", np.int16, "Number of visits in M", -1), Column("nVisit_r", np.int16, "Number of visits in R", -1), diff --git a/tests/test_PfsStarCatalog.py b/tests/test_PfsStarCatalog.py index 5dcff43d..1f9ece6d 100644 --- a/tests/test_PfsStarCatalog.py +++ b/tests/test_PfsStarCatalog.py @@ -8,6 +8,7 @@ from pfs.datamodel import TargetType, Observations from pfs.datamodel import PfsStarCatalog, StarCatalogTable + class PfsStarCatalogTestCase(lsst.utils.tests.TestCase): """ Check the format of example datamodel files are consistent with that specified in the corresponding @@ -49,6 +50,9 @@ def makePfsStarCatalog(self): proposalId=np.array(['A', 'A', 'A']), obCode=np.array(['A', 'B', 'C']), + fiberId=np.array([-1, -1, -1]), + cobraId=np.array([-1, -1, -1]), + nVisit_b=np.array([3, 3, 3]), nVisit_m=np.array([3, 3, 3]), nVisit_r=np.array([3, 3, 3]), From a54f59aa1d3c7c43728e22b3daff65f98a8032a3 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sat, 21 Feb 2026 14:51:34 -0500 Subject: [PATCH 21/24] Removed cobraId which is not available in the observations collection. --- python/pfs/datamodel/ga.py | 1 - tests/test_PfsStarCatalog.py | 1 - 2 files changed, 2 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index e4ccd653..a8cea14d 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -390,7 +390,6 @@ class StarCatalogTable(PfsTable): Column("obCode", str, "Observing Block", ""), Column("fiberId", np.int16, "Fiber identifier", -1), - Column("cobraId", np.int16, "Cobra identifier", -1), Column("nVisit_b", np.int16, "Number of visits in B", -1), Column("nVisit_m", np.int16, "Number of visits in M", -1), diff --git a/tests/test_PfsStarCatalog.py b/tests/test_PfsStarCatalog.py index 1f9ece6d..c0aa04e5 100644 --- a/tests/test_PfsStarCatalog.py +++ b/tests/test_PfsStarCatalog.py @@ -51,7 +51,6 @@ def makePfsStarCatalog(self): obCode=np.array(['A', 'B', 'C']), fiberId=np.array([-1, -1, -1]), - cobraId=np.array([-1, -1, -1]), nVisit_b=np.array([3, 3, 3]), nVisit_m=np.array([3, 3, 3]), From 73ff1b37c0b6930737829e02d52a4addae1bd6e2 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sun, 22 Feb 2026 06:31:56 -0500 Subject: [PATCH 22/24] Add new flags and parameters to PfsStar and StarCatalogTable; update tests accordingly --- python/pfs/datamodel/ga.py | 36 +++++++++++++++++++++--------------- tests/test_PfsStar.py | 20 ++++++++++++-------- tests/test_PfsStarCatalog.py | 12 ++++++++++++ 3 files changed, 45 insertions(+), 23 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index a8cea14d..0659214c 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -45,6 +45,7 @@ class TempFitFlag(IntFlag): UNLIKELYPRIOR = 1 << 7 # "Unlikely prior" RVGUESSMULTIMODAL = 1 << 8 # "log L(RV) during guess is multimodal" RVFITMULTIMODAL = 1 << 9 # "log L(RV) around best fit is multimodal" + NODATA = 1 << 10 # No data to fit, possibly missing 2d output class VelocityCorrections(PfsTable): @@ -285,8 +286,8 @@ class PfsStar(PfsFiberArray): NotesClass = PfsStarNotes FluxTableClass = StarFluxTable - StellarParamsFitsExtName = "STELLARCOVAR" - AbundancesFitsExtName = "ABUNDCOVAR" + StellarParamsCovarFitsExtName = "STELLARCOVAR" + AbundancesCovarFitsExtName = "ABUNDCOVAR" def __init__( self, @@ -329,16 +330,18 @@ def validate(self): def _readImpl(cls, fits): data = super()._readImpl(fits) - # TODO: handle missing extensions - - data["velocityCorrections"] = VelocityCorrections.readHdu(fits) - data["stellarParams"] = StellarParams.readHdu(fits) - if cls.StellarParamsFitsExtName in fits: - data["paramsCovar"] = fits[cls.StellarParamsFitsExtName].data.astype(np.float32) - data["abundances"] = Abundances.readHdu(fits) - if cls.AbundancesFitsExtName in fits: - data["abundCovar"] = fits[cls.AbundancesFitsExtName].data.astype(np.float32) - data["measurementFlags"] = MeasurementFlags.readHdu(fits) + if VelocityCorrections.fitsExtName in fits: + data["velocityCorrections"] = VelocityCorrections.readHdu(fits) + if StellarParams.fitsExtName in fits: + data["stellarParams"] = StellarParams.readHdu(fits) + if cls.StellarParamsCovarFitsExtName in fits: + data["paramsCovar"] = fits[cls.StellarParamsCovarFitsExtName].data.astype(np.float32) + if Abundances.fitsExtName in fits: + data["abundances"] = Abundances.readHdu(fits) + if cls.AbundancesCovarFitsExtName in fits: + data["abundCovar"] = fits[cls.AbundancesCovarFitsExtName].data.astype(np.float32) + if MeasurementFlags.fitsExtName in fits: + data["measurementFlags"] = MeasurementFlags.readHdu(fits) return data @@ -354,13 +357,13 @@ def _writeImpl(self, fits): if self.paramsCovar is not None: fits.append(ImageHDU(self.paramsCovar.astype(np.float32), header=header, - name=self.StellarParamsFitsExtName)) + name=self.StellarParamsCovarFitsExtName)) if self.abundances is not None: self.abundances.writeHdu(fits) if self.abundCovar is not None: fits.append(ImageHDU(self.abundCovar.astype(np.float32), header=header, - name=self.AbundancesFitsExtName)) + name=self.AbundancesCovarFitsExtName)) if self.measurementFlags is not None: self.measurementFlags.writeHdu(fits) @@ -373,7 +376,7 @@ class StarCatalogTable(PfsTable): damdVer = GA_DAMD_VER schema = [ Column("catId", np.int32, "PFS catalog identifier", -1), - Column("objId", np.int32, "Object identifier.", -1), + Column("objId", np.int32, "GA object identifier.", -1), Column("gaiaId", np.int64, "GAIA identifier", -1), Column("ps1Id", np.int64, "PS1 identifier", -1), Column("hscId", np.int64, "HSC identifier", -1), @@ -386,10 +389,13 @@ class StarCatalogTable(PfsTable): Column("pmDec", np.float32, "Proper motion dec [mas/yr]", np.nan), Column("parallax", np.float32, "Parallax [mas]", np.nan), Column("targetType", np.int16, "Target type.", -1), + Column("targetPriority", np.int16, "Target priority.", -1), Column("proposalId", str, "Proposal ID", ""), Column("obCode", str, "Observing Block", ""), Column("fiberId", np.int16, "Fiber identifier", -1), + Column("nVisit", np.int16, "Number of visits", -1), + Column("pfsVisitHash", np.int64, "PFS visit hash", -1), Column("nVisit_b", np.int16, "Number of visits in B", -1), Column("nVisit_m", np.int16, "Number of visits in M", -1), diff --git a/tests/test_PfsStar.py b/tests/test_PfsStar.py index 78da8cbe..2ccc5485 100644 --- a/tests/test_PfsStar.py +++ b/tests/test_PfsStar.py @@ -177,11 +177,15 @@ def testWriteFitsFromFits(self): if os.path.exists(filename): os.unlink(filename) - try: - pfsStar.writeFits(filename) - other = PfsStar.readFits(filename) - self.assertPfsStar(pfsStar, other) - except Exception as e: - raise e - finally: - os.unlink(filename) + pfsStar.writeFits(filename) + other = PfsStar.readFits(filename) + self.assertPfsStar(pfsStar, other) + + # try: + # pfsStar.writeFits(filename) + # other = PfsStar.readFits(filename) + # self.assertPfsStar(pfsStar, other) + # except Exception as e: + # raise e + # finally: + # os.unlink(filename) diff --git a/tests/test_PfsStarCatalog.py b/tests/test_PfsStarCatalog.py index c0aa04e5..c6bc70b4 100644 --- a/tests/test_PfsStarCatalog.py +++ b/tests/test_PfsStarCatalog.py @@ -7,6 +7,7 @@ from pfs.datamodel import TargetType, Observations from pfs.datamodel import PfsStarCatalog, StarCatalogTable +from pfs.datamodel.pfsTable import Column class PfsStarCatalogTestCase(lsst.utils.tests.TestCase): @@ -39,6 +40,7 @@ def makePfsStarCatalog(self): gaiaId=np.array([11, 12, 13]), ps1Id=np.array([21, 22, 23]), hscId=np.array([21, 22, 23]), + sdssId=np.array([21, 22, 23]), miscId=np.array([-1, -1, 33]), ra=np.array([210.01, 210.02, 210.03]), dec=np.array([67.1, 67.2, 67.3]), @@ -47,10 +49,13 @@ def makePfsStarCatalog(self): pmDec=np.array([0.0, 0.0, 0.0]), parallax=np.array([0.0, 0.0, 0.0]), targetType=np.array([TargetType.SCIENCE, TargetType.SCIENCE, TargetType.SCIENCE]), + targetPriority=np.array([1, 1, 1]), proposalId=np.array(['A', 'A', 'A']), obCode=np.array(['A', 'B', 'C']), fiberId=np.array([-1, -1, -1]), + nVisit=np.array([3, 3, 3]), + pfsVisitHash=np.array([81985529216486895, 81985529216486895, 81985529216486895]), nVisit_b=np.array([3, 3, 3]), nVisit_m=np.array([3, 3, 3]), @@ -69,18 +74,25 @@ def makePfsStarCatalog(self): v_los=np.array([0.0, 0.0, 0.0]), v_losErr=np.array([0.0, 0.0, 0.0]), + v_losStatus=np.array(['', '', '']), EBV=np.array([0.0, 0.0, 0.0]), EBVErr=np.array([0.0, 0.0, 0.0]), + EBVStatus=np.array(['', '', '']), T_eff=np.array([0.0, 0.0, 0.0]), T_effErr=np.array([0.0, 0.0, 0.0]), + T_effStatus=np.array(['', '', '']), M_H=np.array([0.0, 0.0, 0.0]), M_HErr=np.array([0.0, 0.0, 0.0]), + M_HStatus=np.array(['', '', '']), a_M=np.array([0.0, 0.0, 0.0]), a_MErr=np.array([0.0, 0.0, 0.0]), + a_MStatus=np.array(['', '', '']), C=np.array([0.0, 0.0, 0.0]), CErr=np.array([0.0, 0.0, 0.0]), + CStatus=np.array(['', '', '']), log_g=np.array([0.0, 0.0, 0.0]), log_gErr=np.array([0.0, 0.0, 0.0]), + log_gStatus=np.array(['', '', '']), flag=np.array([False, False, False]), status=np.array(['', '', '']), From fa8224c7c1a20afbbafec242c916b9c1e1d91f66 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Sun, 22 Feb 2026 12:11:23 -0500 Subject: [PATCH 23/24] Removed survey IDs --- python/pfs/datamodel/ga.py | 5 ----- tests/test_PfsStarCatalog.py | 5 ----- 2 files changed, 10 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 0659214c..59d09594 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -377,11 +377,6 @@ class StarCatalogTable(PfsTable): schema = [ Column("catId", np.int32, "PFS catalog identifier", -1), Column("objId", np.int32, "GA object identifier.", -1), - Column("gaiaId", np.int64, "GAIA identifier", -1), - Column("ps1Id", np.int64, "PS1 identifier", -1), - Column("hscId", np.int64, "HSC identifier", -1), - Column("sdssId", np.int64, "SDSS identifier", -1), - Column("miscId", np.int64, "Miscellaneous identifier", -1), Column("ra", np.float32, "Right ascension ICRS [deg]", np.nan), Column("dec", np.float32, "Declination ICRS [deg]", np.nan), Column("epoch", str, "Coordinate epoch [Jyr]", ''), diff --git a/tests/test_PfsStarCatalog.py b/tests/test_PfsStarCatalog.py index c6bc70b4..bb0ce07d 100644 --- a/tests/test_PfsStarCatalog.py +++ b/tests/test_PfsStarCatalog.py @@ -37,11 +37,6 @@ def makePfsStarCatalog(self): catalog = StarCatalogTable( catId=np.array([10001, 10001, 10001]), objId=np.array([1, 2, 3]), - gaiaId=np.array([11, 12, 13]), - ps1Id=np.array([21, 22, 23]), - hscId=np.array([21, 22, 23]), - sdssId=np.array([21, 22, 23]), - miscId=np.array([-1, -1, 33]), ra=np.array([210.01, 210.02, 210.03]), dec=np.array([67.1, 67.2, 67.3]), epoch=np.array(['J2016.0', 'J2016.0', 'J2016.0']), From 1a97b142c43079fc2df75c4f27a4cffda7d78b58 Mon Sep 17 00:00:00 2001 From: Laszlo Dobos Date: Fri, 6 Mar 2026 14:10:55 -0500 Subject: [PATCH 24/24] Added spectrograph, renamed flags. --- python/pfs/datamodel/ga.py | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/python/pfs/datamodel/ga.py b/python/pfs/datamodel/ga.py index 59d09594..1b37b2aa 100644 --- a/python/pfs/datamodel/ga.py +++ b/python/pfs/datamodel/ga.py @@ -389,6 +389,7 @@ class StarCatalogTable(PfsTable): Column("obCode", str, "Observing Block", ""), Column("fiberId", np.int16, "Fiber identifier", -1), + Column("spectrograph", np.int16, "Spectrograph identifier", -1), Column("nVisit", np.int16, "Number of visits", -1), Column("pfsVisitHash", np.int64, "PFS visit hash", -1), @@ -429,8 +430,8 @@ class StarCatalogTable(PfsTable): Column("log_gErr", np.float32, "log g error", np.nan), Column("log_gStatus", str, "log g flag", ""), - Column("flag", bool, "Measurement flag (true means bad)", False), - Column("status", str, "Measurement flags", ""), + Column("tempfitflag", bool, "Measurement flag (true means bad)", False), + Column("tempfitstatus", str, "Measurement flags", ""), ] fitsExtName = 'STARCATALOG'