Coverage for tests/test_diaCalculationPlugins.py: 13%

1# This file is part of ap_association.

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.

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)

48import lsst.utils.tests

51def run_single_plugin(diaObjectCat,

52 diaObjectId,

53 diaSourceCat,

54 band,

55 plugin):

56 """Wrapper for running single plugins.

58 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`

60 Parameters

61 ----------

62 diaObjectCat : `pandas.DataFrame`

63 Input object catalog to store data into and read from.

64 diaSourcesCat : `pandas.DataFrame`

65 DiaSource catalog to read data from and groupby on.

66 fitlerName : `str`

67 String name of the filter to process.

68 plugin : `lsst.ap.association.DiaCalculationPlugin`

69 Plugin to run.

70 """

71 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)

72 diaSourceCat.set_index(

73 ["diaObjectId", "band", "diaSourceId"],

74 inplace=True,

75 drop=False)

77 objDiaSources = diaSourceCat.loc[diaObjectId]

78 updatingFilterDiaSources = diaSourceCat.loc[

79 (diaObjectId, band), :

80 ]

82 plugin.calculate(diaObjects=diaObjectCat,

83 diaObjectId=diaObjectId,

84 diaSources=objDiaSources,

85 filterDiaSources=updatingFilterDiaSources,

86 band=band)

89def run_multi_plugin(diaObjectCat, diaSourceCat, band, plugin):

90 """Wrapper for running multi plugins.

92 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`

94 Parameters

95 ----------

96 diaObjectCat : `pandas.DataFrame`

97 Input object catalog to store data into and read from.

98 diaSourcesCat : `pandas.DataFrame`

99 DiaSource catalog to read data from and groupby on.

100 filterName : `str`

101 String name of the filter to process.

102 plugin : `lsst.ap.association.DiaCalculationPlugin`

103 Plugin to run.

104 """

105 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)

106 diaSourceCat.set_index(

107 ["diaObjectId", "band", "diaSourceId"],

108 inplace=True,

109 drop=False)

110

111 updatingFilterDiaSources = diaSourceCat.loc[

112 (slice(None), band), :

113 ]

114

115 diaSourcesGB = diaSourceCat.groupby(level=0)

116 filterDiaSourcesGB = updatingFilterDiaSources.groupby(level=0)

117

118 plugin.calculate(diaObjects=diaObjectCat,

119 diaSources=diaSourcesGB,

120 filterDiaSources=filterDiaSourcesGB,

121 band=band)

122

123

124def run_multiband_plugin(diaObjectCat, diaSourceCat, plugin):

125 """Wrapper for running multi plugins.

126

127 Reproduces some of the behavior of `lsst.ap.association.DiaCalcuation.run`

128

129 Parameters

130 ----------

131 diaObjectCat : `pandas.DataFrame`

132 Input object catalog to store data into and read from.

133 diaSourcesCat : `pandas.DataFrame`

134 DiaSource catalog to read data from and groupby on.

135 plugin : `lsst.ap.association.DiaCalculationPlugin`

136 Plugin to run.

137 """

138 diaObjectCat.set_index("diaObjectId", inplace=True, drop=False)

139 diaSourceCat.set_index(

140 ["diaObjectId", "band", "diaSourceId"],

141 inplace=True,

142 drop=False)

143

144 diaSourcesGB = diaSourceCat.groupby(level=0)

145

146 plugin.calculate(diaObjects=diaObjectCat,

147 diaSources=diaSourcesGB,

148 )

149

150

151def make_diaObject_table(objId, plugin, default_value=np.nan, band=None):

152 """Create a minimal diaObject table with columns required for the plugin

153

154 Parameters

155 ----------

156 objId : `int`

157 The diaObjectId

158 plugin : `lsst.ap.association.DiaCalculationPlugin`

159 The plugin that will be run.

160 default_value : `float` or `int`, optional

161 Value to set new columns to.

162 band : `str`, optional

163 Band designation to append to the plugin columns.

164

165 Returns

166 -------

167 diaObjects : `pandas.DataFrame`

168 Output catalog with the required columns for the plugin.

169 """

170 # Add an extra empty diaObject here. This ensures that

171 # we properly test the source/object matching implicit

172 # in the plugin calculations.

173 diaObjects = {"diaObjectId": [objId, objId + 1]}

174 for col in plugin.outputCols:

175 if band is not None:

176 diaObjects[f"{band}_{col}"] = default_value

177 else:

178 diaObjects[col] = default_value

179 return pd.DataFrame(diaObjects)

180

181

182class TestMeanPosition(unittest.TestCase):

183

184 def testCalculate(self):

185 """Test mean position calculation.

186 """

187 n_sources = 10

188 objId = 0

189

190 plug = MeanDiaPosition(MeanDiaPositionConfig(),

191 "ap_meanPosition",

192 None)

193

194 # Test expected means in RA.

195 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

196 diaSources = pd.DataFrame(data={"ra": np.linspace(-1, 1, n_sources),

197 "dec": np.zeros(n_sources),

198 "midpointMjdTai": np.linspace(0, n_sources, n_sources),

199 "diaObjectId": n_sources * [objId],

200 "band": n_sources * ["g"],

201 "diaSourceId": np.arange(n_sources,

202 dtype=int)})

203 run_multi_plugin(diaObjects, diaSources, "g", plug)

204

205 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)

206 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)

207

208 # Test expected means in DEC.

209 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

210 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),

211 "dec": np.linspace(-1, 1, n_sources),

212 "midpointMjdTai": np.linspace(0, n_sources, n_sources),

213 "diaObjectId": n_sources * [objId],

214 "band": n_sources * ["g"],

215 "diaSourceId": np.arange(n_sources,

216 dtype=int)})

217 run_multi_plugin(diaObjects, diaSources, "g", plug)

