Coverage for tests/test_subtractTask.py: 99%
784 statements
« prev ^ index » next coverage.py v7.16.1, created at 2026-09-22 10:45 +0000
1# This file is part of ip_diffim.
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 unittest
24from astropy import units as u
26import lsst.afw.math as afwMath
27import lsst.afw.table as afwTable
28import lsst.geom
29import lsst.meas.algorithms as measAlg
30from lsst.ip.diffim import subtractImages, InsufficientKernelSourcesError
31from lsst.pex.config import FieldValidationError
32from lsst.pipe.base import NoWorkFound
33import lsst.utils.tests
34import numpy as np
35from lsst.ip.diffim.utils import (computeRobustStatistics, computePSFNoiseEquivalentArea,
36 evaluateMeanPsfFwhm, getPsfFwhm)
37from lsst.pex.exceptions import InvalidParameterError
39from utils import makeStats, makeTestImage, CustomCoaddPsf
42class AlardLuptonSubtractTestBase:
43 goodPsfSize = 2.0
44 midPsfSize = 2.4
45 badPsfSize = 2.8
47 def _setup_subtraction(self, fluxField="truth_instFlux", errField="truth_instFluxErr", **kwargs):
48 """Setup and configure the image subtraction PipelineTask.
50 Parameters
51 ----------
52 fluxField : `str`, optional
53 Name of the flux field in the source catalog.
54 errField : `str`, optional
55 Name of the flux error field in the source catalog.
56 **kwargs
57 Any additional config parameters to set.
59 Returns
60 -------
61 `lsst.pipe.base.PipelineTask`
62 The configured Task to use for detection and measurement.
63 """
64 config = self.subtractTask.ConfigClass()
65 config.doSubtractBackground = False
66 config.restrictKernelEdgeSources = False
67 config.sourceSelector.signalToNoise.fluxField = fluxField
68 config.sourceSelector.signalToNoise.errField = errField
69 config.sourceSelector.doUnresolved = True
70 config.sourceSelector.doIsolated = True
71 config.sourceSelector.doRequirePrimary = True
72 config.sourceSelector.doFlags = True
73 config.sourceSelector.doSignalToNoise = True
74 config.sourceSelector.flags.bad = ["base_PsfFlux_flag", ]
75 config.update(**kwargs)
77 return self.subtractTask(config=config)
80class AlardLuptonSubtractTest(AlardLuptonSubtractTestBase, lsst.utils.tests.TestCase):
81 subtractTask = subtractImages.AlardLuptonSubtractTask
83 def test_allowed_config_modes(self):
84 """Verify the allowable modes for convolution.
85 """
86 config = subtractImages.AlardLuptonSubtractTask.ConfigClass()
87 config.mode = 'auto'
88 config.mode = 'convolveScience'
89 config.mode = 'convolveTemplate'
91 with self.assertRaises(FieldValidationError):
92 config.mode = 'aotu'
94 def test_mismatched_template(self):
95 """Test that an error is raised if the template
96 does not fully contain the science image.
97 """
98 xSize = 200
99 ySize = 200
100 science, sources = makeTestImage(psfSize=self.midPsfSize, xSize=xSize + 20, ySize=ySize + 20)
101 template, _ = makeTestImage(psfSize=self.midPsfSize, xSize=xSize, ySize=ySize,
102 doApplyCalibration=True)
103 task = self._setup_subtraction()
104 with self.assertRaises(AssertionError):
105 task.run(template, science, sources)
107 def test_mismatched_filter(self):
108 """Test that an error is raised if the science and template have
109 different bands.
110 """
111 xSize = 200
112 ySize = 200
113 science, sources = makeTestImage(psfSize=self.midPsfSize, xSize=xSize + 20, ySize=ySize + 20,
114 band="g", physicalFilter="g noCamera")
115 template, _ = makeTestImage(psfSize=self.midPsfSize, xSize=xSize, ySize=ySize,
116 doApplyCalibration=True, band="not-g", physicalFilter="not-g noCamera")
117 task = self._setup_subtraction()
118 with self.assertRaises(AssertionError):
119 task.run(template, science, sources)
121 def test_incomplete_template_coverage(self):
122 noiseLevel = 1.
123 border = 20
124 xSize = 400
125 ySize = 400
126 nSources = 80
127 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
128 nSrc=nSources, xSize=xSize, ySize=ySize)
129 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
130 nSrc=nSources, templateBorderSize=border, doApplyCalibration=True,
131 xSize=xSize, ySize=ySize)
133 science_height = science.getBBox().getDimensions().getY()
135 def _run_and_check_coverage(template_coverage,
136 requiredTemplateFraction=0.1,
137 minTemplateFractionForExpectedSuccess=0.2):
138 template_cut = template.clone()
139 template_height = int(science_height*template_coverage + border)
140 template_cut.image.array[:, template_height:] = 0
141 template_cut.mask.array[:, template_height:] = template_cut.mask.getPlaneBitMask('NO_DATA')
142 task = self._setup_subtraction(
143 requiredTemplateFraction=requiredTemplateFraction,
144 minTemplateFractionForExpectedSuccess=minTemplateFractionForExpectedSuccess
145 )
146 if template_coverage < requiredTemplateFraction:
147 doRaise = True
148 elif template_coverage < minTemplateFractionForExpectedSuccess:
149 doRaise = True
150 else:
151 doRaise = False
152 if doRaise:
153 with self.assertRaises(NoWorkFound):
154 task.run(template_cut, science.clone(), sources.copy(deep=True))
155 else:
156 task.run(template_cut, science.clone(), sources.copy(deep=True))
157 _run_and_check_coverage(template_coverage=0.09)
158 _run_and_check_coverage(template_coverage=0.15)
159 _run_and_check_coverage(template_coverage=0.7)
161 def test_clear_template_mask(self):
162 noiseLevel = 1.
163 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
164 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
165 templateBorderSize=20, doApplyCalibration=True)
166 diffimEmptyMaskPlanes = ["DETECTED", "DETECTED_NEGATIVE"]
167 task = self._setup_subtraction(mode="convolveTemplate")
168 # Ensure that each each mask plane is set for some pixels
169 mask = template.mask
170 x0 = 50
171 x1 = 75
172 y0 = 150
173 y1 = 175
174 scienceMaskCheck = {}
175 for maskPlane in mask.getMaskPlaneDict().keys():
176 scienceMaskCheck[maskPlane] = np.sum(science.mask.array & mask.getPlaneBitMask(maskPlane) > 0)
177 mask.array[x0: x1, y0: y1] |= mask.getPlaneBitMask(maskPlane)
178 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) > 0))
180 output = task.run(template, science, sources)
181 # Verify that the template mask has been modified in place
182 for maskPlane in mask.getMaskPlaneDict().keys():
183 if maskPlane in diffimEmptyMaskPlanes:
184 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) == 0))
185 elif maskPlane in task.config.preserveTemplateMask:
186 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) > 0))
187 else:
188 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) == 0))
189 # Mask planes set in the science image should also be set in the difference
190 # Except the "DETECTED" planes should have been cleared
191 diffimMask = output.difference.mask
192 for maskPlane, scienceSum in scienceMaskCheck.items():
193 diffimSum = np.sum(diffimMask.array & mask.getPlaneBitMask(maskPlane) > 0)
194 if maskPlane in diffimEmptyMaskPlanes:
195 self.assertEqual(diffimSum, 0)
196 else:
197 self.assertTrue(diffimSum >= scienceSum)
199 def test_equal_images(self):
200 """Test that running with enough sources produces reasonable output,
201 with the same size psf in the template and science.
202 """
203 noiseLevel = 1.
204 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
205 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
206 templateBorderSize=20, doApplyCalibration=True)
207 task = self._setup_subtraction()
208 output = task.run(template, science, sources)
209 # There shoud be no NaN values in the difference image
210 self.assertTrue(np.all(np.isfinite(output.difference.image.array)))
211 # Mean of difference image should be close to zero.
212 meanError = noiseLevel/np.sqrt(output.difference.image.array.size)
213 # Make sure to include pixels with the DETECTED mask bit set.
214 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA", "DETECTED", "DETECTED_NEGATIVE"))
215 differenceMean = computeRobustStatistics(output.difference.image, output.difference.mask, statsCtrl)
216 self.assertFloatsAlmostEqual(differenceMean, 0, atol=5*meanError)
217 # stddev of difference image should be close to expected value.
218 differenceStd = computeRobustStatistics(output.difference.image, output.difference.mask,
219 makeStats(), statistic=afwMath.STDEV)
220 self.assertFloatsAlmostEqual(differenceStd, np.sqrt(2)*noiseLevel, rtol=0.1)
222 def test_equal_images_missing_mask_planes(self):
223 """Test that running with enough sources produces reasonable output,
224 with the same size psf in the template and science and with missing
225 mask planes.
226 """
227 noiseLevel = 1.
228 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
229 addMaskPlanes=[])
230 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
231 templateBorderSize=20, doApplyCalibration=True, addMaskPlanes=[])
232 task = self._setup_subtraction()
233 output = task.run(template, science, sources)
234 # There shoud be no NaN values in the difference image
235 self.assertTrue(np.all(np.isfinite(output.difference.image.array)))
236 # Mean of difference image should be close to zero.
237 meanError = noiseLevel/np.sqrt(output.difference.image.array.size)
238 # Make sure to include pixels with the DETECTED mask bit set.
239 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA", "DETECTED", "DETECTED_NEGATIVE"))
240 differenceMean = computeRobustStatistics(output.difference.image, output.difference.mask, statsCtrl)
241 self.assertFloatsAlmostEqual(differenceMean, 0, atol=5*meanError)
242 # stddev of difference image should be close to expected value.
243 differenceStd = computeRobustStatistics(output.difference.image, output.difference.mask,
244 makeStats(), statistic=afwMath.STDEV)
245 self.assertFloatsAlmostEqual(differenceStd, np.sqrt(2)*noiseLevel, rtol=0.1)
247 def test_psf_size(self):
248 """Test that the image subtract task runs without failing, if
249 fwhmExposureBuffer and fwhmExposureGrid parameters are set.
250 """
251 noiseLevel = 1.
252 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
253 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
254 templateBorderSize=20, doApplyCalibration=True)
256 schema = afwTable.ExposureTable.makeMinimalSchema()
257 weightKey = schema.addField("weight", type="D", doc="Coadd weight")
258 exposureCatalog = afwTable.ExposureCatalog(schema)
259 kernel = measAlg.DoubleGaussianPsf(7, 7, 2.0).getKernel()
260 psf = measAlg.KernelPsf(kernel, template.getBBox().getCenter())
262 record = exposureCatalog.addNew()
263 record.setPsf(psf)
264 record.setWcs(template.wcs)
265 record.setD(weightKey, 1.0)
266 record.setBBox(template.getBBox())
268 customPsf = CustomCoaddPsf(exposureCatalog, template.wcs)
269 template.setPsf(customPsf)
271 # Test that we get an exception if we simply get the FWHM at center.
272 with self.assertRaises(InvalidParameterError):
273 getPsfFwhm(template.psf, True)
275 with self.assertRaises(InvalidParameterError):
276 getPsfFwhm(template.psf, False)
278 # Test that evaluateMeanPsfFwhm runs successfully on the template.
279 evaluateMeanPsfFwhm(template, fwhmExposureBuffer=0.05, fwhmExposureGrid=10)
281 # Since the PSF is spatially invariant, the FWHM should be the same at
282 # all points in the science image.
283 fwhm1 = getPsfFwhm(science.psf, False)
284 fwhm2 = evaluateMeanPsfFwhm(science, fwhmExposureBuffer=0.05, fwhmExposureGrid=10)
285 self.assertAlmostEqual(fwhm1[0], fwhm2, places=13)
286 self.assertAlmostEqual(fwhm1[1], fwhm2, places=13)
288 self.assertAlmostEqual(evaluateMeanPsfFwhm(science, fwhmExposureBuffer=0.05,
289 fwhmExposureGrid=10),
290 getPsfFwhm(science.psf, True), places=7
291 )
293 # Test that the image subtraction task runs successfully.
294 task = self._setup_subtraction()
296 # Test that the task runs if we take the mean FWHM on a grid.
297 with self.assertLogs(level="INFO") as cm:
298 task.run(template, science, sources)
300 # Check that evaluateMeanPsfFwhm was called.
301 # This tests that getPsfFwhm failed raising InvalidParameterError,
302 # that is caught and handled appropriately.
303 logMessage = ("INFO:lsst.alardLuptonSubtract:Unable to evaluate PSF at the average position. "
304 "Evaluting PSF on a grid of points."
305 )
306 self.assertIn(logMessage, cm.output)
308 def test_auto_convolveTemplate(self):
309 """Test that auto mode gives the same result as convolveTemplate when
310 the template psf is the smaller.
311 """
312 noiseLevel = 1.
313 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
314 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
315 templateBorderSize=20, doApplyCalibration=True)
316 task = self._setup_subtraction(mode="convolveTemplate")
317 output = task.run(template.clone(), science.clone(), sources)
319 task = self._setup_subtraction(mode="auto")
320 outputAuto = task.run(template, science, sources)
321 self.assertMaskedImagesEqual(output.difference.maskedImage, outputAuto.difference.maskedImage)
323 def test_auto_convolveScience(self):
324 """Test that auto mode gives the same result as convolveScience when
325 the science psf is the smaller.
326 """
327 noiseLevel = 1.
328 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
329 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
330 templateBorderSize=20, doApplyCalibration=True)
331 task = self._setup_subtraction(mode="convolveScience")
332 output = task.run(template.clone(), science.clone(), sources)
334 task = self._setup_subtraction(mode="auto")
335 outputAuto = task.run(template, science, sources)
336 self.assertMaskedImagesEqual(output.difference.maskedImage, outputAuto.difference.maskedImage)
338 def test_science_better(self):
339 """Test that running with enough sources produces reasonable output,
340 with the science psf being smaller than the template.
341 """
342 statsCtrl = makeStats()
343 statsCtrlDetect = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
345 def _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel, templateNoiseLevel):
346 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=scienceNoiseLevel,
347 noiseSeed=6)
348 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=templateNoiseLevel, noiseSeed=7,
349 templateBorderSize=20, doApplyCalibration=True)
350 task = self._setup_subtraction(mode="convolveScience")
351 output = task.run(template, science, sources)
352 self.assertFloatsAlmostEqual(task.metadata["scaleScienceVarianceFactor"], 1., atol=.05)
353 # Mean of difference image should be close to zero.
354 nGoodPix = np.sum(np.isfinite(output.difference.image.array))
355 meanError = (scienceNoiseLevel + templateNoiseLevel)/np.sqrt(nGoodPix)
356 diffimMean = computeRobustStatistics(output.difference.image, output.difference.mask,
357 statsCtrlDetect)
359 self.assertFloatsAlmostEqual(diffimMean, 0, atol=5*meanError)
360 # stddev of difference image should be close to expected value.
361 noiseLevel = np.sqrt(scienceNoiseLevel**2 + templateNoiseLevel**2)
362 varianceMean = computeRobustStatistics(output.difference.variance, output.difference.mask,
363 statsCtrl)
364 diffimStd = computeRobustStatistics(output.difference.image, output.difference.mask,
365 statsCtrl, statistic=afwMath.STDEV)
366 self.assertFloatsAlmostEqual(varianceMean, noiseLevel**2, rtol=0.1)
367 self.assertFloatsAlmostEqual(diffimStd, noiseLevel, rtol=0.1)
369 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=1., templateNoiseLevel=1.)
370 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=1., templateNoiseLevel=.1)
371 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=.1, templateNoiseLevel=.1)
373 def test_template_better(self):
374 """Test that running with enough sources produces reasonable output,
375 with the template psf being smaller than the science.
376 """
377 statsCtrl = makeStats()
378 statsCtrlDetect = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
380 def _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel, templateNoiseLevel):
381 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=scienceNoiseLevel,
382 noiseSeed=6)
383 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=templateNoiseLevel, noiseSeed=7,
384 templateBorderSize=20, doApplyCalibration=True)
385 task = self._setup_subtraction()
386 output = task.run(template, science, sources)
387 self.assertFloatsAlmostEqual(task.metadata["scaleScienceVarianceFactor"], 1., atol=.05)
388 # There should be no NaNs in the image if we convolve the template with a buffer
389 self.assertTrue(np.all(np.isfinite(output.difference.image.array)))
390 # Mean of difference image should be close to zero.
391 meanError = (scienceNoiseLevel + templateNoiseLevel)/np.sqrt(output.difference.image.array.size)
393 diffimMean = computeRobustStatistics(output.difference.image, output.difference.mask,
394 statsCtrlDetect)
395 self.assertFloatsAlmostEqual(diffimMean, 0, atol=5*meanError)
396 # stddev of difference image should be close to expected value.
397 noiseLevel = np.sqrt(scienceNoiseLevel**2 + templateNoiseLevel**2)
398 varianceMean = computeRobustStatistics(output.difference.variance, output.difference.mask,
399 statsCtrl)
400 diffimStd = computeRobustStatistics(output.difference.image, output.difference.mask,
401 statsCtrl, statistic=afwMath.STDEV)
402 self.assertFloatsAlmostEqual(varianceMean, noiseLevel**2, rtol=0.1)
403 self.assertFloatsAlmostEqual(diffimStd, noiseLevel, rtol=0.1)
405 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=1., templateNoiseLevel=1.)
406 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=1., templateNoiseLevel=.1)
407 _run_and_check_images(statsCtrl, statsCtrlDetect, scienceNoiseLevel=.1, templateNoiseLevel=.1)
409 def test_symmetry(self):
410 """Test that convolving the science and convolving the template are
411 symmetric: if the psfs are switched between them, the difference image
412 should be nearly the same.
413 """
414 noiseLevel = 1.
415 # Don't include a border for the template, in order to make the results
416 # comparable when we swap which image is treated as the "science" image.
417 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel,
418 noiseSeed=6, templateBorderSize=0, doApplyCalibration=True)
419 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel,
420 noiseSeed=7, templateBorderSize=0, doApplyCalibration=True)
421 task = self._setup_subtraction(mode='auto')
423 # The science image will be modified in place, so use a copy for the second run.
424 science_better = task.run(template.clone(), science.clone(), sources)
425 template_better = task.run(science, template, sources)
427 delta = template_better.difference.clone()
428 delta.image -= science_better.difference.image
429 delta.variance -= science_better.difference.variance
430 delta.mask.array -= science_better.difference.mask.array
432 statsCtrl = makeStats()
433 # Mean of delta should be very close to zero.
434 nGoodPix = np.sum(np.isfinite(delta.image.array))
435 meanError = 2*noiseLevel/np.sqrt(nGoodPix)
436 deltaMean = computeRobustStatistics(delta.image, delta.mask, statsCtrl)
437 deltaStd = computeRobustStatistics(delta.image, delta.mask, statsCtrl, statistic=afwMath.STDEV)
438 self.assertFloatsAlmostEqual(deltaMean, 0, atol=5*meanError)
439 # stddev of difference image should be close to expected value
440 self.assertFloatsAlmostEqual(deltaStd, 2*np.sqrt(2)*noiseLevel, rtol=.1)
442 def test_few_sources(self):
443 """Test with only 1 source, to check that we get a useful error.
444 """
445 xSize = 256
446 ySize = 256
447 science, sources = makeTestImage(psfSize=self.midPsfSize, nSrc=10, xSize=xSize, ySize=ySize)
448 template, _ = makeTestImage(psfSize=self.goodPsfSize, nSrc=10, xSize=xSize, ySize=ySize,
449 doApplyCalibration=True)
450 task = self._setup_subtraction()
451 sources = sources[0:1]
452 with self.assertRaises(InsufficientKernelSourcesError):
453 task.run(template, science, sources)
455 def test_kernel_source_selector(self):
456 """Check that kernel source selection behaves as expected.
457 """
458 xSize = 256
459 ySize = 256
460 nSourcesSimulated = 20
461 sciencePsfSize = self.midPsfSize
462 templatePsfSize = self.goodPsfSize
463 science, sources = makeTestImage(psfSize=sciencePsfSize, nSrc=nSourcesSimulated,
464 xSize=xSize, ySize=ySize)
465 template, _ = makeTestImage(psfSize=templatePsfSize, nSrc=nSourcesSimulated,
466 xSize=xSize, ySize=ySize, doApplyCalibration=True)
468 def _run_and_check_sources(sourcesIn, maxKernelSources=1000, minKernelSources=3):
469 sources = sourcesIn.copy(deep=True)
471 task = self._setup_subtraction(maxKernelSources=maxKernelSources,
472 minKernelSources=minKernelSources,
473 )
474 task.templatePsfSize = templatePsfSize
475 task.sciencePsfSize = sciencePsfSize
476 task.matchedPsfSize = sciencePsfSize
477 # Verify that source flags are not set in the input catalog
478 # Note that this will use the last flag in the list for the rest of
479 # the test.
480 for badSourceFlag in task.sourceSelector.config.flags.bad:
481 self.assertEqual(np.sum(sources[badSourceFlag]), 0)
482 nSources = len(sources)
483 # Flag a third of the sources
484 sources[0:: 3][badSourceFlag] = True
485 rejectRadius = 2*task.config.makeKernel.kernel.active.kernelSize
486 bbox = science.getBBox()
487 bbox.grow(-rejectRadius)
488 edgeSources = ~bbox.contains(sources.getX(), sources.getY())
489 sources[edgeSources][badSourceFlag] = True
490 nBadSources = np.sum(sources[badSourceFlag])
491 if maxKernelSources > 0:
492 nGoodSources = np.minimum(nSources - nBadSources, maxKernelSources)
493 else:
494 nGoodSources = nSources - nBadSources
496 signalToNoise = sources.getPsfInstFlux()/sources.getPsfInstFluxErr()
497 signalToNoise = signalToNoise[~sources[badSourceFlag]]
498 signalToNoise.sort()
499 selectSources = task._sourceSelector(template, science, sources)
500 self.assertEqual(nGoodSources, len(selectSources))
501 signalToNoiseOut = selectSources.getPsfInstFlux()/selectSources.getPsfInstFluxErr()
502 signalToNoiseOut.sort()
503 self.assertFloatsAlmostEqual(signalToNoise[-nGoodSources:], signalToNoiseOut)
505 _run_and_check_sources(sources)
506 _run_and_check_sources(sources, maxKernelSources=len(sources)//3)
507 _run_and_check_sources(sources, maxKernelSources=-1)
508 with self.assertRaises(RuntimeError):
509 _run_and_check_sources(sources, minKernelSources=1000)
511 def test_order_equal_images(self):
512 """Verify that the result is the same regardless of convolution mode
513 if the images are equivalent.
514 """
515 noiseLevel = .1
516 seed1 = 6
517 seed2 = 7
518 for psfSize in [self.midPsfSize, self.goodPsfSize, self.badPsfSize]:
519 science1, sources1 = makeTestImage(psfSize=psfSize, noiseLevel=noiseLevel, noiseSeed=seed1,
520 clearEdgeMask=True)
521 template1, _ = makeTestImage(psfSize=psfSize, noiseLevel=noiseLevel, noiseSeed=seed2,
522 templateBorderSize=0, doApplyCalibration=True,
523 clearEdgeMask=True)
524 task1 = self._setup_subtraction(mode="convolveTemplate")
525 results_convolveTemplate = task1.run(template1, science1, sources1)
527 science2, sources2 = makeTestImage(psfSize=psfSize, noiseLevel=noiseLevel, noiseSeed=seed1,
528 clearEdgeMask=True)
529 template2, _ = makeTestImage(psfSize=psfSize, noiseLevel=noiseLevel, noiseSeed=seed2,
530 templateBorderSize=0, doApplyCalibration=True,
531 clearEdgeMask=True)
532 task2 = self._setup_subtraction(mode="convolveScience")
533 results_convolveScience = task2.run(template2, science2, sources2)
534 bbox = results_convolveTemplate.difference.getBBox().clippedTo(
535 results_convolveScience.difference.getBBox())
536 diff1 = science1.maskedImage.clone()[bbox]
537 diff1 -= template1.maskedImage[bbox]
538 diff2 = science2.maskedImage.clone()[bbox]
539 diff2 -= template2.maskedImage[bbox]
540 self.assertFloatsAlmostEqual(results_convolveTemplate.difference[bbox].image.array,
541 diff1.image.array,
542 atol=noiseLevel*5.)
543 self.assertFloatsAlmostEqual(results_convolveScience.difference[bbox].image.array,
544 diff2.image.array,
545 atol=noiseLevel*5.)
546 diffErr = noiseLevel*2
547 self.assertMaskedImagesAlmostEqual(results_convolveTemplate.difference[bbox].maskedImage,
548 results_convolveScience.difference[bbox].maskedImage,
549 atol=diffErr*5.)
551 def test_background_subtraction(self):
552 """Check that we can recover the background,
553 and that it is subtracted correctly in the difference image.
555 NOTE: Background subtraction is now turned off by default in
556 subtractImages. It is now run in detectAndMeasure instead, but since the
557 code to run background subtraction is not being removed this test should
558 stay to make sure it continues functioning as intended.
559 """
560 noiseLevel = 1.
561 xSize = 512
562 ySize = 512
563 x0 = 123
564 y0 = 456
565 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
566 templateBorderSize=20,
567 xSize=xSize, ySize=ySize, x0=x0, y0=y0,
568 doApplyCalibration=True)
569 params = [2.2, 2.1, 2.0, 1.2, 1.1, 1.0]
571 bbox2D = lsst.geom.Box2D(lsst.geom.Point2D(x0, y0), lsst.geom.Extent2D(xSize, ySize))
572 background_model = afwMath.Chebyshev1Function2D(params, bbox2D)
573 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
574 background=background_model,
575 xSize=xSize, ySize=ySize, x0=x0, y0=y0)
576 # Don't use ``self._setup_subtraction()`` here.
577 # Modifying the config of a subtask is messy.
578 config = subtractImages.AlardLuptonSubtractTask.ConfigClass()
580 config.sourceSelector.signalToNoise.fluxField = "truth_instFlux"
581 config.sourceSelector.signalToNoise.errField = "truth_instFluxErr"
582 config.doSubtractBackground = True
584 config.makeKernel.kernel.name = "AL"
585 config.makeKernel.kernel.active.fitForBackground = True
586 config.makeKernel.kernel.active.spatialKernelOrder = 1
587 config.makeKernel.kernel.active.spatialBgOrder = 2
588 statsCtrl = makeStats()
590 def _run_and_check_images(config, statsCtrl, mode):
591 """Check that the fit background matches the input model.
592 """
593 config.mode = mode
594 task = subtractImages.AlardLuptonSubtractTask(config=config)
595 output = task.run(template.clone(), science.clone(), sources)
597 # We should be fitting the same number of parameters as were in the input model
598 self.assertEqual(output.backgroundModel.getNParameters(), background_model.getNParameters())
600 # The parameters of the background fit should be close to the input model
601 self.assertFloatsAlmostEqual(np.array(output.backgroundModel.getParameters()),
602 np.array(params), rtol=0.3)
604 # stddev of difference image should be close to expected value.
605 # This will fail if we have mis-subtracted the background.
606 stdVal = computeRobustStatistics(output.difference.image, output.difference.mask,
607 statsCtrl, statistic=afwMath.STDEV)
608 self.assertFloatsAlmostEqual(stdVal, np.sqrt(2)*noiseLevel, rtol=0.1)
610 _run_and_check_images(config, statsCtrl, "convolveTemplate")
611 _run_and_check_images(config, statsCtrl, "convolveScience")
613 def test_scale_variance_convolve_template(self):
614 """Check variance scaling of the image difference.
615 """
616 scienceNoiseLevel = 4.
617 templateNoiseLevel = 2.
618 scaleFactor = 1.345
619 # Make sure to include pixels with the DETECTED mask bit set.
620 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
622 def _run_and_check_images(science, template, sources, statsCtrl,
623 doDecorrelation, doScaleVariance, scaleFactor=1.):
624 """Check that the variance plane matches the expected value for
625 different configurations of ``doDecorrelation`` and ``doScaleVariance``.
626 """
628 task = self._setup_subtraction(doDecorrelation=doDecorrelation,
629 doScaleVariance=doScaleVariance,
630 )
631 output = task.run(template.clone(), science.clone(), sources)
632 if doScaleVariance:
633 self.assertFloatsAlmostEqual(task.metadata["scaleScienceVarianceFactor"],
634 scaleFactor, atol=0.05)
636 scienceNoise = computeRobustStatistics(science.variance, science.mask, statsCtrl)
637 if doDecorrelation:
638 templateNoise = computeRobustStatistics(template.variance, template.mask, statsCtrl)
639 else:
640 templateNoise = computeRobustStatistics(output.matchedTemplate.variance,
641 output.matchedTemplate.mask,
642 statsCtrl)
644 if doScaleVariance:
645 # Only the science variance is scaled here. The template
646 # variance is scaled independently in ``GetTemplateTask``.
647 scienceNoise *= scaleFactor
648 varMean = computeRobustStatistics(output.difference.variance, output.difference.mask, statsCtrl)
649 self.assertFloatsAlmostEqual(varMean, scienceNoise + templateNoise, rtol=0.1)
651 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=scienceNoiseLevel, noiseSeed=6)
652 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=templateNoiseLevel, noiseSeed=7,
653 templateBorderSize=20, doApplyCalibration=True)
654 # Verify that the variance plane of the difference image is correct
655 # when the template and science variance planes are correct
656 _run_and_check_images(science, template, sources, statsCtrl,
657 doDecorrelation=True, doScaleVariance=True)
658 _run_and_check_images(science, template, sources, statsCtrl,
659 doDecorrelation=True, doScaleVariance=False)
660 _run_and_check_images(science, template, sources, statsCtrl,
661 doDecorrelation=False, doScaleVariance=True)
662 _run_and_check_images(science, template, sources, statsCtrl,
663 doDecorrelation=False, doScaleVariance=False)
665 # Verify that the variance plane of the difference image is correct
666 # when the input science variance plane is incorrect
667 science.variance.array /= scaleFactor
668 _run_and_check_images(science, template, sources, statsCtrl,
669 doDecorrelation=True, doScaleVariance=True, scaleFactor=scaleFactor)
670 _run_and_check_images(science, template, sources, statsCtrl,
671 doDecorrelation=True, doScaleVariance=False, scaleFactor=scaleFactor)
672 _run_and_check_images(science, template, sources, statsCtrl,
673 doDecorrelation=False, doScaleVariance=True, scaleFactor=scaleFactor)
674 _run_and_check_images(science, template, sources, statsCtrl,
675 doDecorrelation=False, doScaleVariance=False, scaleFactor=scaleFactor)
677 def test_scale_variance_convolve_science(self):
678 """Check variance scaling of the image difference.
679 """
680 scienceNoiseLevel = 4.
681 templateNoiseLevel = 2.
682 scaleFactor = 1.345
683 # Make sure to include pixels with the DETECTED mask bit set.
684 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
686 def _run_and_check_images(science, template, sources, statsCtrl,
687 doDecorrelation, doScaleVariance, scaleFactor=1.):
688 """Check that the variance plane matches the expected value for
689 different configurations of ``doDecorrelation`` and ``doScaleVariance``.
690 """
692 task = self._setup_subtraction(mode="convolveScience",
693 doDecorrelation=doDecorrelation,
694 doScaleVariance=doScaleVariance,
695 restrictKernelEdgeSources=False,
696 )
697 output = task.run(template.clone(), science.clone(), sources)
698 if doScaleVariance:
699 self.assertFloatsAlmostEqual(task.metadata["scaleScienceVarianceFactor"],
700 scaleFactor, atol=0.05)
702 templateNoise = computeRobustStatistics(template.variance, template.mask, statsCtrl)
703 if doDecorrelation:
704 scienceNoise = computeRobustStatistics(science.variance, science.mask, statsCtrl)
705 else:
706 scienceNoise = computeRobustStatistics(output.matchedScience.variance,
707 output.matchedScience.mask,
708 statsCtrl)
710 if doScaleVariance:
711 # Only the science variance is scaled here. The template
712 # variance is scaled independently in ``GetTemplateTask``.
713 scienceNoise *= scaleFactor
715 varMean = computeRobustStatistics(output.difference.variance, output.difference.mask, statsCtrl)
716 self.assertFloatsAlmostEqual(varMean, scienceNoise + templateNoise, rtol=0.1)
718 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=scienceNoiseLevel, noiseSeed=6)
719 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=templateNoiseLevel, noiseSeed=7,
720 templateBorderSize=20, doApplyCalibration=True)
721 # Verify that the variance plane of the difference image is correct
722 # when the template and science variance planes are correct
723 _run_and_check_images(science, template, sources, statsCtrl,
724 doDecorrelation=True, doScaleVariance=True)
725 _run_and_check_images(science, template, sources, statsCtrl,
726 doDecorrelation=True, doScaleVariance=False)
727 _run_and_check_images(science, template, sources, statsCtrl,
728 doDecorrelation=False, doScaleVariance=True)
729 _run_and_check_images(science, template, sources, statsCtrl,
730 doDecorrelation=False, doScaleVariance=False)
732 # Verify that the variance plane of the difference image is correct
733 # when the input science variance plane is incorrect
734 science.variance.array /= scaleFactor
735 _run_and_check_images(science, template, sources, statsCtrl,
736 doDecorrelation=True, doScaleVariance=True, scaleFactor=scaleFactor)
737 _run_and_check_images(science, template, sources, statsCtrl,
738 doDecorrelation=True, doScaleVariance=False, scaleFactor=scaleFactor)
739 _run_and_check_images(science, template, sources, statsCtrl,
740 doDecorrelation=False, doScaleVariance=True, scaleFactor=scaleFactor)
741 _run_and_check_images(science, template, sources, statsCtrl,
742 doDecorrelation=False, doScaleVariance=False, scaleFactor=scaleFactor)
744 def test_exposure_properties_convolve_template(self):
745 """Check that all necessary exposure metadata is included
746 when the template is convolved.
747 """
748 noiseLevel = 1.
749 seed = 37
750 rng = np.random.RandomState(seed)
751 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
752 psf = science.psf
753 psfAvgPos = psf.getAveragePosition()
754 psfSize = getPsfFwhm(science.psf)
755 psfImg = psf.computeKernelImage(psfAvgPos)
756 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
757 templateBorderSize=20, doApplyCalibration=True)
759 # Generate a random aperture correction map
760 apCorrMap = lsst.afw.image.ApCorrMap()
761 for name in ("a", "b", "c"):
762 apCorrMap.set(name, lsst.afw.math.ChebyshevBoundedField(science.getBBox(), rng.randn(3, 3)))
763 science.info.setApCorrMap(apCorrMap)
765 def _run_and_check_images(doDecorrelation):
766 """Check that the metadata is correct with or without decorrelation.
767 """
768 task = self._setup_subtraction(mode="convolveTemplate",
769 doDecorrelation=doDecorrelation,
770 )
771 output = task.run(template.clone(), science.clone(), sources)
772 psfOut = output.difference.psf
773 psfAvgPos = psfOut.getAveragePosition()
774 if doDecorrelation:
775 # Decorrelation requires recalculating the PSF,
776 # so it will not be the same as the input
777 psfOutSize = getPsfFwhm(science.psf)
778 self.assertFloatsAlmostEqual(psfSize, psfOutSize)
779 else:
780 psfOutImg = psfOut.computeKernelImage(psfAvgPos)
781 self.assertImagesAlmostEqual(psfImg, psfOutImg)
783 # check PSF, WCS, bbox, filterLabel, photoCalib, aperture correction
784 self._compare_apCorrMaps(apCorrMap, output.difference.info.getApCorrMap())
785 self.assertWcsAlmostEqualOverBBox(science.wcs, output.difference.wcs, science.getBBox())
786 self.assertEqual(science.filter, output.difference.filter)
787 self.assertEqual(science.photoCalib, output.difference.photoCalib)
788 _run_and_check_images(doDecorrelation=True)
789 _run_and_check_images(doDecorrelation=False)
791 def test_exposure_properties_convolve_science(self):
792 """Check that all necessary exposure metadata is included
793 when the science image is convolved.
794 """
795 noiseLevel = 1.
796 seed = 37
797 rng = np.random.RandomState(seed)
798 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
799 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
800 templateBorderSize=20, doApplyCalibration=True)
801 psf = template.psf
802 psfAvgPos = psf.getAveragePosition()
803 psfSize = getPsfFwhm(template.psf)
804 psfImg = psf.computeKernelImage(psfAvgPos)
806 # Generate a random aperture correction map
807 apCorrMap = lsst.afw.image.ApCorrMap()
808 for name in ("a", "b", "c"):
809 apCorrMap.set(name, lsst.afw.math.ChebyshevBoundedField(science.getBBox(), rng.randn(3, 3)))
810 science.info.setApCorrMap(apCorrMap)
812 def _run_and_check_images(doDecorrelation):
813 """Check that the metadata is correct with or without decorrelation.
814 """
815 task = self._setup_subtraction(mode="convolveScience",
816 doDecorrelation=doDecorrelation,
817 )
818 output = task.run(template.clone(), science.clone(), sources)
819 if doDecorrelation:
820 # Decorrelation requires recalculating the PSF,
821 # so it will not be the same as the input
822 psfOutSize = getPsfFwhm(template.psf)
823 self.assertFloatsAlmostEqual(psfSize, psfOutSize)
824 else:
825 psfOut = output.difference.psf
826 psfAvgPos = psfOut.getAveragePosition()
827 psfOutImg = psfOut.computeKernelImage(psfAvgPos)
828 self.assertImagesAlmostEqual(psfImg, psfOutImg)
830 # check PSF, WCS, bbox, filterLabel, photoCalib, aperture correction
831 self._compare_apCorrMaps(apCorrMap, output.difference.info.getApCorrMap())
832 self.assertWcsAlmostEqualOverBBox(science.wcs, output.difference.wcs, science.getBBox())
833 self.assertEqual(science.filter, output.difference.filter)
834 self.assertEqual(science.photoCalib, output.difference.photoCalib)
836 _run_and_check_images(doDecorrelation=True)
837 _run_and_check_images(doDecorrelation=False)
839 def _compare_apCorrMaps(self, a, b):
840 """Compare two ApCorrMaps for equality, without assuming that their BoundedFields have the
841 same addresses (i.e. so we can compare after serialization).
843 This function is taken from ``ApCorrMapTestCase`` in afw/tests/.
845 Parameters
846 ----------
847 a, b : `lsst.afw.image.ApCorrMap`
848 The two aperture correction maps to compare.
849 """
850 self.assertEqual(len(a), len(b))
851 for name, value in list(a.items()):
852 value2 = b.get(name)
853 self.assertIsNotNone(value2)
854 self.assertEqual(value.getBBox(), value2.getBBox())
855 self.assertFloatsAlmostEqual(
856 value.getCoefficients(), value2.getCoefficients(), rtol=0.0)
858 def test_fake_mask_plane_propagation(self):
859 """Test that we have the mask planes related to fakes in diffim images.
860 This is testing method called updateMasks
861 """
862 xSize = 200
863 ySize = 200
864 science, sources = makeTestImage(psfSize=self.midPsfSize, xSize=xSize, ySize=ySize)
865 science_fake_img, science_fake_sources = makeTestImage(
866 psfSize=self.midPsfSize, xSize=xSize, ySize=ySize, seed=7, nSrc=2, noiseLevel=0.25, fluxRange=1
867 )
868 template, _ = makeTestImage(psfSize=self.midPsfSize, xSize=xSize, ySize=ySize,
869 doApplyCalibration=True)
870 tmplt_fake_img, tmplt_fake_sources = makeTestImage(
871 psfSize=self.midPsfSize, xSize=xSize, ySize=ySize, seed=9, nSrc=2, noiseLevel=0.25, fluxRange=1
872 )
873 # created fakes and added them to the images
874 science.image.array += science_fake_img.image.array
875 template.image.array += tmplt_fake_img.image.array
877 # TODO: DM-40796 update to INJECTED names when source injection gets refactored
878 # adding mask planes to both science and template images
879 science_mask_planes = science.mask.addMaskPlane("FAKE")
880 template_mask_planes = template.mask.addMaskPlane("FAKE")
882 for a_science_source in science_fake_sources:
883 # 3 x 3 masking of the source locations is fine
884 bbox = lsst.geom.Box2I(
885 lsst.geom.Point2I(a_science_source.getX(), a_science_source.getY()), lsst.geom.Extent2I(3, 3)
886 )
887 science[bbox].mask.array |= science_mask_planes
889 for a_template_source in tmplt_fake_sources:
890 # 3 x 3 masking of the source locations is fine
891 bbox = lsst.geom.Box2I(
892 lsst.geom.Point2I(a_template_source.getX(), a_template_source.getY()),
893 lsst.geom.Extent2I(3, 3)
894 )
895 template[bbox].mask.array |= template_mask_planes
897 science_fake_masked = (science.mask.array & science.mask.getPlaneBitMask("FAKE")) > 0
898 template_fake_masked = (template.mask.array & template.mask.getPlaneBitMask("FAKE")) > 0
900 task = self._setup_subtraction()
901 subtraction = task.run(template, science, sources)
903 # check subtraction mask plane is set where we set the previous masks
904 diff_mask = subtraction.difference.mask
906 # science mask should be now in INJECTED
907 inj_masked = (diff_mask.array & diff_mask.getPlaneBitMask("INJECTED")) > 0
909 # template mask should be now in INJECTED_TEMPLATE
910 injTmplt_masked = (diff_mask.array & diff_mask.getPlaneBitMask("INJECTED_TEMPLATE")) > 0
912 self.assertEqual(np.sum(inj_masked.astype(int)-science_fake_masked.astype(int)), 0)
913 self.assertEqual(np.sum(injTmplt_masked.astype(int)-template_fake_masked.astype(int)), 0)
915 def test_metadata_metrics(self):
916 """Verify fields are added to metadata when subtraction is run, and
917 that the difference image limiting magnitude is calculated correctly,
918 both with a "good" and "bad" seeing template.
919 """
920 science, sources = makeTestImage(psfSize=self.badPsfSize, noiseLevel=1)
921 template_good, _ = makeTestImage(psfSize=self.midPsfSize, doApplyCalibration=True, noiseLevel=0.25,
922 templateBorderSize=20)
923 template_bad, _ = makeTestImage(psfSize=9.5, doApplyCalibration=True, noiseLevel=0.25,
924 templateBorderSize=20)
926 # Add a few sky objects; sky footprints are needed for some metrics.
927 config = measAlg.SkyObjectsTask.ConfigClass()
928 config.nSources = 3
929 skyTask = measAlg.SkyObjectsTask(config=config, name="skySources")
930 skyTask.skySourceKey = sources.schema["sky_source"].asKey()
931 skyTask.run(science.mask, 10, catalog=sources)
932 sources = sources.copy(deep=True)
933 # Add centroids, since these sources were added post-measurement.
934 for record in sources[sources["sky_source"]]:
935 record["truth_x"] = record.getFootprint().getPeaks()[0].getFx()
936 record["truth_y"] = record.getFootprint().getPeaks()[0].getFy()
938 # The metadata fields are attached to the subtractTask, so we do
939 # need to run that; run it for both "good" and "bad" seeing templates
941 subtractTask_good = self._setup_subtraction()
942 _ = subtractTask_good.run(template_good.clone(), science.clone(), sources)
943 subtractTask_bad = self._setup_subtraction()
944 _ = subtractTask_bad.run(template_bad.clone(), science.clone(), sources)
946 # Test that the diffim limiting magnitudes are computed correctly
947 maglim_science = subtractTask_good._calculateMagLim(science)
948 fluxlim_science = (maglim_science*u.ABmag).to_value(u.nJy)
949 maglim_template_good = subtractTask_good._calculateMagLim(template_good)
950 fluxlim_template_good = (maglim_template_good*u.ABmag).to_value(u.nJy)
951 maglim_template_bad = subtractTask_bad._calculateMagLim(template_bad)
952 fluxlim_template_bad = (maglim_template_bad*u.ABmag).to_value(u.nJy)
954 maglim_good = (np.sqrt(fluxlim_science**2 + fluxlim_template_good**2)*u.nJy).to(u.ABmag).value
955 maglim_bad = (np.sqrt(fluxlim_science**2 + fluxlim_template_bad**2)*u.nJy).to(u.ABmag).value
957 self.assertFloatsAlmostEqual(subtractTask_good.metadata['diffimLimitingMagnitude'],
958 maglim_good, atol=1e-6)
959 self.assertFloatsAlmostEqual(subtractTask_bad.metadata['diffimLimitingMagnitude'],
960 maglim_bad, atol=1e-6)
962 # Create a template with a PSF that is not defined at the image center.
963 # First, make an exposure catalog so we can force the template to have
964 # a bad (off-image) PSF. It must have a record with a weight field
965 # and a BBox in order to let us set the PSF manually.
966 template_offimage, _ = makeTestImage()
967 schema = afwTable.ExposureTable.makeMinimalSchema()
968 weightKey = schema.addField("weight", type="D", doc="Coadd weight")
969 exposureCatalog = afwTable.ExposureCatalog(schema)
970 record = exposureCatalog.addNew()
971 record.setD(weightKey, 1.0)
972 record.setBBox(template_offimage.getBBox())
973 kernel = measAlg.DoubleGaussianPsf(7, 7, 2.0).getKernel()
974 psf = measAlg.KernelPsf(kernel, template_offimage.getBBox().getCenter())
975 record.setPsf(psf)
976 record.setWcs(template_offimage.wcs)
977 custom_offimage_psf = CustomCoaddPsf(exposureCatalog, template_offimage.wcs)
978 template_offimage.setPsf(custom_offimage_psf)
980 # Test that template PSF size edge cases are handled correctly.
981 subtractTask_offimage = self._setup_subtraction()
982 _ = subtractTask_offimage.run(template_offimage.clone(), science.clone(), sources)
983 # Test that providing no fallbackPsfSize results in a nan template
984 # limiting magnitude.
985 maglim_template_offimage = subtractTask_offimage._calculateMagLim(template_offimage)
986 self.assertTrue(np.isnan(maglim_template_offimage))
987 # Test that given the provided fallbackPsfSize, the diffim limiting
988 # magnitude is calculated correctly.
989 maglim_template_offimage = 28.182284789714952
990 fluxlim_template_offimage = (maglim_template_offimage*u.ABmag).to_value(u.nJy)
991 maglim_offimage = (np.sqrt(fluxlim_science**2 + fluxlim_template_offimage**2)*u.nJy).to(u.ABmag).value
992 self.assertEqual(subtractTask_offimage.metadata['diffimLimitingMagnitude'], maglim_offimage)
994 # Test that several other expected metadata metrics exist
995 self.assertIn('scienceLimitingMagnitude', subtractTask_good.metadata)
996 self.assertIn('templateLimitingMagnitude', subtractTask_good.metadata)
998 # The mean ratio metric should be much worse on the "bad" subtraction.
999 self.assertLess(subtractTask_good.metadata['differenceFootprintRatioMean'], 0.02)
1000 self.assertGreater(subtractTask_bad.metadata['differenceFootprintRatioMean'], 0.12)
1003class AlardLuptonPreconvolveSubtractTest(AlardLuptonSubtractTestBase, lsst.utils.tests.TestCase):
1004 subtractTask = subtractImages.AlardLuptonPreconvolveSubtractTask
1006 def test_mismatched_template(self):
1007 """Test that an error is raised if the template
1008 does not fully contain the science image.
1009 """
1010 xSize = 200
1011 ySize = 200
1012 science, sources = makeTestImage(psfSize=self.midPsfSize, xSize=xSize + 20, ySize=ySize + 20)
1013 template, _ = makeTestImage(psfSize=self.midPsfSize, xSize=xSize, ySize=ySize,
1014 doApplyCalibration=True)
1015 task = self._setup_subtraction()
1016 with self.assertRaises(AssertionError):
1017 task.run(template, science, sources)
1019 def test_equal_images(self):
1020 """Test that running with enough sources produces reasonable output,
1021 with the same size psf in the template and science.
1022 """
1023 noiseLevel = 1.
1024 xSize = 400
1025 ySize = 400
1026 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
1027 xSize=xSize, ySize=ySize)
1028 template, _ = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1029 templateBorderSize=20, doApplyCalibration=True,
1030 xSize=xSize, ySize=ySize)
1031 task = self._setup_subtraction()
1032 output = task.run(template, science, sources)
1033 # There shoud be no NaN values in the Score image
1034 self.assertTrue(np.all(np.isfinite(output.scoreExposure.image.array)))
1035 # Mean of Score image should be close to zero.
1036 meanError = noiseLevel/np.sqrt(output.scoreExposure.image.array.size)
1037 # Make sure to include pixels with the DETECTED mask bit set.
1038 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
1039 scoreMean = computeRobustStatistics(output.scoreExposure.image,
1040 output.scoreExposure.mask,
1041 statsCtrl)
1042 self.assertFloatsAlmostEqual(scoreMean, 0, atol=5*meanError)
1043 nea = computePSFNoiseEquivalentArea(science.psf)
1044 # stddev of Score image should be close to expected value.
1045 scoreStd = computeRobustStatistics(output.scoreExposure.image, output.scoreExposure.mask,
1046 statsCtrl=statsCtrl, statistic=afwMath.STDEV)
1047 self.assertFloatsAlmostEqual(scoreStd, np.sqrt(2)*noiseLevel/np.sqrt(nea), rtol=0.1)
1049 def test_incomplete_template_coverage(self):
1050 noiseLevel = 1.
1051 border = 20
1052 xSize = 400
1053 ySize = 400
1054 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
1055 xSize=xSize, ySize=ySize)
1056 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1057 templateBorderSize=border, doApplyCalibration=True,
1058 xSize=xSize, ySize=ySize)
1060 science_height = science.getBBox().getDimensions().getY()
1062 def _run_and_check_coverage(template_coverage,
1063 requiredTemplateFraction=0.1,
1064 minTemplateFractionForExpectedSuccess=0.2):
1065 template_cut = template.clone()
1066 template_height = int(science_height*template_coverage + border)
1067 template_cut.image.array[:, template_height:] = 0
1068 template_cut.mask.array[:, template_height:] = template_cut.mask.getPlaneBitMask('NO_DATA')
1069 task = self._setup_subtraction(
1070 requiredTemplateFraction=requiredTemplateFraction,
1071 minTemplateFractionForExpectedSuccess=minTemplateFractionForExpectedSuccess
1072 )
1073 if template_coverage < task.config.requiredTemplateFraction:
1074 doRaise = True
1075 elif template_coverage < task.config.minTemplateFractionForExpectedSuccess:
1076 doRaise = True
1077 else:
1078 doRaise = False
1079 if doRaise:
1080 with self.assertRaises(NoWorkFound):
1081 task.run(template_cut, science.clone(), sources.copy(deep=True))
1082 else:
1083 task.run(template_cut, science.clone(), sources.copy(deep=True))
1084 _run_and_check_coverage(template_coverage=0.09)
1085 _run_and_check_coverage(template_coverage=0.15)
1086 _run_and_check_coverage(template_coverage=.7)
1088 def test_clear_template_mask(self):
1089 noiseLevel = 1.
1090 xSize = 400
1091 ySize = 400
1092 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
1093 xSize=xSize, ySize=ySize)
1094 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1095 templateBorderSize=20, doApplyCalibration=True,
1096 xSize=xSize, ySize=ySize)
1097 diffimEmptyMaskPlanes = ["DETECTED", "DETECTED_NEGATIVE"]
1098 task = self._setup_subtraction()
1099 # Ensure that each each mask plane is set for some pixels
1100 mask = template.mask
1101 x0 = 50
1102 x1 = 75
1103 y0 = 150
1104 y1 = 175
1105 scienceMaskCheck = {}
1106 for maskPlane in mask.getMaskPlaneDict().keys():
1107 scienceMaskCheck[maskPlane] = np.sum(science.mask.array & mask.getPlaneBitMask(maskPlane) > 0)
1108 mask.array[x0: x1, y0: y1] |= mask.getPlaneBitMask(maskPlane)
1109 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) > 0))
1111 output = task.run(template, science, sources)
1112 # Verify that the template mask has been modified in place
1113 for maskPlane in mask.getMaskPlaneDict().keys():
1114 if maskPlane in diffimEmptyMaskPlanes:
1115 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) == 0))
1116 elif maskPlane in task.config.preserveTemplateMask:
1117 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) > 0))
1118 else:
1119 self.assertTrue(np.sum(mask.array & mask.getPlaneBitMask(maskPlane) == 0))
1120 # Mask planes set in the science image should also be set in the difference
1121 # Except the "DETECTED" planes should have been cleared
1122 diffimMask = output.scoreExposure.mask
1123 for maskPlane, scienceSum in scienceMaskCheck.items():
1124 diffimSum = np.sum(diffimMask.array & mask.getPlaneBitMask(maskPlane) > 0)
1125 if maskPlane in diffimEmptyMaskPlanes:
1126 self.assertEqual(diffimSum, 0)
1127 else:
1128 self.assertTrue(diffimSum >= scienceSum)
1130 def test_agnostic_template_psf(self):
1131 """Test that the Score image is the same whether the template PSF is
1132 larger or smaller than the science image PSF.
1133 """
1134 noiseLevel = .3
1135 xSize = 400
1136 ySize = 400
1137 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel,
1138 noiseSeed=6, templateBorderSize=0,
1139 xSize=xSize, ySize=ySize)
1140 template1, _ = makeTestImage(psfSize=self.badPsfSize, noiseLevel=noiseLevel,
1141 noiseSeed=7, doApplyCalibration=True,
1142 xSize=xSize, ySize=ySize)
1143 template2, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel,
1144 noiseSeed=8, doApplyCalibration=True,
1145 xSize=xSize, ySize=ySize)
1146 task = self._setup_subtraction()
1148 science_better = task.run(template1, science.clone(), sources)
1149 template_better = task.run(template2, science, sources)
1150 bbox = science_better.scoreExposure.getBBox().clippedTo(template_better.scoreExposure.getBBox())
1152 delta = template_better.scoreExposure[bbox].clone()
1153 delta.image -= science_better.scoreExposure[bbox].image
1154 delta.variance -= science_better.scoreExposure[bbox].variance
1155 delta.mask.array &= science_better.scoreExposure[bbox].mask.array
1157 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
1158 # Mean of delta should be very close to zero.
1159 nGoodPix = np.sum(np.isfinite(delta.image.array))
1160 meanError = 2*noiseLevel/np.sqrt(nGoodPix)
1161 deltaMean = computeRobustStatistics(delta.image, delta.mask, statsCtrl)
1162 deltaStd = computeRobustStatistics(delta.image, delta.mask, statsCtrl,
1163 statistic=afwMath.STDEV)
1164 self.assertFloatsAlmostEqual(deltaMean, 0, atol=5*meanError)
1165 nea = computePSFNoiseEquivalentArea(science.psf)
1166 # stddev of Score image should be close to expected value
1167 self.assertFloatsAlmostEqual(deltaStd, np.sqrt(2)*noiseLevel/np.sqrt(nea), rtol=.1)
1169 def test_few_sources(self):
1170 """Test with only 1 source, to check that we get a useful error.
1171 """
1172 xSize = 256
1173 ySize = 256
1174 science, sources = makeTestImage(psfSize=self.midPsfSize, nSrc=10, xSize=xSize, ySize=ySize)
1175 template, _ = makeTestImage(psfSize=self.goodPsfSize, nSrc=10, xSize=xSize, ySize=ySize,
1176 doApplyCalibration=True)
1177 task = self._setup_subtraction()
1178 sources = sources[0:1]
1179 with self.assertRaises(InsufficientKernelSourcesError):
1180 task.run(template, science, sources)
1182 def test_background_subtraction(self):
1183 """Check that we can recover the background,
1184 and that it is subtracted correctly in the Score image.
1185 """
1186 noiseLevel = 1.
1187 xSize = 512
1188 ySize = 512
1189 x0 = 123
1190 y0 = 456
1191 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1192 templateBorderSize=20,
1193 xSize=xSize, ySize=ySize, x0=x0, y0=y0,
1194 doApplyCalibration=True)
1195 params = [2.2, 2.1, 2.0, 1.2, 1.1, 1.0]
1197 bbox2D = lsst.geom.Box2D(lsst.geom.Point2D(x0, y0), lsst.geom.Extent2D(xSize, ySize))
1198 background_model = afwMath.Chebyshev1Function2D(params, bbox2D)
1199 science, sources = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
1200 background=background_model,
1201 xSize=xSize, ySize=ySize, x0=x0, y0=y0)
1202 # Don't use ``self._setup_subtraction()`` here.
1203 # Modifying the config of a subtask is messy.
1204 config = subtractImages.AlardLuptonPreconvolveSubtractTask.ConfigClass()
1206 config.sourceSelector.signalToNoise.fluxField = "truth_instFlux"
1207 config.sourceSelector.signalToNoise.errField = "truth_instFluxErr"
1208 config.doSubtractBackground = True
1210 config.makeKernel.kernel.name = "AL"
1211 config.makeKernel.kernel.active.fitForBackground = True
1212 config.makeKernel.kernel.active.spatialKernelOrder = 1
1213 config.makeKernel.kernel.active.spatialBgOrder = 2
1214 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
1216 task = subtractImages.AlardLuptonPreconvolveSubtractTask(config=config)
1217 output = task.run(template.clone(), science.clone(), sources)
1219 # We should be fitting the same number of parameters as were in the input model
1220 self.assertEqual(output.backgroundModel.getNParameters(), background_model.getNParameters())
1222 # The parameters of the background fit should be close to the input model
1223 self.assertFloatsAlmostEqual(np.array(output.backgroundModel.getParameters()),
1224 np.array(params), rtol=0.2)
1226 # stddev of Score image should be close to expected value.
1227 # This will fail if we have mis-subtracted the background.
1228 stdVal = computeRobustStatistics(output.scoreExposure.image, output.scoreExposure.mask,
1229 statsCtrl, statistic=afwMath.STDEV)
1230 # get the img psf Noise Equivalent Area value
1231 nea = computePSFNoiseEquivalentArea(science.psf)
1232 self.assertFloatsAlmostEqual(stdVal, np.sqrt(2)*noiseLevel/np.sqrt(nea), rtol=0.12)
1234 def test_scale_variance(self):
1235 """Check variance scaling of the Score image.
1236 """
1237 scienceNoiseLevel = 4.
1238 templateNoiseLevel = 2.
1239 scaleFactor = 1.345
1240 xSize = 400
1241 ySize = 400
1242 # Make sure to include pixels with the DETECTED mask bit set.
1243 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA"))
1245 def _run_and_check_images(science, template, sources, statsCtrl,
1246 doDecorrelation, doScaleVariance, scaleFactor=1.):
1247 """Check that the variance plane matches the expected value for
1248 different configurations of ``doDecorrelation`` and ``doScaleVariance``.
1249 """
1251 task = self._setup_subtraction(doDecorrelation=doDecorrelation,
1252 doScaleVariance=doScaleVariance,
1253 )
1254 output = task.run(template.clone(), science.clone(), sources)
1255 if doScaleVariance:
1256 self.assertFloatsAlmostEqual(task.metadata["scaleScienceVarianceFactor"],
1257 scaleFactor, atol=0.05)
1259 scienceNoise = computeRobustStatistics(science.variance, science.mask, statsCtrl)
1260 # get the img psf Noise Equivalent Area value
1261 nea = computePSFNoiseEquivalentArea(science.psf)
1262 scienceNoise /= nea
1263 if doDecorrelation:
1264 templateNoise = computeRobustStatistics(template.variance, template.mask, statsCtrl)
1265 templateNoise /= nea
1266 else:
1267 # Don't divide by NEA in this case, since the template is convolved
1268 # and in the same units as the Score exposure.
1269 templateNoise = computeRobustStatistics(output.matchedTemplate.variance,
1270 output.matchedTemplate.mask,
1271 statsCtrl)
1272 if doScaleVariance:
1273 # Only the science variance is scaled here. The template
1274 # variance is scaled independently in ``GetTemplateTask``.
1275 scienceNoise *= scaleFactor
1276 varMean = computeRobustStatistics(output.scoreExposure.variance,
1277 output.scoreExposure.mask,
1278 statsCtrl)
1279 self.assertFloatsAlmostEqual(varMean, scienceNoise + templateNoise, rtol=0.1)
1281 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=scienceNoiseLevel, noiseSeed=6,
1282 xSize=xSize, ySize=ySize)
1283 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=templateNoiseLevel, noiseSeed=7,
1284 templateBorderSize=20, doApplyCalibration=True,
1285 xSize=xSize, ySize=ySize)
1286 # Verify that the variance plane of the Score image is correct
1287 # when the template and science variance planes are correct
1288 _run_and_check_images(science, template, sources, statsCtrl,
1289 doDecorrelation=True, doScaleVariance=True)
1290 _run_and_check_images(science, template, sources, statsCtrl,
1291 doDecorrelation=True, doScaleVariance=False)
1292 _run_and_check_images(science, template, sources, statsCtrl,
1293 doDecorrelation=False, doScaleVariance=True)
1294 _run_and_check_images(science, template, sources, statsCtrl,
1295 doDecorrelation=False, doScaleVariance=False)
1297 # Verify that the variance plane of the Score image is correct
1298 # when the input science variance plane is incorrect
1299 science.variance.array /= scaleFactor
1300 _run_and_check_images(science, template, sources, statsCtrl,
1301 doDecorrelation=True, doScaleVariance=True, scaleFactor=scaleFactor)
1302 _run_and_check_images(science, template, sources, statsCtrl,
1303 doDecorrelation=True, doScaleVariance=False, scaleFactor=scaleFactor)
1304 _run_and_check_images(science, template, sources, statsCtrl,
1305 doDecorrelation=False, doScaleVariance=True, scaleFactor=scaleFactor)
1306 _run_and_check_images(science, template, sources, statsCtrl,
1307 doDecorrelation=False, doScaleVariance=False, scaleFactor=scaleFactor)
1309 def test_exposure_properties(self):
1310 """Check that all necessary exposure metadata is included
1311 with the Score image.
1312 """
1313 noiseLevel = 1.
1314 xSize = 400
1315 ySize = 400
1316 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6,
1317 xSize=xSize, ySize=ySize)
1318 psf = science.psf
1319 psfAvgPos = psf.getAveragePosition()
1320 psfSize = getPsfFwhm(science.psf)
1321 psfImg = psf.computeKernelImage(psfAvgPos)
1322 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1323 templateBorderSize=20, doApplyCalibration=True,
1324 xSize=xSize, ySize=ySize)
1326 def _run_and_check_images(doDecorrelation):
1327 """Check that the metadata is correct with or without decorrelation.
1328 """
1329 task = self._setup_subtraction(doDecorrelation=doDecorrelation)
1330 output = task.run(template.clone(), science.clone(), sources)
1331 psfOut = output.scoreExposure.psf
1332 psfAvgPos = psfOut.getAveragePosition()
1333 if doDecorrelation:
1334 # Decorrelation requires recalculating the PSF,
1335 # so it will not be the same as the input
1336 psfOutSize = getPsfFwhm(science.psf)
1337 self.assertFloatsAlmostEqual(psfSize, psfOutSize)
1338 else:
1339 psfOutImg = psfOut.computeKernelImage(psfAvgPos)
1340 self.assertImagesAlmostEqual(psfImg, psfOutImg)
1342 # check PSF, WCS, bbox, filterLabel, photoCalib
1343 self.assertWcsAlmostEqualOverBBox(science.wcs, output.scoreExposure.wcs, science.getBBox())
1344 self.assertEqual(science.filter, output.scoreExposure.filter)
1345 self.assertEqual(science.photoCalib, output.scoreExposure.photoCalib)
1346 _run_and_check_images(doDecorrelation=True)
1347 _run_and_check_images(doDecorrelation=False)
1350class SimplifiedSubtractTest(AlardLuptonSubtractTestBase, lsst.utils.tests.TestCase):
1351 subtractTask = subtractImages.SimplifiedSubtractTask
1353 def test_runSimplifiedTaskWithExistingKernel(self):
1354 """Test that the simplified task produces the same output as
1355 `AlardLuptonSubtractTask` if it uses the AL kernel.
1356 """
1357 noiseLevel = 1.
1358 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
1359 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1360 templateBorderSize=20, doApplyCalibration=True)
1361 alTask = AlardLuptonSubtractTest._setup_subtraction(AlardLuptonSubtractTest())
1362 task = self._setup_subtraction(useExistingKernel=True)
1364 alResults = alTask.run(template.clone(), science.clone(), sources)
1365 results = task.run(template.clone(), science.clone(),
1366 inputPsfMatchingKernel=alResults.psfMatchingKernel)
1368 self.assertMaskedImagesEqual(alResults.difference, results.difference)
1370 def test_runSimplifiedTaskWithSourceDetection(self):
1371 """Test that the simplified task with source detection produces
1372 reasonable output.
1373 """
1374 noiseLevel = 1.
1375 science, sources = makeTestImage(psfSize=self.midPsfSize, noiseLevel=noiseLevel, noiseSeed=6)
1376 template, _ = makeTestImage(psfSize=self.goodPsfSize, noiseLevel=noiseLevel, noiseSeed=7,
1377 templateBorderSize=20, doApplyCalibration=True)
1378 task = self._setup_subtraction(useExistingKernel=False,
1379 fluxField="base_PsfFlux_instFlux",
1380 errField="base_PsfFlux_instFluxErr",
1381 )
1383 output = task.run(template, science)
1385 # There shoud be no NaN values in the difference image
1386 self.assertTrue(np.all(np.isfinite(output.difference.image.array)))
1387 # Mean of difference image should be close to zero.
1388 meanError = noiseLevel/np.sqrt(output.difference.image.array.size)
1389 # Make sure to include pixels with the DETECTED mask bit set.
1390 statsCtrl = makeStats(badMaskPlanes=("EDGE", "BAD", "NO_DATA", "DETECTED", "DETECTED_NEGATIVE"))
1391 differenceMean = computeRobustStatistics(output.difference.image, output.difference.mask, statsCtrl)
1392 self.assertFloatsAlmostEqual(differenceMean, 0, atol=5*meanError)
1393 # stddev of difference image should be close to expected value.
1394 differenceStd = computeRobustStatistics(output.difference.image, output.difference.mask,
1395 makeStats(), statistic=afwMath.STDEV)
1396 self.assertFloatsAlmostEqual(differenceStd, np.sqrt(2)*noiseLevel, rtol=0.1)
1399def setup_module(module):
1400 lsst.utils.tests.init()
1403class MemoryTestCase(lsst.utils.tests.MemoryTestCase):
1404 pass
1407if __name__ == "__main__": 1407 ↛ 1408line 1407 didn't jump to line 1408 because the condition on line 1407 was never true
1408 lsst.utils.tests.init()
1409 unittest.main()