Fix issues pointed out by reviewer

python/lsst/ip/diffim/zogy.py

@@ -228,7 +228,7 @@ def selectPsf(psf, exposure):
         self.im1_psf = selectPsf(psf1, self.template)
         self.im2_psf = selectPsf(psf2, self.science)
 
-        # Make sure PSFs are the same size. Assume they're square...
+        # Make sure PSFs are the same size.
         if self.im1_psf.shape[0] < self.im2_psf.shape[0]:
             self.im1_psf = np.pad(self.im1_psf, (((self.im2_psf.shape[0] - self.im1_psf.shape[0])//2,
                                                   (self.im2_psf.shape[1] - self.im1_psf.shape[1])//2)),
@@ -284,36 +284,10 @@ def _padPsfToSize(psf, size):
         psf : 2D numpy.array
             The padded copy of the input `psf`.
         """
-        padSize0 = size[0]  # im.shape[0]//2 - psf.shape[0]//2
-        padSize1 = size[1]  # im.shape[1]//2 - psf.shape[1]//2
-        # Hastily assume the psf is odd-sized...
-        if padSize0 > 0 or padSize1 > 0:
-            if padSize0 < 0:
-                padSize0 = 0
-            if padSize1 < 0:
-                padSize1 = 0
-            psf = np.pad(psf, ((padSize0, padSize0-1), (padSize1, padSize1-1)), mode='constant',
-                         constant_values=0)
-        return psf
-
-    @staticmethod
-    def _padPsfToImageSize(psf, im):
-        """Zero-pad `psf` to same dimensions as im.
-
-        Parameters
-        ----------
-        psf : 2D numpy.array
-            Input psf to be padded
-        im : lsst.afw.Image, MaskedImage or Exposure
-            Dimensions of this image are used to set the PSF pad dimensions
-
-        Returns
-        -------
-        psf : 2D numpy.array
-            The padded copy of the input `psf`.
-        """
-        return ZogyTask._padPsfToSize(psf, (im.shape[0]//2 - psf.shape[0]//2,
-                                            im.shape[1]//2 - psf.shape[1]//2))
+        new = np.zeros(size)
+        offset = [size[0]//2 - psf.shape[0]//2 - 1, size[1]//2 - psf.shape[1]//2 - 1]
+        new[offset[0]:(psf.shape[0] + offset[0]), offset[1]:(psf.shape[1] + offset[1])] = psf
+        return new
 
     def computePrereqs(self, psf1=None, psf2=None, padSize=0):
         """Compute standard ZOGY quantities used by (nearly) all methods.
@@ -345,8 +319,8 @@ def computePrereqs(self, psf1=None, psf2=None, padSize=0):
         padSize = self.padSize if padSize is None else padSize
         Pr, Pn = psf1, psf2
         if padSize > 0:
-            Pr = ZogyTask._padPsfToSize(psf1, (padSize, padSize))
-            Pn = ZogyTask._padPsfToSize(psf2, (padSize, padSize))
+            Pr = ZogyTask._padPsfToSize(psf1, (psf1.shape[0]+padSize, psf1.shape[0]+padSize))
+            Pn = ZogyTask._padPsfToSize(psf2, (psf2.shape[0]+padSize, psf2.shape[0]+padSize))
 
         sigR, sigN = self.sig1, self.sig2
         Pr_hat = np.fft.fft2(Pr)
@@ -384,8 +358,8 @@ def computeDiffimFourierSpace(self, debug=False, returnMatchedTemplate=False, **
         - D_var : 2D numpy.array, the variance image for `D`
         """
         # Do all in fourier space (needs image-sized PSFs)
-        psf1 = ZogyTask._padPsfToImageSize(self.im1_psf, self.im1)
-        psf2 = ZogyTask._padPsfToImageSize(self.im2_psf, self.im2)
+        psf1 = ZogyTask._padPsfToSize(self.im1_psf, self.im1.shape)
+        psf2 = ZogyTask._padPsfToSize(self.im2_psf, self.im2.shape)
 
         preqs = self.computePrereqs(psf1, psf2, padSize=0)  # already padded the PSFs
 
@@ -696,8 +670,8 @@ def computeScorrFourierSpace(self, xVarAst=0., yVarAst=0., **kwargs):
         - Dpsf : the PSF of the diffim D, likely never to be used.
         """
         # Do all in fourier space (needs image-sized PSFs)
-        psf1 = ZogyTask._padPsfToImageSize(self.im1_psf, self.im1)
-        psf2 = ZogyTask._padPsfToImageSize(self.im2_psf, self.im2)
+        psf1 = ZogyTask._padPsfToSize(self.im1_psf, self.im1.shape)
+        psf2 = ZogyTask._padPsfToSize(self.im2_psf, self.im2.shape)
 
         preqs = self.computePrereqs(psf1, psf2, padSize=0)  # already padded the PSFs
 
@@ -990,7 +964,7 @@ class ZogyImagePsfMatchConfig(ImagePsfMatchConfig):
 
     def setDefaults(self):
         self.zogyMapReduceConfig.gridStepX = self.zogyMapReduceConfig.gridStepY = 19
-        self.zogyMapReduceConfig.gridSizeX = self.zogyMapReduceConfig.gridSizeY = 20
+        self.zogyMapReduceConfig.cellSizeX = self.zogyMapReduceConfig.cellSizeY = 20
         self.zogyMapReduceConfig.borderSizeX = self.zogyMapReduceConfig.borderSizeY = 6
         self.zogyMapReduceConfig.reducerSubtask.reduceOperation = 'average'
 
@@ -1042,11 +1016,15 @@ def subtractExposures(self, templateExposure, scienceExposure,
         if not self._validateWcs(templateExposure, scienceExposure):
             if doWarping:
                 self.log.info("Astrometrically registering template to science image")
-                templatePsf = templateExposure.getPsf()
+                # Also warp the PSF
+                xyTransform = afwImage.XYTransformFromWcsPair(scienceExposure.getWcs(),
+                                                              templateExposure.getWcs())
+                psfWarped = measAlg.WarpedPsf(templateExposure.getPsf(), xyTransform)
                 templateExposure = self._warper.warpExposure(scienceExposure.getWcs(),
                                                              templateExposure,
                                                              destBBox=scienceExposure.getBBox())
-                templateExposure.setPsf(templatePsf)
+
+                templateExposure.setPsf(psfWarped)
                 templateExposure.writeFits('WARPEDTEMPLATE_ZOGY.fits')
             else:
                 self.log.error("ERROR: Input images not registered")

tests/testImageDecorrelation.py

@@ -120,7 +120,7 @@ def makeFakeImages(size=(256, 256), svar=0.04, tvar=0.04, psf1=3.3, psf2=2.2, of
 
     xim = np.arange(-size[0]//2, size[0]//2, 1)
     yim = np.arange(-size[1]//2, size[1]//2, 1)
-    x0im, y0im = np.meshgrid(xim, yim)
+    x0im, y0im = np.meshgrid(yim, xim)
     im1 = np.random.normal(scale=np.sqrt(svar), size=x0im.shape)  # variance of science image
     im2 = np.random.normal(scale=np.sqrt(tvar), size=x0im.shape)  # variance of template
 
@@ -189,7 +189,7 @@ def makeExposure(imgArray, psfArray, imgVariance):
         @return a new exposure containing the image, PSF and desired variance plane
         """
         # All this code to convert the template image array/psf array into an exposure.
-        bbox = afwGeom.Box2I(afwGeom.Point2I(0, 0), afwGeom.Point2I(imgArray.shape[0]-1, imgArray.shape[1]-1))
+        bbox = afwGeom.Box2I(afwGeom.Point2I(0, 0), afwGeom.Point2I(imgArray.shape[1]-1, imgArray.shape[0]-1))
         im1ex = afwImage.ExposureD(bbox)
         im1ex.getMaskedImage().getImage().getArray()[:, :] = imgArray
         im1ex.getMaskedImage().getVariance().getArray()[:, :] = imgVariance

tests/testZogy.py

@@ -65,7 +65,7 @@ def _setUpImages(self, svar=100., tvar=100., varyPsf=0.):
 
         seed = 666
         self.im1ex, self.im2ex \
-            = makeFakeImages(svar=self.svar, tvar=self.tvar,
+            = makeFakeImages(size=(255, 257), svar=self.svar, tvar=self.tvar,
                              psf1=self.psf1_sigma, psf2=self.psf2_sigma,
                              n_sources=10, psf_yvary_factor=varyPsf,
                              seed=seed, verbose=False)