218

219 self.assertAlmostEqual(diaObjects.loc[objId, "ra"], 0.0)

220 self.assertAlmostEqual(diaObjects.loc[objId, "dec"], 0.0)

221

222 # Test failure mode RA is nan.

223 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

224 diaSources = pd.DataFrame(data={"ra": np.full(n_sources, np.nan),

225 "dec": np.zeros(n_sources),

226 "midpointMjdTai": np.linspace(0, n_sources, n_sources),

227 "diaObjectId": n_sources * [objId],

228 "band": n_sources * ["g"],

229 "diaSourceId": np.arange(n_sources,

230 dtype=int)})

231 run_multi_plugin(diaObjects, diaSources, "g", plug)

232

233 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))

234 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))

235

236 # Test failure mode DEC is nan.

237 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

238 diaSources = pd.DataFrame(data={"ra": np.zeros(n_sources),

239 "dec": np.full(n_sources, np.nan),

240 "midpointMjdTai": np.linspace(0, n_sources, n_sources),

241 "diaObjectId": n_sources * [objId],

242 "band": n_sources * ["g"],

243 "diaSourceId": np.arange(n_sources,

244 dtype=int)})

245 run_multi_plugin(diaObjects, diaSources, "g", plug)

246

247 self.assertTrue(np.isnan(diaObjects.loc[objId, "ra"]))

248 self.assertTrue(np.isnan(diaObjects.loc[objId, "dec"]))

249

250

251class TestHTMIndexPosition(unittest.TestCase):

252

253 def testCalculate(self):

254 """Test HTMPixel assignment calculation.

255 """

256 # Test expected pixelId at RA, DEC = 0

257 objId = 0

258 n_sources = 10

259 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

260 diaObjects.loc[objId, "ra"] = 0.

261 diaObjects.loc[objId, "dec"] = 0.

262 diaSources = pd.DataFrame(

263 data={"diaObjectId": n_sources * [objId],

264 "band": n_sources * ["g"],

265 "diaSourceId": np.arange(n_sources, dtype=int)})

266 plug = HTMIndexDiaPosition(HTMIndexDiaPositionConfig(),

267 "ap_HTMIndex",

268 None)

269

270 run_single_plugin(diaObjectCat=diaObjects,

271 diaObjectId=objId,

272 diaSourceCat=diaSources,

273 band="g",

274 plugin=plug)

275 self.assertEqual(diaObjects.at[objId, "pixelId"],

276 17042430230528)

277

278 # Test expected pixelId at some value of RA and DEC.

279 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

280 diaObjects.loc[objId, "ra"] = 45.37

281 diaObjects.loc[objId, "dec"] = 13.67

282 diaSources = pd.DataFrame(

283 data={"diaObjectId": n_sources * [objId],

284 "band": n_sources * ["g"],

285 "diaSourceId": np.arange(n_sources, dtype=int)})

286 run_single_plugin(diaObjectCat=diaObjects,

287 diaObjectId=objId,

288 diaSourceCat=diaSources,

289 band="g",

290 plugin=plug)

291 self.assertEqual(diaObjects.at[objId, "pixelId"],

292 17450571968473)

293

294

295class TestNDiaSourcesDiaPlugin(unittest.TestCase):

296

297 def testCalculate(self):

298 """Test that the number of DiaSources is correct.

299 """

300

301 for n_sources in [1, 8, 10]:

302 # Test expected number of sources per object.

303 objId = 0

304 diaSources = pd.DataFrame(

305 data={"diaObjectId": n_sources * [objId],

306 "band": n_sources * ["g"],

307 "diaSourceId": np.arange(n_sources, dtype=int)})

308 plug = NumDiaSourcesDiaPlugin(NumDiaSourcesDiaPluginConfig(),

309 "ap_nDiaSources",

310 None)

311 diaObjects = make_diaObject_table(objId, plug, default_value=int(-1))

312 run_multi_plugin(diaObjects, diaSources, "g", plug)

313

314 self.assertEqual(n_sources, diaObjects.at[objId, "nDiaSources"])

315 self.assertEqual(diaObjects["nDiaSources"].dtype, np.int64)

316

317

318class TestSimpleSourceFlagDiaPlugin(unittest.TestCase):

319

320 def testCalculate(self):

321 """Test that DiaObject flags are set.

322 """

323 objId = 0

324 n_sources = 10

325

326 # Test expected flags, no flags set.

327 diaSources = pd.DataFrame(

328 data={"diaObjectId": n_sources * [objId],

329 "band": n_sources * ["g"],

330 "diaSourceId": np.arange(n_sources, dtype=int),

331 "flags": np.zeros(n_sources, dtype=np.uint64)})

332 plug = SimpleSourceFlagDiaPlugin(SimpleSourceFlagDiaPluginConfig(),

333 "ap_diaObjectFlag",

334 None)

335

336 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))

337 run_multi_plugin(diaObjects, diaSources, "g", plug)

338 self.assertEqual(diaObjects.at[objId, "flags"], 0)

339 self.assertEqual(diaObjects["flags"].dtype, np.uint64)

340

341 # Test expected flags, all flags set.

342 diaSources = pd.DataFrame(

343 data={"diaObjectId": n_sources * [objId],

344 "band": n_sources * ["g"],

345 "diaSourceId": np.arange(n_sources, dtype=int),

346 "flags": np.ones(n_sources, dtype=np.uint64)})

347 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))

348 run_multi_plugin(diaObjects, diaSources, "g", plug)

349 self.assertEqual(diaObjects.at[objId, "flags"], 1)

350 self.assertEqual(diaObjects["flags"].dtype, np.uint64)

351

352 # Test expected flags, random flags.

