tickets/DM-21531: Update to latest version of scarlet

python/lsst/meas/extensions/scarlet/deblend.py

@@ -42,12 +42,10 @@ def _estimateRMS(exposure, statsMask):
 
 
 def _getTargetPsf(shape, sigma=1/np.sqrt(2)):
-    X = np.arange(shape[2])
-    Y = np.arange(shape[1])
-    X, Y = np.meshgrid(X, Y)
-    coords = np.stack([Y, X])
+    x = np.arange(shape[2])
+    y = np.arange(shape[1])
     y0, x0 = (shape[1]-1) // 2, (shape[2]-1) // 2
-    target_psf = gaussian(coords, y0, x0, 1, sigma)
+    target_psf = gaussian(y, x, y0, x0, 1, sigma).astype(np.float32)
     target_psf /= target_psf.sum()
     return target_psf
 
@@ -74,7 +72,7 @@ def deblend(mExposure, footprint, log, config):
     mask = (mExposure.mask[:, bbox].array & badPixels) | fpMask[None, :]
     weights[mask > 0] = 0
 
-    psfs = mExposure.computePsfImage(footprint.getCentroid()).array
+    psfs = mExposure.computePsfImage(footprint.getCentroid()).array.astype(np.float32)
     target_psf = _getTargetPsf(psfs.shape)
 
     frame = LsstFrame(images.shape, psfs=target_psf[None])
@@ -279,7 +277,7 @@ def _addSchemaKeys(self, schema):
                                                  doc="Deblender skipped this source")
         self.modelCenter = afwTable.Point2DKey.addFields(schema, name="deblend_peak_center",
                                                          doc="Center used to apply constraints in scarlet",
-                                                         unit="Pixel")
+                                                         unit="pixel")
         self.modelCenterFlux = schema.addField('deblend_peak_instFlux', type=float, units='count',
                                                doc="The instFlux at the peak position of deblended mode")
         # self.log.trace('Added keys to schema: %s', ", ".join(str(x) for x in

tests/test_deblend.py

@@ -1,5 +1,3 @@
-
-
 # This file is part of meas_extensions_scarlet.
 #
 # Developed for the LSST Data Management System.
@@ -52,12 +50,11 @@ def test_deblend_task(self):
 
         # Add some noise, otherwise the task will blow up due to
         # zero variance
-        noise = 10*(np.random.rand(*images.shape)-.5)
+        noise = 10*(np.random.rand(*images.shape).astype(np.float32)-.5)
         images += noise
 
         filters = "grizy"
-        _images = afwImage.MultibandMaskedImage.fromArrays(filters, images.astype(np.float32), None,
-                                                           noise.astype(np.float32))
+        _images = afwImage.MultibandMaskedImage.fromArrays(filters, images.astype(np.float32), None, noise)
         coadds = [afwImage.Exposure(img, dtype=img.image.array.dtype) for img in _images]
         coadds = afwImage.MultibandExposure.fromExposures(filters, coadds)
         for b, coadd in enumerate(coadds):
@@ -89,12 +86,12 @@ def test_deblend_task(self):
 
         seds = np.array([heavy.getImage(fill=0).image.array.sum(axis=(1, 2)) for heavy in heavies])
         true_seds = np.array([
-            [[1665.726318359375, 1745.5401611328125, 1525.91796875, 997.3868408203125, 0.0],
-             [767.100341796875, 1057.0374755859375, 1312.89111328125, 1694.7535400390625, 2069.294921875],
-             [8.08012580871582, 879.344970703125, 2246.90087890625, 4212.82470703125, 6987.0849609375]]
+            [1665.726318359375, 1745.5401611328125, 1525.91796875, 997.3868408203125, 0.0],
+            [767.100341796875, 1057.0374755859375, 1312.89111328125, 1694.7535400390625, 2069.294921875],
+            [8.08012580871582, 879.344970703125, 2246.90087890625, 4212.82470703125, 6987.0849609375]
         ])
 
-        self.assertFloatsAlmostEqual(true_seds, seds, rtol=1e-8, atol=1e-8)
+        self.assertFloatsAlmostEqual(true_seds, seds, rtol=1e-4, atol=1e-4)
 
         bbox = parent.getFootprint().getBBox()
         data = coadds[:, bbox]
@@ -106,8 +103,8 @@ def test_deblend_task(self):
             model[:, heavy.getBBox()].array += heavy.getImage(fill=0).image.array
 
         residual = data.image.array - model.array
