24 from itertools
import product
29 import astropy.units
as u
32 from lsst.daf.butler
import DeferredDatasetHandle
36 from .parquetTable
import ParquetTable, MultilevelParquetTable
40 typeKey='functor', name=None):
41 """Initialize an object defined in a dictionary
43 The object needs to be importable as
44 f'{basePath}.{initDict[typeKey]}'
45 The positional and keyword arguments (if any) are contained in
46 "args" and "kwargs" entries in the dictionary, respectively.
47 This is used in `functors.CompositeFunctor.from_yaml` to initialize
48 a composite functor from a specification in a YAML file.
53 Dictionary describing object's initialization. Must contain
54 an entry keyed by ``typeKey`` that is the name of the object,
55 relative to ``basePath``.
57 Path relative to module in which ``initDict[typeKey]`` is defined.
59 Key of ``initDict`` that is the name of the object
60 (relative to `basePath`).
62 initDict = initDict.copy()
64 pythonType = doImport(f
'{basePath}.{initDict.pop(typeKey)}')
66 if 'args' in initDict:
67 args = initDict.pop(
'args')
68 if isinstance(args, str):
71 element = pythonType(*args, **initDict)
72 except Exception
as e:
73 message = f
'Error in constructing functor "{name}" of type {pythonType.__name__} with args: {args}'
74 raise type(e)(message, e.args)
79 """Define and execute a calculation on a ParquetTable
81 The `__call__` method accepts either a `ParquetTable` object or a
82 `DeferredDatasetHandle`, and returns the
83 result of the calculation as a single column. Each functor defines what
84 columns are needed for the calculation, and only these columns are read
85 from the `ParquetTable`.
87 The action of `__call__` consists of two steps: first, loading the
88 necessary columns from disk into memory as a `pandas.DataFrame` object;
89 and second, performing the computation on this dataframe and returning the
93 To define a new `Functor`, a subclass must define a `_func` method,
94 that takes a `pandas.DataFrame` and returns result in a `pandas.Series`.
95 In addition, it must define the following attributes
97 * `_columns`: The columns necessary to perform the calculation
98 * `name`: A name appropriate for a figure axis label
99 * `shortname`: A name appropriate for use as a dictionary key
101 On initialization, a `Functor` should declare what band (`filt` kwarg)
102 and dataset (e.g. `'ref'`, `'meas'`, `'forced_src'`) it is intended to be
103 applied to. This enables the `_get_data` method to extract the proper
104 columns from the parquet file. If not specified, the dataset will fall back
105 on the `_defaultDataset`attribute. If band is not specified and `dataset`
106 is anything other than `'ref'`, then an error will be raised when trying to
107 perform the calculation.
109 Originally, `Functor` was set up to expect
110 datasets formatted like the `deepCoadd_obj` dataset; that is, a
111 dataframe with a multi-level column index, with the levels of the
112 column index being `band`, `dataset`, and `column`.
113 It has since been generalized to apply to dataframes without mutli-level
114 indices and multi-level indices with just `dataset` and `column` levels.
115 In addition, the `_get_data` method that reads
116 the dataframe from the `ParquetTable` will return a dataframe with column
117 index levels defined by the `_dfLevels` attribute; by default, this is
120 The `_dfLevels` attributes should generally not need to
121 be changed, unless `_func` needs columns from multiple filters or datasets
122 to do the calculation.
123 An example of this is the `lsst.pipe.tasks.functors.Color` functor, for
124 which `_dfLevels = ('band', 'column')`, and `_func` expects the dataframe
125 it gets to have those levels in the column index.
130 Filter upon which to do the calculation
133 Dataset upon which to do the calculation
134 (e.g., 'ref', 'meas', 'forced_src').
138 _defaultDataset =
'ref'
139 _dfLevels = (
'column',)
140 _defaultNoDup =
False
142 def __init__(self, filt=None, dataset=None, noDup=None):
149 if self.
_noDup_noDup
is not None:
156 """Columns required to perform calculation
158 if not hasattr(self,
'_columns'):
159 raise NotImplementedError(
'Must define columns property or _columns attribute')
162 def _get_data_columnLevels(self, data, columnIndex=None):
163 """Gets the names of the column index levels
165 This should only be called in the context of a multilevel table.
166 The logic here is to enable this to work both with the gen2 `MultilevelParquetTable`
167 and with the gen3 `DeferredDatasetHandle`.
171 data : `MultilevelParquetTable` or `DeferredDatasetHandle`
173 columnnIndex (optional): pandas `Index` object
174 if not passed, then it is read from the `DeferredDatasetHandle`
176 if isinstance(data, DeferredDatasetHandle):
177 if columnIndex
is None:
178 columnIndex = data.get(component=
"columns")
179 if columnIndex
is not None:
180 return columnIndex.names
181 if isinstance(data, MultilevelParquetTable):
182 return data.columnLevels
184 raise TypeError(f
"Unknown type for data: {type(data)}!")
186 def _get_data_columnLevelNames(self, data, columnIndex=None):
187 """Gets the content of each of the column levels for a multilevel table
189 Similar to `_get_data_columnLevels`, this enables backward compatibility with gen2.
191 Mirrors original gen2 implementation within `pipe.tasks.parquetTable.MultilevelParquetTable`
193 if isinstance(data, DeferredDatasetHandle):
194 if columnIndex
is None:
195 columnIndex = data.get(component=
"columns")
196 if columnIndex
is not None:
197 columnLevels = columnIndex.names
199 level: list(np.unique(np.array([c
for c
in columnIndex])[:, i]))
200 for i, level
in enumerate(columnLevels)
202 return columnLevelNames
203 if isinstance(data, MultilevelParquetTable):
204 return data.columnLevelNames
206 raise TypeError(f
"Unknown type for data: {type(data)}!")
208 def _colsFromDict(self, colDict, columnIndex=None):
209 """Converts dictionary column specficiation to a list of columns
211 This mirrors the original gen2 implementation within `pipe.tasks.parquetTable.MultilevelParquetTable`
216 for i, lev
in enumerate(columnLevels):
218 if isinstance(colDict[lev], str):
219 new_colDict[lev] = [colDict[lev]]
221 new_colDict[lev] = colDict[lev]
223 new_colDict[lev] = columnIndex.levels[i]
225 levelCols = [new_colDict[lev]
for lev
in columnLevels]
226 cols = product(*levelCols)
230 """Returns columns needed by functor from multilevel dataset
232 To access tables with multilevel column structure, the `MultilevelParquetTable`
233 or `DeferredDatasetHandle` need to be passed either a list of tuples or a
238 data : `MultilevelParquetTable` or `DeferredDatasetHandle`
240 columnIndex (optional): pandas `Index` object
241 either passed or read in from `DeferredDatasetHandle`.
244 If true, then return a list of tuples rather than the column dictionary
245 specification. This is set to `True` by `CompositeFunctor` in order to be able to
246 combine columns from the various component functors.
249 if isinstance(data, DeferredDatasetHandle)
and columnIndex
is None:
250 columnIndex = data.get(component=
"columns")
255 columnDict = {
'column': self.
columnscolumns,
256 'dataset': self.
datasetdataset}
257 if self.
filtfilt
is None:
259 if "band" in columnLevels:
260 if self.
datasetdataset ==
"ref":
261 columnDict[
"band"] = columnLevelNames[
"band"][0]
263 raise ValueError(f
"'filt' not set for functor {self.name}"
264 f
"(dataset {self.dataset}) "
266 "contains multiple filters in column index. "
267 "Set 'filt' or set 'dataset' to 'ref'.")
269 columnDict[
'band'] = self.
filtfilt
271 if isinstance(data, MultilevelParquetTable):
272 return data._colsFromDict(columnDict)
273 elif isinstance(data, DeferredDatasetHandle):
275 return self.
_colsFromDict_colsFromDict(columnDict, columnIndex=columnIndex)
279 def _func(self, df, dropna=True):
280 raise NotImplementedError(
'Must define calculation on dataframe')
282 def _get_columnIndex(self, data):
283 """Return columnIndex
286 if isinstance(data, DeferredDatasetHandle):
287 return data.get(component=
"columns")
291 def _get_data(self, data):
292 """Retrieve dataframe necessary for calculation.
294 The data argument can be a DataFrame, a ParquetTable instance, or a gen3 DeferredDatasetHandle
296 Returns dataframe upon which `self._func` can act.
298 N.B. while passing a raw pandas `DataFrame` *should* work here, it has not been tested.
300 if isinstance(data, pd.DataFrame):
305 is_multiLevel = isinstance(data, MultilevelParquetTable)
or isinstance(columnIndex, pd.MultiIndex)
308 if isinstance(data, ParquetTable)
and not is_multiLevel:
310 df = data.toDataFrame(columns=columns)
319 if isinstance(data, MultilevelParquetTable):
321 df = data.toDataFrame(columns=columns, droplevels=
False)
322 elif isinstance(data, DeferredDatasetHandle):
324 df = data.get(parameters={
"columns": columns})
332 def _setLevels(self, df):
333 levelsToDrop = [n
for n
in df.columns.names
if n
not in self.
_dfLevels_dfLevels]
334 df.columns = df.columns.droplevel(levelsToDrop)
337 def _dropna(self, vals):
343 vals = self.
_func_func(df)
345 vals = self.
failfail(df)
347 vals = self.
_dropna_dropna(vals)
352 """Computes difference between functor called on two different ParquetTable objects
354 return self(data1, **kwargs) - self(data2, **kwargs)
357 return pd.Series(np.full(len(df), np.nan), index=df.index)
361 """Full name of functor (suitable for figure labels)
363 return NotImplementedError
367 """Short name of functor (suitable for column name/dict key)
373 """Perform multiple calculations at once on a catalog
375 The role of a `CompositeFunctor` is to group together computations from
376 multiple functors. Instead of returning `pandas.Series` a
377 `CompositeFunctor` returns a `pandas.Dataframe`, with the column names
378 being the keys of `funcDict`.
380 The `columns` attribute of a `CompositeFunctor` is the union of all columns
381 in all the component functors.
383 A `CompositeFunctor` does not use a `_func` method itself; rather,
384 when a `CompositeFunctor` is called, all its columns are loaded
385 at once, and the resulting dataframe is passed to the `_func` method of each component
386 functor. This has the advantage of only doing I/O (reading from parquet file) once,
387 and works because each individual `_func` method of each component functor does not
388 care if there are *extra* columns in the dataframe being passed; only that it must contain
389 *at least* the `columns` it expects.
391 An important and useful class method is `from_yaml`, which takes as argument the path to a YAML
392 file specifying a collection of functors.
396 funcs : `dict` or `list`
397 Dictionary or list of functors. If a list, then it will be converted
398 into a dictonary according to the `.shortname` attribute of each functor.
405 if type(funcs) == dict:
408 self.
funcDictfuncDict = {f.shortname: f
for f
in funcs}
410 self.
_filt_filt =
None
416 return self.
_filt_filt
421 for _, f
in self.
funcDictfuncDict.items():
423 self.
_filt_filt = filt
426 if isinstance(new, dict):
428 elif isinstance(new, CompositeFunctor):
431 raise TypeError(
'Can only update with dictionary or CompositeFunctor.')
439 return list(set([x
for y
in [f.columns
for f
in self.
funcDictfuncDict.values()]
for x
in y]))
448 f.multilevelColumns(data, returnTuple=
True, **kwargs)
for f
in self.
funcDictfuncDict.values()
456 """Apply the functor to the data table
460 data : `lsst.daf.butler.DeferredDatasetHandle`,
461 `lsst.pipe.tasks.parquetTable.MultilevelParquetTable`,
462 `lsst.pipe.tasks.parquetTable.ParquetTable`,
463 or `pandas.DataFrame`.
464 The table or a pointer to a table on disk from which columns can
470 is_multiLevel = isinstance(data, MultilevelParquetTable)
or isinstance(columnIndex, pd.MultiIndex)
476 if isinstance(data, MultilevelParquetTable):
478 df = data.toDataFrame(columns=columns, droplevels=
False)
479 elif isinstance(data, DeferredDatasetHandle):
481 df = data.get(parameters={
"columns": columns})
484 for k, f
in self.
funcDictfuncDict.items():
486 subdf = f._setLevels(
487 df[f.multilevelColumns(data, returnTuple=
True, columnIndex=columnIndex)]
489 valDict[k] = f._func(subdf)
490 except Exception
as e:
492 valDict[k] = f.fail(subdf)
497 if isinstance(data, DeferredDatasetHandle):
500 elif isinstance(data, pd.DataFrame):
507 valDict = {k: f._func(df)
for k, f
in self.
funcDictfuncDict.items()}
510 valDf = pd.concat(valDict, axis=1)
512 print([(k, type(v))
for k, v
in valDict.items()])
515 if kwargs.get(
'dropna',
False):
516 valDf = valDf.dropna(how=
'any')
522 if renameRules
is None:
524 for old, new
in renameRules:
525 if col.startswith(old):
526 col = col.replace(old, new)
532 filename = os.path.expandvars(filename)
533 with open(filename)
as f:
534 translationDefinition = yaml.safe_load(f)
536 return cls.
from_yamlfrom_yaml(translationDefinition, **kwargs)
541 for func, val
in translationDefinition[
'funcs'].items():
544 if 'flag_rename_rules' in translationDefinition:
545 renameRules = translationDefinition[
'flag_rename_rules']
549 if 'calexpFlags' in translationDefinition:
550 for flag
in translationDefinition[
'calexpFlags']:
551 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'calexp')
553 if 'refFlags' in translationDefinition:
554 for flag
in translationDefinition[
'refFlags']:
555 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'ref')
557 if 'forcedFlags' in translationDefinition:
558 for flag
in translationDefinition[
'forcedFlags']:
559 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'forced_src')
561 if 'flags' in translationDefinition:
562 for flag
in translationDefinition[
'flags']:
563 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'meas')
565 return cls(funcs, **kwargs)
569 """Evaluate an expression on a DataFrame, knowing what the 'mag' function means
571 Builds on `pandas.DataFrame.eval`, which parses and executes math on dataframes.
575 df : pandas.DataFrame
576 Dataframe on which to evaluate expression.
582 expr_new = re.sub(
r'mag\((\w+)\)',
r'-2.5*log(\g<1>)/log(10)', expr)
583 val = df.eval(expr_new, truediv=
True)
585 expr_new = re.sub(
r'mag\((\w+)\)',
r'-2.5*log(\g<1>_instFlux)/log(10)', expr)
586 val = df.eval(expr_new, truediv=
True)
591 """Arbitrary computation on a catalog
593 Column names (and thus the columns to be loaded from catalog) are found
594 by finding all words and trying to ignore all "math-y" words.
599 Expression to evaluate, to be parsed and executed by `mag_aware_eval`.
601 _ignore_words = (
'mag',
'sin',
'cos',
'exp',
'log',
'sqrt')
613 flux_cols = re.findall(
r'mag\(\s*(\w+)\s*\)', self.
exprexpr)
615 cols = [c
for c
in re.findall(
r'[a-zA-Z_]+', self.
exprexpr)
if c
not in self.
_ignore_words_ignore_words]
618 if not re.search(
'_instFlux$', c):
619 cols.append(f
'{c}_instFlux')
624 return list(set([c
for c
in cols
if c
not in not_a_col]))
631 """Get column with specified name
647 return df[self.
colcol]
651 """Return the value of the index for each object
654 columns = [
'coord_ra']
655 _defaultDataset =
'ref'
659 return pd.Series(df.index, index=df.index)
664 _allow_difference =
False
668 return pd.Series(df.index, index=df.index)
672 col =
'base_Footprint_nPix'
676 """Base class for coordinate column, in degrees
685 output = df[self.
colcol] * 180 / np.pi
if self.
_radians_radians
else df[self.
colcol]
690 """Right Ascension, in degrees
696 super().
__init__(
'coord_ra', **kwargs)
699 return super().
__call__(catalog, **kwargs)
703 """Declination, in degrees
709 super().
__init__(
'coord_dec', **kwargs)
712 return super().
__call__(catalog, **kwargs)
716 """Compute the level 20 HtmIndex for the catalog.
720 This functor was implemented to satisfy requirements of old APDB interface
721 which required ``pixelId`` column in DiaObject with HTM20 index. APDB
722 interface had migrated to not need that information, but we keep this
723 class in case it may be useful for something else.
738 def computePixel(row):
747 return self.
pixelatorpixelator.index(sphPoint.getVector())
749 return df.apply(computePixel, axis=1)
753 if not col.endswith(
'_instFlux'):
759 if not col.endswith(
'_instFluxErr'):
760 col +=
'_instFluxErr'
765 """Compute calibrated magnitude
767 Takes a `calib` argument, which returns the flux at mag=0
768 as `calib.getFluxMag0()`. If not provided, then the default
769 `fluxMag0` is 63095734448.0194, which is default for HSC.
770 This default should be removed in DM-21955
772 This calculation hides warnings about invalid values and dividing by zero.
774 As for all functors, a `dataset` and `filt` kwarg should be provided upon
775 initialization. Unlike the default `Functor`, however, the default dataset
776 for a `Mag` is `'meas'`, rather than `'ref'`.
781 Name of flux column from which to compute magnitude. Can be parseable
782 by `lsst.pipe.tasks.functors.fluxName` function---that is, you can pass
783 `'modelfit_CModel'` instead of `'modelfit_CModel_instFlux'`) and it will
785 calib : `lsst.afw.image.calib.Calib` (optional)
786 Object that knows zero point.
788 _defaultDataset =
'meas'
793 if calib
is not None:
797 self.
fluxMag0fluxMag0 = 63095734448.0194
806 with np.warnings.catch_warnings():
807 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
808 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
809 return -2.5*np.log10(df[self.
colcol] / self.
fluxMag0fluxMag0)
813 return f
'mag_{self.col}'
817 """Compute calibrated magnitude uncertainty
819 Takes the same `calib` object as `lsst.pipe.tasks.functors.Mag`.
824 calib : `lsst.afw.image.calib.Calib` (optional)
825 Object that knows zero point.
830 if self.
calibcalib
is not None:
837 return [self.
colcol, self.
colcol +
'Err']
840 with np.warnings.catch_warnings():
841 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
842 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
844 x = df[fluxErrCol] / df[fluxCol]
846 magErr = (2.5 / np.log(10.)) * np.sqrt(x*x + y*y)
851 return super().name +
'_err'
859 return (df[self.
colcol] / self.
fluxMag0fluxMag0) * 1e9
863 _defaultDataset =
'meas'
865 """Functor to calculate magnitude difference"""
874 return [self.
col1col1, self.
col2col2]
877 with np.warnings.catch_warnings():
878 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
879 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
880 return -2.5*np.log10(df[self.
col1col1]/df[self.
col2col2])
884 return f
'(mag_{self.col1} - mag_{self.col2})'
888 return f
'magDiff_{self.col1}_{self.col2}'
892 """Compute the color between two filters
894 Computes color by initializing two different `Mag`
895 functors based on the `col` and filters provided, and
896 then returning the difference.
898 This is enabled by the `_func` expecting a dataframe with a
899 multilevel column index, with both `'band'` and `'column'`,
900 instead of just `'column'`, which is the `Functor` default.
901 This is controlled by the `_dfLevels` attribute.
903 Also of note, the default dataset for `Color` is `forced_src'`,
904 whereas for `Mag` it is `'meas'`.
909 Name of flux column from which to compute; same as would be passed to
910 `lsst.pipe.tasks.functors.Mag`.
913 Filters from which to compute magnitude difference.
914 Color computed is `Mag(filt2) - Mag(filt1)`.
916 _defaultDataset =
'forced_src'
917 _dfLevels = (
'band',
'column')
923 raise RuntimeError(
"Cannot compute Color for %s: %s - %s " % (col, filt2, filt1))
927 self.
mag2mag2 =
Mag(col, filt=filt2, **kwargs)
928 self.
mag1mag1 =
Mag(col, filt=filt1, **kwargs)
941 mag2 = self.mag2._func(df[self.filt2])
942 mag1 = self.mag1._func(df[self.filt1])
947 return [self.
mag1mag1.col, self.
mag2mag2.col]
954 return f
'{self.filt2} - {self.filt1} ({self.col})'
958 return f
"{self.col}_{self.filt2.replace('-', '')}m{self.filt1.replace('-', '')}"
962 """Main function of this subclass is to override the dropna=True
965 _allow_difference =
False
970 return super().
__call__(parq, dropna=
False, **kwargs)
974 _columns = [
"base_ClassificationExtendedness_value"]
975 _column =
"base_ClassificationExtendedness_value"
980 test = (x < 0.5).astype(int)
981 test = test.mask(mask, 2)
985 categories = [
'galaxy',
'star', self.
_null_label_null_label]
986 label = pd.Series(pd.Categorical.from_codes(test, categories=categories),
987 index=x.index, name=
'label')
989 label = label.astype(str)
994 _columns = [
'numStarFlags']
995 labels = {
"star": 0,
"maybe": 1,
"notStar": 2}
1001 n = len(x.unique()) - 1
1003 labels = [
'noStar',
'maybe',
'star']
1004 label = pd.Series(pd.cut(x, [-1, 0, n-1, n], labels=labels),
1005 index=x.index, name=
'label')
1008 label = label.astype(str)
1014 name =
'Deconvolved Moments'
1015 shortname =
'deconvolvedMoments'
1016 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
1017 "ext_shapeHSM_HsmSourceMoments_yy",
1018 "base_SdssShape_xx",
"base_SdssShape_yy",
1019 "ext_shapeHSM_HsmPsfMoments_xx",
1020 "ext_shapeHSM_HsmPsfMoments_yy")
1022 def _func(self, df):
1023 """Calculate deconvolved moments"""
1024 if "ext_shapeHSM_HsmSourceMoments_xx" in df.columns:
1025 hsm = df[
"ext_shapeHSM_HsmSourceMoments_xx"] + df[
"ext_shapeHSM_HsmSourceMoments_yy"]
1027 hsm = np.ones(len(df))*np.nan
1028 sdss = df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]
1029 if "ext_shapeHSM_HsmPsfMoments_xx" in df.columns:
1030 psf = df[
"ext_shapeHSM_HsmPsfMoments_xx"] + df[
"ext_shapeHSM_HsmPsfMoments_yy"]
1035 raise RuntimeError(
'No psf shape parameter found in catalog')
1037 return hsm.where(np.isfinite(hsm), sdss) - psf
1041 """Functor to calculate SDSS trace radius size for sources"""
1042 name =
"SDSS Trace Size"
1043 shortname =
'sdssTrace'
1044 _columns = (
"base_SdssShape_xx",
"base_SdssShape_yy")
1046 def _func(self, df):
1047 srcSize = np.sqrt(0.5*(df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]))
1052 """Functor to calculate SDSS trace radius size difference (%) between object and psf model"""
1053 name =
"PSF - SDSS Trace Size"
1054 shortname =
'psf_sdssTrace'
1055 _columns = (
"base_SdssShape_xx",
"base_SdssShape_yy",
1056 "base_SdssShape_psf_xx",
"base_SdssShape_psf_yy")
1058 def _func(self, df):
1059 srcSize = np.sqrt(0.5*(df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]))
1060 psfSize = np.sqrt(0.5*(df[
"base_SdssShape_psf_xx"] + df[
"base_SdssShape_psf_yy"]))
1061 sizeDiff = 100*(srcSize - psfSize)/(0.5*(srcSize + psfSize))
1066 """Functor to calculate HSM trace radius size for sources"""
1067 name =
'HSM Trace Size'
1068 shortname =
'hsmTrace'
1069 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
1070 "ext_shapeHSM_HsmSourceMoments_yy")
1072 def _func(self, df):
1073 srcSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmSourceMoments_xx"]
1074 + df[
"ext_shapeHSM_HsmSourceMoments_yy"]))
1079 """Functor to calculate HSM trace radius size difference (%) between object and psf model"""
1080 name =
'PSF - HSM Trace Size'
1081 shortname =
'psf_HsmTrace'
1082 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
1083 "ext_shapeHSM_HsmSourceMoments_yy",
1084 "ext_shapeHSM_HsmPsfMoments_xx",
1085 "ext_shapeHSM_HsmPsfMoments_yy")
1087 def _func(self, df):
1088 srcSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmSourceMoments_xx"]
1089 + df[
"ext_shapeHSM_HsmSourceMoments_yy"]))
1090 psfSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmPsfMoments_xx"]
1091 + df[
"ext_shapeHSM_HsmPsfMoments_yy"]))
1092 sizeDiff = 100*(srcSize - psfSize)/(0.5*(srcSize + psfSize))
1097 name =
'HSM Psf FWHM'
1098 _columns = (
'ext_shapeHSM_HsmPsfMoments_xx',
'ext_shapeHSM_HsmPsfMoments_yy')
1101 SIGMA2FWHM = 2*np.sqrt(2*np.log(2))
1103 def _func(self, df):
1105 0.5*(df[
'ext_shapeHSM_HsmPsfMoments_xx'] + df[
'ext_shapeHSM_HsmPsfMoments_yy']))
1109 name =
"Distortion Ellipticity (e1)"
1110 shortname =
"Distortion"
1123 def _func(self, df):
1128 name =
"Ellipticity e2"
1140 def _func(self, df):
1141 return 2*df[self.
colXYcolXY] / (df[self.
colXXcolXX] + df[self.
colYYcolYY])
1156 def _func(self, df):
1157 return (df[self.
colXXcolXX]*df[self.
colYYcolYY] - df[self.
colXYcolXY]**2)**0.25
1161 """Computations using the stored localWcs.
1163 name =
"LocalWcsOperations"
1178 """Compute the distance on the sphere from x2, y1 to x1, y1.
1186 cd11 : `pandas.Series`
1187 [1, 1] element of the local Wcs affine transform.
1188 cd11 : `pandas.Series`
1189 [1, 1] element of the local Wcs affine transform.
1190 cd12 : `pandas.Series`
1191 [1, 2] element of the local Wcs affine transform.
1192 cd21 : `pandas.Series`
1193 [2, 1] element of the local Wcs affine transform.
1194 cd22 : `pandas.Series`
1195 [2, 2] element of the local Wcs affine transform.
1200 RA and dec conversion of x and y given the local Wcs. Returned
1201 units are in radians.
1204 return (x * cd11 + y * cd12, x * cd21 + y * cd22)
1207 """Compute the local pixel scale conversion.
1211 ra1 : `pandas.Series`
1212 Ra of the first coordinate in radians.
1213 dec1 : `pandas.Series`
1214 Dec of the first coordinate in radians.
1215 ra2 : `pandas.Series`
1216 Ra of the second coordinate in radians.
1217 dec2 : `pandas.Series`
1218 Dec of the second coordinate in radians.
1222 dist : `pandas.Series`
1223 Distance on the sphere in radians.
1225 deltaDec = dec2 - dec1
1227 return 2 * np.arcsin(
1229 np.sin(deltaDec / 2) ** 2
1230 + np.cos(dec2) * np.cos(dec1) * np.sin(deltaRa / 2) ** 2))
1233 """Compute the distance on the sphere from x2, y1 to x1, y1.
1237 x1 : `pandas.Series`
1239 y1 : `pandas.Series`
1241 x2 : `pandas.Series`
1243 y2 : `pandas.Series`
1245 cd11 : `pandas.Series`
1246 [1, 1] element of the local Wcs affine transform.
1247 cd11 : `pandas.Series`
1248 [1, 1] element of the local Wcs affine transform.
1249 cd12 : `pandas.Series`
1250 [1, 2] element of the local Wcs affine transform.
1251 cd21 : `pandas.Series`
1252 [2, 1] element of the local Wcs affine transform.
1253 cd22 : `pandas.Series`
1254 [2, 2] element of the local Wcs affine transform.
1258 Distance : `pandas.Series`
1259 Arcseconds per pixel at the location of the local WC
1261 ra1, dec1 = self.
computeDeltaRaDeccomputeDeltaRaDec(x1, y1, cd11, cd12, cd21, cd22)
1262 ra2, dec2 = self.
computeDeltaRaDeccomputeDeltaRaDec(x2, y2, cd11, cd12, cd21, cd22)
1268 """Compute the local pixel scale from the stored CDMatrix.
1280 """Compute the local pixel to scale conversion in arcseconds.
1284 cd11 : `pandas.Series`
1285 [1, 1] element of the local Wcs affine transform in radians.
1286 cd11 : `pandas.Series`
1287 [1, 1] element of the local Wcs affine transform in radians.
1288 cd12 : `pandas.Series`
1289 [1, 2] element of the local Wcs affine transform in radians.
1290 cd21 : `pandas.Series`
1291 [2, 1] element of the local Wcs affine transform in radians.
1292 cd22 : `pandas.Series`
1293 [2, 2] element of the local Wcs affine transform in radians.
1297 pixScale : `pandas.Series`
1298 Arcseconds per pixel at the location of the local WC
1300 return 3600 * np.degrees(np.sqrt(np.fabs(cd11 * cd22 - cd12 * cd21)))
1302 def _func(self, df):
1310 """Convert a value in units pixels squared to units arcseconds squared.
1329 return f
"{self.col}_asArcseconds"
1333 return [self.
colcol,
1339 def _func(self, df):
1347 """Convert a value in units pixels to units arcseconds.
1366 return f
"{self.col}_asArcsecondsSq"
1370 return [self.
colcol,
1376 def _func(self, df):
1381 return df[self.
colcol] * pixScale * pixScale
1385 name =
'Reference Band'
1386 shortname =
'refBand'
1390 return [
"merge_measurement_i",
1391 "merge_measurement_r",
1392 "merge_measurement_z",
1393 "merge_measurement_y",
1394 "merge_measurement_g"]
1396 def _func(self, df):
1397 def getFilterAliasName(row):
1399 colName = row.idxmax()
1400 return colName.replace(
'merge_measurement_',
'')
1402 return df[self.
columnscolumnscolumns].apply(getFilterAliasName, axis=1)
1407 AB_FLUX_SCALE = (0 * u.ABmag).to_value(u.nJy)
1408 LOG_AB_FLUX_SCALE = 12.56
1409 FIVE_OVER_2LOG10 = 1.085736204758129569
1413 def __init__(self, colFlux, colFluxErr=None, calib=None, **kwargs):
1419 if calib
is not None:
1429 return [self.
colcol]
1433 return f
'mag_{self.col}'
1437 if np.abs(a) < np.abs(b):
1442 return np.abs(a) * np.sqrt(1. + q*q)
1448 with np.warnings.catch_warnings():
1449 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
1450 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
1451 return -2.5 * np.log10(dn/fluxMag0)
1454 retVal = self.
vhypotvhypot(dn * fluxMag0Err, dnErr * fluxMag0)
1455 retVal *= self.
AB_FLUX_SCALEAB_FLUX_SCALE / fluxMag0 / fluxMag0
1459 retVal = self.
dn2fluxErrdn2fluxErr(dn, dnErr, fluxMag0, fluxMag0Err) / self.
dn2fluxdn2flux(dn, fluxMag0)
1464 def _func(self, df):
1473 def _func(self, df):
1475 return pd.Series(retArr, index=df.index)
1479 def _func(self, df):
1488 def _func(self, df):
1490 return pd.Series(retArr, index=df.index)
1494 """Base class for calibrating the specified instrument flux column using
1495 the local photometric calibration.
1500 Name of the instrument flux column.
1501 instFluxErrCol : `str`
1502 Name of the assocated error columns for ``instFluxCol``.
1503 photoCalibCol : `str`
1504 Name of local calibration column.
1505 photoCalibErrCol : `str`
1506 Error associated with ``photoCalibCol``
1516 logNJanskyToAB = (1 * u.nJy).to_value(u.ABmag)
1531 """Convert instrument flux to nanojanskys.
1535 instFlux : `numpy.ndarray` or `pandas.Series`
1536 Array of instrument flux measurements
1537 localCalib : `numpy.ndarray` or `pandas.Series`
1538 Array of local photometric calibration estimates.
1542 calibFlux : `numpy.ndarray` or `pandas.Series`
1543 Array of calibrated flux measurements.
1545 return instFlux * localCalib
1548 """Convert instrument flux to nanojanskys.
1552 instFlux : `numpy.ndarray` or `pandas.Series`
1553 Array of instrument flux measurements
1554 instFluxErr : `numpy.ndarray` or `pandas.Series`
1555 Errors on associated ``instFlux`` values
1556 localCalib : `numpy.ndarray` or `pandas.Series`
1557 Array of local photometric calibration estimates.
1558 localCalibErr : `numpy.ndarray` or `pandas.Series`
1559 Errors on associated ``localCalib`` values
1563 calibFluxErr : `numpy.ndarray` or `pandas.Series`
1564 Errors on calibrated flux measurements.
1566 return np.hypot(instFluxErr * localCalib, instFlux * localCalibErr)
1569 """Convert instrument flux to nanojanskys.
1573 instFlux : `numpy.ndarray` or `pandas.Series`
1574 Array of instrument flux measurements
1575 localCalib : `numpy.ndarray` or `pandas.Series`
1576 Array of local photometric calibration estimates.
1580 calibMag : `numpy.ndarray` or `pandas.Series`
1581 Array of calibrated AB magnitudes.
1586 """Convert instrument flux err to nanojanskys.
1590 instFlux : `numpy.ndarray` or `pandas.Series`
1591 Array of instrument flux measurements
1592 instFluxErr : `numpy.ndarray` or `pandas.Series`
1593 Errors on associated ``instFlux`` values
1594 localCalib : `numpy.ndarray` or `pandas.Series`
1595 Array of local photometric calibration estimates.
1596 localCalibErr : `numpy.ndarray` or `pandas.Series`
1597 Errors on associated ``localCalib`` values
1601 calibMagErr: `numpy.ndarray` or `pandas.Series`
1602 Error on calibrated AB magnitudes.
1605 return 2.5 / np.log(10) * err / self.
instFluxToNanojanskyinstFluxToNanojansky(instFlux, instFluxErr)
1609 """Compute calibrated fluxes using the local calibration value.
1625 return f
'flux_{self.instFluxCol}'
1627 def _func(self, df):
1632 """Compute calibrated flux errors using the local calibration value.
1649 return f
'fluxErr_{self.instFluxCol}'
1651 def _func(self, df):
1657 """Compute calibrated AB magnitudes using the local calibration value.
1673 return f
'mag_{self.instFluxCol}'
1675 def _func(self, df):
1681 """Compute calibrated AB magnitude errors using the local calibration value.
1698 return f
'magErr_{self.instFluxCol}'
1700 def _func(self, df):
1708 """Compute absolute mean of dipole fluxes.
1717 LocalDipoleMeanFluxErr
1719 LocalDipoleDiffFluxErr
1749 return f
'dipMeanFlux_{self.instFluxPosCol}_{self.instFluxNegCol}'
1751 def _func(self, df):
1757 """Compute the error on the absolute mean of dipole fluxes.
1766 LocalDipoleMeanFluxErr
1768 LocalDipoleDiffFluxErr
1782 return f
'dipMeanFluxErr_{self.instFluxPosCol}_{self.instFluxNegCol}'
1784 def _func(self, df):
1793 """Compute the absolute difference of dipole fluxes.
1795 Value is (abs(pos) - abs(neg))
1804 LocalDipoleMeanFluxErr
1806 LocalDipoleDiffFluxErr
1817 return f
'dipDiffFlux_{self.instFluxPosCol}_{self.instFluxNegCol}'
1819 def _func(self, df):
1825 """Compute the error on the absolute difference of dipole fluxes.
1834 LocalDipoleMeanFluxErr
1836 LocalDipoleDiffFluxErr
1850 return f
'dipDiffFluxErr_{self.instFluxPosCol}_{self.instFluxNegCol}'
1852 def _func(self, df):
1861 """Base class for returning the ratio of 2 columns.
1863 Can be used to compute a Signal to Noise ratio for any input flux.
1868 Name of the column to use at the numerator in the ratio
1870 Name of the column to use as the denominator in the ratio.
1886 return f
'ratio_{self.numerator}_{self.denominator}'
1888 def _func(self, df):
1889 with np.warnings.catch_warnings():
1890 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
1891 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
def multilevelColumns(self, parq, **kwargs)
def __init__(self, col, filt2, filt1, **kwargs)
def __init__(self, col, **kwargs)
def __init__(self, funcs, **kwargs)
def __call__(self, data, **kwargs)
def from_file(cls, filename, **kwargs)
def from_yaml(cls, translationDefinition, **kwargs)
def renameCol(cls, col, renameRules)
def multilevelColumns(self, data, **kwargs)
def pixelScaleArcseconds(self, cd11, cd12, cd21, cd22)
def __init__(self, col, colCD_1_1, colCD_1_2, colCD_2_1, colCD_2_2, **kwargs)
def __init__(self, col, colCD_1_1, colCD_1_2, colCD_2_1, colCD_2_2, **kwargs)
def __init__(self, col, **kwargs)
def __init__(self, expr, **kwargs)
def __init__(self, **kwargs)
def __call__(self, catalog, **kwargs)
def __init__(self, colXX, colXY, colYY, **kwargs)
def __init__(self, colXX, colXY, colYY, **kwargs)
def __call__(self, data, dropna=False)
def _get_data(self, data)
def _func(self, df, dropna=True)
def multilevelColumns(self, data, columnIndex=None, returnTuple=False)
def _get_data_columnLevelNames(self, data, columnIndex=None)
def difference(self, data1, data2, **kwargs)
def __init__(self, filt=None, dataset=None, noDup=None)
def _get_columnIndex(self, data)
def _colsFromDict(self, colDict, columnIndex=None)
def _get_data_columnLevels(self, data, columnIndex=None)
def __init__(self, ra, decl, **kwargs)
def __call__(self, parq, dropna=False, **kwargs)
def __init__(self, instFluxPosCol, instFluxNegCol, instFluxPosErrCol, instFluxNegErrCol, photoCalibCol, photoCalibErrCol, **kwargs)
def instFluxToNanojansky(self, instFlux, localCalib)
def instFluxErrToMagnitudeErr(self, instFlux, instFluxErr, localCalib, localCalibErr)
def __init__(self, instFluxCol, instFluxErrCol, photoCalibCol, photoCalibErrCol, **kwargs)
def instFluxErrToNanojanskyErr(self, instFlux, instFluxErr, localCalib, localCalibErr)
def instFluxToMagnitude(self, instFlux, localCalib)
def __init__(self, colCD_1_1, colCD_1_2, colCD_2_1, colCD_2_2, **kwargs)
def computeDeltaRaDec(self, x, y, cd11, cd12, cd21, cd22)
def computeSkySeperation(self, ra1, dec1, ra2, dec2)
def getSkySeperationFromPixel(self, x1, y1, x2, y2, cd11, cd12, cd21, cd22)
def __init__(self, col1, col2, **kwargs)
def __init__(self, *args, **kwargs)
def __init__(self, col, calib=None, **kwargs)
def dn2mag(self, dn, fluxMag0)
def dn2flux(self, dn, fluxMag0)
def dn2fluxErr(self, dn, dnErr, fluxMag0, fluxMag0Err)
def dn2MagErr(self, dn, dnErr, fluxMag0, fluxMag0Err)
def __init__(self, colFlux, colFluxErr=None, calib=None, **kwargs)
def __call__(self, catalog, **kwargs)
def __init__(self, **kwargs)
def __init__(self, colXX, colXY, colYY, **kwargs)
def __init__(self, numerator, denominator, **kwargs)
def mag_aware_eval(df, expr)
def init_fromDict(initDict, basePath='lsst.pipe.tasks.functors', typeKey='functor', name=None)