353 diaSources = pd.DataFrame(

354 data={"diaObjectId": n_sources * [objId],

355 "band": n_sources * ["g"],

356 "diaSourceId": np.arange(n_sources, dtype=int),

357 "flags": np.random.randint(0, 2 ** 16, size=n_sources)})

358

359 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))

360 run_multi_plugin(diaObjects, diaSources, "g", plug)

361 self.assertEqual(diaObjects.at[objId, "flags"], 1)

362 self.assertEqual(diaObjects["flags"].dtype, np.uint64)

363

364 # Test expected flags, one flag set.

365 flag_array = np.zeros(n_sources, dtype=np.uint64)

366 flag_array[4] = 256

367 diaSources = pd.DataFrame(

368 data={"diaObjectId": n_sources * [objId],

369 "band": n_sources * ["g"],

370 "diaSourceId": np.arange(n_sources, dtype=int),

371 "flags": flag_array})

372 diaObjects = make_diaObject_table(objId, plug, default_value=np.uint64(0))

373 run_multi_plugin(diaObjects, diaSources, "g", plug)

374 self.assertEqual(diaObjects.at[objId, "flags"], 1)

375 self.assertEqual(diaObjects["flags"].dtype, np.uint64)

376

377

378class TestWeightedMeanDiaPsfFlux(unittest.TestCase):

379

380 def testCalculate(self):

381 """Test mean value calculation.

382 """

383 n_sources = 10

384 objId = 0

385

386 # Test expected mean.

387 # In the first test, we have only one object.

388 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

389 diaSources = pd.DataFrame(

390 data={"diaObjectId": n_sources * [objId],

391 "band": n_sources * ["u"],

392 "diaSourceId": np.arange(n_sources, dtype=int),

393 "psfFlux": np.linspace(-1, 1, n_sources),

394 "psfFluxErr": np.ones(n_sources)})

395

396 plug = WeightedMeanDiaPsfFlux(WeightedMeanDiaPsfFluxConfig(),

397 "ap_meanFlux",

398 None)

399 run_multi_plugin(diaObjects, diaSources, "u", plug)

400

401 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMean"], 0.0)

402 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxMeanErr"],

403 np.sqrt(1 / n_sources))

404 self.assertEqual(diaObjects.loc[objId, "u_psfFluxNdata"], n_sources)

405 # We expect this to be converted to float.

406 # TODO DM-53254: This should be an integer (and should be checked

407 # to be an integer).

408 self.assertEqual(diaObjects["u_psfFluxNdata"].dtype, np.float64)

409

410 # Test expected mean with a nan value.

411 # In the second test, we have two objects (one empty).

412 diaObjects = pd.DataFrame({"diaObjectId": [objId, objId + 1]})

413 fluxes = np.linspace(-1, 1, n_sources)

414 fluxes[4] = np.nan

415 diaSources = pd.DataFrame(

416 data={"diaObjectId": n_sources * [objId],

417 "band": n_sources * ["r"],

418 "diaSourceId": np.arange(n_sources, dtype=int),

419 "psfFlux": fluxes,

420 "psfFluxErr": np.ones(n_sources)})

421 run_multi_plugin(diaObjects, diaSources, "r", plug)

422

423 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMean"],

424 np.nanmean(fluxes))

425 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMeanErr"],

426 np.sqrt(1 / (n_sources - 1)))

427 self.assertEqual(diaObjects.loc[objId, "r_psfFluxNdata"], n_sources - 1)

428 # We expect this to be converted to float.

429 # TODO DM-53254: This should be an integer (and should be checked

430 # to be an integer).

431 self.assertEqual(diaObjects["r_psfFluxNdata"].dtype, np.float64)

432

433

434class TestPercentileDiaPsfFlux(unittest.TestCase):

435

436 def testCalculate(self):

437 """Test flux percentile calculation.

438 """

439 n_sources = 10

440 objId = 0

441

442 # Test expected percentile values.

443 fluxes = np.linspace(-1, 1, n_sources)

444 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

445 diaSources = pd.DataFrame(

446 data={"diaObjectId": n_sources * [objId],

447 "band": n_sources * ["u"],

448 "diaSourceId": np.arange(n_sources, dtype=int),

449 "psfFlux": fluxes,

450 "psfFluxErr": np.ones(n_sources)})

451

452 plug = PercentileDiaPsfFlux(PercentileDiaPsfFluxConfig(),

453 "ap_percentileFlux",

454 None)

455 run_multi_plugin(diaObjects, diaSources, "u", plug)

456 for pTile, testVal in zip(plug.config.percentiles,

457 np.nanpercentile(

458 fluxes,

459 plug.config.percentiles)):

460 self.assertAlmostEqual(

461 diaObjects.at[objId, "u_psfFluxPercentile{:02d}".format(pTile)],

462 testVal)

463

464 # Test expected percentile values with a nan value.

465 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

466 fluxes[4] = np.nan

467 diaSources = pd.DataFrame(

468 data={"diaObjectId": n_sources * [objId],

469 "band": n_sources * ["r"],

470 "diaSourceId": np.arange(n_sources, dtype=int),

471 "psfFlux": fluxes,

472 "psfFluxErr": np.ones(n_sources)})

473 run_multi_plugin(diaObjects, diaSources, "r", plug)

474 for pTile, testVal in zip(plug.config.percentiles,

475 np.nanpercentile(

476 fluxes,

477 plug.config.percentiles)):

478 self.assertAlmostEqual(

479 diaObjects.at[objId, "r_psfFluxPercentile{:02d}".format(pTile)],

480 testVal)

481

482

483class TestSigmaDiaPsfFlux(unittest.TestCase):

484

485 def testCalculate(self):

486 """Test flux scatter calculation.

487 """

488 n_sources = 10

489 objId = 0

490

491 # Test expected sigma scatter of fluxes.