-        self.assertFloatsAlmostEqual(np.abs(residual).sum(), 11601.3867187500)
-        self.assertFloatsAlmostEqual(np.max(np.abs(residual)), 56.1048278809, rtol=1e-8, atol=1e-8)
+        self.assertFloatsAlmostEqual(np.abs(residual).sum(), 11601.3867187500, rtol=1e-5, atol=1e-5)
+        self.assertFloatsAlmostEqual(np.max(np.abs(residual)), 56.1048278809, rtol=1e-5, atol=1e-5)
 
 
 class MemoryTester(lsst.utils.tests.MemoryTestCase):

tests/test_scarlet.py

@@ -72,7 +72,7 @@ def test_init(self):
         self.assertEqual(src.delay_thresh, 10)
 
         self.assertFloatsAlmostEqual(src.sed/3, trueSed)
-        self.assertFloatsAlmostEqual(src.morph*3, trueMorph)
+        self.assertFloatsAlmostEqual(src.morph*3, trueMorph, rtol=1e-7)
         self.assertEqual(src.detectedPeak, peak)
         self.assertEqual(foot.getBBox(), bbox)
 
@@ -142,7 +142,7 @@ def test_to_heavy(self):
         hModel = hFoot.getImage(fill=0).image.array
 
         self.assertEqual(bbox, hFoot.getBBox())
-        self.assertFloatsAlmostEqual(hModel, src.get_model(observation=observation), rtol=1e-6, atol=1e-6)
+        self.assertFloatsAlmostEqual(hModel, src.get_model(observation=observation), rtol=1e-4, atol=1e-4)
 
 
 class TestLsstBlend(lsst.utils.tests.TestCase):
@@ -157,7 +157,7 @@ def test_fit_pointSource(self):
         targetPsfImage, psfImages, images, channels, seds, morphs, targetPsf, psfs = result
         B, Ny, Nx = shape
 
-        frame = lmeScarlet.LsstFrame(shape, psfs=targetPsfImage[None])
+        frame = lmeScarlet.LsstFrame(shape, psfs=targetPsfImage[None], dtype=np.float)
         observation = lmeScarlet.LsstObservation(images, psfs=psfImages).match(frame)
         bgRms = np.ones((B, )) * 1e-3
         sources = []
@@ -171,10 +171,8 @@ def test_fit_pointSource(self):
         blend.fit(10)
 
         self.assertEqual(blend.it, 2)
-        self.assertFloatsAlmostEqual(blend.L_sed, 2.5481250470053265, rtol=1e-10, atol=1e-10)
-        self.assertFloatsAlmostEqual(blend.L_morph, 9024.538938935855, rtol=1e-10, atol=1e-10)
-        self.assertFloatsAlmostEqual(np.array(blend.mse),
-                                     np.array([3.875628098330452e-15, 3.875598349723412e-15]))
+        self.assertFloatsAlmostEqual(blend.L_sed, 2.5481250470053265, rtol=1e-5, atol=1e-5)
+        self.assertFloatsAlmostEqual(blend.L_morph, 9024.538938935855, rtol=1e-5, atol=1e-5)
         self.assertTrue(blend.mse[0] > blend.mse[1])
 
     def test_get_model(self):
@@ -198,7 +196,7 @@ def test_get_model(self):
         self.assertEqual(blend.observations[0], observation)
         self.assertEqual(blend.mse, [])
         model = blend.get_model(observation=observation)
-        self.assertFloatsAlmostEqual(model, images, rtol=1e-6, atol=1e-6)
+        self.assertFloatsAlmostEqual(model, images, rtol=1e-5, atol=1e-5)
 
 
 class MemoryTester(lsst.utils.tests.MemoryTestCase):

tests/utils.py

@@ -57,9 +57,9 @@ def initData(shape, coords, amplitudes=None, convolve=True):
         np.arange(B, dtype=float)[::-1],
         np.ones((B,), dtype=float)
     ]
-    seds = np.array([_seds[n % 3]*amplitudes[n] for n in range(K)])
+    seds = np.array([_seds[n % 3]*amplitudes[n] for n in range(K)], dtype=np.float32)
 
-    morphs = np.zeros((K, Ny, Nx))
+    morphs = np.zeros((K, Ny, Nx), dtype=np.float32)
     for k, coord in enumerate(coords):
         morphs[k, coord[0], coord[1]] = 1
     images = seds.T.dot(morphs.reshape(K, -1)).reshape(shape)