Coverage for tests/test_diaCalculationPlugins.py: 99%
640 statements
« prev ^ index » next coverage.py v7.16.2, created at 2026-11-06 01:51 -0800
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/>.
22import warnings
24from astropy.stats import median_absolute_deviation
25import numpy as np
26import pandas as pd
27from scipy.stats import skew
28import unittest
30from lsst.meas.base import (
31 MeanDiaPosition, MeanDiaPositionConfig,
32 FirstLastDiaSourceTime, FirstLastDiaSourceTimeConfig,
33 HTMIndexDiaPosition, HTMIndexDiaPositionConfig,
34 NumDiaSourcesDiaPlugin, NumDiaSourcesDiaPluginConfig,
35 SimpleSourceFlagDiaPlugin, SimpleSourceFlagDiaPluginConfig,
36 WeightedMeanDiaPsfFlux, WeightedMeanDiaPsfFluxConfig,
37 PercentileDiaPsfFlux, PercentileDiaPsfFluxConfig,
38 SigmaDiaPsfFlux, SigmaDiaPsfFluxConfig,
39 Chi2DiaPsfFlux, Chi2DiaPsfFluxConfig,
40 MadDiaPsfFlux, MadDiaPsfFluxConfig,
41 SkewDiaPsfFlux, SkewDiaPsfFluxConfig,
42 MinMaxDiaPsfFlux, MinMaxDiaPsfFluxConfig,
43 MaxSlopeDiaPsfFlux, MaxSlopeDiaPsfFluxConfig,
44 ErrMeanDiaPsfFlux, ErrMeanDiaPsfFluxConfig,
45 LinearFitDiaPsfFlux, LinearFitDiaPsfFluxConfig,
46 StetsonJDiaPsfFlux, StetsonJDiaPsfFluxConfig,
47 WeightedMeanDiaTotFlux, WeightedMeanDiaTotFluxConfig,
48 SigmaDiaTotFlux, SigmaDiaTotFluxConfig,
49 LombScarglePeriodogram, LombScarglePeriodogramConfig,
50 LombScarglePeriodogramMulti, LombScarglePeriodogramMultiConfig,
51 UnphysicalDiaSourceSeparation)
52import lsst.utils.tests
55def run_single_plugin(diaObjectCat,
56 diaObjectId,
57 diaSourceCat,
58 band,
59 plugin):
60 """Wrapper for running single plugins.
62 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
64 Parameters
65 ----------
66 diaObjectCat : `pandas.DataFrame`
67 Input object catalog to store data into and read from.
68 diaSourcesCat : `pandas.DataFrame`
69 DiaSource catalog to read data from and groupby on.
70 fitlerName : `str`
71 String name of the filter to process.
72 plugin : `lsst.ap.association.DiaCalculationPlugin`
73 Plugin to run.
74 """
75 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
76 diaSourceCat.set_index(
77 ["diaObjectId", "band", "diaSourceId"],
78 inplace=True,
79 drop=False)
81 objDiaSources = diaSourceCat.loc[diaObjectId]
82 updatingFilterDiaSources = diaSourceCat.loc[
83 (diaObjectId, band), :
84 ]
86 plugin.calculate(diaObjects=diaObjectCat,
87 diaObjectId=diaObjectId,
88 diaSources=objDiaSources,
89 filterDiaSources=updatingFilterDiaSources,
90 band=band)
93def run_multi_plugin(diaObjectCat, diaSourceCat, band, plugin):
94 """Wrapper for running multi plugins.
96 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
98 Parameters
99 ----------
100 diaObjectCat : `pandas.DataFrame`
101 Input object catalog to store data into and read from.
102 diaSourcesCat : `pandas.DataFrame`
103 DiaSource catalog to read data from and groupby on.
104 filterName : `str`
105 String name of the filter to process.
106 plugin : `lsst.ap.association.DiaCalculationPlugin`
107 Plugin to run.
108 """
109 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
110 diaSourceCat.set_index(
111 ["diaObjectId", "band", "diaSourceId"],
112 inplace=True,
113 drop=False)
115 updatingFilterDiaSources = diaSourceCat.loc[
116 (slice(None), band), :
117 ]
119 diaSourcesGB = diaSourceCat.groupby(level=0)
120 filterDiaSourcesGB = updatingFilterDiaSources.groupby(level=0)
122 plugin.calculate(diaObjects=diaObjectCat,
123 diaSources=diaSourcesGB,
124 filterDiaSources=filterDiaSourcesGB,
125 band=band)
128def run_multiband_plugin(diaObjectCat, diaSourceCat, plugin):
129 """Wrapper for running multi plugins.
131 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`
133 Parameters
134 ----------
135 diaObjectCat : `pandas.DataFrame`
136 Input object catalog to store data into and read from.
137 diaSourcesCat : `pandas.DataFrame`
138 DiaSource catalog to read data from and groupby on.
139 plugin : `lsst.ap.association.DiaCalculationPlugin`
140 Plugin to run.
141 """
142 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)
143 diaSourceCat.set_index(
144 ["diaObjectId", "band", "diaSourceId"],
145 inplace=True,
146 drop=False)
148 diaSourcesGB = diaSourceCat.groupby(level=0)
150 plugin.calculate(diaObjects=diaObjectCat,
151 diaSources=diaSourcesGB,
152 )
155def make_diaObject_table(objId, plugin, default_value=np.nan, band=None):
156 """Create a minimal diaObject table with columns required for the plugin
158 Parameters
159 ----------
160 objId : `int`
161 The diaObjectId
162 plugin : `lsst.ap.association.DiaCalculationPlugin`
163 The plugin that will be run.
164 default_value : `float` or `int`, optional
165 Value to set new columns to.
166 band : `str`, optional
167 Band designation to append to the plugin columns.
169 Returns
170 -------
171 diaObjects : `pandas.DataFrame`
172 Output catalog with the required columns for the plugin.
173 """
174 # Add an extra empty diaObject here. This ensures that
175 # we properly test the source/object matching implicit
176 # in the plugin calculations.
177 diaObjects = {"diaObjectId": [objId, objId + 1]}
178 for col in plugin.outputCols:
179 if band is not None:
180 diaObjects[f"{band}_{col}"] = default_value
181 else:
182 diaObjects[col] = default_value
183 return pd.DataFrame(diaObjects)
186class TestMeanPosition(unittest.TestCase):
188 def testCalculate(self):
189 """Test mean position calculation.
191 DiaSources here are constructed without per-source coordinate
192 errors, so each ``run_multi_plugin`` call legitimately triggers
193 the no-errors warning from the plugin. Suppress it here so the
194 signal of an unexpected warning elsewhere is not drowned out.
195 """
196 n_sources = 10
197 objId = 0
199 # configure a 2 degree max separation
200 plug = MeanDiaPosition(MeanDiaPositionConfig(MaxAllowedDiaSourceSeparation=7200.0),
201 "ap_meanPosition",
202 None)
204 warnings.filterwarnings("ignore",
205 message="No DiaSources with finite coordinate errors",
206 category=UserWarning)
207 self.addCleanup(warnings.resetwarnings)
209 # Test expected means in RA.
210 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
211 diaSources = pd.DataFrame(data={"ra": np.linspace(-1, 1, n_sources),
212 "dec": np.zeros(n_sources),
213 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
214 "diaObjectId": n_sources * [objId],
215 "band": n_sources * ["g"],
216 "diaSourceId": np.arange(n_sources,
217 dtype=int)})
218 run_multi_plugin(diaObjects, diaSources, "g", plug)
220 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
221 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
223 # Test expected means in DEC.
224 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
225 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),
226 "dec": np.linspace(-1, 1, n_sources),
227 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
228 "diaObjectId": n_sources * [objId],
229 "band": n_sources * ["g"],
230 "diaSourceId": np.arange(n_sources,
231 dtype=int)})
232 run_multi_plugin(diaObjects, diaSources, "g", plug)
234 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
235 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
237 # Test failure mode RA is nan.
238 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
239 diaSources = pd.DataFrame(data={"ra": np.full(n_sources, np.nan),
240 "dec": np.zeros(n_sources),
241 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
242 "diaObjectId": n_sources * [objId],
243 "band": n_sources * ["g"],
244 "diaSourceId": np.arange(n_sources,
245 dtype=int)})
246 run_multi_plugin(diaObjects, diaSources, "g", plug)
248 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))
249 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))
251 # Test failure mode DEC is nan.
252 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
253 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),
254 "dec": np.full(n_sources, np.nan),
255 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
256 "diaObjectId": n_sources * [objId],
257 "band": n_sources * ["g"],
258 "diaSourceId": np.arange(n_sources,
259 dtype=int)})
260 run_multi_plugin(diaObjects, diaSources, "g", plug)
262 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))
263 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))
265 # configure the default 3 arcsecond separation
266 plug = MeanDiaPosition(MeanDiaPositionConfig(MaxAllowedDiaSourceSeparation=3.0),
267 "ap_meanPosition",
268 None)
270 # These 1 degree separations should raise
271 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
272 diaSources = pd.DataFrame(data={"ra": np.linspace(-1, 1, n_sources),
273 "dec": np.zeros(n_sources),
274 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
275 "diaObjectId": n_sources * [objId],
276 "band": n_sources * ["g"],
277 "diaSourceId": np.arange(n_sources,
278 dtype=int)})
279 with self.assertRaises(UnphysicalDiaSourceSeparation):
280 run_multi_plugin(diaObjects, diaSources, "g", plug)
282 # 1 arcsecond separations should not raise
283 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
284 diaSources = pd.DataFrame(data={"ra": np.linspace(-1/3600., 1/3600., n_sources),
285 "dec": np.zeros(n_sources),
286 "midpointMjdTai": np.linspace(0, n_sources, n_sources),
287 "diaObjectId": n_sources * [objId],
288 "band": n_sources * ["g"],
289 "diaSourceId": np.arange(n_sources,
290 dtype=int)})
291 run_multi_plugin(diaObjects, diaSources, "g", plug)
293 def _makeDiaSourcesWithUncertainties(self, raErr, decErr, ra_dec_Cov, ra=None, dec=None, objId=0):
294 """Build a tiny DiaSource DataFrame with per-source uncertainties.
295 """
296 n = len(raErr)
297 if ra is None:
298 ra = np.zeros(n)
299 if dec is None:
300 dec = np.zeros(n)
301 return pd.DataFrame(data={
302 "ra": ra,
303 "dec": dec,
304 "raErr": raErr,
305 "decErr": decErr,
306 "ra_dec_Cov": ra_dec_Cov,
307 "midpointMjdTai": np.arange(n, dtype=float),
308 "diaObjectId": n * [objId],
309 "band": n * ["g"],
310 "diaSourceId": np.arange(n, dtype=int),
311 })
313 def testUncertaintyDiagonalOnly(self):
314 """Two coincident sources, no off-diagonal: chi^2 = 0, scale
315 factor is 1, so the output is the diagonal inverse-variance
316 weighted-mean covariance.
317 """
318 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
319 objId = 0
320 raErr = np.array([1e-6, 2e-6])
321 decErr = np.array([1e-6, 2e-6])
322 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
323 diaSources = self._makeDiaSourcesWithUncertainties(
324 raErr=raErr,
325 decErr=decErr,
326 ra_dec_Cov=np.array([np.nan, np.nan]),
327 )
328 run_multi_plugin(diaObjects, diaSources, "g", plug)
330 # Var(weighted mean) = 1 / sum(1/sigma^2) per axis.
331 expectedRaErr = 1.0/np.sqrt(np.sum(1.0/raErr**2))
332 expectedDecErr = 1.0/np.sqrt(np.sum(1.0/decErr**2))
333 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], expectedRaErr)
334 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], expectedDecErr)
335 # No per-source ra_dec_Cov, so output ra_dec_Cov is NaN.
336 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra_dec_Cov"]))
338 def testUncertaintyFullCovariance(self):
339 """Two coincident sources with off-diagonal covariance: chi^2 = 0,
340 so the output is C_formal = (sum_i inv(C_i))^-1.
341 """
342 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
343 objId = 0
344 raErr = np.array([1e-6, 2e-6])
345 decErr = np.array([1.5e-6, 1.0e-6])
346 rho = np.array([0.3, -0.2]) # correlation coefficient
347 raDecCov = rho * raErr * decErr
349 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
350 diaSources = self._makeDiaSourcesWithUncertainties(
351 raErr=raErr, decErr=decErr, ra_dec_Cov=raDecCov)
352 run_multi_plugin(diaObjects, diaSources, "g", plug)
354 n = len(raErr)
355 cov = np.zeros((n, 2, 2))
356 cov[:, 0, 0] = raErr**2
357 cov[:, 1, 1] = decErr**2
358 cov[:, 0, 1] = raDecCov
359 cov[:, 1, 0] = raDecCov
360 covObj = np.linalg.inv(np.linalg.inv(cov).sum(axis=0))
362 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], np.sqrt(covObj[0, 0]))
363 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], np.sqrt(covObj[1, 1]))
364 self.assertAlmostEqual(diaObjects.loc[objId, "ra_dec_Cov"], covObj[0, 1])
366 def testUncertaintySingleSourceCopiesThrough(self):
367 """A single-source group: outputs equal the source's uncertainties.
368 """
369 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
370 objId = 0
371 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
372 diaSources = self._makeDiaSourcesWithUncertainties(
373 raErr=np.array([1.5e-6]),
374 decErr=np.array([0.8e-6]),
375 ra_dec_Cov=np.array([2e-13]),
376 )
377 run_multi_plugin(diaObjects, diaSources, "g", plug)
379 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], 1.5e-6)
380 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], 0.8e-6)
381 self.assertAlmostEqual(diaObjects.loc[objId, "ra_dec_Cov"], 2e-13)
383 def testUncertaintyMissingColumnsFallsBackToScatter(self):
384 """No raErr/decErr columns + N>=2 spread-out sources -> the
385 weighted-mean path is unusable, so the DiaObject falls back to
386 the unweighted mean position with the scatter-only SEM, and a
387 warning is emitted.
388 """
389 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
390 objId = 0
391 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
392 n = 3
393 ra = np.linspace(-1e-4, 1e-4, n)
394 dec = np.zeros(n)
395 diaSources = pd.DataFrame(data={
396 "ra": ra,
397 "dec": dec,
398 "midpointMjdTai": np.arange(n, dtype=float),
399 "diaObjectId": n * [objId],
400 "band": n * ["g"],
401 "diaSourceId": np.arange(n, dtype=int),
402 })
403 with self.assertWarnsRegex(UserWarning, "No DiaSources with finite coordinate errors"):
404 run_multi_plugin(diaObjects, diaSources, "g", plug)
406 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
407 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
409 # Expected scatter term: at dec=0, the tangent-plane east offsets
410 # equal the RA deltas (in degrees) to high precision, so the
411 # standard error of the mean in RA is sample_std(ra, ddof=1) /
412 # sqrt(N). Dec is identically zero so its scatter is zero.
413 expectedRaErr = np.std(ra, ddof=1)/np.sqrt(n)
414 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], expectedRaErr)
415 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], 0.0)
416 self.assertAlmostEqual(diaObjects.loc[objId, "ra_dec_Cov"], 0.0)
418 def testUncertaintyMissingColumnsSingleSourceEmitsNaN(self):
419 """With a single source and no raErr/decErr the uncertainty is
420 undefined but the mean position is still computed; a warning is
421 emitted that the per-source errors are unusable.
422 """
423 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
424 objId = 0
425 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
426 diaSources = pd.DataFrame(data={
427 "ra": [0.0],
428 "dec": [0.0],
429 "midpointMjdTai": [0.0],
430 "diaObjectId": [objId],
431 "band": ["g"],
432 "diaSourceId": [0],
433 })
434 with self.assertWarnsRegex(UserWarning, "No DiaSources with finite coordinate errors"):
435 run_multi_plugin(diaObjects, diaSources, "g", plug)
437 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)
438 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
439 self.assertTrue(np.isnan(diaObjects.loc[objId, "raErr"]))
440 self.assertTrue(np.isnan(diaObjects.loc[objId, "decErr"]))
441 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra_dec_Cov"]))
443 def testUncertaintyPartialCovarianceFallsBackToDiagonal(self):
444 """Some sources have NaN ra_dec_Cov but valid raErr/decErr.
446 When at least one included source lacks a finite ra_dec_Cov,
447 diagonal weights are used for all included sources, the off-
448 diagonal output is NaN, and the chi-squared falls back to its
449 diagonal form. Sources are coincident here, so chi^2 = 0 and
450 there is no rescaling.
451 """
452 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
453 objId = 0
454 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
455 raErr = np.array([1e-6, 2e-6, 3e-6])
456 decErr = np.array([1e-6, 2e-6, 3e-6])
457 # Only one source has a finite covariance
458 raDecCov = np.array([1e-13, np.nan, np.nan])
459 diaSources = self._makeDiaSourcesWithUncertainties(raErr=raErr, decErr=decErr, ra_dec_Cov=raDecCov)
460 run_multi_plugin(diaObjects, diaSources, "g", plug)
462 # Diagonal inverse-variance weighted mean covariance.
463 expectedRaErr = 1.0/np.sqrt(np.sum(1.0/raErr**2))
464 expectedDecErr = 1.0/np.sqrt(np.sum(1.0/decErr**2))
465 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], expectedRaErr)
466 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], expectedDecErr)
467 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra_dec_Cov"]))
469 def testUncertaintyScaleFactorInflatesWhenScatterExceedsErrors(self):
470 """Two sources whose positional separation is much larger than
471 their per-source errors: chi^2 >> dof, so the chi-squared scale
472 factor inflates the formal covariance by chi^2 / dof.
473 """
474 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
475 objId = 0
476 # Two sources symmetric around RA = 0, Dec = 0. Spread of
477 # 0.36 arcsec on each axis (within MaxAllowedDiaSourceSeparation
478 # of 3 arcsec), still 100x the per-source sigma of 1e-6 deg = 3.6 mas,
479 # so chi^2 / dof ~ 4e4 / 2 ~ 2e4.
480 delta = 1e-4 # degrees
481 ra = np.array([-delta, delta])
482 dec = np.array([-delta, delta])
483 raErr = np.array([1e-6, 1e-6])
484 decErr = np.array([1e-6, 1e-6])
485 raDecCov = np.array([np.nan, np.nan])
487 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
488 diaSources = self._makeDiaSourcesWithUncertainties(
489 raErr=raErr, decErr=decErr, ra_dec_Cov=raDecCov, ra=ra, dec=dec)
490 run_multi_plugin(diaObjects, diaSources, "g", plug)
492 # With equal per-source errors the weighted mean equals the
493 # unweighted mean = (0, 0), so the tangent-plane residuals are
494 # essentially (ra, dec).
495 # chi^2 = sum_i (r_east^2/sigma_ra^2 + r_north^2/sigma_dec^2)
496 # dof = 2 * (N - 1) = 2.
497 chi2 = np.sum(ra**2/raErr**2 + dec**2/decErr**2)
498 dof = 2*(len(raErr) - 1)
499 scale = max(1.0, chi2/dof)
500 # Diagonal weighted-mean variance: 1/sum(1/sigma^2) per axis.
501 varRaFormal = 1.0/np.sum(1.0/raErr**2)
502 varDecFormal = 1.0/np.sum(1.0/decErr**2)
503 expectedRaErr = np.sqrt(scale*varRaFormal)
504 expectedDecErr = np.sqrt(scale*varDecFormal)
506 # The inflation should be large (scale ~ 1e6):
507 self.assertGreater(scale, 1e3)
509 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"] / expectedRaErr, 1.0, places=5)
510 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"] / expectedDecErr, 1.0, places=5)
511 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra_dec_Cov"]))
513 def testUncertaintyScaleFactorNoInflationWhenConsistent(self):
514 """When per-source positions are coincident, chi^2 = 0 and the
515 chi-squared scale factor is exactly 1: the output equals
516 C_formal = (sum_i inv(C_i))^-1.
517 """
518 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
519 objId = 0
520 raErr = np.array([1e-6, 2e-6, 1.5e-6])
521 decErr = np.array([1e-6, 1.5e-6, 2e-6])
522 raDecCov = np.array([2e-13, -1e-13, 3e-14])
524 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
525 diaSources = self._makeDiaSourcesWithUncertainties(raErr=raErr, decErr=decErr, ra_dec_Cov=raDecCov)
526 run_multi_plugin(diaObjects, diaSources, "g", plug)
528 n = len(raErr)
529 cov = np.zeros((n, 2, 2))
530 cov[:, 0, 0] = raErr**2
531 cov[:, 1, 1] = decErr**2
532 cov[:, 0, 1] = raDecCov
533 cov[:, 1, 0] = raDecCov
534 covObj = np.linalg.inv(np.linalg.inv(cov).sum(axis=0))
536 self.assertAlmostEqual(diaObjects.loc[objId, "raErr"], np.sqrt(covObj[0, 0]))
537 self.assertAlmostEqual(diaObjects.loc[objId, "decErr"], np.sqrt(covObj[1, 1]))
538 self.assertAlmostEqual(diaObjects.loc[objId, "ra_dec_Cov"], covObj[0, 1])
540 def testWeightedMeanPositionDiffersFromUnweighted(self):
541 """Two sources at asymmetric (ra) positions with very unequal
542 per-source errors: the reported position is the inverse-variance
543 weighted mean, which is much closer to the more-precise source
544 than to the unweighted midpoint.
545 """
546 plug = MeanDiaPosition(MeanDiaPositionConfig(), "ap_meanPosition", None)
547 objId = 0
548 d = 1e-5 # degrees; well inside MaxAllowedDiaSourceSeparation.
549 ra = np.array([-d, d])
550 dec = np.array([0.0, 0.0])
551 # Source 0 is 100x more precise than source 1.
552 raErr = np.array([1e-6, 1e-4])
553 decErr = np.array([1e-6, 1e-4])
554 raDecCov = np.array([np.nan, np.nan])
556 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
557 diaSources = self._makeDiaSourcesWithUncertainties(
558 raErr=raErr, decErr=decErr, ra_dec_Cov=raDecCov, ra=ra, dec=dec)
559 run_multi_plugin(diaObjects, diaSources, "g", plug)
561 # Diagonal weighted mean of (ra_0, ra_1) with weights w_i = 1/raErr_i^2.
562 w = 1.0/raErr**2
563 expectedRa = np.sum(w*ra)/np.sum(w)
564 # Output RA is wrapped to [0, 360); normalize back to a signed
565 # offset near zero before comparing.
566 outRa = ((diaObjects.loc[objId, "ra"] + 180.0) % 360.0) - 180.0
567 # The unweighted mean would be 0; the weighted mean should be
568 # very close to source 0 at ra = -d.
569 self.assertAlmostEqual(outRa, expectedRa)
570 self.assertLess(abs(outRa - (-d)), 0.01*d)
571 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)
574class TestFirstLastDiaSourceTime(unittest.TestCase):
576 def _makeSources(self, objId, times):
577 n = len(times)
578 return pd.DataFrame(
579 data={"diaObjectId": n * [objId],
580 "band": n * ["g"],
581 "diaSourceId": np.arange(n, dtype=int),
582 "midpointMjdTai": np.asarray(times, dtype=float)})
584 def testNewObject(self):
585 """First/last times are the min/max source times for a new object.
586 """
587 objId = 0
588 times = [60010.0, 60005.5, 60020.25]
589 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
590 "ap_firstLastDiaSourceTime",
591 None)
592 diaObjects = make_diaObject_table(objId, plug)
593 # Include an extra object with no DiaSources
594 diaObjects.loc[1, "firstDiaSourceMjdTai"] = 59000.0
595 diaObjects.loc[1, "lastDiaSourceMjdTai"] = 59100.0
596 diaSources = self._makeSources(objId, times)
597 run_multiband_plugin(diaObjects, diaSources, plug)
599 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 60005.5)
600 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 60020.25)
601 # The object with no sources must be left untouched.
602 self.assertEqual(diaObjects.loc[objId + 1, "firstDiaSourceMjdTai"], 59000.0)
603 self.assertEqual(diaObjects.loc[objId + 1, "lastDiaSourceMjdTai"], 59100.0)
605 def testStoredValuesOutsideRangeAreKept(self):
606 """Stored first/last times outside the loaded range are retained.
607 """
608 objId = 0
609 times = [60010.0, 60005.5, 60020.25]
610 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
611 "ap_firstLastDiaSourceTime",
612 None)
613 diaObjects = make_diaObject_table(objId, plug)
614 diaObjects.loc[0, "firstDiaSourceMjdTai"] = 59000.0
615 diaObjects.loc[0, "lastDiaSourceMjdTai"] = 61000.0
616 diaSources = self._makeSources(objId, times)
617 run_multiband_plugin(diaObjects, diaSources, plug)
619 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 59000.0)
620 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 61000.0)
622 def testStoredValuesInsideRangeAreExtended(self):
623 """Stored first/last times inside the loaded range are extended.
624 """
625 objId = 0
626 times = [60010.0, 60005.5, 60020.25]
627 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
628 "ap_firstLastDiaSourceTime",
629 None)
630 diaObjects = make_diaObject_table(objId, plug)
631 diaObjects.loc[0, "firstDiaSourceMjdTai"] = 60008.0
632 diaObjects.loc[0, "lastDiaSourceMjdTai"] = 60012.0
633 diaSources = self._makeSources(objId, times)
634 run_multiband_plugin(diaObjects, diaSources, plug)
636 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 60005.5)
637 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 60020.25)
639 def testOnlyOneStoredValue(self):
640 """A stored first time and an unset last time are handled
641 independently.
642 """
643 objId = 0
644 times = [60010.0, 60005.5, 60020.25]
645 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
646 "ap_firstLastDiaSourceTime",
647 None)
648 diaObjects = make_diaObject_table(objId, plug)
649 diaObjects.loc[0, "firstDiaSourceMjdTai"] = 59000.0
650 diaSources = self._makeSources(objId, times)
651 run_multiband_plugin(diaObjects, diaSources, plug)
653 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 59000.0)
654 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 60020.25)
656 def testNanSourceTimesIgnored(self):
657 """DiaSources with NaN midpointMjdTai do not affect the result.
658 """
659 objId = 0
660 times = [60010.0, np.nan, 60020.25]
661 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
662 "ap_firstLastDiaSourceTime",
663 None)
664 diaObjects = make_diaObject_table(objId, plug)
665 diaSources = self._makeSources(objId, times)
666 run_multiband_plugin(diaObjects, diaSources, plug)
668 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 60010.0)
669 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 60020.25)
671 def testMissingOutputColumns(self):
672 """Output columns are created when absent from the DiaObject table.
673 """
674 objId = 0
675 times = [60010.0, 60005.5]
676 plug = FirstLastDiaSourceTime(FirstLastDiaSourceTimeConfig(),
677 "ap_firstLastDiaSourceTime",
678 None)
679 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
680 diaSources = self._makeSources(objId, times)
681 run_multiband_plugin(diaObjects, diaSources, plug)
683 self.assertEqual(diaObjects.loc[objId, "firstDiaSourceMjdTai"], 60005.5)
684 self.assertEqual(diaObjects.loc[objId, "lastDiaSourceMjdTai"], 60010.0)
687class TestHTMIndexPosition(unittest.TestCase):
689 def testCalculate(self):
690 """Test HTMPixel assignment calculation.
691 """
692 # Test expected pixelId at RA, DEC = 0
693 objId = 0
694 n_sources = 10
695 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
696 diaObjects.loc[objId, "ra"] = 0.
697 diaObjects.loc[objId, "dec"] = 0.
698 diaSources = pd.DataFrame(
699 data={"diaObjectId": n_sources * [objId],
700 "band": n_sources * ["g"],
701 "diaSourceId": np.arange(n_sources, dtype=int)})
702 plug = HTMIndexDiaPosition(HTMIndexDiaPositionConfig(),
703 "ap_HTMIndex",
704 None)
706 run_single_plugin(diaObjectCat=diaObjects,
707 diaObjectId=objId,
708 diaSourceCat=diaSources,
709 band="g",
710 plugin=plug)
711 self.assertEqual(diaObjects.at[objId, "pixelId"],
712 17042430230528)
714 # Test expected pixelId at some value of RA and DEC.
715 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
716 diaObjects.loc[objId, "ra"] = 45.37
717 diaObjects.loc[objId, "dec"] = 13.67
718 diaSources = pd.DataFrame(
719 data={"diaObjectId": n_sources * [objId],
720 "band": n_sources * ["g"],
721 "diaSourceId": np.arange(n_sources, dtype=int)})
722 run_single_plugin(diaObjectCat=diaObjects,
723 diaObjectId=objId,
724 diaSourceCat=diaSources,
725 band="g",
726 plugin=plug)
727 self.assertEqual(diaObjects.at[objId, "pixelId"],
728 17450571968473)
731class TestNDiaSourcesDiaPlugin(unittest.TestCase):
733 def testCalculate(self):
734 """Test that the number of DiaSources is correct.
735 """
737 for n_sources in [1, 8, 10]:
738 # Test expected number of sources per object.
739 objId = 0
740 diaSources = pd.DataFrame(
741 data={"diaObjectId": n_sources * [objId],
742 "band": n_sources * ["g"],
743 "diaSourceId": np.arange(n_sources, dtype=int)})
744 plug = NumDiaSourcesDiaPlugin(NumDiaSourcesDiaPluginConfig(),
745 "ap_nDiaSources",
746 None)
747 diaObjects = make_diaObject_table(objId, plug, default_value=int(-1))
748 run_multi_plugin(diaObjects, diaSources, "g", plug)
750 self.assertEqual(n_sources, diaObjects.at[objId, "nDiaSources"])
751 self.assertEqual(diaObjects["nDiaSources"].dtype, np.int64)
754class TestSimpleSourceFlagDiaPlugin(unittest.TestCase):
756 def testCalculate(self):
757 """Test that DiaObject flags are set.
758 """
759 objId = 0
760 n_sources = 10
762 # Test expected flags, no flags set.
763 diaSources = pd.DataFrame(
764 data={"diaObjectId": n_sources * [objId],
765 "band": n_sources * ["g"],
766 "diaSourceId": np.arange(n_sources, dtype=int),
767 "flags": np.zeros(n_sources, dtype=np.uint64)})
768 plug = SimpleSourceFlagDiaPlugin(SimpleSourceFlagDiaPluginConfig(),
769 "ap_diaObjectFlag",
770 None)
772 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
773 run_multi_plugin(diaObjects, diaSources, "g", plug)
774 self.assertEqual(diaObjects.at[objId, "flags"], 0)
775 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
777 # Test expected flags, all flags set.
778 diaSources = pd.DataFrame(
779 data={"diaObjectId": n_sources * [objId],
780 "band": n_sources * ["g"],
781 "diaSourceId": np.arange(n_sources, dtype=int),
782 "flags": np.ones(n_sources, dtype=np.uint64)})
783 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
784 run_multi_plugin(diaObjects, diaSources, "g", plug)
785 self.assertEqual(diaObjects.at[objId, "flags"], 1)
786 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
788 # Test expected flags, random flags.
789 diaSources = pd.DataFrame(
790 data={"diaObjectId": n_sources * [objId],
791 "band": n_sources * ["g"],
792 "diaSourceId": np.arange(n_sources, dtype=int),
793 "flags": np.random.randint(0, 2 ** 16, size=n_sources)})
795 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
796 run_multi_plugin(diaObjects, diaSources, "g", plug)
797 self.assertEqual(diaObjects.at[objId, "flags"], 1)
798 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
800 # Test expected flags, one flag set.
801 flag_array = np.zeros(n_sources, dtype=np.uint64)
802 flag_array[4] = 256
803 diaSources = pd.DataFrame(
804 data={"diaObjectId": n_sources * [objId],
805 "band": n_sources * ["g"],
806 "diaSourceId": np.arange(n_sources, dtype=int),
807 "flags": flag_array})
808 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))
809 run_multi_plugin(diaObjects, diaSources, "g", plug)
810 self.assertEqual(diaObjects.at[objId, "flags"], 1)
811 self.assertEqual(diaObjects["flags"].dtype, np.uint64)
814class TestWeightedMeanDiaPsfFlux(unittest.TestCase):
816 def testCalculate(self):
817 """Test mean value calculation.
818 """
819 n_sources = 10
820 objId = 0
822 # Test expected mean.
823 # In the first test, we have only one object.
824 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
825 diaSources = pd.DataFrame(
826 data={"diaObjectId": n_sources * [objId],
827 "band": n_sources * ["u"],
828 "diaSourceId": np.arange(n_sources, dtype=int),
829 "psfFlux": np.linspace(-1, 1, n_sources),
830 "psfFluxErr": np.ones(n_sources)})
832 plug = WeightedMeanDiaPsfFlux(WeightedMeanDiaPsfFluxConfig(),
833 "ap_meanFlux",
834 None)
835 run_multi_plugin(diaObjects, diaSources, "u", plug)
837 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMean"], 0.0)
838 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMeanErr"],
839 np.sqrt(1 / n_sources))
840 self.assertEqual(diaObjects.loc[objId, "u_psfFluxNdata"], n_sources)
841 # We expect this to be converted to float.
842 # TODO DM-53254: This should be an integer (and should be checked
843 # to be an integer).
844 self.assertEqual(diaObjects["u_psfFluxNdata"].dtype, np.float64)
846 # Test expected mean with a nan value.
847 # In the second test, we have two objects (one empty).
848 diaObjects = pd.DataFrame({"diaObjectId": [objId, objId + 1]})
849 fluxes = np.linspace(-1, 1, n_sources)
850 fluxes[4] = np.nan
851 diaSources = pd.DataFrame(
852 data={"diaObjectId": n_sources * [objId],
853 "band": n_sources * ["r"],
854 "diaSourceId": np.arange(n_sources, dtype=int),
855 "psfFlux": fluxes,
856 "psfFluxErr": np.ones(n_sources)})
857 run_multi_plugin(diaObjects, diaSources, "r", plug)
859 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMean"],
860 np.nanmean(fluxes))
861 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMeanErr"],
862 np.sqrt(1 / (n_sources - 1)))
863 self.assertEqual(diaObjects.loc[objId, "r_psfFluxNdata"], n_sources - 1)
864 # We expect this to be converted to float.
865 # TODO DM-53254: This should be an integer (and should be checked
866 # to be an integer).
867 self.assertEqual(diaObjects["r_psfFluxNdata"].dtype, np.float64)
870class TestPercentileDiaPsfFlux(unittest.TestCase):
872 def testCalculate(self):
873 """Test flux percentile calculation.
874 """
875 n_sources = 10
876 objId = 0
878 # Test expected percentile values.
879 fluxes = np.linspace(-1, 1, n_sources)
880 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
881 diaSources = pd.DataFrame(
882 data={"diaObjectId": n_sources * [objId],
883 "band": n_sources * ["u"],
884 "diaSourceId": np.arange(n_sources, dtype=int),
885 "psfFlux": fluxes,
886 "psfFluxErr": np.ones(n_sources)})
888 plug = PercentileDiaPsfFlux(PercentileDiaPsfFluxConfig(),
889 "ap_percentileFlux",
890 None)
891 run_multi_plugin(diaObjects, diaSources, "u", plug)
892 for pTile, testVal in zip(plug.config.percentiles,
893 np.nanpercentile(
894 fluxes,
895 plug.config.percentiles)):
896 self.assertAlmostEqual(
897 diaObjects.at[objId, "u_psfFluxPercentile{:02d}".format(pTile)],
898 testVal)
900 # Test expected percentile values with a nan value.
901 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
902 fluxes[4] = np.nan
903 diaSources = pd.DataFrame(
904 data={"diaObjectId": n_sources * [objId],
905 "band": n_sources * ["r"],
906 "diaSourceId": np.arange(n_sources, dtype=int),
907 "psfFlux": fluxes,
908 "psfFluxErr": np.ones(n_sources)})
909 run_multi_plugin(diaObjects, diaSources, "r", plug)
910 for pTile, testVal in zip(plug.config.percentiles,
911 np.nanpercentile(
912 fluxes,
913 plug.config.percentiles)):
914 self.assertAlmostEqual(
915 diaObjects.at[objId, "r_psfFluxPercentile{:02d}".format(pTile)],
916 testVal)
919class TestSigmaDiaPsfFlux(unittest.TestCase):
921 def testCalculate(self):
922 """Test flux scatter calculation.
923 """
924 n_sources = 10
925 objId = 0
927 # Test expected sigma scatter of fluxes.
928 fluxes = np.linspace(-1, 1, n_sources)
929 diaSources = pd.DataFrame(
930 data={"diaObjectId": n_sources * [objId],
931 "band": n_sources * ["u"],
932 "diaSourceId": np.arange(n_sources, dtype=int),
933 "psfFlux": fluxes,
934 "psfFluxErr": np.ones(n_sources)})
936 plug = SigmaDiaPsfFlux(SigmaDiaPsfFluxConfig(),
937 "ap_sigmaFlux",
938 None)
939 diaObjects = make_diaObject_table(objId, plug, band='u')
940 run_multi_plugin(diaObjects, diaSources, "u", plug)
941 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxSigma"],
942 np.nanstd(fluxes, ddof=1))
944 # test one input, returns nan.
945 diaSources = pd.DataFrame(
946 data={"diaObjectId": 1 * [objId],
947 "band": 1 * ["g"],
948 "diaSourceId": [0],
949 "psfFlux": [fluxes[0]],
950 "psfFluxErr": [1.]})
952 diaObjects = make_diaObject_table(objId, plug, band='g')
953 run_multi_plugin(diaObjects, diaSources, "g", plug)
954 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxSigma"]))
956 # Test expected sigma scatter of fluxes with a nan value.
957 fluxes[4] = np.nan
958 diaSources = pd.DataFrame(
959 data={"diaObjectId": n_sources * [objId],
960 "band": n_sources * ["r"],
961 "diaSourceId": np.arange(n_sources, dtype=int),
962 "psfFlux": fluxes,
963 "psfFluxErr": np.ones(n_sources)})
965 diaObjects = make_diaObject_table(objId, plug, band='r')
966 run_multi_plugin(diaObjects, diaSources, "r", plug)
967 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxSigma"],
968 np.nanstd(fluxes, ddof=1))
971class TestChi2DiaPsfFlux(unittest.TestCase):
973 def testCalculate(self):
974 """Test flux chi2 calculation.
975 """
976 n_sources = 10
977 objId = 0
979 # Test expected chi^2 value.
980 fluxes = np.linspace(-1, 1, n_sources)
981 diaObjects = pd.DataFrame({"diaObjectId": [objId],
982 "u_psfFluxMean": [0.0]})
983 diaSources = pd.DataFrame(
984 data={"diaObjectId": n_sources * [objId],
985 "band": n_sources * ["u"],
986 "diaSourceId": np.arange(n_sources, dtype=int),
987 "psfFlux": fluxes,
988 "psfFluxErr": np.ones(n_sources)})
990 plug = Chi2DiaPsfFlux(Chi2DiaPsfFluxConfig(),
991 "ap_chi2Flux",
992 None)
993 run_multi_plugin(diaObjects, diaSources, "u", plug)
994 self.assertAlmostEqual(
995 diaObjects.loc[objId, "u_psfFluxChi2"],
996 np.nansum(((diaSources["psfFlux"]
997 - np.nanmean(diaSources["psfFlux"]))
998 / diaSources["psfFluxErr"]) ** 2))
1000 # Test expected chi^2 value with a nan value set.
1001 fluxes[4] = np.nan
1002 diaObjects = pd.DataFrame({"diaObjectId": [objId],
1003 "r_psfFluxMean": [np.nanmean(fluxes)]})
1004 diaSources = pd.DataFrame(
1005 data={"diaObjectId": n_sources * [objId],
1006 "band": n_sources * ["r"],
1007 "diaSourceId": np.arange(n_sources, dtype=int),
1008 "psfFlux": fluxes,
1009 "psfFluxErr": np.ones(n_sources)})
1010 run_multi_plugin(diaObjects, diaSources, "r", plug)
1011 self.assertAlmostEqual(
1012 diaObjects.loc[objId, "r_psfFluxChi2"],
1013 np.nansum(((diaSources["psfFlux"]
1014 - np.nanmean(diaSources["psfFlux"]))
1015 / diaSources["psfFluxErr"]) ** 2))
1018class TestMadDiaPsfFlux(unittest.TestCase):
1020 def testCalculate(self):
1021 """Test flux median absolute deviation calculation.
1022 """
1023 n_sources = 10
1024 objId = 0
1026 # Test expected MAD value.
1027 fluxes = np.linspace(-1, 1, n_sources)
1028 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1029 diaSources = pd.DataFrame(
1030 data={"diaObjectId": n_sources * [objId],
1031 "band": n_sources * ["u"],
1032 "diaSourceId": np.arange(n_sources, dtype=int),
1033 "psfFlux": fluxes,
1034 "psfFluxErr": np.ones(n_sources)})
1036 plug = MadDiaPsfFlux(MadDiaPsfFluxConfig(),
1037 "ap_madFlux",
1038 None)
1039 run_multi_plugin(diaObjects, diaSources, "u", plug)
1040 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMAD"],
1041 median_absolute_deviation(fluxes,
1042 ignore_nan=True))
1044 # Test expected MAD value with a nan set.
1045 fluxes[4] = np.nan
1046 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1047 diaSources = pd.DataFrame(
1048 data={"diaObjectId": n_sources * [objId],
1049 "band": n_sources * ["r"],
1050 "diaSourceId": np.arange(n_sources, dtype=int),
1051 "psfFlux": fluxes,
1052 "psfFluxErr": np.ones(n_sources)})
1053 run_multi_plugin(diaObjects, diaSources, "r", plug)
1054 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMAD"],
1055 median_absolute_deviation(fluxes,
1056 ignore_nan=True))
1059class TestSkewDiaPsfFlux(unittest.TestCase):
1061 def testCalculate(self):
1062 """Test flux skew calculation.
1063 """
1064 n_sources = 10
1065 objId = 0
1067 # Test expected skew value.
1068 fluxes = np.linspace(-1, 1, n_sources)
1069 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1070 diaSources = pd.DataFrame(
1071 data={"diaObjectId": n_sources * [objId],
1072 "band": n_sources * ["u"],
1073 "diaSourceId": np.arange(n_sources, dtype=int),
1074 "psfFlux": fluxes,
1075 "psfFluxErr": np.ones(n_sources)})
1077 plug = SkewDiaPsfFlux(SkewDiaPsfFluxConfig(),
1078 "ap_skewFlux",
1079 None)
1080 run_multi_plugin(diaObjects, diaSources, "u", plug)
1081 self.assertAlmostEqual(
1082 diaObjects.loc[objId, "u_psfFluxSkew"],
1083 skew_wrapper(fluxes))
1085 # Test expected skew value with a nan set.
1086 fluxes[4] = np.nan
1087 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1088 diaSources = pd.DataFrame(
1089 data={"diaObjectId": n_sources * [objId],
1090 "band": n_sources * ["r"],
1091 "diaSourceId": np.arange(n_sources, dtype=int),
1092 "psfFlux": fluxes,
1093 "psfFluxErr": np.ones(n_sources)})
1094 run_multi_plugin(diaObjects, diaSources, "r", plug)
1096 self.assertAlmostEqual(
1097 diaObjects.at[objId, "r_psfFluxSkew"],
1098 skew_wrapper(fluxes))
1101class TestMinMaxDiaPsfFlux(unittest.TestCase):
1103 def testCalculate(self):
1104 """Test flux min/max calculation.
1105 """
1106 n_sources = 10
1107 objId = 0
1109 # Test expected MinMax fluxes.
1110 fluxes = np.linspace(-1, 1, n_sources)
1111 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1112 diaSources = pd.DataFrame(
1113 data={"diaObjectId": n_sources * [objId],
1114 "band": n_sources * ["u"],
1115 "diaSourceId": np.arange(n_sources, dtype=int),
1116 "psfFlux": fluxes,
1117 "psfFluxErr": np.ones(n_sources)})
1119 plug = MinMaxDiaPsfFlux(MinMaxDiaPsfFluxConfig(),
1120 "ap_minMaxFlux",
1121 None)
1122 run_multi_plugin(diaObjects, diaSources, "u", plug)
1123 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMin"], -1)
1124 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMax"], 1)
1126 # Test expected MinMax fluxes with a nan set.
1127 fluxes[4] = np.nan
1128 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1129 diaSources = pd.DataFrame(
1130 data={"diaObjectId": n_sources * [objId],
1131 "band": n_sources * ["r"],
1132 "diaSourceId": np.arange(n_sources, dtype=int),
1133 "psfFlux": fluxes,
1134 "psfFluxErr": np.ones(n_sources)})
1135 run_multi_plugin(diaObjects, diaSources, "r", plug)
1136 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMin"], -1)
1137 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMax"], 1)
1140class TestMaxSlopeDiaPsfFlux(unittest.TestCase):
1142 def testCalculate(self):
1143 """Test flux maximum slope.
1144 """
1145 n_sources = 10
1146 objId = 0
1148 # Test max slope value.
1149 fluxes = np.linspace(-1, 1, n_sources)
1150 times = np.concatenate([np.linspace(0, 1, n_sources)[:-1], [1 - 1/90]])
1151 diaSources = pd.DataFrame(
1152 data={"diaObjectId": n_sources * [objId],
1153 "band": n_sources * ["u"],
1154 "diaSourceId": np.arange(n_sources, dtype=int),
1155 "psfFlux": fluxes,
1156 "psfFluxErr": np.ones(n_sources),
1157 "midpointMjdTai": times})
1159 plug = MaxSlopeDiaPsfFlux(MaxSlopeDiaPsfFluxConfig(),
1160 "ap_maxSlopeFlux",
1161 None)
1162 diaObjects = make_diaObject_table(objId, plug, band='u')
1163 run_multi_plugin(diaObjects, diaSources, "u", plug)
1164 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMaxSlope"], 2 + 2/9)
1166 # Test max slope value returns nan on 1 input.
1167 diaSources = pd.DataFrame(
1168 data={"diaObjectId": 1 * [objId],
1169 "band": 1 * ["g"],
1170 "diaSourceId": np.arange(1, dtype=int),
1171 "psfFlux": fluxes[0],
1172 "psfFluxErr": np.ones(1),
1173 "midpointMjdTai": times[0]})
1174 diaObjects = make_diaObject_table(objId, plug, band='g')
1175 run_multi_plugin(diaObjects, diaSources, "g", plug)
1176 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxMaxSlope"]))
1178 # Test max slope value inputing nan values.
1179 fluxes[4] = np.nan
1180 times[7] = np.nan
1181 diaSources = pd.DataFrame(
1182 data={"diaObjectId": n_sources * [objId],
1183 "band": n_sources * ["r"],
1184 "diaSourceId": np.arange(n_sources, dtype=int),
1185 "psfFlux": fluxes,
1186 "psfFluxErr": np.ones(n_sources),
1187 "midpointMjdTai": times})
1188 diaObjects = make_diaObject_table(objId, plug, band='r')
1189 run_multi_plugin(diaObjects, diaSources, "r", plug)
1190 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMaxSlope"], 2 + 2 / 9)
1193class TestErrMeanDiaPsfFlux(unittest.TestCase):
1195 def testCalculate(self):
1196 """Test error mean calculation.
1197 """
1198 n_sources = 10
1199 objId = 0
1201 # Test mean of the errors.
1202 fluxes = np.linspace(-1, 1, n_sources)
1203 errors = np.linspace(1, 2, n_sources)
1204 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1205 diaSources = pd.DataFrame(
1206 data={"diaObjectId": n_sources * [objId],
1207 "band": n_sources * ["u"],
1208 "diaSourceId": np.arange(n_sources, dtype=int),
1209 "psfFlux": fluxes,
1210 "psfFluxErr": errors})
1212 plug = ErrMeanDiaPsfFlux(ErrMeanDiaPsfFluxConfig(),
1213 "ap_errMeanFlux",
1214 None)
1215 run_multi_plugin(diaObjects, diaSources, "u", plug)
1216 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxErrMean"],
1217 np.nanmean(errors).astype(np.float32))
1219 # Test mean of the errors with input nan value.
1220 errors[4] = np.nan
1221 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1222 diaSources = pd.DataFrame(
1223 data={"diaObjectId": n_sources * [objId],
1224 "band": n_sources * ["r"],
1225 "diaSourceId": np.arange(n_sources, dtype=int),
1226 "psfFlux": fluxes,
1227 "psfFluxErr": errors})
1228 run_multi_plugin(diaObjects, diaSources, "r", plug)
1229 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxErrMean"],
1230 np.nanmean(errors).astype(np.float32))
1233class TestLinearFitDiaPsfFlux(unittest.TestCase):
1235 def testCalculate(self):
1236 """Test a linear fit to flux vs time.
1237 """
1238 n_sources = 10
1239 objId = 0
1241 # Test best fit linear model.
1242 fluxes = np.linspace(-1, 1, n_sources)
1243 errors = np.linspace(1, 2, n_sources)
1244 times = np.linspace(0, 1, n_sources)
1245 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1246 diaSources = pd.DataFrame(
1247 data={"diaObjectId": n_sources * [objId],
1248 "band": n_sources * ["u"],
1249 "diaSourceId": np.arange(n_sources, dtype=int),
1250 "psfFlux": fluxes,
1251 "psfFluxErr": errors,
1252 "midpointMjdTai": times})
1254 plug = LinearFitDiaPsfFlux(LinearFitDiaPsfFluxConfig(),
1255 "ap_LinearFit",
1256 None)
1257 run_multi_plugin(diaObjects, diaSources, "u", plug)
1258 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearSlope"],
1259 2.)
1260 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearIntercept"],
1261 -1.)
1263 # Test best fit linear model with input nans.
1264 fluxes[7] = np.nan
1265 errors[4] = np.nan
1266 times[2] = np.nan
1267 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1268 diaSources = pd.DataFrame(
1269 data={"diaObjectId": n_sources * [objId],
1270 "band": n_sources * ["r"],
1271 "diaSourceId": np.arange(n_sources, dtype=int),
1272 "psfFlux": fluxes,
1273 "psfFluxErr": errors,
1274 "midpointMjdTai": times})
1275 run_multi_plugin(diaObjects, diaSources, "r", plug)
1276 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearSlope"], 2.)
1277 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearIntercept"],
1278 -1.)
1281class TestStetsonJDiaPsfFlux(unittest.TestCase):
1283 def testCalculate(self):
1284 """Test the stetsonJ statistic.
1285 """
1286 n_sources = 10
1287 objId = 0
1289 # Test stetsonJ calculation.
1290 fluxes = np.linspace(-1, 1, n_sources)
1291 errors = np.ones(n_sources)
1292 diaObjects = pd.DataFrame({"diaObjectId": [objId],
1293 "u_psfFluxMean": [np.nanmean(fluxes)]})
1294 diaSources = pd.DataFrame(
1295 data={"diaObjectId": n_sources * [objId],
1296 "band": n_sources * ["u"],
1297 "diaSourceId": np.arange(n_sources, dtype=int),
1298 "psfFlux": fluxes,
1299 "psfFluxErr": errors})
1301 plug = StetsonJDiaPsfFlux(StetsonJDiaPsfFluxConfig(),
1302 "ap_StetsonJ",
1303 None)
1304 run_multi_plugin(diaObjects, diaSources, "u", plug)
1305 # Expected StetsonJ for the values created. Confirmed using Cesimum's
1306 # implementation. http://github.com/cesium-ml/cesium
1307 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxStetsonJ"],
1308 -0.5958393936080928)
1310 # Test stetsonJ calculation returns nan on single input.
1311 diaObjects = pd.DataFrame({"diaObjectId": [objId],
1312 "g_psfFluxMean": [np.nanmean(fluxes)]})
1313 diaSources = pd.DataFrame(
1314 data={"diaObjectId": 1 * [objId],
1315 "band": 1 * ["g"],
1316 "diaSourceId": np.arange(1, dtype=int),
1317 "psfFlux": fluxes[0],
1318 "psfFluxErr": errors[0]})
1319 run_multi_plugin(diaObjects, diaSources, "g", plug)
1320 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxStetsonJ"]))
1322 # Test stetsonJ calculation returns when nans are input.
1323 fluxes[7] = np.nan
1324 errors[4] = np.nan
1325 nonNanMask = np.logical_and(~np.isnan(fluxes),
1326 ~np.isnan(errors))
1327 diaObjects = pd.DataFrame(
1328 {"diaObjectId": [objId],
1329 "r_psfFluxMean": [np.average(fluxes[nonNanMask],
1330 weights=errors[nonNanMask])]})
1331 diaSources = pd.DataFrame(
1332 data={"diaObjectId": n_sources * [objId],
1333 "band": n_sources * ["r"],
1334 "diaSourceId": np.arange(n_sources, dtype=int),
1335 "psfFlux": fluxes,
1336 "psfFluxErr": errors})
1337 run_multi_plugin(diaObjects, diaSources, "r", plug)
1338 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxStetsonJ"],
1339 -0.5412797916187173)
1342class TestWeightedMeanDiaTotFlux(unittest.TestCase):
1344 def testCalculate(self):
1345 """Test mean value calculation.
1346 """
1347 n_sources = 10
1348 objId = 0
1350 # Test test mean on scienceFlux.
1351 fluxes = np.linspace(-1, 1, n_sources)
1352 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1353 diaSources = pd.DataFrame(
1354 data={"diaObjectId": n_sources * [objId],
1355 "band": n_sources * ["u"],
1356 "diaSourceId": np.arange(n_sources, dtype=int),
1357 "scienceFlux": fluxes,
1358 "scienceFluxErr": np.ones(n_sources)})
1360 plug = WeightedMeanDiaTotFlux(WeightedMeanDiaTotFluxConfig(),
1361 "ap_meanTotFlux",
1362 None)
1363 run_multi_plugin(diaObjects, diaSources, "u", plug)
1365 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMean"], 0.0)
1366 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMeanErr"],
1367 np.sqrt(1 / n_sources))
1369 # Test test mean on scienceFlux with input nans
1370 fluxes[4] = np.nan
1371 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1372 diaSources = pd.DataFrame(
1373 data={"diaObjectId": n_sources * [objId],
1374 "band": n_sources * ["r"],
1375 "diaSourceId": np.arange(n_sources, dtype=int),
1376 "scienceFlux": fluxes,
1377 "scienceFluxErr": np.ones(n_sources)})
1378 run_multi_plugin(diaObjects, diaSources, "r", plug)
1380 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMean"],
1381 np.nanmean(fluxes))
1382 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMeanErr"],
1383 np.sqrt(1 / (n_sources - 1)))
1386def generatePeriodicData(n=10, period=10):
1387 """Generate noisy, sinusoidally-varying periodic data for testing Lomb-
1388 Scargle Periodogram.
1390 The returned fluxes will have, within the errors, the passed-in period and
1391 a power close to 1, because the fluxes are purely sinusoidal.
1393 Parameters
1394 ----------
1395 n : int
1396 Number of data points to generate.
1397 period : float
1398 Period of the periodic signal.
1400 Returns
1401 -------
1402 t : np.ndarray
1403 Time values.
1404 y_obs : np.ndarray
1405 Observed flux values.
1406 """
1407 np.random.seed(42)
1409 t = np.linspace(-2*np.pi, 2*np.pi, n) + 100*np.random.random(n)
1410 y = 10 + np.sin(2 * np.pi * t / period)
1411 y_obs = np.random.normal(y, 0.001)
1413 return t, y_obs
1416class TestMultiLombScarglePeriodogram(lsst.utils.tests.TestCase):
1418 def testCalculate(self):
1419 """Test Mulitband Lomb Scargle Periodogram."""
1420 n_sources = 10
1421 objId = 0
1423 # Create synthetic multi-band data
1424 times, fluxes = generatePeriodicData(n_sources, period=10)
1425 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1426 diaSources = pd.DataFrame(
1427 data={"diaObjectId": n_sources * [objId],
1428 "band": n_sources//2 * ["u"] + n_sources//2 * ["g"],
1429 "diaSourceId": np.arange(n_sources, dtype=int),
1430 "midpointMjdTai": times,
1431 "psfFlux": fluxes,
1432 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1434 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),
1435 "ap_lombScarglePeriodogramMulti",
1436 None)
1438 run_multiband_plugin(diaObjects, diaSources, plugin)
1439 self.assertAlmostEqual(diaObjects.at[objId, "multiPeriod"], 10, delta=0.04)
1440 self.assertAlmostEqual(diaObjects.at[objId, "multiPower"], 1, delta=1e-2)
1441 # This implementation of LS returns a normalized power < 1.
1442 self.assertLess(diaObjects.at[objId, "multiPower"], 1)
1443 self.assertAlmostEqual(diaObjects.at[objId, "multiFap"], 0, delta=0.04)
1444 # Note: The below values are empirical, but seem reasonable, and
1445 # test that we get values for each band.
1446 self.assertAlmostEqual(diaObjects.at[objId, "u_multiAmp"], 0.029, delta=0.01)
1447 self.assertAlmostEqual(diaObjects.at[objId, "g_multiAmp"], 0.029, delta=0.01)
1448 self.assertAlmostEqual(diaObjects.at[objId, "u_multiPhase"], -2.0, delta=0.2)
1449 self.assertAlmostEqual(diaObjects.at[objId, "g_multiPhase"], 1.0, delta=0.1)
1451 def testCalculateTwoSources(self):
1452 """Test Mulitband Lomb Scargle Periodogram with 2 sources (minimum
1453 detections = 5), which will result in NaN output."""
1454 objId = 0
1455 n_sources = 2
1456 times, fluxes = generatePeriodicData(n_sources, period=10)
1457 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1458 diaSources = pd.DataFrame(
1459 data={"diaObjectId": n_sources * [objId],
1460 "band": n_sources * ["u"],
1461 "diaSourceId": np.arange(n_sources, dtype=int),
1462 "midpointMjdTai": times,
1463 "psfFlux": fluxes,
1464 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1466 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),
1467 "ap_lombScarglePeriodogramMulti",
1468 None)
1470 run_multi_plugin(diaObjects, diaSources, "u", plugin)
1471 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPeriod"]))
1472 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPower"]))
1473 self.assertTrue(np.isnan(diaObjects.at[objId, "multiFap"]))
1476class TestLombScarglePeriodogram(lsst.utils.tests.TestCase):
1478 def testCalculate(self):
1479 """Test Lomb Scargle Periodogram."""
1480 n_sources = 10
1481 objId = 0
1483 # Test period calculation.
1484 times, fluxes = generatePeriodicData(n_sources, period=10)
1485 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1486 diaSources = pd.DataFrame(
1487 data={"diaObjectId": n_sources * [objId],
1488 "band": n_sources * ["u"],
1489 "diaSourceId": np.arange(n_sources, dtype=int),
1490 "midpointMjdTai": times,
1491 "psfFlux": fluxes,
1492 "psfFluxErr": 1e-3+np.zeros(n_sources)})
1494 plugin = LombScarglePeriodogram(LombScarglePeriodogramConfig(),
1495 "ap_lombScarglePeriodogram",
1496 None)
1498 run_multi_plugin(diaObjects, diaSources, "u", plugin)
1499 self.assertAlmostEqual(diaObjects.at[objId, "u_period"], 10, delta=0.04)
1500 # This implementation of LS returns a normalized power < 1.
1501 self.assertAlmostEqual(diaObjects.at[objId, "u_power"], 1, delta=1e-2)
1502 self.assertLess(diaObjects.at[objId, "u_power"], 1)
1504 # Test that we get the same result with a NaN flux.
1505 diaSources.loc[4, "psfFlux"] = np.nan
1506 diaObjects = pd.DataFrame({"diaObjectId": [objId]})
1507 diaSources = pd.DataFrame(
1508 data={"diaObjectId": n_sources * [objId],
1509 "band": n_sources * ["r"],
1510 "diaSourceId": np.arange(n_sources, dtype=int),
1511 "midpointMjdTai": times,
1512 "psfFlux": fluxes,
1513 "psfFluxErr": np.ones(n_sources)})
1514 run_multi_plugin(diaObjects, diaSources, "r", plugin)
1515 self.assertAlmostEqual(diaObjects.at[objId, "r_period"], 10, delta=0.04)
1516 self.assertAlmostEqual(diaObjects.at[objId, "r_power"], 1, delta=1e-2)
1517 # This implementation of LS returns a normalized power < 1.
1518 self.assertLess(diaObjects.at[objId, "r_power"], 1)
1521class TestSigmaDiaTotFlux(unittest.TestCase):
1523 def testCalculate(self):
1524 """Test flux scatter calculation.
1525 """
1526 n_sources = 10
1527 objId = 0
1529 # Test test scatter on scienceFlux.
1530 fluxes = np.linspace(-1, 1, n_sources)
1531 diaSources = pd.DataFrame(
1532 data={"diaObjectId": n_sources * [objId],
1533 "band": n_sources * ["u"],
1534 "diaSourceId": np.arange(n_sources, dtype=int),
1535 "scienceFlux": fluxes,
1536 "scienceFluxErr": np.ones(n_sources)})
1538 plug = SigmaDiaTotFlux(SigmaDiaTotFluxConfig(),
1539 "ap_sigmaTotFlux",
1540 None)
1541 diaObjects = make_diaObject_table(objId, plug, band='u')
1542 run_multi_plugin(diaObjects, diaSources, "u", plug)
1543 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxSigma"],
1544 np.nanstd(fluxes, ddof=1))
1546 # Test test scatter on scienceFlux returns nan on 1 input.
1547 diaSources = pd.DataFrame(
1548 data={"diaObjectId": 1 * [objId],
1549 "band": 1 * ["g"],
1550 "diaSourceId": np.arange(1, dtype=int),
1551 "scienceFlux": fluxes[0],
1552 "scienceFluxErr": np.ones(1)})
1553 diaObjects = make_diaObject_table(objId, plug, band='g')
1554 run_multi_plugin(diaObjects, diaSources, "g", plug)
1555 self.assertTrue(np.isnan(diaObjects.at[objId, "g_scienceFluxSigma"]))
1557 # Test test scatter on scienceFlux takes input nans.
1558 fluxes[4] = np.nan
1559 diaSources = pd.DataFrame(
1560 data={"diaObjectId": n_sources * [objId],
1561 "band": n_sources * ["r"],
1562 "diaSourceId": np.arange(n_sources, dtype=int),
1563 "scienceFlux": fluxes,
1564 "scienceFluxErr": np.ones(n_sources)})
1565 diaObjects = make_diaObject_table(objId, plug, band='r')
1566 run_multi_plugin(diaObjects, diaSources, "r", plug)
1567 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxSigma"],
1568 np.nanstd(fluxes, ddof=1))
1571def skew_wrapper(values):
1572 """Compute scipy skew, omitting nans.
1574 This version works with both scipy<1.9 (where it erroneously returns a
1575 masked array) and scipy>=1.9 (where it correctly returns a float).
1577 Parameters
1578 ----------
1579 values : `np.ndarray`
1581 Returns
1582 -------
1583 skew_value : `float`
1584 """
1585 value = skew(values, bias=False, nan_policy="omit")
1586 if isinstance(value, np.ma.masked_array): 1586 ↛ 1587line 1586 didn't jump to line 1587 because the condition on line 1586 was never true
1587 return value.data
1588 else:
1589 return value
1592class MemoryTester(lsst.utils.tests.MemoryTestCase):
1593 pass
1596def setup_module(module):
1597 lsst.utils.tests.init()
1600if __name__ == "__main__": 1600 ↛ 1601line 1600 didn't jump to line 1601 because the condition on line 1600 was never true
1601 lsst.utils.tests.init()
1602 unittest.main()