492 fluxes = np.linspace(-1, 1, n_sources)

493 diaSources = pd.DataFrame(

494 data={"diaObjectId": n_sources * [objId],

495 "band": n_sources * ["u"],

496 "diaSourceId": np.arange(n_sources, dtype=int),

497 "psfFlux": fluxes,

498 "psfFluxErr": np.ones(n_sources)})

499

500 plug = SigmaDiaPsfFlux(SigmaDiaPsfFluxConfig(),

501 "ap_sigmaFlux",

502 None)

503 diaObjects = make_diaObject_table(objId, plug, band='u')

504 run_multi_plugin(diaObjects, diaSources, "u", plug)

505 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxSigma"],

506 np.nanstd(fluxes, ddof=1))

507

508 # test one input, returns nan.

509 diaSources = pd.DataFrame(

510 data={"diaObjectId": 1 * [objId],

511 "band": 1 * ["g"],

512 "diaSourceId": [0],

513 "psfFlux": [fluxes[0]],

514 "psfFluxErr": [1.]})

515

516 diaObjects = make_diaObject_table(objId, plug, band='g')

517 run_multi_plugin(diaObjects, diaSources, "g", plug)

518 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxSigma"]))

519

520 # Test expected sigma scatter of fluxes with a nan value.

521 fluxes[4] = np.nan

522 diaSources = pd.DataFrame(

523 data={"diaObjectId": n_sources * [objId],

524 "band": n_sources * ["r"],

525 "diaSourceId": np.arange(n_sources, dtype=int),

526 "psfFlux": fluxes,

527 "psfFluxErr": np.ones(n_sources)})

528

529 diaObjects = make_diaObject_table(objId, plug, band='r')

530 run_multi_plugin(diaObjects, diaSources, "r", plug)

531 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxSigma"],

532 np.nanstd(fluxes, ddof=1))

533

534

535class TestChi2DiaPsfFlux(unittest.TestCase):

536

537 def testCalculate(self):

538 """Test flux chi2 calculation.

539 """

540 n_sources = 10

541 objId = 0

542

543 # Test expected chi^2 value.

544 fluxes = np.linspace(-1, 1, n_sources)

545 diaObjects = pd.DataFrame({"diaObjectId": [objId],

546 "u_psfFluxMean": [0.0]})

547 diaSources = pd.DataFrame(

548 data={"diaObjectId": n_sources * [objId],

549 "band": n_sources * ["u"],

550 "diaSourceId": np.arange(n_sources, dtype=int),

551 "psfFlux": fluxes,

552 "psfFluxErr": np.ones(n_sources)})

553

554 plug = Chi2DiaPsfFlux(Chi2DiaPsfFluxConfig(),

555 "ap_chi2Flux",

556 None)

557 run_multi_plugin(diaObjects, diaSources, "u", plug)

558 self.assertAlmostEqual(

559 diaObjects.loc[objId, "u_psfFluxChi2"],

560 np.nansum(((diaSources["psfFlux"]

561 - np.nanmean(diaSources["psfFlux"]))

562 / diaSources["psfFluxErr"]) ** 2))

563

564 # Test expected chi^2 value with a nan value set.

565 fluxes[4] = np.nan

566 diaObjects = pd.DataFrame({"diaObjectId": [objId],

567 "r_psfFluxMean": [np.nanmean(fluxes)]})

568 diaSources = pd.DataFrame(

569 data={"diaObjectId": n_sources * [objId],

570 "band": n_sources * ["r"],

571 "diaSourceId": np.arange(n_sources, dtype=int),

572 "psfFlux": fluxes,

573 "psfFluxErr": np.ones(n_sources)})

574 run_multi_plugin(diaObjects, diaSources, "r", plug)

575 self.assertAlmostEqual(

576 diaObjects.loc[objId, "r_psfFluxChi2"],

577 np.nansum(((diaSources["psfFlux"]

578 - np.nanmean(diaSources["psfFlux"]))

579 / diaSources["psfFluxErr"]) ** 2))

580

581

582class TestMadDiaPsfFlux(unittest.TestCase):

583

584 def testCalculate(self):

585 """Test flux median absolute deviation calculation.

586 """

587 n_sources = 10

588 objId = 0

589

590 # Test expected MAD value.

591 fluxes = np.linspace(-1, 1, n_sources)

592 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

593 diaSources = pd.DataFrame(

594 data={"diaObjectId": n_sources * [objId],

595 "band": n_sources * ["u"],

596 "diaSourceId": np.arange(n_sources, dtype=int),

597 "psfFlux": fluxes,

598 "psfFluxErr": np.ones(n_sources)})

599

600 plug = MadDiaPsfFlux(MadDiaPsfFluxConfig(),

601 "ap_madFlux",

602 None)

603 run_multi_plugin(diaObjects, diaSources, "u", plug)

604 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMAD"],

605 median_absolute_deviation(fluxes,

606 ignore_nan=True))

607

608 # Test expected MAD value with a nan set.

609 fluxes[4] = np.nan

610 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

611 diaSources = pd.DataFrame(

612 data={"diaObjectId": n_sources * [objId],

613 "band": n_sources * ["r"],

614 "diaSourceId": np.arange(n_sources, dtype=int),

615 "psfFlux": fluxes,

616 "psfFluxErr": np.ones(n_sources)})

617 run_multi_plugin(diaObjects, diaSources, "r", plug)

618 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMAD"],

619 median_absolute_deviation(fluxes,

620 ignore_nan=True))

621

622

623class TestSkewDiaPsfFlux(unittest.TestCase):

624

625 def testCalculate(self):

626 """Test flux skew calculation.

627 """

628 n_sources = 10

629 objId = 0

630

631 # Test expected skew value.

632 fluxes = np.linspace(-1, 1, n_sources)

