Loading [MathJax]/extensions/tex2jax.js
LSST Applications g0fba68d861+330b683386,g1ec0fe41b4+3ea9d11450,g1fd858c14a+9be2b0f3b9,g2440f9efcc+8c5ae1fdc5,g33b6eb7922+23bc9e47ac,g35bb328faa+8c5ae1fdc5,g4a4af6cd76+d25431c27e,g4d2262a081+da6e3a53c4,g53246c7159+8c5ae1fdc5,g55585698de+be1c65ba71,g56a49b3a55+92a7603e7a,g60b5630c4e+be1c65ba71,g67b6fd64d1+3fc8cb0b9e,g78460c75b0+7e33a9eb6d,g786e29fd12+668abc6043,g8352419a5c+8c5ae1fdc5,g8852436030+60e38ee5ff,g89139ef638+3fc8cb0b9e,g8de97aba08+a886b35a30,g94187f82dc+be1c65ba71,g989de1cb63+3fc8cb0b9e,g9d31334357+be1c65ba71,g9f33ca652e+69d6bbdd4b,gabe3b4be73+8856018cbb,gabf8522325+977d9fabaf,gb1101e3267+b0077987df,gb89ab40317+3fc8cb0b9e,gc91f06edcd+2e2ca305f6,gcf25f946ba+60e38ee5ff,gd6cbbdb0b4+1cc2750d2e,gdb1c4ca869+be65c9c1d7,gde0f65d7ad+05f0c6b5f8,ge278dab8ac+6b863515ed,ge410e46f29+3fc8cb0b9e,gf35d7ec915+97dd712d81,gf5e32f922b+8c5ae1fdc5,gf618743f1b+3164b09b60,gf67bdafdda+3fc8cb0b9e,w.2025.18
LSST Data Management Base Package
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
diaCalculationPlugins.py
Go to the documentation of this file.
1# This file is part of ap_association.
2#
3# Developed for the LSST Data Management System.
4# This product includes software developed by the LSST Project
5# (https://www.lsst.org).
6# See the COPYRIGHT file at the top-level directory of this distribution
7# for details of code ownership.
8#
9# This program is free software: you can redistribute it and/or modify
10# it under the terms of the GNU General Public License as published by
11# the Free Software Foundation, either version 3 of the License, or
12# (at your option) any later version.
13#
14# This program is distributed in the hope that it will be useful,
15# but WITHOUT ANY WARRANTY; without even the implied warranty of
16# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17# GNU General Public License for more details.
18#
19# You should have received a copy of the GNU General Public License
20# along with this program. If not, see <https://www.gnu.org/licenses/>.
21
22"""Plugins for use in DiaSource summary statistics.
23
24Output columns must be
25as defined in the schema of the Apdb both in name and units.
26"""
27
28import functools
29import warnings
30
31from astropy.stats import median_absolute_deviation
32import numpy as np
33import pandas as pd
34from scipy.optimize import lsq_linear
35
36import lsst.geom as geom
37import lsst.pex.config as pexConfig
38import lsst.sphgeom as sphgeom
39from astropy.timeseries import LombScargle
40from astropy.timeseries import LombScargleMultiband
41import math
42import statistics
43
44from .diaCalculation import (
45 DiaObjectCalculationPluginConfig,
46 DiaObjectCalculationPlugin)
47from .pluginRegistry import register
48
49
50__all__ = ("MeanDiaPositionConfig", "MeanDiaPosition",
51 "HTMIndexDiaPosition", "HTMIndexDiaPositionConfig",
52 "NumDiaSourcesDiaPlugin", "NumDiaSourcesDiaPluginConfig",
53 "SimpleSourceFlagDiaPlugin", "SimpleSourceFlagDiaPluginConfig",
54 "WeightedMeanDiaPsfFluxConfig", "WeightedMeanDiaPsfFlux",
55 "PercentileDiaPsfFlux", "PercentileDiaPsfFluxConfig",
56 "SigmaDiaPsfFlux", "SigmaDiaPsfFluxConfig",
57 "Chi2DiaPsfFlux", "Chi2DiaPsfFluxConfig",
58 "MadDiaPsfFlux", "MadDiaPsfFluxConfig",
59 "SkewDiaPsfFlux", "SkewDiaPsfFluxConfig",
60 "MinMaxDiaPsfFlux", "MinMaxDiaPsfFluxConfig",
61 "MaxSlopeDiaPsfFlux", "MaxSlopeDiaPsfFluxConfig",
62 "ErrMeanDiaPsfFlux", "ErrMeanDiaPsfFluxConfig",
63 "LinearFitDiaPsfFlux", "LinearFitDiaPsfFluxConfig",
64 "StetsonJDiaPsfFlux", "StetsonJDiaPsfFluxConfig",
65 "WeightedMeanDiaTotFlux", "WeightedMeanDiaTotFluxConfig",
66 "SigmaDiaTotFlux", "SigmaDiaTotFluxConfig",
67 "LombScarglePeriodogram", "LombScarglePeriodogramConfig",
68 "LombScarglePeriodogramMulti", "LombScarglePeriodogramMultiConfig")
69
70
71def catchWarnings(_func=None, *, warns=[]):
72 """Decorator for generically catching numpy warnings.
73 """
74 def decoratorCatchWarnings(func):
75 @functools.wraps(func)
76 def wrapperCatchWarnings(*args, **kwargs):
77 with warnings.catch_warnings():
78 for val in warns:
79 warnings.filterwarnings("ignore", val)
80 return func(*args, **kwargs)
81 return wrapperCatchWarnings
82
83 if _func is None:
84 return decoratorCatchWarnings
85 else:
86 return decoratorCatchWarnings(_func)
87
88
89def compute_optimized_periodogram_grid(x0, oversampling_factor=5, nyquist_factor=100):
90 """
91 Computes an optimized periodogram frequency grid for a given time series.
92
93 Parameters
94 ----------
95 x0 : `array`
96 The input time axis.
97 oversampling_factor : `int`, optional
98 The oversampling factor for frequency grid.
99 nyquist_factor : `int`, optional
100 The Nyquist factor for frequency grid.
101
102 Returns
103 -------
104 frequencies : `array`
105 The computed optimized periodogram frequency grid.
106 """
107
108 num_points = len(x0)
109 baseline = np.max(x0) - np.min(x0)
110
111 # Calculate the frequency resolution based on oversampling factor and baseline
112 frequency_resolution = 1. / baseline / oversampling_factor
113
114 num_frequencies = int(
115 0.5 * oversampling_factor * nyquist_factor * num_points)
116 frequencies = frequency_resolution + \
117 frequency_resolution * np.arange(num_frequencies)
118
119 return frequencies
120
121
123 pass
124
125
126@register("ap_lombScarglePeriodogram")
128 """Compute the single-band period of a DiaObject given a set of DiaSources.
129 """
130 ConfigClass = LombScarglePeriodogramConfig
131
132 plugType = "multi"
133 outputCols = ["period", "power"]
134 needsFilter = True
135
136 @classmethod
139
140 @catchWarnings(warns=["All-NaN slice encountered"])
141 def calculate(self,
142 diaObjects,
143 diaSources,
144 filterDiaSources,
145 band):
146 """Compute the periodogram.
147
148 Parameters
149 ----------
150 diaObjects : `pandas.DataFrame`
151 Summary objects to store values in.
152 diaSources : `pandas.DataFrame` or `pandas.DataFrameGroupBy`
153 Catalog of DiaSources summarized by this DiaObject.
154 """
155
156 # Check and initialize output columns in diaObjects.
157 if (periodCol := f"{band}_period") not in diaObjects.columns:
158 diaObjects[periodCol] = np.nan
159 if (powerCol := f"{band}_power") not in diaObjects.columns:
160 diaObjects[powerCol] = np.nan
161
162 def _calculate_period(df, min_detections=5, nterms=1, oversampling_factor=5, nyquist_factor=100):
163 """Compute the Lomb-Scargle periodogram given a set of DiaSources.
164
165 Parameters
166 ----------
167 df : `pandas.DataFrame`
168 The input DataFrame.
169 min_detections : `int`, optional
170 The minimum number of detections.
171 nterms : `int`, optional
172 The number of terms in the Lomb-Scargle model.
173 oversampling_factor : `int`, optional
174 The oversampling factor for frequency grid.
175 nyquist_factor : `int`, optional
176 The Nyquist factor for frequency grid.
177
178 Returns
179 -------
180 pd_tab : `pandas.Series`
181 The output DataFrame with the Lomb-Scargle parameters.
182 """
183 tmpDf = df[~np.logical_or(np.isnan(df["psfFlux"]),
184 np.isnan(df["midpointMjdTai"]))]
185
186 if len(tmpDf) < min_detections:
187 return pd.Series({periodCol: np.nan, powerCol: np.nan})
188
189 time = tmpDf["midpointMjdTai"].to_numpy()
190 flux = tmpDf["psfFlux"].to_numpy()
191 flux_err = tmpDf["psfFluxErr"].to_numpy()
192
193 lsp = LombScargle(time, flux, dy=flux_err, nterms=nterms)
195 time, oversampling_factor=oversampling_factor, nyquist_factor=nyquist_factor)
196 period = 1/f_grid
197 power = lsp.power(f_grid)
198
199 return pd.Series({periodCol: period[np.argmax(power)],
200 powerCol: np.max(power)})
201
202 diaObjects.loc[:, [periodCol, powerCol]
203 ] = filterDiaSources.apply(_calculate_period)
204
205
207 pass
208
209
210@register("ap_lombScarglePeriodogramMulti")
212 """Compute the multi-band LombScargle periodogram of a DiaObject given a set of DiaSources.
213 """
214 ConfigClass = LombScarglePeriodogramMultiConfig
215
216 plugType = "multi"
217 outputCols = ["multiPeriod", "multiPower",
218 "multiFap", "multiAmp", "multiPhase"]
219 needsFilter = True
220
221 @classmethod
224
225 @staticmethod
226 def calculate_baluev_fap(time, n, maxPeriod, zmax):
227 """Calculate the False-Alarm probability using the Baluev approximation.
228
229 Parameters
230 ----------
231 time : `array`
232 The input time axis.
233 n : `int`
234 The number of detections.
235 maxPeriod : `float`
236 The maximum period in the grid.
237 zmax : `float`
238 The maximum power in the grid.
239
240 Returns
241 -------
242 fap_estimate : `float`
243 The False-Alarm probability Baluev approximation.
244
245 Notes
246 ----------
247 .. [1] Baluev, R. V. 2008, MNRAS, 385, 1279
248 .. [2] Süveges, M., Guy, L.P., Eyer, L., et al. 2015, MNRAS, 450, 2052
249 """
250 if n <= 2:
251 return np.nan
252
253 gam_ratio = math.factorial(int((n - 1)/2)) / math.factorial(int((n - 2)/2))
254 fu = 1/maxPeriod
255 return gam_ratio * np.sqrt(
256 4*np.pi*statistics.variance(time)
257 ) * fu * (1-zmax)**((n-4)/2) * np.sqrt(zmax)
258
259 @staticmethod
260 def generate_lsp_params(lsp_model, fbest, bands):
261 """Generate the Lomb-Scargle parameters.
262 Parameters
263 ----------
264 lsp_model : `astropy.timeseries.LombScargleMultiband`
265 The Lomb-Scargle model.
266 fbest : `float`
267 The best period.
268 bands : `array`
269 The bands of the time series.
270
271 Returns
272 -------
273 Amp : `array`
274 The amplitude of the time series.
275 Ph : `array`
276 The phase of the time series.
277
278 Notes
279 ----------
280 .. [1] VanderPlas, J. T., & Ivezic, Z. 2015, ApJ, 812, 18
281 """
282 best_params = lsp_model.model_parameters(fbest, units=True)
283
284 name_params = [f"theta_base_{i}" for i in range(3)]
285 name_params += [f"theta_band_{band}_{i}" for band in np.unique(bands) for i in range(3)]
286
287 df_params = pd.DataFrame([best_params], columns=name_params)
288
289 unique_bands = np.unique(bands)
290
291 amplitude_band = [np.sqrt(df_params[f"theta_band_{band}_1"]**2
292 + df_params[f"theta_band_{band}_2"]**2)
293 for band in unique_bands]
294 phase_bands = [np.arctan2(df_params[f"theta_band_{band}_2"],
295 df_params[f"theta_band_{band}_1"]) for band in unique_bands]
296
297 amp = [a[0] for a in amplitude_band]
298 ph = [p[0] for p in phase_bands]
299
300 return amp, ph
301
302 @catchWarnings(warns=["All-NaN slice encountered"])
303 def calculate(self,
304 diaObjects,
305 diaSources,
306 **kwargs):
307 """Compute the multi-band LombScargle periodogram of a DiaObject given
308 a set of DiaSources.
309
310 Parameters
311 ----------
312 diaObjects : `pandas.DataFrame`
313 Summary objects to store values in.
314 diaSources : `pandas.DataFrame` or `pandas.DataFrameGroupBy`
315 Catalog of DiaSources summarized by this DiaObject.
316 **kwargs : `dict`
317 Unused kwargs that are always passed to a plugin.
318 """
319
320 bands_arr = diaSources['band'].unique().values
321 unique_bands = np.unique(np.concatenate(bands_arr))
322 # Check and initialize output columns in diaObjects.
323 if (periodCol := "multiPeriod") not in diaObjects.columns:
324 diaObjects[periodCol] = np.nan
325 if (powerCol := "multiPower") not in diaObjects.columns:
326 diaObjects[powerCol] = np.nan
327 if (fapCol := "multiFap") not in diaObjects.columns:
328 diaObjects[fapCol] = np.nan
329 ampCol = "multiAmp"
330 phaseCol = "multiPhase"
331 for i in range(len(unique_bands)):
332 ampCol_band = f"{unique_bands[i]}_{ampCol}"
333 if ampCol_band not in diaObjects.columns:
334 diaObjects[ampCol_band] = np.nan
335 phaseCol_band = f"{unique_bands[i]}_{phaseCol}"
336 if phaseCol_band not in diaObjects.columns:
337 diaObjects[phaseCol_band] = np.nan
338
339 def _calculate_period_multi(df, all_unique_bands,
340 min_detections=9, oversampling_factor=5, nyquist_factor=100):
341 """Calculate the multi-band Lomb-Scargle periodogram.
342
343 Parameters
344 ----------
345 df : `pandas.DataFrame`
346 The input DataFrame.
347 all_unique_bands : `list` of `str`
348 List of all bands present in the diaSource table that is being worked on.
349 min_detections : `int`, optional
350 The minimum number of detections, including all bands.
351 oversampling_factor : `int`, optional
352 The oversampling factor for frequency grid.
353 nyquist_factor : `int`, optional
354 The Nyquist factor for frequency grid.
355
356 Returns
357 -------
358 pd_tab : `pandas.Series`
359 The output DataFrame with the Lomb-Scargle parameters.
360 """
361 tmpDf = df[~np.logical_or(np.isnan(df["psfFlux"]),
362 np.isnan(df["midpointMjdTai"]))]
363
364 if (len(tmpDf)) < min_detections:
365 pd_tab_nodet = pd.Series({periodCol: np.nan,
366 powerCol: np.nan,
367 fapCol: np.nan})
368 for band in all_unique_bands:
369 pd_tab_nodet[f"{band}_{ampCol}"] = np.nan
370 pd_tab_nodet[f"{band}_{phaseCol}"] = np.nan
371
372 return pd_tab_nodet
373
374 time = tmpDf["midpointMjdTai"].to_numpy()
375 flux = tmpDf["psfFlux"].to_numpy()
376 flux_err = tmpDf["psfFluxErr"].to_numpy()
377 bands = tmpDf["band"].to_numpy()
378
379 lsp = LombScargleMultiband(time, flux, bands, dy=flux_err,
380 nterms_base=1, nterms_band=1)
381
383 time, oversampling_factor=oversampling_factor, nyquist_factor=nyquist_factor)
384 period = 1/f_grid
385 power = lsp.power(f_grid)
386
387 fap_estimate = self.calculate_baluev_fap(
388 time, len(time), period[np.argmax(power)], np.max(power))
389
390 params_table_new = self.generate_lsp_params(lsp, f_grid[np.argmax(power)], bands)
391
392 pd_tab = pd.Series({periodCol: period[np.argmax(power)],
393 powerCol: np.max(power),
394 fapCol: fap_estimate
395 })
396
397 # Initialize the per-band amplitude/phase columns as NaNs
398 for band in all_unique_bands:
399 pd_tab[f"{band}_{ampCol}"] = np.nan
400 pd_tab[f"{band}_{phaseCol}"] = np.nan
401
402 # Populate the values of only the bands that have data for this diaSource
403 unique_bands = np.unique(bands)
404 for i in range(len(unique_bands)):
405 pd_tab[f"{unique_bands[i]}_{ampCol}"] = params_table_new[0][i]
406 pd_tab[f"{unique_bands[i]}_{phaseCol}"] = params_table_new[1][i]
407
408 return pd_tab
409
410 columns_list = [periodCol, powerCol, fapCol]
411 for i in range(len(unique_bands)):
412 columns_list.append(f"{unique_bands[i]}_{ampCol}")
413 columns_list.append(f"{unique_bands[i]}_{phaseCol}")
414
415 diaObjects.loc[:, columns_list
416 ] = diaSources.apply(_calculate_period_multi, unique_bands)
417
418
420 pass
421
422
423@register("ap_meanPosition")
425 """Compute the mean position of a DiaObject given a set of DiaSources.
426 """
427
428 ConfigClass = MeanDiaPositionConfig
429
430 plugType = 'multi'
431
432 outputCols = ["ra", "dec", "radecMjdTai"]
433 needsFilter = False
434
435 @classmethod
438
439 def calculate(self, diaObjects, diaSources, **kwargs):
440 """Compute the mean ra/dec position of the diaObject given the
441 diaSource locations.
442
443 Parameters
444 ----------
445 diaObjects : `pandas.DataFrame`
446 Summary objects to store values in.
447 diaSources : `pandas.DataFrame` or `pandas.DataFrameGroupBy`
448 Catalog of DiaSources summarized by this DiaObject.
449 **kwargs
450 Any additional keyword arguments that may be passed to the plugin.
451 """
452 for outCol in self.outputCols:
453 if outCol not in diaObjects.columns:
454 diaObjects[outCol] = np.nan
455
456 def _computeMeanPos(df):
457 aveCoord = geom.averageSpherePoint(
458 list(geom.SpherePoint(src["ra"], src["dec"], geom.degrees)
459 for idx, src in df.iterrows()))
460 ra = aveCoord.getRa().asDegrees()
461 dec = aveCoord.getDec().asDegrees()
462 if np.isnan(ra) or np.isnan(dec):
463 radecMjdTai = np.nan
464 else:
465 radecMjdTai = df["midpointMjdTai"].max()
466
467 return pd.Series({"ra": aveCoord.getRa().asDegrees(),
468 "dec": aveCoord.getDec().asDegrees(),
469 "radecMjdTai": radecMjdTai})
470
471 ans = diaSources.apply(_computeMeanPos)
472 diaObjects.loc[:, ["ra", "dec", "radecMjdTai"]] = ans
473
474
476
477 htmLevel = pexConfig.Field(
478 dtype=int,
479 doc="Level of the HTM pixelization.",
480 default=20,
481 )
482
483
484@register("ap_HTMIndex")
486 """Compute the mean position of a DiaObject given a set of DiaSources.
487
488 Notes
489 -----
490 This plugin was implemented to satisfy requirements of old APDB interface
491 which required ``pixelId`` column in DiaObject with HTM20 index. APDB
492 interface had migrated to not need that information, but we keep this
493 plugin in case it may be useful for something else.
494 """
495 ConfigClass = HTMIndexDiaPositionConfig
496
497 plugType = 'single'
498
499 inputCols = ["ra", "dec"]
500 outputCols = ["pixelId"]
501 needsFilter = False
502
503 def __init__(self, config, name, metadata):
504 DiaObjectCalculationPlugin.__init__(self, config, name, metadata)
506
507 @classmethod
509 return cls.FLUX_MOMENTS_CALCULATED
510
511 def calculate(self, diaObjects, diaObjectId, **kwargs):
512 """Compute the mean position of a DiaObject given a set of DiaSources
513
514 Parameters
515 ----------
516 diaObjects : `pandas.dataFrame`
517 Summary objects to store values in and read ra/dec from.
518 diaObjectId : `int`
519 Id of the diaObject to update.
520 **kwargs
521 Any additional keyword arguments that may be passed to the plugin.
522 """
523 sphPoint = geom.SpherePoint(
524 diaObjects.at[diaObjectId, "ra"] * geom.degrees,
525 diaObjects.at[diaObjectId, "dec"] * geom.degrees)
526 diaObjects.at[diaObjectId, "pixelId"] = self.pixelator.index(
527 sphPoint.getVector())
528
529
531 pass
532
533
534@register("ap_nDiaSources")
536 """Compute the total number of DiaSources associated with this DiaObject.
537 """
538
539 ConfigClass = NumDiaSourcesDiaPluginConfig
540 outputCols = ["nDiaSources"]
541 plugType = "multi"
542 needsFilter = False
543
544 @classmethod
547
548 def calculate(self, diaObjects, diaSources, **kwargs):
549 """Compute the total number of DiaSources associated with this DiaObject.
550
551 Parameters
552 ----------
553 diaObject : `dict`
554 Summary object to store values in and read ra/dec from.
555 **kwargs
556 Any additional keyword arguments that may be passed to the plugin.
557 """
558 dtype = diaObjects["nDiaSources"].dtype
559 diaObjects.loc[:, "nDiaSources"] = diaSources.diaObjectId.count().astype(dtype)
560
561
563 pass
564
565
566@register("ap_diaObjectFlag")
568 """Find if any DiaSource is flagged.
569
570 Set the DiaObject flag if any DiaSource is flagged.
571 """
572
573 ConfigClass = NumDiaSourcesDiaPluginConfig
574 outputCols = ["flags"]
575 plugType = "multi"
576 needsFilter = False
577
578 @classmethod
581
582 def calculate(self, diaObjects, diaSources, **kwargs):
583 """Find if any DiaSource is flagged.
584
585 Set the DiaObject flag if any DiaSource is flagged.
586
587 Parameters
588 ----------
589 diaObject : `dict`
590 Summary object to store values in and read ra/dec from.
591 **kwargs
592 Any additional keyword arguments that may be passed to the plugin.
593 """
594 dtype = diaObjects["flags"].dtype
595 diaObjects.loc[:, "flags"] = diaSources.flags.any().astype(dtype)
596
597
604 """Compute the weighted mean and mean error on the point source fluxes
605 of the DiaSource measured on the difference image.
606
607 Additionally store number of usable data points.
608 """
609
610 ConfigClass = WeightedMeanDiaPsfFluxConfig
611 outputCols = ["psfFluxMean", "psfFluxMeanErr", "psfFluxNdata"]
612 plugType = "multi"
613 needsFilter = True
614
615 @classmethod
618
619 @catchWarnings(warns=["invalid value encountered",
620 "divide by zero"])
621 def calculate(self,
622 diaObjects,
623 diaSources,
624 filterDiaSources,
625 band,
626 **kwargs):
627 """Compute the weighted mean and mean error of the point source flux.
628
629 Parameters
630 ----------
631 diaObject : `dict`
632 Summary object to store values in.
633 diaSources : `pandas.DataFrame`
634 DataFrame representing all diaSources associated with this
635 diaObject.
636 filterDiaSources : `pandas.DataFrame`
637 DataFrame representing diaSources associated with this
638 diaObject that are observed in the band pass ``band``.
639 band : `str`
640 Simple, string name of the filter for the flux being calculated.
641 **kwargs
642 Any additional keyword arguments that may be passed to the plugin.
643 """
644 meanName = "{}_psfFluxMean".format(band)
645 errName = "{}_psfFluxMeanErr".format(band)
646 nDataName = "{}_psfFluxNdata".format(band)
647 if meanName not in diaObjects.columns:
648 diaObjects[meanName] = np.nan
649 if errName not in diaObjects.columns:
650 diaObjects[errName] = np.nan
651 if nDataName not in diaObjects.columns:
652 diaObjects[nDataName] = 0
653
654 def _weightedMean(df):
655 tmpDf = df[~np.logical_or(np.isnan(df["psfFlux"]),
656 np.isnan(df["psfFluxErr"]))]
657 tot_weight = np.nansum(1 / tmpDf["psfFluxErr"] ** 2)
658 fluxMean = np.nansum(tmpDf["psfFlux"]
659 / tmpDf["psfFluxErr"] ** 2)
660 fluxMean /= tot_weight
661 if tot_weight > 0:
662 fluxMeanErr = np.sqrt(1 / tot_weight)
663 else:
664 fluxMeanErr = np.nan
665 nFluxData = len(tmpDf)
666
667 return pd.Series({meanName: fluxMean,
668 errName: fluxMeanErr,
669 nDataName: nFluxData},
670 dtype="object")
671 df = filterDiaSources.apply(_weightedMean).astype(diaObjects.dtypes[[meanName, errName, nDataName]])
672
673 diaObjects.loc[:, [meanName, errName, nDataName]] = df
674
675
677 percentiles = pexConfig.ListField(
678 dtype=int,
679 default=[5, 25, 50, 75, 95],
680 doc="Percentiles to calculate to compute values for. Should be "
681 "integer values."
682 )
683
684
685@register("ap_percentileFlux")
687 """Compute percentiles of diaSource fluxes.
688 """
689
690 ConfigClass = PercentileDiaPsfFluxConfig
691 # Output columns are created upon instantiation of the class.
692 outputCols = []
693 plugType = "multi"
694 needsFilter = True
695
696 def __init__(self, config, name, metadata, **kwargs):
697 DiaObjectCalculationPlugin.__init__(self,
698 config,
699 name,
700 metadata,
701 **kwargs)
702 self.outputCols = ["psfFluxPercentile{:02d}".format(percent)
703 for percent in self.config.percentiles]
704
705 @classmethod
708
709 @catchWarnings(warns=["All-NaN slice encountered"])
710 def calculate(self,
711 diaObjects,
712 diaSources,
713 filterDiaSources,
714 band,
715 **kwargs):
716 """Compute the percentile fluxes of the point source flux.
717
718 Parameters
719 ----------
720 diaObject : `dict`
721 Summary object to store values in.
722 diaSources : `pandas.DataFrame`
723 DataFrame representing all diaSources associated with this
724 diaObject.
725 filterDiaSources : `pandas.DataFrame`
726 DataFrame representing diaSources associated with this
727 diaObject that are observed in the band pass ``band``.
728 band : `str`
729 Simple, string name of the filter for the flux being calculated.
730 **kwargs
731 Any additional keyword arguments that may be passed to the plugin.
732 """
733 pTileNames = []
734 dtype = None
735 for tilePercent in self.config.percentiles:
736 pTileName = "{}_psfFluxPercentile{:02d}".format(band,
737 tilePercent)
738 pTileNames.append(pTileName)
739 if pTileName not in diaObjects.columns:
740 diaObjects[pTileName] = np.nan
741 elif dtype is None:
742 dtype = diaObjects[pTileName].dtype
743
744 def _fluxPercentiles(df):
745 pTiles = np.nanpercentile(df["psfFlux"], self.config.percentiles)
746 return pd.Series(
747 dict((tileName, pTile)
748 for tileName, pTile in zip(pTileNames, pTiles)))
749 if dtype is None:
750 diaObjects.loc[:, pTileNames] = filterDiaSources.apply(_fluxPercentiles)
751 else:
752 diaObjects.loc[:, pTileNames] = filterDiaSources.apply(_fluxPercentiles).astype(dtype)
753
754
756 pass
757
758
759@register("ap_sigmaFlux")
761 """Compute scatter of diaSource fluxes.
762 """
763
764 ConfigClass = SigmaDiaPsfFluxConfig
765 # Output columns are created upon instantiation of the class.
766 outputCols = ["psfFluxSigma"]
767 plugType = "multi"
768 needsFilter = True
769
770 @classmethod
773
774 def calculate(self,
775 diaObjects,
776 diaSources,
777 filterDiaSources,
778 band,
779 **kwargs):
780 """Compute the sigma fluxes of the point source flux.
781
782 Parameters
783 ----------
784 diaObject : `dict`
785 Summary object to store values in.
786 diaSources : `pandas.DataFrame`
787 DataFrame representing all diaSources associated with this
788 diaObject.
789 filterDiaSources : `pandas.DataFrame`
790 DataFrame representing diaSources associated with this
791 diaObject that are observed in the band pass ``band``.
792 band : `str`
793 Simple, string name of the filter for the flux being calculated.
794 **kwargs
795 Any additional keyword arguments that may be passed to the plugin.
796 """
797 # Set "delta degrees of freedom (ddf)" to 1 to calculate the unbiased
798 # estimator of scatter (i.e. 'N - 1' instead of 'N').
799 column = "{}_psfFluxSigma".format(band)
800 if column in diaObjects:
801 dtype = diaObjects[column].dtype
802 diaObjects.loc[:, column] = filterDiaSources.psfFlux.std().astype(dtype)
803 else:
804 diaObjects.loc[:, column] = filterDiaSources.psfFlux.std()
805
806
808 pass
809
810
811@register("ap_chi2Flux")
813 """Compute chi2 of diaSource fluxes.
814 """
815
816 ConfigClass = Chi2DiaPsfFluxConfig
817
818 # Required input Cols
819 inputCols = ["psfFluxMean"]
820 # Output columns are created upon instantiation of the class.
821 outputCols = ["psfFluxChi2"]
822 plugType = "multi"
823 needsFilter = True
824
825 @classmethod
827 return cls.FLUX_MOMENTS_CALCULATED
828
829 @catchWarnings(warns=["All-NaN slice encountered"])
830 def calculate(self,
831 diaObjects,
832 diaSources,
833 filterDiaSources,
834 band,
835 **kwargs):
836 """Compute the chi2 of the point source fluxes.
837
838 Parameters
839 ----------
840 diaObject : `dict`
841 Summary object to store values in.
842 diaSources : `pandas.DataFrame`
843 DataFrame representing all diaSources associated with this
844 diaObject.
845 filterDiaSources : `pandas.DataFrame`
846 DataFrame representing diaSources associated with this
847 diaObject that are observed in the band pass ``band``.
848 band : `str`
849 Simple, string name of the filter for the flux being calculated.
850 **kwargs
851 Any additional keyword arguments that may be passed to the plugin.
852 """
853 meanName = "{}_psfFluxMean".format(band)
854 column = "{}_psfFluxChi2".format(band)
855
856 def _chi2(df):
857 delta = (df["psfFlux"]
858 - diaObjects.at[df.diaObjectId.iat[0], meanName])
859 return np.nansum((delta / df["psfFluxErr"]) ** 2)
860
861 if column in diaObjects:
862 dtype = diaObjects[column].dtype
863 diaObjects.loc[:, column] = filterDiaSources.apply(_chi2).astype(dtype)
864 else:
865 diaObjects.loc[:, column] = filterDiaSources.apply(_chi2)
866
867
869 pass
870
871
872@register("ap_madFlux")
874 """Compute median absolute deviation of diaSource fluxes.
875 """
876
877 ConfigClass = MadDiaPsfFluxConfig
878
879 # Required input Cols
880 # Output columns are created upon instantiation of the class.
881 outputCols = ["psfFluxMAD"]
882 plugType = "multi"
883 needsFilter = True
884
885 @classmethod
888
889 @catchWarnings(warns=["All-NaN slice encountered"])
890 def calculate(self,
891 diaObjects,
892 diaSources,
893 filterDiaSources,
894 band,
895 **kwargs):
896 """Compute the median absolute deviation of the point source fluxes.
897
898 Parameters
899 ----------
900 diaObject : `dict`
901 Summary object to store values in.
902 diaSources : `pandas.DataFrame`
903 DataFrame representing all diaSources associated with this
904 diaObject.
905 filterDiaSources : `pandas.DataFrame`
906 DataFrame representing diaSources associated with this
907 diaObject that are observed in the band pass ``band``.
908 band : `str`
909 Simple, string name of the filter for the flux being calculated.
910 **kwargs
911 Any additional keyword arguments that may be passed to the plugin.
912 """
913 column = "{}_psfFluxMAD".format(band)
914 if column in diaObjects:
915 dtype = diaObjects[column].dtype
916 diaObjects.loc[:, column] = \
917 filterDiaSources.psfFlux.apply(median_absolute_deviation,
918 ignore_nan=True).astype(dtype)
919 else:
920 diaObjects.loc[:, column] = \
921 filterDiaSources.psfFlux.apply(median_absolute_deviation,
922 ignore_nan=True)
923
924
926 pass
927
928
929@register("ap_skewFlux")
931 """Compute the skew of diaSource fluxes.
932 """
933
934 ConfigClass = SkewDiaPsfFluxConfig
935
936 # Required input Cols
937 # Output columns are created upon instantiation of the class.
938 outputCols = ["psfFluxSkew"]
939 plugType = "multi"
940 needsFilter = True
941
942 @classmethod
945
946 def calculate(self,
947 diaObjects,
948 diaSources,
949 filterDiaSources,
950 band,
951 **kwargs):
952 """Compute the skew of the point source fluxes.
953
954 Parameters
955 ----------
956 diaObject : `dict`
957 Summary object to store values in.
958 diaSources : `pandas.DataFrame`
959 DataFrame representing all diaSources associated with this
960 diaObject.
961 filterDiaSources : `pandas.DataFrame`
962 DataFrame representing diaSources associated with this
963 diaObject that are observed in the band pass ``band``.
964 band : `str`
965 Simple, string name of the filter for the flux being calculated.
966 **kwargs
967 Any additional keyword arguments that may be passed to the plugin.
968 """
969 column = "{}_psfFluxSkew".format(band)
970 if column in diaObjects:
971 dtype = diaObjects[column].dtype
972 diaObjects.loc[:, column] = filterDiaSources.psfFlux.skew().astype(dtype)
973 else:
974 diaObjects.loc[:, column] = filterDiaSources.psfFlux.skew()
975
976
978 pass
979
980
981@register("ap_minMaxFlux")
983 """Compute min/max of diaSource fluxes.
984 """
985
986 ConfigClass = MinMaxDiaPsfFluxConfig
987
988 # Required input Cols
989 # Output columns are created upon instantiation of the class.
990 outputCols = ["psfFluxMin", "psfFluxMax"]
991 plugType = "multi"
992 needsFilter = True
993
994 @classmethod
997
998 def calculate(self,
999 diaObjects,
1000 diaSources,
1001 filterDiaSources,
1002 band,
1003 **kwargs):
1004 """Compute min/max of the point source fluxes.
1005
1006 Parameters
1007 ----------
1008 diaObject : `dict`
1009 Summary object to store values in.
1010 diaSources : `pandas.DataFrame`
1011 DataFrame representing all diaSources associated with this
1012 diaObject.
1013 filterDiaSources : `pandas.DataFrame`
1014 DataFrame representing diaSources associated with this
1015 diaObject that are observed in the band pass ``band``.
1016 band : `str`
1017 Simple, string name of the filter for the flux being calculated.
1018 **kwargs
1019 Any additional keyword arguments that may be passed to the plugin.
1020 """
1021 minName = "{}_psfFluxMin".format(band)
1022 if minName not in diaObjects.columns:
1023 diaObjects[minName] = np.nan
1024 maxName = "{}_psfFluxMax".format(band)
1025 if maxName not in diaObjects.columns:
1026 diaObjects[maxName] = np.nan
1027
1028 dtype = diaObjects[minName].dtype
1029 diaObjects.loc[:, minName] = filterDiaSources.psfFlux.min().astype(dtype)
1030 diaObjects.loc[:, maxName] = filterDiaSources.psfFlux.max().astype(dtype)
1031
1032
1034 pass
1035
1036
1037@register("ap_maxSlopeFlux")
1039 """Compute the maximum ratio time ordered deltaFlux / deltaTime.
1040 """
1041
1042 ConfigClass = MinMaxDiaPsfFluxConfig
1043
1044 # Required input Cols
1045 # Output columns are created upon instantiation of the class.
1046 outputCols = ["psfFluxMaxSlope"]
1047 plugType = "multi"
1048 needsFilter = True
1049
1050 @classmethod
1053
1054 def calculate(self,
1055 diaObjects,
1056 diaSources,
1057 filterDiaSources,
1058 band,
1059 **kwargs):
1060 """Compute the maximum ratio time ordered deltaFlux / deltaTime.
1061
1062 Parameters
1063 ----------
1064 diaObject : `dict`
1065 Summary object to store values in.
1066 diaSources : `pandas.DataFrame`
1067 DataFrame representing all diaSources associated with this
1068 diaObject.
1069 filterDiaSources : `pandas.DataFrame`
1070 DataFrame representing diaSources associated with this
1071 diaObject that are observed in the band pass ``band``.
1072 band : `str`
1073 Simple, string name of the filter for the flux being calculated.
1074 **kwargs
1075 Any additional keyword arguments that may be passed to the plugin.
1076 """
1077
1078 def _maxSlope(df):
1079 tmpDf = df[~np.logical_or(np.isnan(df["psfFlux"]),
1080 np.isnan(df["midpointMjdTai"]))]
1081 if len(tmpDf) < 2:
1082 return np.nan
1083 times = tmpDf["midpointMjdTai"].to_numpy()
1084 timeArgs = times.argsort()
1085 times = times[timeArgs]
1086 fluxes = tmpDf["psfFlux"].to_numpy()[timeArgs]
1087 return (np.diff(fluxes) / np.diff(times)).max()
1088
1089 column = "{}_psfFluxMaxSlope".format(band)
1090 if column in diaObjects:
1091 dtype = diaObjects[column].dtype
1092 diaObjects.loc[:, column] = filterDiaSources.apply(_maxSlope).astype(dtype)
1093 else:
1094 diaObjects.loc[:, column] = filterDiaSources.apply(_maxSlope)
1095
1096
1098 pass
1099
1100
1101@register("ap_meanErrFlux")
1103 """Compute the mean of the dia source errors.
1104 """
1105
1106 ConfigClass = ErrMeanDiaPsfFluxConfig
1107
1108 # Required input Cols
1109 # Output columns are created upon instantiation of the class.
1110 outputCols = ["psfFluxErrMean"]
1111 plugType = "multi"
1112 needsFilter = True
1113
1114 @classmethod
1117
1118 def calculate(self,
1119 diaObjects,
1120 diaSources,
1121 filterDiaSources,
1122 band,
1123 **kwargs):
1124 """Compute the mean of the dia source errors.
1125
1126 Parameters
1127 ----------
1128 diaObject : `dict`
1129 Summary object to store values in.
1130 diaSources : `pandas.DataFrame`
1131 DataFrame representing all diaSources associated with this
1132 diaObject.
1133 filterDiaSources : `pandas.DataFrame`
1134 DataFrame representing diaSources associated with this
1135 diaObject that are observed in the band pass ``band``.
1136 band : `str`
1137 Simple, string name of the filter for the flux being calculated.
1138 **kwargs
1139 Any additional keyword arguments that may be passed to the plugin.
1140 """
1141 column = "{}_psfFluxErrMean".format(band)
1142 if column in diaObjects:
1143 dtype = diaObjects[column].dtype
1144 diaObjects.loc[:, column] = filterDiaSources.psfFluxErr.mean().astype(dtype)
1145 else:
1146 diaObjects.loc[:, column] = filterDiaSources.psfFluxErr.mean()
1147
1148
1150 pass
1151
1152
1153@register("ap_linearFit")
1155 """Compute fit a linear model to flux vs time.
1156 """
1157
1158 ConfigClass = LinearFitDiaPsfFluxConfig
1159
1160 # Required input Cols
1161 # Output columns are created upon instantiation of the class.
1162 outputCols = ["psfFluxLinearSlope", "psfFluxLinearIntercept"]
1163 plugType = "multi"
1164 needsFilter = True
1165
1166 @classmethod
1169
1170 def calculate(self,
1171 diaObjects,
1172 diaSources,
1173 filterDiaSources,
1174 band,
1175 **kwargs):
1176 """Compute fit a linear model to flux vs time.
1177
1178 Parameters
1179 ----------
1180 diaObject : `dict`
1181 Summary object to store values in.
1182 diaSources : `pandas.DataFrame`
1183 DataFrame representing all diaSources associated with this
1184 diaObject.
1185 filterDiaSources : `pandas.DataFrame`
1186 DataFrame representing diaSources associated with this
1187 diaObject that are observed in the band pass ``band``.
1188 band : `str`
1189 Simple, string name of the filter for the flux being calculated.
1190 **kwargs
1191 Any additional keyword arguments that may be passed to the plugin.
1192 """
1193
1194 mName = "{}_psfFluxLinearSlope".format(band)
1195 if mName not in diaObjects.columns:
1196 diaObjects[mName] = np.nan
1197 bName = "{}_psfFluxLinearIntercept".format(band)
1198 if bName not in diaObjects.columns:
1199 diaObjects[bName] = np.nan
1200 dtype = diaObjects[mName].dtype
1201
1202 def _linearFit(df):
1203 tmpDf = df[~np.logical_or(
1204 np.isnan(df["psfFlux"]),
1205 np.logical_or(np.isnan(df["psfFluxErr"]),
1206 np.isnan(df["midpointMjdTai"])))]
1207 if len(tmpDf) < 2:
1208 return pd.Series({mName: np.nan, bName: np.nan})
1209 fluxes = tmpDf["psfFlux"].to_numpy()
1210 errors = tmpDf["psfFluxErr"].to_numpy()
1211 times = tmpDf["midpointMjdTai"].to_numpy()
1212 A = np.array([times / errors, 1 / errors]).transpose()
1213 m, b = lsq_linear(A, fluxes / errors).x
1214 return pd.Series({mName: m, bName: b}, dtype=dtype)
1215
1216 diaObjects.loc[:, [mName, bName]] = filterDiaSources.apply(_linearFit)
1217
1218
1220 pass
1221
1222
1223@register("ap_stetsonJ")
1225 """Compute the StetsonJ statistic on the DIA point source fluxes.
1226 """
1227
1228 ConfigClass = LinearFitDiaPsfFluxConfig
1229
1230 # Required input Cols
1231 inputCols = ["psfFluxMean"]
1232 # Output columns are created upon instantiation of the class.
1233 outputCols = ["psfFluxStetsonJ"]
1234 plugType = "multi"
1235 needsFilter = True
1236
1237 @classmethod
1239 return cls.FLUX_MOMENTS_CALCULATED
1240
1241 def calculate(self,
1242 diaObjects,
1243 diaSources,
1244 filterDiaSources,
1245 band,
1246 **kwargs):
1247 """Compute the StetsonJ statistic on the DIA point source fluxes.
1248
1249 Parameters
1250 ----------
1251 diaObject : `dict`
1252 Summary object to store values in.
1253 diaSources : `pandas.DataFrame`
1254 DataFrame representing all diaSources associated with this
1255 diaObject.
1256 filterDiaSources : `pandas.DataFrame`
1257 DataFrame representing diaSources associated with this
1258 diaObject that are observed in the band pass ``band``.
1259 band : `str`
1260 Simple, string name of the filter for the flux being calculated.
1261 **kwargs
1262 Any additional keyword arguments that may be passed to the plugin.
1263 """
1264 meanName = "{}_psfFluxMean".format(band)
1265
1266 def _stetsonJ(df):
1267 tmpDf = df[~np.logical_or(np.isnan(df["psfFlux"]),
1268 np.isnan(df["psfFluxErr"]))]
1269 if len(tmpDf) < 2:
1270 return np.nan
1271 fluxes = tmpDf["psfFlux"].to_numpy()
1272 errors = tmpDf["psfFluxErr"].to_numpy()
1273
1274 return self._stetson_J(
1275 fluxes,
1276 errors,
1277 diaObjects.at[tmpDf.diaObjectId.iat[0], meanName])
1278
1279 column = "{}_psfFluxStetsonJ".format(band)
1280 if column in diaObjects:
1281 dtype = diaObjects[column].dtype
1282 diaObjects.loc[:, column] = filterDiaSources.apply(_stetsonJ).astype(dtype)
1283 else:
1284 diaObjects.loc[:, column] = filterDiaSources.apply(_stetsonJ)
1285
1286 def _stetson_J(self, fluxes, errors, mean=None):
1287 """Compute the single band stetsonJ statistic.
1288
1289 Parameters
1290 ----------
1291 fluxes : `numpy.ndarray` (N,)
1292 Calibrated lightcurve flux values.
1293 errors : `numpy.ndarray` (N,)
1294 Errors on the calibrated lightcurve fluxes.
1295 mean : `float`
1296 Starting mean from previous plugin.
1297
1298 Returns
1299 -------
1300 stetsonJ : `float`
1301 stetsonJ statistic for the input fluxes and errors.
1302
1303 References
1304 ----------
1305 .. [1] Stetson, P. B., "On the Automatic Determination of Light-Curve
1306 Parameters for Cepheid Variables", PASP, 108, 851S, 1996
1307 """
1308 n_points = len(fluxes)
1309 flux_mean = self._stetson_mean(fluxes, errors, mean)
1310 delta_val = (
1311 np.sqrt(n_points / (n_points - 1)) * (fluxes - flux_mean) / errors)
1312 p_k = delta_val ** 2 - 1
1313
1314 return np.mean(np.sign(p_k) * np.sqrt(np.fabs(p_k)))
1315
1317 values,
1318 errors,
1319 mean=None,
1320 alpha=2.,
1321 beta=2.,
1322 n_iter=20,
1323 tol=1e-6):
1324 """Compute the stetson mean of the fluxes which down-weights outliers.
1325
1326 Weighted biased on an error weighted difference scaled by a constant
1327 (1/``a``) and raised to the power beta. Higher betas more harshly
1328 penalize outliers and ``a`` sets the number of sigma where a weighted
1329 difference of 1 occurs.
1330
1331 Parameters
1332 ----------
1333 values : `numpy.dnarray`, (N,)
1334 Input values to compute the mean of.
1335 errors : `numpy.ndarray`, (N,)
1336 Errors on the input values.
1337 mean : `float`
1338 Starting mean value or None.
1339 alpha : `float`
1340 Scalar down-weighting of the fractional difference. lower->more
1341 clipping. (Default value is 2.)
1342 beta : `float`
1343 Power law slope of the used to down-weight outliers. higher->more
1344 clipping. (Default value is 2.)
1345 n_iter : `int`
1346 Number of iterations of clipping.
1347 tol : `float`
1348 Fractional and absolute tolerance goal on the change in the mean
1349 before exiting early. (Default value is 1e-6)
1350
1351 Returns
1352 -------
1353 mean : `float`
1354 Weighted stetson mean result.
1355
1356 References
1357 ----------
1358 .. [1] Stetson, P. B., "On the Automatic Determination of Light-Curve
1359 Parameters for Cepheid Variables", PASP, 108, 851S, 1996
1360 """
1361 n_points = len(values)
1362 n_factor = np.sqrt(n_points / (n_points - 1))
1363 inv_var = 1 / errors ** 2
1364
1365 if mean is None:
1366 mean = np.average(values, weights=inv_var)
1367 for iter_idx in range(n_iter):
1368 chi = np.fabs(n_factor * (values - mean) / errors)
1369 tmp_mean = np.average(
1370 values,
1371 weights=inv_var / (1 + (chi / alpha) ** beta))
1372 diff = np.fabs(tmp_mean - mean)
1373 mean = tmp_mean
1374 if diff / mean < tol and diff < tol:
1375 break
1376 return mean
1377
1378
1380 pass
1381
1382
1383@register("ap_meanTotFlux")
1385 """Compute the weighted mean and mean error on the point source fluxes
1386 forced photometered at the DiaSource location in the calibrated image.
1387
1388 Additionally store number of usable data points.
1389 """
1390
1391 ConfigClass = WeightedMeanDiaPsfFluxConfig
1392 outputCols = ["scienceFluxMean", "scienceFluxMeanErr"]
1393 plugType = "multi"
1394 needsFilter = True
1395
1396 @classmethod
1399
1400 @catchWarnings(warns=["invalid value encountered",
1401 "divide by zero"])
1402 def calculate(self,
1403 diaObjects,
1404 diaSources,
1405 filterDiaSources,
1406 band,
1407 **kwargs):
1408 """Compute the weighted mean and mean error of the point source flux.
1409
1410 Parameters
1411 ----------
1412 diaObject : `dict`
1413 Summary object to store values in.
1414 diaSources : `pandas.DataFrame`
1415 DataFrame representing all diaSources associated with this
1416 diaObject.
1417 filterDiaSources : `pandas.DataFrame`
1418 DataFrame representing diaSources associated with this
1419 diaObject that are observed in the band pass ``band``.
1420 band : `str`
1421 Simple, string name of the filter for the flux being calculated.
1422 **kwargs
1423 Any additional keyword arguments that may be passed to the plugin.
1424 """
1425 totMeanName = "{}_scienceFluxMean".format(band)
1426 if totMeanName not in diaObjects.columns:
1427 diaObjects[totMeanName] = np.nan
1428 totErrName = "{}_scienceFluxMeanErr".format(band)
1429 if totErrName not in diaObjects.columns:
1430 diaObjects[totErrName] = np.nan
1431
1432 def _meanFlux(df):
1433 tmpDf = df[~np.logical_or(np.isnan(df["scienceFlux"]),
1434 np.isnan(df["scienceFluxErr"]))]
1435 tot_weight = np.nansum(1 / tmpDf["scienceFluxErr"] ** 2)
1436 fluxMean = np.nansum(tmpDf["scienceFlux"]
1437 / tmpDf["scienceFluxErr"] ** 2)
1438 fluxMean /= tot_weight
1439 fluxMeanErr = np.sqrt(1 / tot_weight)
1440
1441 return pd.Series({totMeanName: fluxMean,
1442 totErrName: fluxMeanErr})
1443
1444 df = filterDiaSources.apply(_meanFlux).astype(diaObjects.dtypes[[totMeanName, totErrName]])
1445 diaObjects.loc[:, [totMeanName, totErrName]] = df
1446
1447
1449 pass
1450
1451
1452@register("ap_sigmaTotFlux")
1454 """Compute scatter of diaSource fluxes.
1455 """
1456
1457 ConfigClass = SigmaDiaPsfFluxConfig
1458 # Output columns are created upon instantiation of the class.
1459 outputCols = ["scienceFluxSigma"]
1460 plugType = "multi"
1461 needsFilter = True
1462
1463 @classmethod
1466
1467 def calculate(self,
1468 diaObjects,
1469 diaSources,
1470 filterDiaSources,
1471 band,
1472 **kwargs):
1473 """Compute the sigma fluxes of the point source flux measured on the
1474 calibrated image.
1475
1476 Parameters
1477 ----------
1478 diaObject : `dict`
1479 Summary object to store values in.
1480 diaSources : `pandas.DataFrame`
1481 DataFrame representing all diaSources associated with this
1482 diaObject.
1483 filterDiaSources : `pandas.DataFrame`
1484 DataFrame representing diaSources associated with this
1485 diaObject that are observed in the band pass ``band``.
1486 band : `str`
1487 Simple, string name of the filter for the flux being calculated.
1488 **kwargs
1489 Any additional keyword arguments that may be passed to the plugin.
1490 """
1491 # Set "delta degrees of freedom (ddf)" to 1 to calculate the unbiased
1492 # estimator of scatter (i.e. 'N - 1' instead of 'N').
1493 column = "{}_scienceFluxSigma".format(band)
1494 if column in diaObjects:
1495 dtype = diaObjects[column].dtype
1496 diaObjects.loc[:, column] = filterDiaSources.scienceFlux.std().astype(dtype)
1497 else:
1498
1499 diaObjects.loc[:, column] = filterDiaSources.scienceFlux.std()
Point in an unspecified spherical coordinate system.
Definition SpherePoint.h:57
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
_stetson_mean(self, values, errors, mean=None, alpha=2., beta=2., n_iter=20, tol=1e-6)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
calculate(self, diaObjects, diaSources, filterDiaSources, band, **kwargs)
HtmPixelization provides HTM indexing of points and regions.
SpherePoint averageSpherePoint(std::vector< SpherePoint > const &coords)
Return the average of a list of coordinates.
compute_optimized_periodogram_grid(x0, oversampling_factor=5, nyquist_factor=100)
register(name, shouldApCorr=False, apCorrList=())