Coverage for tests / test_diaCalculationPlugins.py: 13%
419 statements
« prev ^ index » next coverage.py v7.13.5, created at 2026-05-06 08:36 +0000
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/>.
22from astropy.stats import median_absolute_deviation
23import numpy as np
24import pandas as pd
25from scipy.stats import skew
26import unittest
28from lsst.meas.base import (
29 MeanDiaPosition, MeanDiaPositionConfig,
30 HTMIndexDiaPosition, HTMIndexDiaPositionConfig,
31 NumDiaSourcesDiaPlugin, NumDiaSourcesDiaPluginConfig,
32 SimpleSourceFlagDiaPlugin, SimpleSourceFlagDiaPluginConfig,
33 WeightedMeanDiaPsfFlux, WeightedMeanDiaPsfFluxConfig,
34 PercentileDiaPsfFlux, PercentileDiaPsfFluxConfig,
35 SigmaDiaPsfFlux, SigmaDiaPsfFluxConfig,
36 Chi2DiaPsfFlux, Chi2DiaPsfFluxConfig,
37 MadDiaPsfFlux, MadDiaPsfFluxConfig,
38 SkewDiaPsfFlux, SkewDiaPsfFluxConfig,
39 MinMaxDiaPsfFlux, MinMaxDiaPsfFluxConfig,
40 MaxSlopeDiaPsfFlux, MaxSlopeDiaPsfFluxConfig,
41 ErrMeanDiaPsfFlux, ErrMeanDiaPsfFluxConfig,
42 LinearFitDiaPsfFlux, LinearFitDiaPsfFluxConfig,
43 StetsonJDiaPsfFlux, StetsonJDiaPsfFluxConfig,
44 WeightedMeanDiaTotFlux, WeightedMeanDiaTotFluxConfig,
45 SigmaDiaTotFlux, SigmaDiaTotFluxConfig,
46 LombScarglePeriodogram, LombScarglePeriodogramConfig,
47 LombScarglePeriodogramMulti, LombScarglePeriodogramMultiConfig,
48 UnphysicalDiaSourceSeparation)
49import lsst.utils.tests
52def run_single_plugin(diaObjectCat,
53 diaObjectId,
54 diaSourceCat,
55 band,
56 plugin):
57 """Wrapper for running single plugins.
59 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
61 Parameters
62 ----------
63 diaObjectCat : `pandas.DataFrame`
64 Input object catalog to store data into and read from.
65 diaSourcesCat : `pandas.DataFrame`
66 DiaSource catalog to read data from and groupby on.
67 fitlerName : `str`
68 String name of the filter to process.
69 plugin : `lsst.ap.association.DiaCalculationPlugin`
70 Plugin to run.
71 """
72 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
73 diaSourceCat.set_index(
74 ["diaObjectId", "band", "diaSourceId"],
75 inplace=True,
76 drop=False)
78 objDiaSources = diaSourceCat.loc[diaObjectId]
79 updatingFilterDiaSources = diaSourceCat.loc[
80 (diaObjectId, band), :
81 ]
83 plugin.calculate(diaObjects=diaObjectCat,
84 diaObjectId=diaObjectId,
85 diaSources=objDiaSources,
86 filterDiaSources=updatingFilterDiaSources,
87 band=band)
90def run_multi_plugin(diaObjectCat, diaSourceCat, band, plugin):
91 """Wrapper for running multi plugins.
93 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
95 Parameters
96 ----------
97 diaObjectCat : `pandas.DataFrame`
98 Input object catalog to store data into and read from.
99 diaSourcesCat : `pandas.DataFrame`
100 DiaSource catalog to read data from and groupby on.
101 filterName : `str`
102 String name of the filter to process.
103 plugin : `lsst.ap.association.DiaCalculationPlugin`
104 Plugin to run.
105 """
106 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
107 diaSourceCat.set_index(
108 ["diaObjectId", "band", "diaSourceId"],
109 inplace=True,
110 drop=False)
112 updatingFilterDiaSources = diaSourceCat.loc[
113 (slice(None), band), :
114 ]
116 diaSourcesGB = diaSourceCat.groupby(level=0)
117 filterDiaSourcesGB = updatingFilterDiaSources.groupby(level=0)
119 plugin.calculate(diaObjects=diaObjectCat,
120 diaSources=diaSourcesGB,
121 filterDiaSources=filterDiaSourcesGB,
122 band=band)
125def run_multiband_plugin(diaObjectCat, diaSourceCat, plugin):
126 """Wrapper for running multi plugins.
128 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
130 Parameters
131 ----------
132 diaObjectCat : `pandas.DataFrame`
133 Input object catalog to store data into and read from.
134 diaSourcesCat : `pandas.DataFrame`
135 DiaSource catalog to read data from and groupby on.
136 plugin : `lsst.ap.association.DiaCalculationPlugin`
137 Plugin to run.
138 """
139 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
140 diaSourceCat.set_index(
141 ["diaObjectId", "band", "diaSourceId"],
142 inplace=True,
143 drop=False)
145 diaSourcesGB = diaSourceCat.groupby(level=0)
147 plugin.calculate(diaObjects=diaObjectCat,
148 diaSources=diaSourcesGB,
149 )
152def make_diaObject_table(objId, plugin, default_value=np.nan, band=None):
153 """Create a minimal diaObject table with columns required for the plugin
155 Parameters
156 ----------
157 objId : `int`
158 The diaObjectId
159 plugin : `lsst.ap.association.DiaCalculationPlugin`
160 The plugin that will be run.
161 default_value : `float` or `int`, optional
162 Value to set new columns to.
163 band : `str`, optional
164 Band designation to append to the plugin columns.
166 Returns
167 -------
168 diaObjects : `pandas.DataFrame`
169 Output catalog with the required columns for the plugin.
170 """
171 # Add an extra empty diaObject here. This ensures that
172 # we properly test the source/object matching implicit
173 # in the plugin calculations.
174 diaObjects = {"diaObjectId": [objId, objId + 1]}
175 for col in plugin.outputCols:
176 if band is not None:
177 diaObjects[f"{band}_{col}"] = default_value
178 else:
179 diaObjects[col] = default_value
180 return pd.DataFrame(diaObjects)
183class TestMeanPosition(unittest.TestCase):
185 def testCalculate(self):
186 """Test mean position calculation.
187 """
188 n_sources = 10
189 objId = 0
191 # configure a 2 degree max separation
192 plug = MeanDiaPosition(MeanDiaPositionConfig(MaxAllowedDiaSourceSeparation=7200.0),
193 "ap_meanPosition",
194 None)
196 # Test expected means in RA.
197 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
198 diaSources = pd.DataFrame(data={"ra": np.linspace(-1, 1, n_sources),
199 "dec": np.zeros(n_sources),
200 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
201 "diaObjectId": n_sources * [objId],
202 "band": n_sources * ["g"],
203 "diaSourceId": np.arange(n_sources,
204 dtype=int)})
205 run_multi_plugin(diaObjects, diaSources, "g", plug)
207 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
208 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
210 # Test expected means in DEC.
211 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
212 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),
213 "dec": np.linspace(-1, 1, n_sources),
214 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
215 "diaObjectId": n_sources * [objId],
216 "band": n_sources * ["g"],
217 "diaSourceId": np.arange(n_sources,
218 dtype=int)})
219 run_multi_plugin(diaObjects, diaSources, "g", plug)
221 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
222 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
224 # Test failure mode RA is nan.
225 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
226 diaSources = pd.DataFrame(data={"ra": np.full(n_sources, np.nan),
227 "dec": np.zeros(n_sources),
228 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
229 "diaObjectId": n_sources * [objId],
230 "band": n_sources * ["g"],
231 "diaSourceId": np.arange(n_sources,
232 dtype=int)})
233 run_multi_plugin(diaObjects, diaSources, "g", plug)
235 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))
236 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))
238 # Test failure mode DEC is nan.
239 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
240 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),
241 "dec": np.full(n_sources, np.nan),
242 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
243 "diaObjectId": n_sources * [objId],
244 "band": n_sources * ["g"],
245 "diaSourceId": np.arange(n_sources,
246 dtype=int)})
247 run_multi_plugin(diaObjects, diaSources, "g", plug)
249 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))
250 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))
252 # configure the default 3 arcsecond separation
253 plug = MeanDiaPosition(MeanDiaPositionConfig(MaxAllowedDiaSourceSeparation=3.0),
254 "ap_meanPosition",
255 None)
257 # These 1 degree separations should raise
258 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
259 diaSources = pd.DataFrame(data={"ra": np.linspace(-1, 1, n_sources),
260 "dec": np.zeros(n_sources),
261 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
262 "diaObjectId": n_sources * [objId],
263 "band": n_sources * ["g"],
264 "diaSourceId": np.arange(n_sources,
265 dtype=int)})
266 with self.assertRaises(UnphysicalDiaSourceSeparation):
267 run_multi_plugin(diaObjects, diaSources, "g", plug)
269 # 1 arcsecond separations should not raise
270 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
271 diaSources = pd.DataFrame(data={"ra": np.linspace(-1/3600., 1/3600., n_sources),
272 "dec": np.zeros(n_sources),
273 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
274 "diaObjectId": n_sources * [objId],
275 "band": n_sources * ["g"],
276 "diaSourceId": np.arange(n_sources,
277 dtype=int)})
278 run_multi_plugin(diaObjects, diaSources, "g", plug)
281class TestHTMIndexPosition(unittest.TestCase):
283 def testCalculate(self):
284 """Test HTMPixel assignment calculation.
285 """
286 # Test expected pixelId at RA, DEC = 0
287 objId = 0
288 n_sources = 10
289 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
290 diaObjects.loc[objId, "ra"] = 0.
291 diaObjects.loc[objId, "dec"] = 0.
292 diaSources = pd.DataFrame(
293 data={"diaObjectId": n_sources * [objId],
294 "band": n_sources * ["g"],
295 "diaSourceId": np.arange(n_sources, dtype=int)})
296 plug = HTMIndexDiaPosition(HTMIndexDiaPositionConfig(),
297 "ap_HTMIndex",
298 None)
300 run_single_plugin(diaObjectCat=diaObjects,
301 diaObjectId=objId,
302 diaSourceCat=diaSources,
303 band="g",
304 plugin=plug)
305 self.assertEqual(diaObjects.at[objId, "pixelId"],
306 17042430230528)
308 # Test expected pixelId at some value of RA and DEC.
309 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
310 diaObjects.loc[objId, "ra"] = 45.37
311 diaObjects.loc[objId, "dec"] = 13.67
312 diaSources = pd.DataFrame(
313 data={"diaObjectId": n_sources * [objId],
314 "band": n_sources * ["g"],
315 "diaSourceId": np.arange(n_sources, dtype=int)})
316 run_single_plugin(diaObjectCat=diaObjects,
317 diaObjectId=objId,
318 diaSourceCat=diaSources,
319 band="g",
320 plugin=plug)
321 self.assertEqual(diaObjects.at[objId, "pixelId"],
322 17450571968473)
325class TestNDiaSourcesDiaPlugin(unittest.TestCase):
327 def testCalculate(self):
328 """Test that the number of DiaSources is correct.
329 """
331 for n_sources in [1, 8, 10]:
332 # Test expected number of sources per object.
333 objId = 0
334 diaSources = pd.DataFrame(
335 data={"diaObjectId": n_sources * [objId],
336 "band": n_sources * ["g"],
337 "diaSourceId": np.arange(n_sources, dtype=int)})
338 plug = NumDiaSourcesDiaPlugin(NumDiaSourcesDiaPluginConfig(),
339 "ap_nDiaSources",
340 None)
341 diaObjects = make_diaObject_table(objId, plug, default_value=int(-1))
342 run_multi_plugin(diaObjects, diaSources, "g", plug)
344 self.assertEqual(n_sources, diaObjects.at[objId, "nDiaSources"])
345 self.assertEqual(diaObjects["nDiaSources"].dtype, np.int64)
348class TestSimpleSourceFlagDiaPlugin(unittest.TestCase):
350 def testCalculate(self):
351 """Test that DiaObject flags are set.
352 """
353 objId = 0
354 n_sources = 10
356 # Test expected flags, no flags set.
357 diaSources = pd.DataFrame(
358 data={"diaObjectId": n_sources * [objId],
359 "band": n_sources * ["g"],
360 "diaSourceId": np.arange(n_sources, dtype=int),
361 "flags": np.zeros(n_sources, dtype=np.uint64)})
362 plug = SimpleSourceFlagDiaPlugin(SimpleSourceFlagDiaPluginConfig(),
363 "ap_diaObjectFlag",
364 None)
366 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
367 run_multi_plugin(diaObjects, diaSources, "g", plug)
368 self.assertEqual(diaObjects.at[objId, "flags"], 0)
369 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
371 # Test expected flags, all flags set.
372 diaSources = pd.DataFrame(
373 data={"diaObjectId": n_sources * [objId],
374 "band": n_sources * ["g"],
375 "diaSourceId": np.arange(n_sources, dtype=int),
376 "flags": np.ones(n_sources, dtype=np.uint64)})
377 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
378 run_multi_plugin(diaObjects, diaSources, "g", plug)
379 self.assertEqual(diaObjects.at[objId, "flags"], 1)
380 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
382 # Test expected flags, random flags.
383 diaSources = pd.DataFrame(
384 data={"diaObjectId": n_sources * [objId],
385 "band": n_sources * ["g"],
386 "diaSourceId": np.arange(n_sources, dtype=int),
387 "flags": np.random.randint(0, 2 ** 16, size=n_sources)})
389 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
390 run_multi_plugin(diaObjects, diaSources, "g", plug)
391 self.assertEqual(diaObjects.at[objId, "flags"], 1)
392 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
394 # Test expected flags, one flag set.
395 flag_array = np.zeros(n_sources, dtype=np.uint64)
396 flag_array[4] = 256
397 diaSources = pd.DataFrame(
398 data={"diaObjectId": n_sources * [objId],
399 "band": n_sources * ["g"],
400 "diaSourceId": np.arange(n_sources, dtype=int),
401 "flags": flag_array})
402 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
403 run_multi_plugin(diaObjects, diaSources, "g", plug)
404 self.assertEqual(diaObjects.at[objId, "flags"], 1)
405 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
408class TestWeightedMeanDiaPsfFlux(unittest.TestCase):
410 def testCalculate(self):
411 """Test mean value calculation.
412 """
413 n_sources = 10
414 objId = 0
416 # Test expected mean.
417 # In the first test, we have only one object.
418 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
419 diaSources = pd.DataFrame(
420 data={"diaObjectId": n_sources * [objId],
421 "band": n_sources * ["u"],
422 "diaSourceId": np.arange(n_sources, dtype=int),
423 "psfFlux": np.linspace(-1, 1, n_sources),
424 "psfFluxErr": np.ones(n_sources)})
426 plug = WeightedMeanDiaPsfFlux(WeightedMeanDiaPsfFluxConfig(),
427 "ap_meanFlux",
428 None)
429 run_multi_plugin(diaObjects, diaSources, "u", plug)
431 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMean"], 0.0)
432 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMeanErr"],
433 np.sqrt(1 / n_sources))
434 self.assertEqual(diaObjects.loc[objId, "u_psfFluxNdata"], n_sources)
435 # We expect this to be converted to float.
436 # TODO DM-53254: This should be an integer (and should be checked
437 # to be an integer).
438 self.assertEqual(diaObjects["u_psfFluxNdata"].dtype, np.float64)
440 # Test expected mean with a nan value.
441 # In the second test, we have two objects (one empty).
442 diaObjects = pd.DataFrame({"diaObjectId": [objId, objId + 1]})
443 fluxes = np.linspace(-1, 1, n_sources)
444 fluxes[4] = np.nan
445 diaSources = pd.DataFrame(
446 data={"diaObjectId": n_sources * [objId],
447 "band": n_sources * ["r"],
448 "diaSourceId": np.arange(n_sources, dtype=int),
449 "psfFlux": fluxes,
450 "psfFluxErr": np.ones(n_sources)})
451 run_multi_plugin(diaObjects, diaSources, "r", plug)
453 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMean"],
454 np.nanmean(fluxes))
455 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMeanErr"],
456 np.sqrt(1 / (n_sources - 1)))
457 self.assertEqual(diaObjects.loc[objId, "r_psfFluxNdata"], n_sources - 1)
458 # We expect this to be converted to float.
459 # TODO DM-53254: This should be an integer (and should be checked
460 # to be an integer).
461 self.assertEqual(diaObjects["r_psfFluxNdata"].dtype, np.float64)
464class TestPercentileDiaPsfFlux(unittest.TestCase):
466 def testCalculate(self):
467 """Test flux percentile calculation.
468 """
469 n_sources = 10
470 objId = 0
472 # Test expected percentile values.
473 fluxes = np.linspace(-1, 1, n_sources)
474 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
475 diaSources = pd.DataFrame(
476 data={"diaObjectId": n_sources * [objId],
477 "band": n_sources * ["u"],
478 "diaSourceId": np.arange(n_sources, dtype=int),
479 "psfFlux": fluxes,
480 "psfFluxErr": np.ones(n_sources)})
482 plug = PercentileDiaPsfFlux(PercentileDiaPsfFluxConfig(),
483 "ap_percentileFlux",
484 None)
485 run_multi_plugin(diaObjects, diaSources, "u", plug)
486 for pTile, testVal in zip(plug.config.percentiles,
487 np.nanpercentile(
488 fluxes,
489 plug.config.percentiles)):
490 self.assertAlmostEqual(
491 diaObjects.at[objId, "u_psfFluxPercentile{:02d}".format(pTile)],
492 testVal)
494 # Test expected percentile values with a nan value.
495 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
496 fluxes[4] = np.nan
497 diaSources = pd.DataFrame(
498 data={"diaObjectId": n_sources * [objId],
499 "band": n_sources * ["r"],
500 "diaSourceId": np.arange(n_sources, dtype=int),
501 "psfFlux": fluxes,
502 "psfFluxErr": np.ones(n_sources)})
503 run_multi_plugin(diaObjects, diaSources, "r", plug)
504 for pTile, testVal in zip(plug.config.percentiles,
505 np.nanpercentile(
506 fluxes,
507 plug.config.percentiles)):
508 self.assertAlmostEqual(
509 diaObjects.at[objId, "r_psfFluxPercentile{:02d}".format(pTile)],
510 testVal)
513class TestSigmaDiaPsfFlux(unittest.TestCase):
515 def testCalculate(self):
516 """Test flux scatter calculation.
517 """
518 n_sources = 10
519 objId = 0
521 # Test expected sigma scatter of fluxes.
522 fluxes = np.linspace(-1, 1, n_sources)
523 diaSources = pd.DataFrame(
524 data={"diaObjectId": n_sources * [objId],
525 "band": n_sources * ["u"],
526 "diaSourceId": np.arange(n_sources, dtype=int),
527 "psfFlux": fluxes,
528 "psfFluxErr": np.ones(n_sources)})
530 plug = SigmaDiaPsfFlux(SigmaDiaPsfFluxConfig(),
531 "ap_sigmaFlux",
532 None)
533 diaObjects = make_diaObject_table(objId, plug, band='u')
534 run_multi_plugin(diaObjects, diaSources, "u", plug)
535 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxSigma"],
536 np.nanstd(fluxes, ddof=1))
538 # test one input, returns nan.
539 diaSources = pd.DataFrame(
540 data={"diaObjectId": 1 * [objId],
541 "band": 1 * ["g"],
542 "diaSourceId": [0],
543 "psfFlux": [fluxes[0]],
544 "psfFluxErr": [1.]})
546 diaObjects = make_diaObject_table(objId, plug, band='g')
547 run_multi_plugin(diaObjects, diaSources, "g", plug)
548 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxSigma"]))
550 # Test expected sigma scatter of fluxes with a nan value.
551 fluxes[4] = np.nan
552 diaSources = pd.DataFrame(
553 data={"diaObjectId": n_sources * [objId],
554 "band": n_sources * ["r"],
555 "diaSourceId": np.arange(n_sources, dtype=int),
556 "psfFlux": fluxes,
557 "psfFluxErr": np.ones(n_sources)})
559 diaObjects = make_diaObject_table(objId, plug, band='r')
560 run_multi_plugin(diaObjects, diaSources, "r", plug)
561 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxSigma"],
562 np.nanstd(fluxes, ddof=1))
565class TestChi2DiaPsfFlux(unittest.TestCase):
567 def testCalculate(self):
568 """Test flux chi2 calculation.
569 """
570 n_sources = 10
571 objId = 0
573 # Test expected chi^2 value.
574 fluxes = np.linspace(-1, 1, n_sources)
575 diaObjects = pd.DataFrame({"diaObjectId": [objId],
576 "u_psfFluxMean": [0.0]})
577 diaSources = pd.DataFrame(
578 data={"diaObjectId": n_sources * [objId],
579 "band": n_sources * ["u"],
580 "diaSourceId": np.arange(n_sources, dtype=int),
581 "psfFlux": fluxes,
582 "psfFluxErr": np.ones(n_sources)})
584 plug = Chi2DiaPsfFlux(Chi2DiaPsfFluxConfig(),
585 "ap_chi2Flux",
586 None)
587 run_multi_plugin(diaObjects, diaSources, "u", plug)
588 self.assertAlmostEqual(
589 diaObjects.loc[objId, "u_psfFluxChi2"],
590 np.nansum(((diaSources["psfFlux"]
591 - np.nanmean(diaSources["psfFlux"]))
592 / diaSources["psfFluxErr"]) ** 2))
594 # Test expected chi^2 value with a nan value set.
595 fluxes[4] = np.nan
596 diaObjects = pd.DataFrame({"diaObjectId": [objId],
597 "r_psfFluxMean": [np.nanmean(fluxes)]})
598 diaSources = pd.DataFrame(
599 data={"diaObjectId": n_sources * [objId],
600 "band": n_sources * ["r"],
601 "diaSourceId": np.arange(n_sources, dtype=int),
602 "psfFlux": fluxes,
603 "psfFluxErr": np.ones(n_sources)})
604 run_multi_plugin(diaObjects, diaSources, "r", plug)
605 self.assertAlmostEqual(
606 diaObjects.loc[objId, "r_psfFluxChi2"],
607 np.nansum(((diaSources["psfFlux"]
608 - np.nanmean(diaSources["psfFlux"]))
609 / diaSources["psfFluxErr"]) ** 2))
612class TestMadDiaPsfFlux(unittest.TestCase):
614 def testCalculate(self):
615 """Test flux median absolute deviation calculation.
616 """
617 n_sources = 10
618 objId = 0
620 # Test expected MAD value.
621 fluxes = np.linspace(-1, 1, n_sources)
622 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
623 diaSources = pd.DataFrame(
624 data={"diaObjectId": n_sources * [objId],
625 "band": n_sources * ["u"],
626 "diaSourceId": np.arange(n_sources, dtype=int),
627 "psfFlux": fluxes,
628 "psfFluxErr": np.ones(n_sources)})
630 plug = MadDiaPsfFlux(MadDiaPsfFluxConfig(),
631 "ap_madFlux",
632 None)
633 run_multi_plugin(diaObjects, diaSources, "u", plug)
634 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMAD"],
635 median_absolute_deviation(fluxes,
636 ignore_nan=True))
638 # Test expected MAD value with a nan set.
639 fluxes[4] = np.nan
640 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
641 diaSources = pd.DataFrame(
642 data={"diaObjectId": n_sources * [objId],
643 "band": n_sources * ["r"],
644 "diaSourceId": np.arange(n_sources, dtype=int),
645 "psfFlux": fluxes,
646 "psfFluxErr": np.ones(n_sources)})
647 run_multi_plugin(diaObjects, diaSources, "r", plug)
648 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMAD"],
649 median_absolute_deviation(fluxes,
650 ignore_nan=True))
653class TestSkewDiaPsfFlux(unittest.TestCase):
655 def testCalculate(self):
656 """Test flux skew calculation.
657 """
658 n_sources = 10
659 objId = 0
661 # Test expected skew value.
662 fluxes = np.linspace(-1, 1, n_sources)
663 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
664 diaSources = pd.DataFrame(
665 data={"diaObjectId": n_sources * [objId],
666 "band": n_sources * ["u"],
667 "diaSourceId": np.arange(n_sources, dtype=int),
668 "psfFlux": fluxes,
669 "psfFluxErr": np.ones(n_sources)})
671 plug = SkewDiaPsfFlux(SkewDiaPsfFluxConfig(),
672 "ap_skewFlux",
673 None)
674 run_multi_plugin(diaObjects, diaSources, "u", plug)
675 self.assertAlmostEqual(
676 diaObjects.loc[objId, "u_psfFluxSkew"],
677 skew_wrapper(fluxes))
679 # Test expected skew value with a nan set.
680 fluxes[4] = np.nan
681 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
682 diaSources = pd.DataFrame(
683 data={"diaObjectId": n_sources * [objId],
684 "band": n_sources * ["r"],
685 "diaSourceId": np.arange(n_sources, dtype=int),
686 "psfFlux": fluxes,
687 "psfFluxErr": np.ones(n_sources)})
688 run_multi_plugin(diaObjects, diaSources, "r", plug)
690 self.assertAlmostEqual(
691 diaObjects.at[objId, "r_psfFluxSkew"],
692 skew_wrapper(fluxes))
695class TestMinMaxDiaPsfFlux(unittest.TestCase):
697 def testCalculate(self):
698 """Test flux min/max calculation.
699 """
700 n_sources = 10
701 objId = 0
703 # Test expected MinMax fluxes.
704 fluxes = np.linspace(-1, 1, n_sources)
705 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
706 diaSources = pd.DataFrame(
707 data={"diaObjectId": n_sources * [objId],
708 "band": n_sources * ["u"],
709 "diaSourceId": np.arange(n_sources, dtype=int),
710 "psfFlux": fluxes,
711 "psfFluxErr": np.ones(n_sources)})
713 plug = MinMaxDiaPsfFlux(MinMaxDiaPsfFluxConfig(),
714 "ap_minMaxFlux",
715 None)
716 run_multi_plugin(diaObjects, diaSources, "u", plug)
717 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMin"], -1)
718 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMax"], 1)
720 # Test expected MinMax fluxes with a nan set.
721 fluxes[4] = np.nan
722 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
723 diaSources = pd.DataFrame(
724 data={"diaObjectId": n_sources * [objId],
725 "band": n_sources * ["r"],
726 "diaSourceId": np.arange(n_sources, dtype=int),
727 "psfFlux": fluxes,
728 "psfFluxErr": np.ones(n_sources)})
729 run_multi_plugin(diaObjects, diaSources, "r", plug)
730 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMin"], -1)
731 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMax"], 1)
734class TestMaxSlopeDiaPsfFlux(unittest.TestCase):
736 def testCalculate(self):
737 """Test flux maximum slope.
738 """
739 n_sources = 10
740 objId = 0
742 # Test max slope value.
743 fluxes = np.linspace(-1, 1, n_sources)
744 times = np.concatenate([np.linspace(0, 1, n_sources)[:-1], [1 - 1/90]])
745 diaSources = pd.DataFrame(
746 data={"diaObjectId": n_sources * [objId],
747 "band": n_sources * ["u"],
748 "diaSourceId": np.arange(n_sources, dtype=int),
749 "psfFlux": fluxes,
750 "psfFluxErr": np.ones(n_sources),
751 "midpointMjdTai": times})
753 plug = MaxSlopeDiaPsfFlux(MaxSlopeDiaPsfFluxConfig(),
754 "ap_maxSlopeFlux",
755 None)
756 diaObjects = make_diaObject_table(objId, plug, band='u')
757 run_multi_plugin(diaObjects, diaSources, "u", plug)
758 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMaxSlope"], 2 + 2/9)
760 # Test max slope value returns nan on 1 input.
761 diaSources = pd.DataFrame(
762 data={"diaObjectId": 1 * [objId],
763 "band": 1 * ["g"],
764 "diaSourceId": np.arange(1, dtype=int),
765 "psfFlux": fluxes[0],
766 "psfFluxErr": np.ones(1),
767 "midpointMjdTai": times[0]})
768 diaObjects = make_diaObject_table(objId, plug, band='g')
769 run_multi_plugin(diaObjects, diaSources, "g", plug)
770 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxMaxSlope"]))
772 # Test max slope value inputing nan values.
773 fluxes[4] = np.nan
774 times[7] = np.nan
775 diaSources = pd.DataFrame(
776 data={"diaObjectId": n_sources * [objId],
777 "band": n_sources * ["r"],
778 "diaSourceId": np.arange(n_sources, dtype=int),
779 "psfFlux": fluxes,
780 "psfFluxErr": np.ones(n_sources),
781 "midpointMjdTai": times})
782 diaObjects = make_diaObject_table(objId, plug, band='r')
783 run_multi_plugin(diaObjects, diaSources, "r", plug)
784 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMaxSlope"], 2 + 2 / 9)
787class TestErrMeanDiaPsfFlux(unittest.TestCase):
789 def testCalculate(self):
790 """Test error mean calculation.
791 """
792 n_sources = 10
793 objId = 0
795 # Test mean of the errors.
796 fluxes = np.linspace(-1, 1, n_sources)
797 errors = np.linspace(1, 2, n_sources)
798 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
799 diaSources = pd.DataFrame(
800 data={"diaObjectId": n_sources * [objId],
801 "band": n_sources * ["u"],
802 "diaSourceId": np.arange(n_sources, dtype=int),
803 "psfFlux": fluxes,
804 "psfFluxErr": errors})
806 plug = ErrMeanDiaPsfFlux(ErrMeanDiaPsfFluxConfig(),
807 "ap_errMeanFlux",
808 None)
809 run_multi_plugin(diaObjects, diaSources, "u", plug)
810 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxErrMean"],
811 np.nanmean(errors).astype(np.float32))
813 # Test mean of the errors with input nan value.
814 errors[4] = np.nan
815 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
816 diaSources = pd.DataFrame(
817 data={"diaObjectId": n_sources * [objId],
818 "band": n_sources * ["r"],
819 "diaSourceId": np.arange(n_sources, dtype=int),
820 "psfFlux": fluxes,
821 "psfFluxErr": errors})
822 run_multi_plugin(diaObjects, diaSources, "r", plug)
823 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxErrMean"],
824 np.nanmean(errors).astype(np.float32))
827class TestLinearFitDiaPsfFlux(unittest.TestCase):
829 def testCalculate(self):
830 """Test a linear fit to flux vs time.
831 """
832 n_sources = 10
833 objId = 0
835 # Test best fit linear model.
836 fluxes = np.linspace(-1, 1, n_sources)
837 errors = np.linspace(1, 2, n_sources)
838 times = np.linspace(0, 1, n_sources)
839 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
840 diaSources = pd.DataFrame(
841 data={"diaObjectId": n_sources * [objId],
842 "band": n_sources * ["u"],
843 "diaSourceId": np.arange(n_sources, dtype=int),
844 "psfFlux": fluxes,
845 "psfFluxErr": errors,
846 "midpointMjdTai": times})
848 plug = LinearFitDiaPsfFlux(LinearFitDiaPsfFluxConfig(),
849 "ap_LinearFit",
850 None)
851 run_multi_plugin(diaObjects, diaSources, "u", plug)
852 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearSlope"],
853 2.)
854 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearIntercept"],
855 -1.)
857 # Test best fit linear model with input nans.
858 fluxes[7] = np.nan
859 errors[4] = np.nan
860 times[2] = np.nan
861 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
862 diaSources = pd.DataFrame(
863 data={"diaObjectId": n_sources * [objId],
864 "band": n_sources * ["r"],
865 "diaSourceId": np.arange(n_sources, dtype=int),
866 "psfFlux": fluxes,
867 "psfFluxErr": errors,
868 "midpointMjdTai": times})
869 run_multi_plugin(diaObjects, diaSources, "r", plug)
870 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearSlope"], 2.)
871 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearIntercept"],
872 -1.)
875class TestStetsonJDiaPsfFlux(unittest.TestCase):
877 def testCalculate(self):
878 """Test the stetsonJ statistic.
879 """
880 n_sources = 10
881 objId = 0
883 # Test stetsonJ calculation.
884 fluxes = np.linspace(-1, 1, n_sources)
885 errors = np.ones(n_sources)
886 diaObjects = pd.DataFrame({"diaObjectId": [objId],
887 "u_psfFluxMean": [np.nanmean(fluxes)]})
888 diaSources = pd.DataFrame(
889 data={"diaObjectId": n_sources * [objId],
890 "band": n_sources * ["u"],
891 "diaSourceId": np.arange(n_sources, dtype=int),
892 "psfFlux": fluxes,
893 "psfFluxErr": errors})
895 plug = StetsonJDiaPsfFlux(StetsonJDiaPsfFluxConfig(),
896 "ap_StetsonJ",
897 None)
898 run_multi_plugin(diaObjects, diaSources, "u", plug)
899 # Expected StetsonJ for the values created. Confirmed using Cesimum's
900 # implementation. http://github.com/cesium-ml/cesium
901 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxStetsonJ"],
902 -0.5958393936080928)
904 # Test stetsonJ calculation returns nan on single input.
905 diaObjects = pd.DataFrame({"diaObjectId": [objId],
906 "g_psfFluxMean": [np.nanmean(fluxes)]})
907 diaSources = pd.DataFrame(
908 data={"diaObjectId": 1 * [objId],
909 "band": 1 * ["g"],
910 "diaSourceId": np.arange(1, dtype=int),
911 "psfFlux": fluxes[0],
912 "psfFluxErr": errors[0]})
913 run_multi_plugin(diaObjects, diaSources, "g", plug)
914 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxStetsonJ"]))
916 # Test stetsonJ calculation returns when nans are input.
917 fluxes[7] = np.nan
918 errors[4] = np.nan
919 nonNanMask = np.logical_and(~np.isnan(fluxes),
920 ~np.isnan(errors))
921 diaObjects = pd.DataFrame(
922 {"diaObjectId": [objId],
923 "r_psfFluxMean": [np.average(fluxes[nonNanMask],
924 weights=errors[nonNanMask])]})
925 diaSources = pd.DataFrame(
926 data={"diaObjectId": n_sources * [objId],
927 "band": n_sources * ["r"],
928 "diaSourceId": np.arange(n_sources, dtype=int),
929 "psfFlux": fluxes,
930 "psfFluxErr": errors})
931 run_multi_plugin(diaObjects, diaSources, "r", plug)
932 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxStetsonJ"],
933 -0.5412797916187173)
936class TestWeightedMeanDiaTotFlux(unittest.TestCase):
938 def testCalculate(self):
939 """Test mean value calculation.
940 """
941 n_sources = 10
942 objId = 0
944 # Test test mean on scienceFlux.
945 fluxes = np.linspace(-1, 1, n_sources)
946 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
947 diaSources = pd.DataFrame(
948 data={"diaObjectId": n_sources * [objId],
949 "band": n_sources * ["u"],
950 "diaSourceId": np.arange(n_sources, dtype=int),
951 "scienceFlux": fluxes,
952 "scienceFluxErr": np.ones(n_sources)})
954 plug = WeightedMeanDiaTotFlux(WeightedMeanDiaTotFluxConfig(),
955 "ap_meanTotFlux",
956 None)
957 run_multi_plugin(diaObjects, diaSources, "u", plug)
959 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMean"], 0.0)
960 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMeanErr"],
961 np.sqrt(1 / n_sources))
963 # Test test mean on scienceFlux with input nans
964 fluxes[4] = np.nan
965 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
966 diaSources = pd.DataFrame(
967 data={"diaObjectId": n_sources * [objId],
968 "band": n_sources * ["r"],
969 "diaSourceId": np.arange(n_sources, dtype=int),
970 "scienceFlux": fluxes,
971 "scienceFluxErr": np.ones(n_sources)})
972 run_multi_plugin(diaObjects, diaSources, "r", plug)
974 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMean"],
975 np.nanmean(fluxes))
976 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMeanErr"],
977 np.sqrt(1 / (n_sources - 1)))
980def generatePeriodicData(n=10, period=10):
981 """Generate noisy, sinusoidally-varying periodic data for testing Lomb-
982 Scargle Periodogram.
984 The returned fluxes will have, within the errors, the passed-in period and
985 a power close to 1, because the fluxes are purely sinusoidal.
987 Parameters
988 ----------
989 n : int
990 Number of data points to generate.
991 period : float
992 Period of the periodic signal.
994 Returns
995 -------
996 t : np.ndarray
997 Time values.
998 y_obs : np.ndarray
999 Observed flux values.
1000 """
1001 np.random.seed(42)
1003 t = np.linspace(-2*np.pi, 2*np.pi, n) + 100*np.random.random(n)
1004 y = 10 + np.sin(2 * np.pi * t / period)
1005 y_obs = np.random.normal(y, 0.001)
1007 return t, y_obs
1010class TestMultiLombScarglePeriodogram(lsst.utils.tests.TestCase):
1012 def testCalculate(self):
1013 """Test Mulitband Lomb Scargle Periodogram."""
1014 n_sources = 10
1015 objId = 0
1017 # Create synthetic multi-band data
1018 times, fluxes = generatePeriodicData(n_sources, period=10)
1019 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1020 diaSources = pd.DataFrame(
1021 data={"diaObjectId": n_sources * [objId],
1022 "band": n_sources//2 * ["u"] + n_sources//2 * ["g"],
1023 "diaSourceId": np.arange(n_sources, dtype=int),
1024 "midpointMjdTai": times,
1025 "psfFlux": fluxes,
1026 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1028 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),
1029 "ap_lombScarglePeriodogramMulti",
1030 None)
1032 run_multiband_plugin(diaObjects, diaSources, plugin)
1033 self.assertAlmostEqual(diaObjects.at[objId, "multiPeriod"], 10, delta=0.04)
1034 self.assertAlmostEqual(diaObjects.at[objId, "multiPower"], 1, delta=1e-2)
1035 # This implementation of LS returns a normalized power < 1.
1036 self.assertLess(diaObjects.at[objId, "multiPower"], 1)
1037 self.assertAlmostEqual(diaObjects.at[objId, "multiFap"], 0, delta=0.04)
1038 # Note: The below values are empirical, but seem reasonable, and
1039 # test that we get values for each band.
1040 self.assertAlmostEqual(diaObjects.at[objId, "u_multiAmp"], 0.029, delta=0.01)
1041 self.assertAlmostEqual(diaObjects.at[objId, "g_multiAmp"], 0.029, delta=0.01)
1042 self.assertAlmostEqual(diaObjects.at[objId, "u_multiPhase"], -2.0, delta=0.2)
1043 self.assertAlmostEqual(diaObjects.at[objId, "g_multiPhase"], 1.0, delta=0.1)
1045 def testCalculateTwoSources(self):
1046 """Test Mulitband Lomb Scargle Periodogram with 2 sources (minimum
1047 detections = 5), which will result in NaN output."""
1048 objId = 0
1049 n_sources = 2
1050 times, fluxes = generatePeriodicData(n_sources, period=10)
1051 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1052 diaSources = pd.DataFrame(
1053 data={"diaObjectId": n_sources * [objId],
1054 "band": n_sources * ["u"],
1055 "diaSourceId": np.arange(n_sources, dtype=int),
1056 "midpointMjdTai": times,
1057 "psfFlux": fluxes,
1058 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1060 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),
1061 "ap_lombScarglePeriodogramMulti",
1062 None)
1064 run_multi_plugin(diaObjects, diaSources, "u", plugin)
1065 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPeriod"]))
1066 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPower"]))
1067 self.assertTrue(np.isnan(diaObjects.at[objId, "multiFap"]))
1070class TestLombScarglePeriodogram(lsst.utils.tests.TestCase):
1072 def testCalculate(self):
1073 """Test Lomb Scargle Periodogram."""
1074 n_sources = 10
1075 objId = 0
1077 # Test period calculation.
1078 times, fluxes = generatePeriodicData(n_sources, period=10)
1079 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1080 diaSources = pd.DataFrame(
1081 data={"diaObjectId": n_sources * [objId],
1082 "band": n_sources * ["u"],
1083 "diaSourceId": np.arange(n_sources, dtype=int),
1084 "midpointMjdTai": times,
1085 "psfFlux": fluxes,
1086 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1088 plugin = LombScarglePeriodogram(LombScarglePeriodogramConfig(),
1089 "ap_lombScarglePeriodogram",
1090 None)
1092 run_multi_plugin(diaObjects, diaSources, "u", plugin)
1093 self.assertAlmostEqual(diaObjects.at[objId, "u_period"], 10, delta=0.04)
1094 # This implementation of LS returns a normalized power < 1.
1095 self.assertAlmostEqual(diaObjects.at[objId, "u_power"], 1, delta=1e-2)
1096 self.assertLess(diaObjects.at[objId, "u_power"], 1)
1098 # Test that we get the same result with a NaN flux.
1099 diaSources.loc[4, "psfFlux"] = np.nan
1100 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1101 diaSources = pd.DataFrame(
1102 data={"diaObjectId": n_sources * [objId],
1103 "band": n_sources * ["r"],
1104 "diaSourceId": np.arange(n_sources, dtype=int),
1105 "midpointMjdTai": times,
1106 "psfFlux": fluxes,
1107 "psfFluxErr": np.ones(n_sources)})
1108 run_multi_plugin(diaObjects, diaSources, "r", plugin)
1109 self.assertAlmostEqual(diaObjects.at[objId, "r_period"], 10, delta=0.04)
1110 self.assertAlmostEqual(diaObjects.at[objId, "r_power"], 1, delta=1e-2)
1111 # This implementation of LS returns a normalized power < 1.
1112 self.assertLess(diaObjects.at[objId, "r_power"], 1)
1115class TestSigmaDiaTotFlux(unittest.TestCase):
1117 def testCalculate(self):
1118 """Test flux scatter calculation.
1119 """
1120 n_sources = 10
1121 objId = 0
1123 # Test test scatter on scienceFlux.
1124 fluxes = np.linspace(-1, 1, n_sources)
1125 diaSources = pd.DataFrame(
1126 data={"diaObjectId": n_sources * [objId],
1127 "band": n_sources * ["u"],
1128 "diaSourceId": np.arange(n_sources, dtype=int),
1129 "scienceFlux": fluxes,
1130 "scienceFluxErr": np.ones(n_sources)})
1132 plug = SigmaDiaTotFlux(SigmaDiaTotFluxConfig(),
1133 "ap_sigmaTotFlux",
1134 None)
1135 diaObjects = make_diaObject_table(objId, plug, band='u')
1136 run_multi_plugin(diaObjects, diaSources, "u", plug)
1137 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxSigma"],
1138 np.nanstd(fluxes, ddof=1))
1140 # Test test scatter on scienceFlux returns nan on 1 input.
1141 diaSources = pd.DataFrame(
1142 data={"diaObjectId": 1 * [objId],
1143 "band": 1 * ["g"],
1144 "diaSourceId": np.arange(1, dtype=int),
1145 "scienceFlux": fluxes[0],
1146 "scienceFluxErr": np.ones(1)})
1147 diaObjects = make_diaObject_table(objId, plug, band='g')
1148 run_multi_plugin(diaObjects, diaSources, "g", plug)
1149 self.assertTrue(np.isnan(diaObjects.at[objId, "g_scienceFluxSigma"]))
1151 # Test test scatter on scienceFlux takes input nans.
1152 fluxes[4] = np.nan
1153 diaSources = pd.DataFrame(
1154 data={"diaObjectId": n_sources * [objId],
1155 "band": n_sources * ["r"],
1156 "diaSourceId": np.arange(n_sources, dtype=int),
1157 "scienceFlux": fluxes,
1158 "scienceFluxErr": np.ones(n_sources)})
1159 diaObjects = make_diaObject_table(objId, plug, band='r')
1160 run_multi_plugin(diaObjects, diaSources, "r", plug)
1161 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxSigma"],
1162 np.nanstd(fluxes, ddof=1))
1165def skew_wrapper(values):
1166 """Compute scipy skew, omitting nans.
1168 This version works with both scipy<1.9 (where it erroneously returns a
1169 masked array) and scipy>=1.9 (where it correctly returns a float).
1171 Parameters
1172 ----------
1173 values : `np.ndarray`
1175 Returns
1176 -------
1177 skew_value : `float`
1178 """
1179 value = skew(values, bias=False, nan_policy="omit")
1180 if isinstance(value, np.ma.masked_array):
1181 return value.data
1182 else:
1183 return value
1186class MemoryTester(lsst.utils.tests.MemoryTestCase):
1187 pass
1190def setup_module(module):
1191 lsst.utils.tests.init()
1194if __name__ == "__main__": 1194 ↛ 1195line 1194 didn't jump to line 1195 because the condition on line 1194 was never true
1195 lsst.utils.tests.init()
1196 unittest.main()