633 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

634 diaSources = pd.DataFrame(

635 data={"diaObjectId": n_sources * [objId],

636 "band": n_sources * ["u"],

637 "diaSourceId": np.arange(n_sources, dtype=int),

638 "psfFlux": fluxes,

639 "psfFluxErr": np.ones(n_sources)})

640

641 plug = SkewDiaPsfFlux(SkewDiaPsfFluxConfig(),

642 "ap_skewFlux",

643 None)

644 run_multi_plugin(diaObjects, diaSources, "u", plug)

645 self.assertAlmostEqual(

646 diaObjects.loc[objId, "u_psfFluxSkew"],

647 skew_wrapper(fluxes))

648

649 # Test expected skew value with a nan set.

650 fluxes[4] = np.nan

651 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

652 diaSources = pd.DataFrame(

653 data={"diaObjectId": n_sources * [objId],

654 "band": n_sources * ["r"],

655 "diaSourceId": np.arange(n_sources, dtype=int),

656 "psfFlux": fluxes,

657 "psfFluxErr": np.ones(n_sources)})

658 run_multi_plugin(diaObjects, diaSources, "r", plug)

659

660 self.assertAlmostEqual(

661 diaObjects.at[objId, "r_psfFluxSkew"],

662 skew_wrapper(fluxes))

663

664

665class TestMinMaxDiaPsfFlux(unittest.TestCase):

666

667 def testCalculate(self):

668 """Test flux min/max calculation.

669 """

670 n_sources = 10

671 objId = 0

672

673 # Test expected MinMax fluxes.

674 fluxes = np.linspace(-1, 1, n_sources)

675 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

676 diaSources = pd.DataFrame(

677 data={"diaObjectId": n_sources * [objId],

678 "band": n_sources * ["u"],

679 "diaSourceId": np.arange(n_sources, dtype=int),

680 "psfFlux": fluxes,

681 "psfFluxErr": np.ones(n_sources)})

682

683 plug = MinMaxDiaPsfFlux(MinMaxDiaPsfFluxConfig(),

684 "ap_minMaxFlux",

685 None)

686 run_multi_plugin(diaObjects, diaSources, "u", plug)

687 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMin"], -1)

688 self.assertEqual(diaObjects.loc[objId, "u_psfFluxMax"], 1)

689

690 # Test expected MinMax fluxes with a nan set.

691 fluxes[4] = np.nan

692 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

693 diaSources = pd.DataFrame(

694 data={"diaObjectId": n_sources * [objId],

695 "band": n_sources * ["r"],

696 "diaSourceId": np.arange(n_sources, dtype=int),

697 "psfFlux": fluxes,

698 "psfFluxErr": np.ones(n_sources)})

699 run_multi_plugin(diaObjects, diaSources, "r", plug)

700 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMin"], -1)

701 self.assertEqual(diaObjects.loc[objId, "r_psfFluxMax"], 1)

702

703

704class TestMaxSlopeDiaPsfFlux(unittest.TestCase):

705

706 def testCalculate(self):

707 """Test flux maximum slope.

708 """

709 n_sources = 10

710 objId = 0

711

712 # Test max slope value.

713 fluxes = np.linspace(-1, 1, n_sources)

714 times = np.concatenate([np.linspace(0, 1, n_sources)[:-1], [1 - 1/90]])

715 diaSources = pd.DataFrame(

716 data={"diaObjectId": n_sources * [objId],

717 "band": n_sources * ["u"],

718 "diaSourceId": np.arange(n_sources, dtype=int),

719 "psfFlux": fluxes,

720 "psfFluxErr": np.ones(n_sources),

721 "midpointMjdTai": times})

722

723 plug = MaxSlopeDiaPsfFlux(MaxSlopeDiaPsfFluxConfig(),

724 "ap_maxSlopeFlux",

725 None)

726 diaObjects = make_diaObject_table(objId, plug, band='u')

727 run_multi_plugin(diaObjects, diaSources, "u", plug)

728 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxMaxSlope"], 2 + 2/9)

729

730 # Test max slope value returns nan on 1 input.

731 diaSources = pd.DataFrame(

732 data={"diaObjectId": 1 * [objId],

733 "band": 1 * ["g"],

734 "diaSourceId": np.arange(1, dtype=int),

735 "psfFlux": fluxes[0],

736 "psfFluxErr": np.ones(1),

737 "midpointMjdTai": times[0]})

738 diaObjects = make_diaObject_table(objId, plug, band='g')

739 run_multi_plugin(diaObjects, diaSources, "g", plug)

740 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxMaxSlope"]))

741

742 # Test max slope value inputing nan values.

743 fluxes[4] = np.nan

744 times[7] = np.nan

745 diaSources = pd.DataFrame(

746 data={"diaObjectId": n_sources * [objId],

747 "band": n_sources * ["r"],

748 "diaSourceId": np.arange(n_sources, dtype=int),

749 "psfFlux": fluxes,

750 "psfFluxErr": np.ones(n_sources),

751 "midpointMjdTai": times})

752 diaObjects = make_diaObject_table(objId, plug, band='r')

753 run_multi_plugin(diaObjects, diaSources, "r", plug)

754 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxMaxSlope"], 2 + 2 / 9)

755

756

757class TestErrMeanDiaPsfFlux(unittest.TestCase):

758

759 def testCalculate(self):

760 """Test error mean calculation.

761 """

762 n_sources = 10

763 objId = 0

764

765 # Test mean of the errors.

766 fluxes = np.linspace(-1, 1, n_sources)

767 errors = np.linspace(1, 2, n_sources)

768 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

769 diaSources = pd.DataFrame(

770 data={"diaObjectId": n_sources * [objId],

771 "band": n_sources * ["u"],

772 "diaSourceId": np.arange(n_sources, dtype=int),

773 "psfFlux": fluxes,

774 "psfFluxErr": errors})

