3 from itertools
import product
7 import astropy.units
as u
10 from lsst.daf.butler
import DeferredDatasetHandle
11 from .parquetTable
import ParquetTable, MultilevelParquetTable
15 typeKey='functor', name=None):
16 """Initialize an object defined in a dictionary
18 The object needs to be importable as
19 f'{basePath}.{initDict[typeKey]}'
20 The positional and keyword arguments (if any) are contained in
21 "args" and "kwargs" entries in the dictionary, respectively.
22 This is used in `functors.CompositeFunctor.from_yaml` to initialize
23 a composite functor from a specification in a YAML file.
28 Dictionary describing object's initialization. Must contain
29 an entry keyed by ``typeKey`` that is the name of the object,
30 relative to ``basePath``.
32 Path relative to module in which ``initDict[typeKey]`` is defined.
34 Key of ``initDict`` that is the name of the object
35 (relative to `basePath`).
37 initDict = initDict.copy()
39 pythonType = doImport(f
'{basePath}.{initDict.pop(typeKey)}')
41 if 'args' in initDict:
42 args = initDict.pop(
'args')
43 if isinstance(args, str):
46 element = pythonType(*args, **initDict)
47 except Exception
as e:
48 message = f
'Error in constructing functor "{name}" of type {pythonType.__name__} with args: {args}'
49 raise type(e)(message, e.args)
54 """Define and execute a calculation on a ParquetTable
56 The `__call__` method accepts either a `ParquetTable` object or a
57 `DeferredDatasetHandle`, and returns the
58 result of the calculation as a single column. Each functor defines what
59 columns are needed for the calculation, and only these columns are read
60 from the `ParquetTable`.
62 The action of `__call__` consists of two steps: first, loading the
63 necessary columns from disk into memory as a `pandas.DataFrame` object;
64 and second, performing the computation on this dataframe and returning the
68 To define a new `Functor`, a subclass must define a `_func` method,
69 that takes a `pandas.DataFrame` and returns result in a `pandas.Series`.
70 In addition, it must define the following attributes
72 * `_columns`: The columns necessary to perform the calculation
73 * `name`: A name appropriate for a figure axis label
74 * `shortname`: A name appropriate for use as a dictionary key
76 On initialization, a `Functor` should declare what band (`filt` kwarg)
77 and dataset (e.g. `'ref'`, `'meas'`, `'forced_src'`) it is intended to be
78 applied to. This enables the `_get_data` method to extract the proper
79 columns from the parquet file. If not specified, the dataset will fall back
80 on the `_defaultDataset`attribute. If band is not specified and `dataset`
81 is anything other than `'ref'`, then an error will be raised when trying to
82 perform the calculation.
84 As currently implemented, `Functor` is only set up to expect a
85 dataset of the format of the `deepCoadd_obj` dataset; that is, a
86 dataframe with a multi-level column index,
87 with the levels of the column index being `band`,
88 `dataset`, and `column`. This is defined in the `_columnLevels` attribute,
89 as well as being implicit in the role of the `filt` and `dataset` attributes
90 defined at initialization. In addition, the `_get_data` method that reads
91 the dataframe from the `ParquetTable` will return a dataframe with column
92 index levels defined by the `_dfLevels` attribute; by default, this is
95 The `_columnLevels` and `_dfLevels` attributes should generally not need to
96 be changed, unless `_func` needs columns from multiple filters or datasets
97 to do the calculation.
98 An example of this is the `lsst.pipe.tasks.functors.Color` functor, for
99 which `_dfLevels = ('band', 'column')`, and `_func` expects the dataframe
100 it gets to have those levels in the column index.
105 Filter upon which to do the calculation
108 Dataset upon which to do the calculation
109 (e.g., 'ref', 'meas', 'forced_src').
113 _defaultDataset =
'ref'
114 _columnLevels = (
'band',
'dataset',
'column')
115 _dfLevels = (
'column',)
116 _defaultNoDup =
False
118 def __init__(self, filt=None, dataset=None, noDup=None):
125 if self.
_noDup_noDup
is not None:
132 """Columns required to perform calculation
134 if not hasattr(self,
'_columns'):
135 raise NotImplementedError(
'Must define columns property or _columns attribute')
138 def _get_data_columnLevels(self, data, columnIndex=None):
139 """Gets the names of the column index levels
141 This should only be called in the context of a multilevel table.
142 The logic here is to enable this to work both with the gen2 `MultilevelParquetTable`
143 and with the gen3 `DeferredDatasetHandle`.
147 data : `MultilevelParquetTable` or `DeferredDatasetHandle`
149 columnnIndex (optional): pandas `Index` object
150 if not passed, then it is read from the `DeferredDatasetHandle`
152 if isinstance(data, DeferredDatasetHandle):
153 if columnIndex
is None:
154 columnIndex = data.get(component=
"columns")
155 if columnIndex
is not None:
156 return columnIndex.names
157 if isinstance(data, MultilevelParquetTable):
158 return data.columnLevels
160 raise TypeError(f
"Unknown type for data: {type(data)}!")
162 def _get_data_columnLevelNames(self, data, columnIndex=None):
163 """Gets the content of each of the column levels for a multilevel table
165 Similar to `_get_data_columnLevels`, this enables backward compatibility with gen2.
167 Mirrors original gen2 implementation within `pipe.tasks.parquetTable.MultilevelParquetTable`
169 if isinstance(data, DeferredDatasetHandle):
170 if columnIndex
is None:
171 columnIndex = data.get(component=
"columns")
172 if columnIndex
is not None:
173 columnLevels = columnIndex.names
175 level: list(np.unique(np.array([c
for c
in columnIndex])[:, i]))
176 for i, level
in enumerate(columnLevels)
178 return columnLevelNames
179 if isinstance(data, MultilevelParquetTable):
180 return data.columnLevelNames
182 raise TypeError(f
"Unknown type for data: {type(data)}!")
184 def _colsFromDict(self, colDict, columnIndex=None):
185 """Converts dictionary column specficiation to a list of columns
187 This mirrors the original gen2 implementation within `pipe.tasks.parquetTable.MultilevelParquetTable`
192 for i, lev
in enumerate(columnLevels):
194 if isinstance(colDict[lev], str):
195 new_colDict[lev] = [colDict[lev]]
197 new_colDict[lev] = colDict[lev]
199 new_colDict[lev] = columnIndex.levels[i]
201 levelCols = [new_colDict[lev]
for lev
in columnLevels]
202 cols = product(*levelCols)
206 """Returns columns needed by functor from multilevel dataset
208 To access tables with multilevel column structure, the `MultilevelParquetTable`
209 or `DeferredDatasetHandle` need to be passed either a list of tuples or a
214 data : `MultilevelParquetTable` or `DeferredDatasetHandle`
216 columnIndex (optional): pandas `Index` object
217 either passed or read in from `DeferredDatasetHandle`.
220 If true, then return a list of tuples rather than the column dictionary
221 specification. This is set to `True` by `CompositeFunctor` in order to be able to
222 combine columns from the various component functors.
225 if isinstance(data, DeferredDatasetHandle)
and columnIndex
is None:
226 columnIndex = data.get(component=
"columns")
231 if not set(columnLevels) == set(self.
_columnLevels_columnLevels):
233 "ParquetTable does not have the expected column levels. "
234 f
"Got {columnLevels}; expected {self._columnLevels}."
237 columnDict = {
'column': self.
columnscolumns,
238 'dataset': self.
datasetdataset}
239 if self.
filtfilt
is None:
241 if "band" in columnLevels:
242 if self.
datasetdataset ==
"ref":
243 columnDict[
"band"] = columnLevelNames[
"band"][0]
245 raise ValueError(f
"'filt' not set for functor {self.name}"
246 f
"(dataset {self.dataset}) "
248 "contains multiple filters in column index. "
249 "Set 'filt' or set 'dataset' to 'ref'.")
251 columnDict[
'band'] = self.
filtfilt
253 if isinstance(data, MultilevelParquetTable):
254 return data._colsFromDict(columnDict)
255 elif isinstance(data, DeferredDatasetHandle):
257 return self.
_colsFromDict_colsFromDict(columnDict, columnIndex=columnIndex)
261 def _func(self, df, dropna=True):
262 raise NotImplementedError(
'Must define calculation on dataframe')
264 def _get_columnIndex(self, data):
265 """Return columnIndex
268 if isinstance(data, DeferredDatasetHandle):
269 return data.get(component=
"columns")
273 def _get_data(self, data):
274 """Retrieve dataframe necessary for calculation.
276 The data argument can be a DataFrame, a ParquetTable instance, or a gen3 DeferredDatasetHandle
278 Returns dataframe upon which `self._func` can act.
280 N.B. while passing a raw pandas `DataFrame` *should* work here, it has not been tested.
282 if isinstance(data, pd.DataFrame):
287 is_multiLevel = isinstance(data, MultilevelParquetTable)
or isinstance(columnIndex, pd.MultiIndex)
290 if isinstance(data, ParquetTable)
and not is_multiLevel:
292 df = data.toDataFrame(columns=columns)
299 if isinstance(data, MultilevelParquetTable):
301 df = data.toDataFrame(columns=columns, droplevels=
False)
302 elif isinstance(data, DeferredDatasetHandle):
304 df = data.get(parameters={
"columns": columns})
310 def _setLevels(self, df):
311 levelsToDrop = [n
for n
in df.columns.names
if n
not in self.
_dfLevels_dfLevels]
312 df.columns = df.columns.droplevel(levelsToDrop)
315 def _dropna(self, vals):
321 vals = self.
_func_func(df)
323 vals = self.
failfail(df)
325 vals = self.
_dropna_dropna(vals)
330 """Computes difference between functor called on two different ParquetTable objects
332 return self(data1, **kwargs) - self(data2, **kwargs)
335 return pd.Series(np.full(len(df), np.nan), index=df.index)
339 """Full name of functor (suitable for figure labels)
341 return NotImplementedError
345 """Short name of functor (suitable for column name/dict key)
351 """Perform multiple calculations at once on a catalog
353 The role of a `CompositeFunctor` is to group together computations from
354 multiple functors. Instead of returning `pandas.Series` a
355 `CompositeFunctor` returns a `pandas.Dataframe`, with the column names
356 being the keys of `funcDict`.
358 The `columns` attribute of a `CompositeFunctor` is the union of all columns
359 in all the component functors.
361 A `CompositeFunctor` does not use a `_func` method itself; rather,
362 when a `CompositeFunctor` is called, all its columns are loaded
363 at once, and the resulting dataframe is passed to the `_func` method of each component
364 functor. This has the advantage of only doing I/O (reading from parquet file) once,
365 and works because each individual `_func` method of each component functor does not
366 care if there are *extra* columns in the dataframe being passed; only that it must contain
367 *at least* the `columns` it expects.
369 An important and useful class method is `from_yaml`, which takes as argument the path to a YAML
370 file specifying a collection of functors.
374 funcs : `dict` or `list`
375 Dictionary or list of functors. If a list, then it will be converted
376 into a dictonary according to the `.shortname` attribute of each functor.
383 if type(funcs) == dict:
386 self.
funcDictfuncDict = {f.shortname: f
for f
in funcs}
388 self.
_filt_filt =
None
394 return self.
_filt_filt
399 for _, f
in self.
funcDictfuncDict.items():
401 self.
_filt_filt = filt
404 if isinstance(new, dict):
406 elif isinstance(new, CompositeFunctor):
409 raise TypeError(
'Can only update with dictionary or CompositeFunctor.')
417 return list(set([x
for y
in [f.columns
for f
in self.
funcDictfuncDict.values()]
for x
in y]))
426 f.multilevelColumns(data, returnTuple=
True, **kwargs)
for f
in self.
funcDictfuncDict.values()
434 """Apply the functor to the data table
438 data : `lsst.daf.butler.DeferredDatasetHandle`,
439 `lsst.pipe.tasks.parquetTable.MultilevelParquetTable`,
440 `lsst.pipe.tasks.parquetTable.ParquetTable`,
441 or `pandas.DataFrame`.
442 The table or a pointer to a table on disk from which columns can
448 is_multiLevel = isinstance(data, MultilevelParquetTable)
or isinstance(columnIndex, pd.MultiIndex)
454 if isinstance(data, MultilevelParquetTable):
456 df = data.toDataFrame(columns=columns, droplevels=
False)
457 elif isinstance(data, DeferredDatasetHandle):
459 df = data.get(parameters={
"columns": columns})
462 for k, f
in self.
funcDictfuncDict.items():
464 subdf = f._setLevels(
465 df[f.multilevelColumns(data, returnTuple=
True, columnIndex=columnIndex)]
467 valDict[k] = f._func(subdf)
469 valDict[k] = f.fail(subdf)
472 if isinstance(data, DeferredDatasetHandle):
475 elif isinstance(data, pd.DataFrame):
482 valDict = {k: f._func(df)
for k, f
in self.
funcDictfuncDict.items()}
485 valDf = pd.concat(valDict, axis=1)
487 print([(k, type(v))
for k, v
in valDict.items()])
490 if kwargs.get(
'dropna',
False):
491 valDf = valDf.dropna(how=
'any')
497 if renameRules
is None:
499 for old, new
in renameRules:
500 if col.startswith(old):
501 col = col.replace(old, new)
506 with open(filename)
as f:
507 translationDefinition = yaml.safe_load(f)
509 return cls.
from_yamlfrom_yaml(translationDefinition, **kwargs)
514 for func, val
in translationDefinition[
'funcs'].items():
517 if 'flag_rename_rules' in translationDefinition:
518 renameRules = translationDefinition[
'flag_rename_rules']
522 if 'refFlags' in translationDefinition:
523 for flag
in translationDefinition[
'refFlags']:
524 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'ref')
526 if 'flags' in translationDefinition:
527 for flag
in translationDefinition[
'flags']:
528 funcs[cls.
renameColrenameCol(flag, renameRules)] =
Column(flag, dataset=
'meas')
530 return cls(funcs, **kwargs)
534 """Evaluate an expression on a DataFrame, knowing what the 'mag' function means
536 Builds on `pandas.DataFrame.eval`, which parses and executes math on dataframes.
540 df : pandas.DataFrame
541 Dataframe on which to evaluate expression.
547 expr_new = re.sub(
r'mag\((\w+)\)',
r'-2.5*log(\g<1>)/log(10)', expr)
548 val = df.eval(expr_new, truediv=
True)
550 expr_new = re.sub(
r'mag\((\w+)\)',
r'-2.5*log(\g<1>_instFlux)/log(10)', expr)
551 val = df.eval(expr_new, truediv=
True)
556 """Arbitrary computation on a catalog
558 Column names (and thus the columns to be loaded from catalog) are found
559 by finding all words and trying to ignore all "math-y" words.
564 Expression to evaluate, to be parsed and executed by `mag_aware_eval`.
566 _ignore_words = (
'mag',
'sin',
'cos',
'exp',
'log',
'sqrt')
578 flux_cols = re.findall(
r'mag\(\s*(\w+)\s*\)', self.
exprexpr)
580 cols = [c
for c
in re.findall(
r'[a-zA-Z_]+', self.
exprexpr)
if c
not in self.
_ignore_words_ignore_words]
583 if not re.search(
'_instFlux$', c):
584 cols.append(f
'{c}_instFlux')
589 return list(set([c
for c
in cols
if c
not in not_a_col]))
596 """Get column with specified name
612 return df[self.
colcol]
616 """Return the value of the index for each object
619 columns = [
'coord_ra']
620 _defaultDataset =
'ref'
624 return pd.Series(df.index, index=df.index)
629 _allow_difference =
False
633 return pd.Series(df.index, index=df.index)
637 col =
'base_Footprint_nPix'
641 """Base class for coordinate column, in degrees
650 output = df[self.
colcol] * 180 / np.pi
if self.
_radians_radians
else df[self.
colcol]
655 """Right Ascension, in degrees
661 super().
__init__(
'coord_ra', **kwargs)
664 return super().
__call__(catalog, **kwargs)
668 """Declination, in degrees
674 super().
__init__(
'coord_dec', **kwargs)
677 return super().
__call__(catalog, **kwargs)
681 if not col.endswith(
'_instFlux'):
687 if not col.endswith(
'_instFluxErr'):
688 col +=
'_instFluxErr'
693 """Compute calibrated magnitude
695 Takes a `calib` argument, which returns the flux at mag=0
696 as `calib.getFluxMag0()`. If not provided, then the default
697 `fluxMag0` is 63095734448.0194, which is default for HSC.
698 This default should be removed in DM-21955
700 This calculation hides warnings about invalid values and dividing by zero.
702 As for all functors, a `dataset` and `filt` kwarg should be provided upon
703 initialization. Unlike the default `Functor`, however, the default dataset
704 for a `Mag` is `'meas'`, rather than `'ref'`.
709 Name of flux column from which to compute magnitude. Can be parseable
710 by `lsst.pipe.tasks.functors.fluxName` function---that is, you can pass
711 `'modelfit_CModel'` instead of `'modelfit_CModel_instFlux'`) and it will
713 calib : `lsst.afw.image.calib.Calib` (optional)
714 Object that knows zero point.
716 _defaultDataset =
'meas'
721 if calib
is not None:
725 self.
fluxMag0fluxMag0 = 63095734448.0194
734 with np.warnings.catch_warnings():
735 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
736 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
737 return -2.5*np.log10(df[self.
colcol] / self.
fluxMag0fluxMag0)
741 return f
'mag_{self.col}'
745 """Compute calibrated magnitude uncertainty
747 Takes the same `calib` object as `lsst.pipe.tasks.functors.Mag`.
752 calib : `lsst.afw.image.calib.Calib` (optional)
753 Object that knows zero point.
758 if self.
calibcalib
is not None:
765 return [self.
colcol, self.
colcol +
'Err']
768 with np.warnings.catch_warnings():
769 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
770 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
772 x = df[fluxErrCol] / df[fluxCol]
774 magErr = (2.5 / np.log(10.)) * np.sqrt(x*x + y*y)
779 return super().name +
'_err'
787 return (df[self.
colcol] / self.
fluxMag0fluxMag0) * 1e9
791 _defaultDataset =
'meas'
793 """Functor to calculate magnitude difference"""
802 return [self.
col1col1, self.
col2col2]
805 with np.warnings.catch_warnings():
806 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
807 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
808 return -2.5*np.log10(df[self.
col1col1]/df[self.
col2col2])
812 return f
'(mag_{self.col1} - mag_{self.col2})'
816 return f
'magDiff_{self.col1}_{self.col2}'
820 """Compute the color between two filters
822 Computes color by initializing two different `Mag`
823 functors based on the `col` and filters provided, and
824 then returning the difference.
826 This is enabled by the `_func` expecting a dataframe with a
827 multilevel column index, with both `'band'` and `'column'`,
828 instead of just `'column'`, which is the `Functor` default.
829 This is controlled by the `_dfLevels` attribute.
831 Also of note, the default dataset for `Color` is `forced_src'`,
832 whereas for `Mag` it is `'meas'`.
837 Name of flux column from which to compute; same as would be passed to
838 `lsst.pipe.tasks.functors.Mag`.
841 Filters from which to compute magnitude difference.
842 Color computed is `Mag(filt2) - Mag(filt1)`.
844 _defaultDataset =
'forced_src'
845 _dfLevels = (
'band',
'column')
851 raise RuntimeError(
"Cannot compute Color for %s: %s - %s " % (col, filt2, filt1))
855 self.
mag2mag2 =
Mag(col, filt=filt2, **kwargs)
856 self.
mag1mag1 =
Mag(col, filt=filt1, **kwargs)
869 mag2 = self.mag2._func(df[self.filt2])
870 mag1 = self.mag1._func(df[self.filt1])
875 return [self.
mag1mag1.col, self.
mag2mag2.col]
882 return f
'{self.filt2} - {self.filt1} ({self.col})'
886 return f
"{self.col}_{self.filt2.replace('-', '')}m{self.filt1.replace('-', '')}"
890 """Main function of this subclass is to override the dropna=True
893 _allow_difference =
False
898 return super().
__call__(parq, dropna=
False, **kwargs)
902 _columns = [
"base_ClassificationExtendedness_value"]
903 _column =
"base_ClassificationExtendedness_value"
908 test = (x < 0.5).astype(int)
909 test = test.mask(mask, 2)
913 categories = [
'galaxy',
'star', self.
_null_label_null_label]
914 label = pd.Series(pd.Categorical.from_codes(test, categories=categories),
915 index=x.index, name=
'label')
917 label = label.astype(str)
922 _columns = [
'numStarFlags']
923 labels = {
"star": 0,
"maybe": 1,
"notStar": 2}
929 n = len(x.unique()) - 1
931 labels = [
'noStar',
'maybe',
'star']
932 label = pd.Series(pd.cut(x, [-1, 0, n-1, n], labels=labels),
933 index=x.index, name=
'label')
936 label = label.astype(str)
942 name =
'Deconvolved Moments'
943 shortname =
'deconvolvedMoments'
944 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
945 "ext_shapeHSM_HsmSourceMoments_yy",
946 "base_SdssShape_xx",
"base_SdssShape_yy",
947 "ext_shapeHSM_HsmPsfMoments_xx",
948 "ext_shapeHSM_HsmPsfMoments_yy")
951 """Calculate deconvolved moments"""
952 if "ext_shapeHSM_HsmSourceMoments_xx" in df.columns:
953 hsm = df[
"ext_shapeHSM_HsmSourceMoments_xx"] + df[
"ext_shapeHSM_HsmSourceMoments_yy"]
955 hsm = np.ones(len(df))*np.nan
956 sdss = df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]
957 if "ext_shapeHSM_HsmPsfMoments_xx" in df.columns:
958 psf = df[
"ext_shapeHSM_HsmPsfMoments_xx"] + df[
"ext_shapeHSM_HsmPsfMoments_yy"]
963 raise RuntimeError(
'No psf shape parameter found in catalog')
965 return hsm.where(np.isfinite(hsm), sdss) - psf
969 """Functor to calculate SDSS trace radius size for sources"""
970 name =
"SDSS Trace Size"
971 shortname =
'sdssTrace'
972 _columns = (
"base_SdssShape_xx",
"base_SdssShape_yy")
975 srcSize = np.sqrt(0.5*(df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]))
980 """Functor to calculate SDSS trace radius size difference (%) between object and psf model"""
981 name =
"PSF - SDSS Trace Size"
982 shortname =
'psf_sdssTrace'
983 _columns = (
"base_SdssShape_xx",
"base_SdssShape_yy",
984 "base_SdssShape_psf_xx",
"base_SdssShape_psf_yy")
987 srcSize = np.sqrt(0.5*(df[
"base_SdssShape_xx"] + df[
"base_SdssShape_yy"]))
988 psfSize = np.sqrt(0.5*(df[
"base_SdssShape_psf_xx"] + df[
"base_SdssShape_psf_yy"]))
989 sizeDiff = 100*(srcSize - psfSize)/(0.5*(srcSize + psfSize))
994 """Functor to calculate HSM trace radius size for sources"""
995 name =
'HSM Trace Size'
996 shortname =
'hsmTrace'
997 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
998 "ext_shapeHSM_HsmSourceMoments_yy")
1000 def _func(self, df):
1001 srcSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmSourceMoments_xx"]
1002 + df[
"ext_shapeHSM_HsmSourceMoments_yy"]))
1007 """Functor to calculate HSM trace radius size difference (%) between object and psf model"""
1008 name =
'PSF - HSM Trace Size'
1009 shortname =
'psf_HsmTrace'
1010 _columns = (
"ext_shapeHSM_HsmSourceMoments_xx",
1011 "ext_shapeHSM_HsmSourceMoments_yy",
1012 "ext_shapeHSM_HsmPsfMoments_xx",
1013 "ext_shapeHSM_HsmPsfMoments_yy")
1015 def _func(self, df):
1016 srcSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmSourceMoments_xx"]
1017 + df[
"ext_shapeHSM_HsmSourceMoments_yy"]))
1018 psfSize = np.sqrt(0.5*(df[
"ext_shapeHSM_HsmPsfMoments_xx"]
1019 + df[
"ext_shapeHSM_HsmPsfMoments_yy"]))
1020 sizeDiff = 100*(srcSize - psfSize)/(0.5*(srcSize + psfSize))
1025 name =
'HSM Psf FWHM'
1026 _columns = (
'ext_shapeHSM_HsmPsfMoments_xx',
'ext_shapeHSM_HsmPsfMoments_yy')
1029 SIGMA2FWHM = 2*np.sqrt(2*np.log(2))
1031 def _func(self, df):
1033 0.5*(df[
'ext_shapeHSM_HsmPsfMoments_xx'] + df[
'ext_shapeHSM_HsmPsfMoments_yy']))
1037 name =
"Distortion Ellipticity (e1)"
1038 shortname =
"Distortion"
1051 def _func(self, df):
1056 name =
"Ellipticity e2"
1068 def _func(self, df):
1069 return 2*df[self.
colXYcolXY] / (df[self.
colXXcolXX] + df[self.
colYYcolYY])
1084 def _func(self, df):
1085 return (df[self.
colXXcolXX]*df[self.
colYYcolYY] - df[self.
colXYcolXY]**2)**0.25
1089 """Computations using the stored localWcs.
1091 name =
"LocalWcsOperations"
1106 """Compute the distance on the sphere from x2, y1 to x1, y1.
1114 cd11 : `pandas.Series`
1115 [1, 1] element of the local Wcs affine transform.
1116 cd11 : `pandas.Series`
1117 [1, 1] element of the local Wcs affine transform.
1118 cd12 : `pandas.Series`
1119 [1, 2] element of the local Wcs affine transform.
1120 cd21 : `pandas.Series`
1121 [2, 1] element of the local Wcs affine transform.
1122 cd22 : `pandas.Series`
1123 [2, 2] element of the local Wcs affine transform.
1128 RA and dec conversion of x and y given the local Wcs. Returned
1129 units are in radians.
1132 return (x * cd11 + y * cd12, x * cd21 + y * cd22)
1135 """Compute the local pixel scale conversion.
1139 ra1 : `pandas.Series`
1140 Ra of the first coordinate in radians.
1141 dec1 : `pandas.Series`
1142 Dec of the first coordinate in radians.
1143 ra2 : `pandas.Series`
1144 Ra of the second coordinate in radians.
1145 dec2 : `pandas.Series`
1146 Dec of the second coordinate in radians.
1150 dist : `pandas.Series`
1151 Distance on the sphere in radians.
1153 deltaDec = dec2 - dec1
1155 return 2 * np.arcsin(
1157 np.sin(deltaDec / 2) ** 2
1158 + np.cos(dec2) * np.cos(dec1) * np.sin(deltaRa / 2) ** 2))
1161 """Compute the distance on the sphere from x2, y1 to x1, y1.
1165 x1 : `pandas.Series`
1167 y1 : `pandas.Series`
1169 x2 : `pandas.Series`
1171 y2 : `pandas.Series`
1173 cd11 : `pandas.Series`
1174 [1, 1] element of the local Wcs affine transform.
1175 cd11 : `pandas.Series`
1176 [1, 1] element of the local Wcs affine transform.
1177 cd12 : `pandas.Series`
1178 [1, 2] element of the local Wcs affine transform.
1179 cd21 : `pandas.Series`
1180 [2, 1] element of the local Wcs affine transform.
1181 cd22 : `pandas.Series`
1182 [2, 2] element of the local Wcs affine transform.
1186 Distance : `pandas.Series`
1187 Arcseconds per pixel at the location of the local WC
1189 ra1, dec1 = self.
computeDeltaRaDeccomputeDeltaRaDec(x1, y1, cd11, cd12, cd21, cd22)
1190 ra2, dec2 = self.
computeDeltaRaDeccomputeDeltaRaDec(x2, y2, cd11, cd12, cd21, cd22)
1196 """Compute the local pixel scale from the stored CDMatrix.
1208 """Compute the local pixel to scale conversion in arcseconds.
1212 cd11 : `pandas.Series`
1213 [1, 1] element of the local Wcs affine transform in radians.
1214 cd11 : `pandas.Series`
1215 [1, 1] element of the local Wcs affine transform in radians.
1216 cd12 : `pandas.Series`
1217 [1, 2] element of the local Wcs affine transform in radians.
1218 cd21 : `pandas.Series`
1219 [2, 1] element of the local Wcs affine transform in radians.
1220 cd22 : `pandas.Series`
1221 [2, 2] element of the local Wcs affine transform in radians.
1225 pixScale : `pandas.Series`
1226 Arcseconds per pixel at the location of the local WC
1228 return 3600 * np.degrees(np.sqrt(np.fabs(cd11 * cd22 - cd12 * cd21)))
1230 def _func(self, df):
1238 """Convert a value in units pixels to units arcseconds.
1257 return f
"{self.col}_asArcseconds"
1261 return [self.
colcol,
1267 def _func(self, df):
1275 name =
'Reference Band'
1276 shortname =
'refBand'
1280 return [
"merge_measurement_i",
1281 "merge_measurement_r",
1282 "merge_measurement_z",
1283 "merge_measurement_y",
1284 "merge_measurement_g"]
1286 def _func(self, df):
1287 def getFilterAliasName(row):
1289 colName = row.idxmax()
1290 return colName.replace(
'merge_measurement_',
'')
1292 return df[self.
columnscolumnscolumns].apply(getFilterAliasName, axis=1)
1297 AB_FLUX_SCALE = (0 * u.ABmag).to_value(u.nJy)
1298 LOG_AB_FLUX_SCALE = 12.56
1299 FIVE_OVER_2LOG10 = 1.085736204758129569
1303 def __init__(self, colFlux, colFluxErr=None, calib=None, **kwargs):
1309 if calib
is not None:
1319 return [self.
colcol]
1323 return f
'mag_{self.col}'
1327 if np.abs(a) < np.abs(b):
1332 return np.abs(a) * np.sqrt(1. + q*q)
1338 with np.warnings.catch_warnings():
1339 np.warnings.filterwarnings(
'ignore',
r'invalid value encountered')
1340 np.warnings.filterwarnings(
'ignore',
r'divide by zero')
1341 return -2.5 * np.log10(dn/fluxMag0)
1344 retVal = self.
vhypotvhypot(dn * fluxMag0Err, dnErr * fluxMag0)
1345 retVal *= self.
AB_FLUX_SCALEAB_FLUX_SCALE / fluxMag0 / fluxMag0
1349 retVal = self.
dn2fluxErrdn2fluxErr(dn, dnErr, fluxMag0, fluxMag0Err) / self.
dn2fluxdn2flux(dn, fluxMag0)
1354 def _func(self, df):
1363 def _func(self, df):
1365 return pd.Series(retArr, index=df.index)
1369 def _func(self, df):
1378 def _func(self, df):
1380 return pd.Series(retArr, index=df.index)
1384 """Base class for calibrating the specified instrument flux column using
1385 the local photometric calibration.
1390 Name of the instrument flux column.
1391 instFluxErrCol : `str`
1392 Name of the assocated error columns for ``instFluxCol``.
1393 photoCalibCol : `str`
1394 Name of local calibration column.
1395 photoCalibErrCol : `str`
1396 Error associated with ``photoCalibCol``
1406 logNJanskyToAB = (1 * u.nJy).to_value(u.ABmag)
1421 """Convert instrument flux to nanojanskys.
1425 instFlux : `numpy.ndarray` or `pandas.Series`
1426 Array of instrument flux measurements
1427 localCalib : `numpy.ndarray` or `pandas.Series`
1428 Array of local photometric calibration estimates.
1432 calibFlux : `numpy.ndarray` or `pandas.Series`
1433 Array of calibrated flux measurements.
1435 return instFlux * localCalib
1438 """Convert instrument flux to nanojanskys.
1442 instFlux : `numpy.ndarray` or `pandas.Series`
1443 Array of instrument flux measurements
1444 instFluxErr : `numpy.ndarray` or `pandas.Series`
1445 Errors on associated ``instFlux`` values
1446 localCalib : `numpy.ndarray` or `pandas.Series`
1447 Array of local photometric calibration estimates.
1448 localCalibErr : `numpy.ndarray` or `pandas.Series`
1449 Errors on associated ``localCalib`` values
1453 calibFluxErr : `numpy.ndarray` or `pandas.Series`
1454 Errors on calibrated flux measurements.
1456 return np.hypot(instFluxErr * localCalib, instFlux * localCalibErr)
1459 """Convert instrument flux to nanojanskys.
1463 instFlux : `numpy.ndarray` or `pandas.Series`
1464 Array of instrument flux measurements
1465 localCalib : `numpy.ndarray` or `pandas.Series`
1466 Array of local photometric calibration estimates.
1470 calibMag : `numpy.ndarray` or `pandas.Series`
1471 Array of calibrated AB magnitudes.
1476 """Convert instrument flux err to nanojanskys.
1480 instFlux : `numpy.ndarray` or `pandas.Series`
1481 Array of instrument flux measurements
1482 instFluxErr : `numpy.ndarray` or `pandas.Series`
1483 Errors on associated ``instFlux`` values
1484 localCalib : `numpy.ndarray` or `pandas.Series`
1485 Array of local photometric calibration estimates.
1486 localCalibErr : `numpy.ndarray` or `pandas.Series`
1487 Errors on associated ``localCalib`` values
1491 calibMagErr: `numpy.ndarray` or `pandas.Series`
1492 Error on calibrated AB magnitudes.
1495 return 2.5 / np.log(10) * err / self.
instFluxToNanojanskyinstFluxToNanojansky(instFlux, instFluxErr)
1499 """Compute calibrated fluxes using the local calibration value.
1515 return f
'flux_{self.instFluxCol}'
1517 def _func(self, df):
1522 """Compute calibrated flux errors using the local calibration value.
1539 return f
'fluxErr_{self.instFluxCol}'
1541 def _func(self, df):
1547 """Compute calibrated AB magnitudes using the local calibration value.
1563 return f
'mag_{self.instFluxCol}'
1565 def _func(self, df):
1571 """Compute calibrated AB magnitude errors using the local calibration value.
1588 return f
'magErr_{self.instFluxCol}'
1590 def _func(self, df):
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, **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 __call__(self, parq, dropna=False, **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 mag_aware_eval(df, expr)
def init_fromDict(initDict, basePath='lsst.pipe.tasks.functors', typeKey='functor', name=None)