775

776 plug = ErrMeanDiaPsfFlux(ErrMeanDiaPsfFluxConfig(),

777 "ap_errMeanFlux",

778 None)

779 run_multi_plugin(diaObjects, diaSources, "u", plug)

780 self.assertAlmostEqual(diaObjects.at[objId, "u_psfFluxErrMean"],

781 np.nanmean(errors).astype(np.float32))

782

783 # Test mean of the errors with input nan value.

784 errors[4] = np.nan

785 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

786 diaSources = pd.DataFrame(

787 data={"diaObjectId": n_sources * [objId],

788 "band": n_sources * ["r"],

789 "diaSourceId": np.arange(n_sources, dtype=int),

790 "psfFlux": fluxes,

791 "psfFluxErr": errors})

792 run_multi_plugin(diaObjects, diaSources, "r", plug)

793 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxErrMean"],

794 np.nanmean(errors).astype(np.float32))

795

796

797class TestLinearFitDiaPsfFlux(unittest.TestCase):

798

799 def testCalculate(self):

800 """Test a linear fit to flux vs time.

801 """

802 n_sources = 10

803 objId = 0

804

805 # Test best fit linear model.

806 fluxes = np.linspace(-1, 1, n_sources)

807 errors = np.linspace(1, 2, n_sources)

808 times = np.linspace(0, 1, n_sources)

809 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

810 diaSources = pd.DataFrame(

811 data={"diaObjectId": n_sources * [objId],

812 "band": n_sources * ["u"],

813 "diaSourceId": np.arange(n_sources, dtype=int),

814 "psfFlux": fluxes,

815 "psfFluxErr": errors,

816 "midpointMjdTai": times})

817

818 plug = LinearFitDiaPsfFlux(LinearFitDiaPsfFluxConfig(),

819 "ap_LinearFit",

820 None)

821 run_multi_plugin(diaObjects, diaSources, "u", plug)

822 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearSlope"],

823 2.)

824 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxLinearIntercept"],

825 -1.)

826

827 # Test best fit linear model with input nans.

828 fluxes[7] = np.nan

829 errors[4] = np.nan

830 times[2] = np.nan

831 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

832 diaSources = pd.DataFrame(

833 data={"diaObjectId": n_sources * [objId],

834 "band": n_sources * ["r"],

835 "diaSourceId": np.arange(n_sources, dtype=int),

836 "psfFlux": fluxes,

837 "psfFluxErr": errors,

838 "midpointMjdTai": times})

839 run_multi_plugin(diaObjects, diaSources, "r", plug)

840 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearSlope"], 2.)

841 self.assertAlmostEqual(diaObjects.loc[objId, "r_psfFluxLinearIntercept"],

842 -1.)

843

844

845class TestStetsonJDiaPsfFlux(unittest.TestCase):

846

847 def testCalculate(self):

848 """Test the stetsonJ statistic.

849 """

850 n_sources = 10

851 objId = 0

852

853 # Test stetsonJ calculation.

854 fluxes = np.linspace(-1, 1, n_sources)

855 errors = np.ones(n_sources)

856 diaObjects = pd.DataFrame({"diaObjectId": [objId],

857 "u_psfFluxMean": [np.nanmean(fluxes)]})

858 diaSources = pd.DataFrame(

859 data={"diaObjectId": n_sources * [objId],

860 "band": n_sources * ["u"],

861 "diaSourceId": np.arange(n_sources, dtype=int),

862 "psfFlux": fluxes,

863 "psfFluxErr": errors})

864

865 plug = StetsonJDiaPsfFlux(StetsonJDiaPsfFluxConfig(),

866 "ap_StetsonJ",

867 None)

868 run_multi_plugin(diaObjects, diaSources, "u", plug)

869 # Expected StetsonJ for the values created. Confirmed using Cesimum's

870 # implementation. http://github.com/cesium-ml/cesium

871 self.assertAlmostEqual(diaObjects.loc[objId, "u_psfFluxStetsonJ"],

872 -0.5958393936080928)

873

874 # Test stetsonJ calculation returns nan on single input.

875 diaObjects = pd.DataFrame({"diaObjectId": [objId],

876 "g_psfFluxMean": [np.nanmean(fluxes)]})

877 diaSources = pd.DataFrame(

878 data={"diaObjectId": 1 * [objId],

879 "band": 1 * ["g"],

880 "diaSourceId": np.arange(1, dtype=int),

881 "psfFlux": fluxes[0],

882 "psfFluxErr": errors[0]})

883 run_multi_plugin(diaObjects, diaSources, "g", plug)

884 self.assertTrue(np.isnan(diaObjects.at[objId, "g_psfFluxStetsonJ"]))

885

886 # Test stetsonJ calculation returns when nans are input.

887 fluxes[7] = np.nan

888 errors[4] = np.nan

889 nonNanMask = np.logical_and(~np.isnan(fluxes),

890 ~np.isnan(errors))

891 diaObjects = pd.DataFrame(

892 {"diaObjectId": [objId],

893 "r_psfFluxMean": [np.average(fluxes[nonNanMask],

894 weights=errors[nonNanMask])]})

895 diaSources = pd.DataFrame(

896 data={"diaObjectId": n_sources * [objId],

897 "band": n_sources * ["r"],

898 "diaSourceId": np.arange(n_sources, dtype=int),

899 "psfFlux": fluxes,

900 "psfFluxErr": errors})

901 run_multi_plugin(diaObjects, diaSources, "r", plug)

902 self.assertAlmostEqual(diaObjects.at[objId, "r_psfFluxStetsonJ"],

903 -0.5412797916187173)

904

905

906class TestWeightedMeanDiaTotFlux(unittest.TestCase):

907

908 def testCalculate(self):

909 """Test mean value calculation.

910 """

911 n_sources = 10

912 objId = 0

913

914 # Test test mean on scienceFlux.

915 fluxes = np.linspace(-1, 1, n_sources)

916 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

917 diaSources = pd.DataFrame(

918 data={"diaObjectId": n_sources * [objId],

919 "band": n_sources * ["u"],

920 "diaSourceId": np.arange(n_sources, dtype=int),

921 "scienceFlux": fluxes,

922 "scienceFluxErr": np.ones(n_sources)})

923

924 plug = WeightedMeanDiaTotFlux(WeightedMeanDiaTotFluxConfig(),

925 "ap_meanTotFlux",

926 None)

927 run_multi_plugin(diaObjects, diaSources, "u", plug)

928

929 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMean"], 0.0)

930 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxMeanErr"],

931 np.sqrt(1 / n_sources))

932

933 # Test test mean on scienceFlux with input nans

934 fluxes[4] = np.nan

935 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

936 diaSources = pd.DataFrame(

937 data={"diaObjectId": n_sources * [objId],

938 "band": n_sources * ["r"],

939 "diaSourceId": np.arange(n_sources, dtype=int),

940 "scienceFlux": fluxes,

941 "scienceFluxErr": np.ones(n_sources)})

942 run_multi_plugin(diaObjects, diaSources, "r", plug)

943

944 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMean"],

945 np.nanmean(fluxes))

946 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxMeanErr"],

947 np.sqrt(1 / (n_sources - 1)))

948

949

950def generatePeriodicData(n=10, period=10):

951 """Generate noisy, sinusoidally-varying periodic data for testing Lomb-

952 Scargle Periodogram.

953

954 The returned fluxes will have, within the errors, the passed-in period and

955 a power close to 1, because the fluxes are purely sinusoidal.

956

957 Parameters

958 ----------

959 n : int

960 Number of data points to generate.

961 period : float

962 Period of the periodic signal.

963

964 Returns

965 -------

966 t : np.ndarray

967 Time values.

968 y_obs : np.ndarray

969 Observed flux values.

970 """

971 np.random.seed(42)

972

973 t = np.linspace(-2*np.pi, 2*np.pi, n) + 100*np.random.random(n)

974 y = 10 + np.sin(2 * np.pi * t / period)

975 y_obs = np.random.normal(y, 0.001)

976

977 return t, y_obs

978

979

980class TestMultiLombScarglePeriodogram(lsst.utils.tests.TestCase):

981

982 def testCalculate(self):

983 """Test Mulitband Lomb Scargle Periodogram."""

984 n_sources = 10

985 objId = 0

986

987 # Create synthetic multi-band data

988 times, fluxes = generatePeriodicData(n_sources, period=10)

989 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

990 diaSources = pd.DataFrame(

991 data={"diaObjectId": n_sources * [objId],

992 "band": n_sources//2 * ["u"] + n_sources//2 * ["g"],

993 "diaSourceId": np.arange(n_sources, dtype=int),

994 "midpointMjdTai": times,

995 "psfFlux": fluxes,

996 "psfFluxErr": 1e-3+np.zeros(n_sources)})

997

998 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),

999 "ap_lombScarglePeriodogramMulti",

1000 None)

1001

1002 run_multiband_plugin(diaObjects, diaSources, plugin)

1003 self.assertAlmostEqual(diaObjects.at[objId, "multiPeriod"], 10, delta=0.04)

1004 self.assertAlmostEqual(diaObjects.at[objId, "multiPower"], 1, delta=1e-2)

1005 # This implementation of LS returns a normalized power < 1.

1006 self.assertLess(diaObjects.at[objId, "multiPower"], 1)

1007 self.assertAlmostEqual(diaObjects.at[objId, "multiFap"], 0, delta=0.04)

1008 # Note: The below values are empirical, but seem reasonable, and

1009 # test that we get values for each band.

1010 self.assertAlmostEqual(diaObjects.at[objId, "u_multiAmp"], 0.029, delta=0.01)

1011 self.assertAlmostEqual(diaObjects.at[objId, "g_multiAmp"], 0.029, delta=0.01)

1012 self.assertAlmostEqual(diaObjects.at[objId, "u_multiPhase"], -2.0, delta=0.2)

1013 self.assertAlmostEqual(diaObjects.at[objId, "g_multiPhase"], 1.0, delta=0.1)

1014

1015 def testCalculateTwoSources(self):

1016 """Test Mulitband Lomb Scargle Periodogram with 2 sources (minimum

1017 detections = 5), which will result in NaN output."""

1018 objId = 0

1019 n_sources = 2

1020 times, fluxes = generatePeriodicData(n_sources, period=10)

1021 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

1022 diaSources = pd.DataFrame(

1023 data={"diaObjectId": n_sources * [objId],

1024 "band": n_sources * ["u"],

1025 "diaSourceId": np.arange(n_sources, dtype=int),

1026 "midpointMjdTai": times,

1027 "psfFlux": fluxes,

1028 "psfFluxErr": 1e-3+np.zeros(n_sources)})

1029

1030 plugin = LombScarglePeriodogramMulti(LombScarglePeriodogramMultiConfig(),

1031 "ap_lombScarglePeriodogramMulti",

1032 None)

1033

1034 run_multi_plugin(diaObjects, diaSources, "u", plugin)

1035 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPeriod"]))

1036 self.assertTrue(np.isnan(diaObjects.at[objId, "multiPower"]))

1037 self.assertTrue(np.isnan(diaObjects.at[objId, "multiFap"]))

1038

1039

1040class TestLombScarglePeriodogram(lsst.utils.tests.TestCase):

1041

1042 def testCalculate(self):

1043 """Test Lomb Scargle Periodogram."""

1044 n_sources = 10

1045 objId = 0

1046

1047 # Test period calculation.

1048 times, fluxes = generatePeriodicData(n_sources, period=10)

1049 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

1050 diaSources = pd.DataFrame(

1051 data={"diaObjectId": n_sources * [objId],

1052 "band": n_sources * ["u"],

1053 "diaSourceId": np.arange(n_sources, dtype=int),

1054 "midpointMjdTai": times,

1055 "psfFlux": fluxes,

1056 "psfFluxErr": 1e-3+np.zeros(n_sources)})

1057

1058 plugin = LombScarglePeriodogram(LombScarglePeriodogramConfig(),

1059 "ap_lombScarglePeriodogram",

1060 None)

1061

1062 run_multi_plugin(diaObjects, diaSources, "u", plugin)

1063 self.assertAlmostEqual(diaObjects.at[objId, "u_period"], 10, delta=0.04)

1064 # This implementation of LS returns a normalized power < 1.

1065 self.assertAlmostEqual(diaObjects.at[objId, "u_power"], 1, delta=1e-2)

1066 self.assertLess(diaObjects.at[objId, "u_power"], 1)

1067

1068 # Test that we get the same result with a NaN flux.

1069 diaSources.loc[4, "psfFlux"] = np.nan

1070 diaObjects = pd.DataFrame({"diaObjectId": [objId]})

1071 diaSources = pd.DataFrame(

1072 data={"diaObjectId": n_sources * [objId],

1073 "band": n_sources * ["r"],

1074 "diaSourceId": np.arange(n_sources, dtype=int),

1075 "midpointMjdTai": times,

1076 "psfFlux": fluxes,

1077 "psfFluxErr": np.ones(n_sources)})

1078 run_multi_plugin(diaObjects, diaSources, "r", plugin)

1079 self.assertAlmostEqual(diaObjects.at[objId, "r_period"], 10, delta=0.04)

1080 self.assertAlmostEqual(diaObjects.at[objId, "r_power"], 1, delta=1e-2)

1081 # This implementation of LS returns a normalized power < 1.

1082 self.assertLess(diaObjects.at[objId, "r_power"], 1)

1083

1084

1085class TestSigmaDiaTotFlux(unittest.TestCase):

1086

1087 def testCalculate(self):

1088 """Test flux scatter calculation.

1089 """

1090 n_sources = 10

1091 objId = 0

1092

1093 # Test test scatter on scienceFlux.

1094 fluxes = np.linspace(-1, 1, n_sources)

1095 diaSources = pd.DataFrame(

1096 data={"diaObjectId": n_sources * [objId],

1097 "band": n_sources * ["u"],

1098 "diaSourceId": np.arange(n_sources, dtype=int),

1099 "scienceFlux": fluxes,

1100 "scienceFluxErr": np.ones(n_sources)})

1101

1102 plug = SigmaDiaTotFlux(SigmaDiaTotFluxConfig(),

1103 "ap_sigmaTotFlux",

1104 None)

1105 diaObjects = make_diaObject_table(objId, plug, band='u')

1106 run_multi_plugin(diaObjects, diaSources, "u", plug)

1107 self.assertAlmostEqual(diaObjects.at[objId, "u_scienceFluxSigma"],

1108 np.nanstd(fluxes, ddof=1))

1109

1110 # Test test scatter on scienceFlux returns nan on 1 input.

1111 diaSources = pd.DataFrame(

1112 data={"diaObjectId": 1 * [objId],

1113 "band": 1 * ["g"],

1114 "diaSourceId": np.arange(1, dtype=int),

1115 "scienceFlux": fluxes[0],

1116 "scienceFluxErr": np.ones(1)})

1117 diaObjects = make_diaObject_table(objId, plug, band='g')

1118 run_multi_plugin(diaObjects, diaSources, "g", plug)

1119 self.assertTrue(np.isnan(diaObjects.at[objId, "g_scienceFluxSigma"]))

1120

1121 # Test test scatter on scienceFlux takes input nans.

1122 fluxes[4] = np.nan

1123 diaSources = pd.DataFrame(

1124 data={"diaObjectId": n_sources * [objId],

1125 "band": n_sources * ["r"],

1126 "diaSourceId": np.arange(n_sources, dtype=int),

1127 "scienceFlux": fluxes,

1128 "scienceFluxErr": np.ones(n_sources)})

1129 diaObjects = make_diaObject_table(objId, plug, band='r')

1130 run_multi_plugin(diaObjects, diaSources, "r", plug)

1131 self.assertAlmostEqual(diaObjects.at[objId, "r_scienceFluxSigma"],

1132 np.nanstd(fluxes, ddof=1))

1133

1134

1135def skew_wrapper(values):

1136 """Compute scipy skew, omitting nans.

1137

1138 This version works with both scipy<1.9 (where it erroneously returns a

1139 masked array) and scipy>=1.9 (where it correctly returns a float).

1140

1141 Parameters

1142 ----------

1143 values : `np.ndarray`

1144

1145 Returns

1146 -------

1147 skew_value : `float`

1148 """

1149 value = skew(values, bias=False, nan_policy="omit")

1150 if isinstance(value, np.ma.masked_array):

1151 return value.data

1152 else:

1153 return value

1154

1155

1156class MemoryTester(lsst.utils.tests.MemoryTestCase):

1157 pass

1158

1159

1160def setup_module(module):

1161 lsst.utils.tests.init()

1162

1163

1164if __name__ == "__main__": 1164 ↛ 1165line 1164 didn't jump to line 1165 because the condition on line 1164 was never true

1165 lsst.utils.tests.init()

1166 unittest.main()

Coverage for tests / test_diaCalculationPlugins.py: 13%

411 statements