From 685ab81c2ce752cd8fe77c8edf7237152ca4e912 Mon Sep 17 00:00:00 2001 From: beckermr Date: Fri, 1 Mar 2024 06:25:12 -0600 Subject: [PATCH 01/27] ENH add fixed apodization for DES --- piff/psf.py | 57 ++++++++++++++++++++++++++++++++++++++++++++++++++++- 1 file changed, 56 insertions(+), 1 deletion(-) diff --git a/piff/psf.py b/piff/psf.py index aa6f12ca..1bead0b9 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -20,10 +20,31 @@ import fitsio import galsim import sys +from numba import vectorize from .star import Star, StarData from .util import write_kwargs, read_kwargs + +@vectorize +def _ap_kern_kern(x, m, h): + # cumulative triweight kernel + y = (x - m) / h + 3 + if y < -3: + return 0 + elif y > 3: + return 1 + else: + val = ( + -5 * y ** 7 / 69984 + + 7 * y ** 5 / 2592 + - 35 * y ** 3 / 864 + + 35 * y / 96 + + 1 / 2 + ) + return val + + class PSF(object): """The base class for describing a PSF model across a field of view. @@ -99,7 +120,7 @@ def parseKwargs(cls, config_psf, logger=None): raise NotImplementedError("Derived classes must define the parseKwargs function") def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_size=48, - image=None, logger=None, **kwargs): + image=None, logger=None, apodize=True, **kwargs): r"""Draws an image of the PSF at a given location. The normal usage would be to specify (chipnum, x, y), in which case Piff will use the @@ -161,6 +182,10 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz [default: 48] :param image: An existing image on which to draw, if desired. [default: None] :param logger: A logger object for logging debug info. [default: None] + :param apodize: Optional float giving the number of half light radii after + which the profile is smoothy apodized to zero at the image edge. + [default: True which uses a default transition width of 1 pixel and + a default apodization radius of 4.25 pixels] :param \**kwargs: Any additional properties required for the interpolation. :returns: A GalSim Image of the PSF @@ -201,6 +226,36 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz prof.drawImage(image, method=method, center=center) + if apodize: + xpix, ypix = image.get_pixel_centers() + dx = xpix - center[0] + dy = ypix - center[1] + r2 = dx**2 + dy**2 + + # this algorithm is adaptive but for DES we'll fix things + # to a constant value for stability + # _im = image.copy() + # _im.wcs = None + # _im.scale = 1.0 + # hlr = _im.calculateHLR(center=galsim.PositionD(center)) + # aprad = apodize * hlr + # msk_nonzero = _im.array != 0 + # max_r = min( + # np.abs(dx[(dx < 0) & msk_nonzero].min()), + # np.abs(dx[(dx > 0) & msk_nonzero].max()), + # np.abs(dy[(dy < 0) & msk_nonzero].min()), + # np.abs(dy[(dy > 0) & msk_nonzero].max()), + # ) + # apwidth = np.abs(hlr) / 2.355 + # apwidth = min(max(apwidth, 0.5), 5.0) + # aprad = max(min(aprad, max_r - 6 * apwidth - 1), 2 * apwidth) + + apwidth = 1.0 + aprad = 4.25 + + apim = image._array * _ap_kern_kern(aprad, np.sqrt(r2), apwidth) + image._array = apim / np.sum(apim) * np.sum(image._array) + return image def get_profile(self, x, y, chipnum=None, flux=1.0, logger=None, **kwargs): From 54fba85bc50e3ebb0d9a9e5c607a3917095dbcac Mon Sep 17 00:00:00 2001 From: beckermr Date: Fri, 1 Mar 2024 06:42:19 -0600 Subject: [PATCH 02/27] REF redo options a bit --- piff/psf.py | 52 +++++++++++++++++++++++++++------------------------- 1 file changed, 27 insertions(+), 25 deletions(-) diff --git a/piff/psf.py b/piff/psf.py index 1bead0b9..5c3ff0c1 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -120,7 +120,7 @@ def parseKwargs(cls, config_psf, logger=None): raise NotImplementedError("Derived classes must define the parseKwargs function") def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_size=48, - image=None, logger=None, apodize=True, **kwargs): + image=None, logger=None, apodize=(1.0, 4.25), **kwargs): r"""Draws an image of the PSF at a given location. The normal usage would be to specify (chipnum, x, y), in which case Piff will use the @@ -182,10 +182,11 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz [default: 48] :param image: An existing image on which to draw, if desired. [default: None] :param logger: A logger object for logging debug info. [default: None] - :param apodize: Optional float giving the number of half light radii after - which the profile is smoothy apodized to zero at the image edge. - [default: True which uses a default transition width of 1 pixel and - a default apodization radius of 4.25 pixels] + :param apodize: Optional parameter to set apodizatoon. If a float/int, gives the + number of half light radii after which the profile is smoothy apodized + to zero a width of ~2.55 half light radii. If a tuple/list, gives + the apodization width and the apodization radius in pixels. + [default: None, which means no apodization.] :param \**kwargs: Any additional properties required for the interpolation. :returns: A GalSim Image of the PSF @@ -232,26 +233,27 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz dy = ypix - center[1] r2 = dx**2 + dy**2 - # this algorithm is adaptive but for DES we'll fix things - # to a constant value for stability - # _im = image.copy() - # _im.wcs = None - # _im.scale = 1.0 - # hlr = _im.calculateHLR(center=galsim.PositionD(center)) - # aprad = apodize * hlr - # msk_nonzero = _im.array != 0 - # max_r = min( - # np.abs(dx[(dx < 0) & msk_nonzero].min()), - # np.abs(dx[(dx > 0) & msk_nonzero].max()), - # np.abs(dy[(dy < 0) & msk_nonzero].min()), - # np.abs(dy[(dy > 0) & msk_nonzero].max()), - # ) - # apwidth = np.abs(hlr) / 2.355 - # apwidth = min(max(apwidth, 0.5), 5.0) - # aprad = max(min(aprad, max_r - 6 * apwidth - 1), 2 * apwidth) - - apwidth = 1.0 - aprad = 4.25 + if isinstance(apodize, (tuple, list)): + apwidth, aprad = apodize + else: + wcs = image.wcs + try: + image.wcs = None + image.scale = 1.0 + hlr = image.calculateHLR(center=galsim.PositionD(center)) + finally: + image.wcs = wcs + aprad = apodize * hlr + msk_nonzero = image.array != 0 + max_r = min( + np.abs(dx[(dx < 0) & msk_nonzero].min()), + np.abs(dx[(dx > 0) & msk_nonzero].max()), + np.abs(dy[(dy < 0) & msk_nonzero].min()), + np.abs(dy[(dy > 0) & msk_nonzero].max()), + ) + apwidth = np.abs(hlr) / 2.355 + apwidth = min(max(apwidth, 0.5), 5.0) + aprad = max(min(aprad, max_r - 6 * apwidth - 1), 2 * apwidth) apim = image._array * _ap_kern_kern(aprad, np.sqrt(r2), apwidth) image._array = apim / np.sum(apim) * np.sum(image._array) From 139066ab8ca1ca23349bd0a3cc475c0c8182172a Mon Sep 17 00:00:00 2001 From: beckermr Date: Fri, 1 Mar 2024 06:45:46 -0600 Subject: [PATCH 03/27] REF no numba --- piff/psf.py | 27 ++++++++++++--------------- 1 file changed, 12 insertions(+), 15 deletions(-) diff --git a/piff/psf.py b/piff/psf.py index 5c3ff0c1..4a6496ad 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -20,29 +20,26 @@ import fitsio import galsim import sys -from numba import vectorize from .star import Star, StarData from .util import write_kwargs, read_kwargs -@vectorize def _ap_kern_kern(x, m, h): # cumulative triweight kernel y = (x - m) / h + 3 - if y < -3: - return 0 - elif y > 3: - return 1 - else: - val = ( - -5 * y ** 7 / 69984 - + 7 * y ** 5 / 2592 - - 35 * y ** 3 / 864 - + 35 * y / 96 - + 1 / 2 - ) - return val + apval = np.zeros_like(m) + msk = y > 3 + apval[msk] = 1 + msk = (y > -3) & (~msk) + apval[msk] = ( + -5 * y[msk] ** 7 / 69984 + + 7 * y[msk] ** 5 / 2592 + - 35 * y[msk] ** 3 / 864 + + 35 * y[msk] / 96 + + 1 / 2 + ) + return apval class PSF(object): From aba23c760f9e2f0d4cc87ef2941d8961cf438d25 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Sat, 2 Mar 2024 06:22:38 -0600 Subject: [PATCH 04/27] Update piff/psf.py --- piff/psf.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/piff/psf.py b/piff/psf.py index 4a6496ad..3850602a 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -117,7 +117,7 @@ def parseKwargs(cls, config_psf, logger=None): raise NotImplementedError("Derived classes must define the parseKwargs function") def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_size=48, - image=None, logger=None, apodize=(1.0, 4.25), **kwargs): + image=None, logger=None, apodize=None, **kwargs): r"""Draws an image of the PSF at a given location. The normal usage would be to specify (chipnum, x, y), in which case Piff will use the From 8e645842068ee49c3d00facd6c869436cd56eb96 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 14:48:17 -0500 Subject: [PATCH 05/27] Update piff/psf.py --- piff/psf.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/piff/psf.py b/piff/psf.py index 3850602a..4a6496ad 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -117,7 +117,7 @@ def parseKwargs(cls, config_psf, logger=None): raise NotImplementedError("Derived classes must define the parseKwargs function") def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_size=48, - image=None, logger=None, apodize=None, **kwargs): + image=None, logger=None, apodize=(1.0, 4.25), **kwargs): r"""Draws an image of the PSF at a given location. The normal usage would be to specify (chipnum, x, y), in which case Piff will use the From c76660a2db16fac4b998b3f8d45fde0e12f9e4c2 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 14:48:59 -0500 Subject: [PATCH 06/27] Update piff/psf.py --- piff/psf.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/piff/psf.py b/piff/psf.py index 4a6496ad..a6e935cf 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -183,7 +183,7 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz number of half light radii after which the profile is smoothy apodized to zero a width of ~2.55 half light radii. If a tuple/list, gives the apodization width and the apodization radius in pixels. - [default: None, which means no apodization.] + [default: (1.0, 4.25), which means a width of 1 pixel and radius of 4.25 pixels.] :param \**kwargs: Any additional properties required for the interpolation. :returns: A GalSim Image of the PSF From cec733c7979d5854571cb91b3b7fde7807f0494e Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 14:51:31 -0500 Subject: [PATCH 07/27] Update _version.py --- piff/_version.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/piff/_version.py b/piff/_version.py index 841c3f9e..c1ba101f 100644 --- a/piff/_version.py +++ b/piff/_version.py @@ -12,5 +12,5 @@ # this list of conditions and the disclaimer given in the documentation # and/or other materials provided with the distribution. -__version__ = '1.3' +__version__ = '1.3.3.1' __version_info__ = tuple(map(int, __version__.split('.'))) From 88755ec6c19cbf7b099ee1d7fe6631b669d41568 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 15:07:30 -0500 Subject: [PATCH 08/27] Update ci.yml --- .github/workflows/ci.yml | 1 + 1 file changed, 1 insertion(+) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 822c8447..5eb0a428 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -16,6 +16,7 @@ jobs: runs-on: ${{ matrix.os }} strategy: + fail-fast: false matrix: # First all python versions in basic linux os: [ ubuntu-latest ] From 0448559af37d3c20b334c034af9ea1b31e2bc5fe Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 16:47:59 -0500 Subject: [PATCH 09/27] Update requirements.txt --- requirements.txt | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/requirements.txt b/requirements.txt index d03b7b19..cfd6693b 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,10 +1,10 @@ -numpy>=1.17 +numpy>=1.17,<2 scipy>=1.2 pyyaml>=5.1 lsstdesc.coord>=1.2 treecorr>=4.3.1 fitsio>=1.0 matplotlib>=3.3 -galsim>=2.3 +galsim>=2.3,<2.5 treegp>=0.6 threadpoolctl>=3.1 From c0161b4634f902d2ac28ee91e5f51f7222f1da30 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 16:50:17 -0500 Subject: [PATCH 10/27] Update ci.yml --- .github/workflows/ci.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index 5eb0a428..db78dbfc 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -82,7 +82,7 @@ jobs: python -m pip install -U pip # Do these first to clarify potential conflicts - pip install -U setuptools numpy + pip install -U setuptools "numpy<2" # Standard dependencies pip install -U -r requirements.txt From 1e8280eaa0b0246fad4bef8ed606c5a38ea303ea Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:26:57 -0500 Subject: [PATCH 11/27] Update test_pixel.py --- tests/test_pixel.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_pixel.py b/tests/test_pixel.py index 2cd33e73..d0b4fc22 100644 --- a/tests/test_pixel.py +++ b/tests/test_pixel.py @@ -1531,7 +1531,7 @@ def test_color(): # Check the convenience function that an end user would typically use offset = s.center_to_offset(s.fit.center) image = psf.draw(x=s['x'], y=s['y'], color=s['color'], - stamp_size=32, flux=s.fit.flux, offset=offset) + stamp_size=32, flux=s.fit.flux, offset=offset, apodize=None) # They may be up to 1 pixel off in either direction, so find the common bounds. b = image.bounds & fit_stamp.bounds np.testing.assert_allclose(image[b].array, fit_stamp[b].array, rtol=1.e-6, atol=1.e-4) From 767f46a4d57ce48d4a4f173ee3727a27a4547f35 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:27:23 -0500 Subject: [PATCH 12/27] Update test_simple.py --- tests/test_simple.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index 57e23910..2678451c 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -309,7 +309,7 @@ def test_single_image(): # test that draw works test_image = psf.draw(x=target['x'], y=target['y'], stamp_size=config['input']['stamp_size'], - flux=target.fit.flux, offset=target.fit.center) + flux=target.fit.flux, offset=target.fit.center, apodize=None) # this image should be the same values as test_star assert test_image == test_star.image # test that draw does not copy the image From 78494ce0e48d709f069fad7345c4602bfa9806ac Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:27:53 -0500 Subject: [PATCH 13/27] Update test_simple.py --- tests/test_simple.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index 2678451c..5f211b8b 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -743,7 +743,7 @@ def test_draw(): # Now use the regular PSF.draw() command. This version is equivalent to the above. # (It's not equal all the way to machine precision, but pretty close.) - im1 = psf.draw(x, y, chipnum, stamp_size=48) + im1 = psf.draw(x, y, chipnum, stamp_size=48, , apodize=None) np.testing.assert_allclose(im1.array, star.data.image.array, rtol=1.e-14, atol=1.e-14) # The wcs in the image is the wcs of the original image From 030911546d10ed302bcbe80db462582571dd98fe Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:28:33 -0500 Subject: [PATCH 14/27] Update test_wcs.py --- tests/test_wcs.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index 2a373671..cbf36e93 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -366,7 +366,7 @@ def test_single(): # This is the more user-friendly way to do this. # Equivalent to ~machine precision. - im = psf.draw(x, y, chipnum=chipnum) + im = psf.draw(x, y, chipnum=chipnum, apodize=None) print('im = ',im) print('star im = ',star.data.image) print('max diff = ',np.max(np.abs(im.array - star.data.image.array))) From 60b4c72235376a3dc49fb3a122e6b2bdc82c2d07 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:29:19 -0500 Subject: [PATCH 15/27] Update test_wcs.py --- tests/test_wcs.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index cbf36e93..9546ed7a 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -612,7 +612,7 @@ def test_olddes(): print('area at 0,0 = ',psf.wcs[0].pixelArea(galsim.PositionD(0,0)),' = %f**2'%( psf.wcs[0].pixelArea(galsim.PositionD(0,0))**0.5)) assert np.isclose(psf.wcs[0].pixelArea(galsim.PositionD(0,0)), 0.2628**2, rtol=1.e-3) - image = psf.draw(x=103.3, y=592.0, logger=logger) + image = psf.draw(x=103.3, y=592.0, logger=logger, apodize=None) print('image shape = ',image.array.shape) print('image near center = ',image.array[23:26,23:26]) print('image sum = ',image.array.sum()) @@ -668,7 +668,7 @@ def test_newdes(): print('area at 0,0 = ',psf.wcs[0].pixelArea(galsim.PositionD(0,0)),' = %f**2'%( psf.wcs[0].pixelArea(galsim.PositionD(0,0))**0.5)) assert np.isclose(psf.wcs[0].pixelArea(galsim.PositionD(0,0)), 0.263021**2, rtol=1.e-3) - image = psf.draw(x=103.3, y=592.0, logger=logger) + image = psf.draw(x=103.3, y=592.0, logger=logger, apodize=None) print('image shape = ',image.array.shape) print('image near center = ',image.array[23:26,23:26]) print('image sum = ',image.array.sum()) From 1a81a807ebf78693cc0f05ca2a37d46b0c8abcae Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:32:36 -0500 Subject: [PATCH 16/27] Update tests/test_simple.py --- tests/test_simple.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index 5f211b8b..c9865790 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -743,7 +743,7 @@ def test_draw(): # Now use the regular PSF.draw() command. This version is equivalent to the above. # (It's not equal all the way to machine precision, but pretty close.) - im1 = psf.draw(x, y, chipnum, stamp_size=48, , apodize=None) + im1 = psf.draw(x, y, chipnum, stamp_size=48, apodize=None) np.testing.assert_allclose(im1.array, star.data.image.array, rtol=1.e-14, atol=1.e-14) # The wcs in the image is the wcs of the original image From b4b5ab3a2edc4e076903c5b328eca11e43e23d4e Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 17:37:59 -0500 Subject: [PATCH 17/27] Update ci.yml --- .github/workflows/ci.yml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml index db78dbfc..c3e15890 100644 --- a/.github/workflows/ci.yml +++ b/.github/workflows/ci.yml @@ -88,7 +88,7 @@ jobs: pip install -U -r requirements.txt # Extra packages needed for testing - pip install -U nose coverage pytest nbval ipykernel + pip install -U nose coverage "pytest<8" nbval ipykernel - name: Install Pixmappy (not on pip) run: | From 5105bf3808ac06a3913d3bb2ef1afd54ac0dd3f2 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 22:14:35 -0500 Subject: [PATCH 18/27] Update test_simple.py --- tests/test_simple.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index c9865790..8c858619 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -314,7 +314,7 @@ def test_single_image(): assert test_image == test_star.image # test that draw does not copy the image image_ref = psf.draw(x=target['x'], y=target['y'], stamp_size=config['input']['stamp_size'], - flux=target.fit.flux, offset=target.fit.center, image=test_image) + flux=target.fit.flux, offset=target.fit.center, image=test_image, apodize=None) image_ref.array[0,0] = 123456789 assert test_image.array[0,0] == image_ref.array[0,0] assert test_star.image.array[0,0] != test_image.array[0,0] From acb61985dde3042e2de2fa16e65db45b797e6317 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 22:17:24 -0500 Subject: [PATCH 19/27] Update test_simple.py --- tests/test_simple.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index 8c858619..09c19a77 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -768,13 +768,13 @@ def test_draw(): # We can center the star at an arbitrary location on the image. # The default is equivalent to center=(x,y). So check that this is equivalent. # Also, 48 is the default stamp size, so that can be omitted here. - im2 = psf.draw(x, y, chipnum, center=(x,y)) + im2 = psf.draw(x, y, chipnum, center=(x,y), apodize=None) assert im2.bounds == im1.bounds np.testing.assert_allclose(im2.array, im1.array, rtol=1.e-14, atol=1.e-14) # Moving by an integer number of pixels should be very close to the same image # over a different slice of the array. - im3 = psf.draw(x, y, chipnum, center=(x+1, y+3)) + im3 = psf.draw(x, y, chipnum, center=(x+1, y+3), apodize=None) assert im3.bounds == im1.bounds # (Remember -- numpy indexing is y,x!) # Also, the FFTs will be different in detail, so only match to 1.e-6. @@ -814,7 +814,7 @@ def test_draw(): # then center=True works fine to draw in the center of that image. im6 = im5.copy() im6.setCenter(0,0) - psf.draw(x, y, chipnum, center=True, image=im6) + psf.draw(x, y, chipnum, center=True, image=im6, apodize=None) assert im6.bounds.center == galsim.PositionI(0,0) np.testing.assert_allclose(im6.array.sum(), 1., rtol=1.e-3) hsm = im6.FindAdaptiveMom() @@ -845,7 +845,7 @@ def test_draw(): # The offset parameter can add an additional to whatever center is used. # Here center=None, so this is equivalent to im4 above. - im9 = psf.draw(x, y, chipnum, offset=(1.3,-0.8)) + im9 = psf.draw(x, y, chipnum, offset=(1.3,-0.8), apodize=None) assert im9.bounds == im1.bounds hsm = im9.FindAdaptiveMom() np.testing.assert_allclose(im9.array, im4.array, rtol=1.e-14, atol=1.e-14) @@ -854,7 +854,7 @@ def test_draw(): # use for this, but it's allowed. (The above with default center is used in unit # tests a number of times, so that version at least is useful if only for us. # I'm hard pressed to imaging end users wanting to specify things this way though.) - im10 = psf.draw(x, y, chipnum, center=(x+0.8, y-0.3), offset=(0.5,-0.5)) + im10 = psf.draw(x, y, chipnum, center=(x+0.8, y-0.3), offset=(0.5,-0.5), apodize=None) assert im10.bounds == im1.bounds np.testing.assert_allclose(im10.array, im4.array, rtol=1.e-14, atol=1.e-14) From 99f5106f4f19dc971015c4c3f941263e1cf965e1 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 22:18:49 -0500 Subject: [PATCH 20/27] Update test_wcs.py --- tests/test_wcs.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index 9546ed7a..eaed388e 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -628,7 +628,7 @@ def test_olddes(): # Also check that it is picklable. psf2 = copy.deepcopy(psf) - image2 = psf2.draw(x=103.3, y=592.0) + image2 = psf2.draw(x=103.3, y=592.0, apodize=None) np.testing.assert_equal(image2.array, image.array) @timer @@ -684,7 +684,7 @@ def test_newdes(): # Also check that it is picklable. psf2 = copy.deepcopy(psf) - image2 = psf2.draw(x=103.3, y=592.0) + image2 = psf2.draw(x=103.3, y=592.0, apodize=None) np.testing.assert_equal(image2.array, image.array) @timer From 2f9d6b42d7c213ce0c19660cc7942b34742119cb Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Mon, 15 Jul 2024 22:25:06 -0500 Subject: [PATCH 21/27] Update test_simple.py --- tests/test_simple.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/tests/test_simple.py b/tests/test_simple.py index 09c19a77..51c85be8 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -788,14 +788,14 @@ def test_draw(): # Can center at other locations, and the hsm centroids should come out centered pretty # close to that location. # (Of course the array will be different here, so can't test that.) - im4 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8)) + im4 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), apodize=None) assert im4.bounds == im1.bounds hsm = im4.FindAdaptiveMom() np.testing.assert_allclose(hsm.moments_centroid.x, x+1.3, atol=0.01) np.testing.assert_allclose(hsm.moments_centroid.y, y-0.8, atol=0.01) # Also allowed is center=True to place in the center of the image. - im5 = psf.draw(x, y, chipnum, center=True) + im5 = psf.draw(x, y, chipnum, center=True, apodize=None) assert im5.bounds == im1.bounds assert im5.array.shape == (48,48) np.testing.assert_allclose(im5.bounds.true_center.x, x, atol=0.5) @@ -824,7 +824,7 @@ def test_draw(): np.testing.assert_allclose(im6.array, im5.array, rtol=1.e-14, atol=1.e-14) # Check non-even stamp size. Also, not unit flux while we're at it. - im7 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), stamp_size=43, flux=23.7) + im7 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), stamp_size=43, flux=23.7, apodize=None) assert im7.array.shape == (43,43) np.testing.assert_allclose(im7.bounds.true_center.x, x, atol=0.5) np.testing.assert_allclose(im7.bounds.true_center.y, y, atol=0.5) @@ -836,7 +836,7 @@ def test_draw(): # Can't do mixed even/odd shape with stamp_size, but it will respect a provided image. im8 = galsim.Image(43,44) im8.setCenter(x,y) # It will respect the given bounds, so put it near the right place. - psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), image=im8, flux=23.7) + psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), image=im8, flux=23.7, apodize=None) assert im8.array.shape == (44,43) np.testing.assert_allclose(im8.array.sum(), 23.7, rtol=1.e-3) hsm = im8.FindAdaptiveMom() From 4351177b758194be22c74792b2d3573dd81239bc Mon Sep 17 00:00:00 2001 From: beckermr Date: Tue, 16 Jul 2024 11:45:32 -0500 Subject: [PATCH 22/27] REF move apodization to the pixelgrid --- piff/pixelgrid.py | 35 ++++++++++++++++++++++++++++ piff/psf.py | 55 +------------------------------------------- tests/test_pixel.py | 6 ++++- tests/test_simple.py | 32 ++++++++++++++++---------- tests/test_wcs.py | 23 ++++++++++++++---- 5 files changed, 79 insertions(+), 72 deletions(-) diff --git a/piff/pixelgrid.py b/piff/pixelgrid.py index 4d113c77..213c174a 100644 --- a/piff/pixelgrid.py +++ b/piff/pixelgrid.py @@ -25,6 +25,26 @@ from .model import Model from .star import Star, StarData, StarFit +APODIZE_PARAMS = (1.0 * 0.263, 4.25 * 0.263) + + +def _ap_kern_kern(x, m, h): + # cumulative triweight kernel + y = (x - m) / h + 3 + apval = np.zeros_like(m) + msk = y > 3 + apval[msk] = 1 + msk = (y > -3) & (~msk) + apval[msk] = ( + -5 * y[msk] ** 7 / 69984 + + 7 * y[msk] ** 5 / 2592 + - 35 * y[msk] ** 3 / 864 + + 35 * y[msk] / 96 + + 1 / 2 + ) + return apval + + class PixelGrid(Model): """A PSF modeled as interpolation between a grid of points. @@ -445,6 +465,21 @@ def getProfile(self, params): :returns: a galsim.GSObject instance """ im = galsim.Image(params.reshape(self.size,self.size), scale=self.scale) + + if APODIZE_PARAMS is not None: + xpix, ypix = im.get_pixel_centers() + # use_true_center = False below + dx = xpix - im.center.x + dy = ypix - im.center.y + r2 = dx**2 + dy**2 + + apwidth, aprad = APODIZE_PARAMS # in arcsec + apwidth /= self.scale # convert to pixels + aprad /= self.scale # convert to pixels + + apim = im._array * _ap_kern_kern(aprad, np.sqrt(r2), apwidth) + im._array = apim / np.sum(apim) * np.sum(im._array) + return galsim.InterpolatedImage(im, x_interpolant=self.interp, use_true_center=False, flux=1.) diff --git a/piff/psf.py b/piff/psf.py index d1a2d9ef..b8f877f1 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -25,23 +25,6 @@ from .util import write_kwargs, read_kwargs -def _ap_kern_kern(x, m, h): - # cumulative triweight kernel - y = (x - m) / h + 3 - apval = np.zeros_like(m) - msk = y > 3 - apval[msk] = 1 - msk = (y > -3) & (~msk) - apval[msk] = ( - -5 * y[msk] ** 7 / 69984 - + 7 * y[msk] ** 5 / 2592 - - 35 * y[msk] ** 3 / 864 - + 35 * y[msk] / 96 - + 1 / 2 - ) - return apval - - class PSF(object): """The base class for describing a PSF model across a field of view. @@ -117,7 +100,7 @@ def parseKwargs(cls, config_psf, logger=None): raise NotImplementedError("Derived classes must define the parseKwargs function") def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_size=48, - image=None, logger=None, apodize=(1.0, 4.25), **kwargs): + image=None, logger=None, **kwargs): r"""Draws an image of the PSF at a given location. The normal usage would be to specify (chipnum, x, y), in which case Piff will use the @@ -179,11 +162,6 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz [default: 48] :param image: An existing image on which to draw, if desired. [default: None] :param logger: A logger object for logging debug info. [default: None] - :param apodize: Optional parameter to set apodizatoon. If a float/int, gives the - number of half light radii after which the profile is smoothy apodized - to zero a width of ~2.55 half light radii. If a tuple/list, gives - the apodization width and the apodization radius in pixels. - [default: (1.0, 4.25), which means a width of 1 pixel and radius of 4.25 pixels.] :param \**kwargs: Any additional properties required for the interpolation. :returns: A GalSim Image of the PSF @@ -224,37 +202,6 @@ def draw(self, x, y, chipnum=None, flux=1.0, center=None, offset=None, stamp_siz prof.drawImage(image, method=method, center=center) - if apodize: - xpix, ypix = image.get_pixel_centers() - dx = xpix - center[0] - dy = ypix - center[1] - r2 = dx**2 + dy**2 - - if isinstance(apodize, (tuple, list)): - apwidth, aprad = apodize - else: - wcs = image.wcs - try: - image.wcs = None - image.scale = 1.0 - hlr = image.calculateHLR(center=galsim.PositionD(center)) - finally: - image.wcs = wcs - aprad = apodize * hlr - msk_nonzero = image.array != 0 - max_r = min( - np.abs(dx[(dx < 0) & msk_nonzero].min()), - np.abs(dx[(dx > 0) & msk_nonzero].max()), - np.abs(dy[(dy < 0) & msk_nonzero].min()), - np.abs(dy[(dy > 0) & msk_nonzero].max()), - ) - apwidth = np.abs(hlr) / 2.355 - apwidth = min(max(apwidth, 0.5), 5.0) - aprad = max(min(aprad, max_r - 6 * apwidth - 1), 2 * apwidth) - - apim = image._array * _ap_kern_kern(aprad, np.sqrt(r2), apwidth) - image._array = apim / np.sum(apim) * np.sum(image._array) - return image def get_profile(self, x, y, chipnum=None, flux=1.0, logger=None, **kwargs): diff --git a/tests/test_pixel.py b/tests/test_pixel.py index d0b4fc22..43799979 100644 --- a/tests/test_pixel.py +++ b/tests/test_pixel.py @@ -1436,6 +1436,10 @@ def test_var(): def test_color(): """Test having a color-dependent size. """ + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + pixel_scale = 0.3 image = galsim.ImageF(1024,1024, scale=pixel_scale) @@ -1531,7 +1535,7 @@ def test_color(): # Check the convenience function that an end user would typically use offset = s.center_to_offset(s.fit.center) image = psf.draw(x=s['x'], y=s['y'], color=s['color'], - stamp_size=32, flux=s.fit.flux, offset=offset, apodize=None) + stamp_size=32, flux=s.fit.flux, offset=offset) # They may be up to 1 pixel off in either direction, so find the common bounds. b = image.bounds & fit_stamp.bounds np.testing.assert_allclose(image[b].array, fit_stamp[b].array, rtol=1.e-6, atol=1.e-4) diff --git a/tests/test_simple.py b/tests/test_simple.py index 51c85be8..9170d818 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -139,6 +139,10 @@ def test_single_image(): else: logger = piff.config.setup_logger(log_file='output/test_single_image.log') + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + # Make the image image = galsim.Image(2048, 2048, scale=0.26) @@ -309,12 +313,12 @@ def test_single_image(): # test that draw works test_image = psf.draw(x=target['x'], y=target['y'], stamp_size=config['input']['stamp_size'], - flux=target.fit.flux, offset=target.fit.center, apodize=None) + flux=target.fit.flux, offset=target.fit.center) # this image should be the same values as test_star assert test_image == test_star.image # test that draw does not copy the image image_ref = psf.draw(x=target['x'], y=target['y'], stamp_size=config['input']['stamp_size'], - flux=target.fit.flux, offset=target.fit.center, image=test_image, apodize=None) + flux=target.fit.flux, offset=target.fit.center, image=test_image) image_ref.array[0,0] = 123456789 assert test_image.array[0,0] == image_ref.array[0,0] assert test_star.image.array[0,0] != test_image.array[0,0] @@ -708,6 +712,10 @@ def test_draw(): else: logger = piff.config.setup_logger(log_file='output/test_draw.log') + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + # Use an existing Piff solution to match as closely as possible how users would actually # use this function. psf = piff.read('input/test_single_py27.piff', logger=logger) @@ -743,7 +751,7 @@ def test_draw(): # Now use the regular PSF.draw() command. This version is equivalent to the above. # (It's not equal all the way to machine precision, but pretty close.) - im1 = psf.draw(x, y, chipnum, stamp_size=48, apodize=None) + im1 = psf.draw(x, y, chipnum, stamp_size=48) np.testing.assert_allclose(im1.array, star.data.image.array, rtol=1.e-14, atol=1.e-14) # The wcs in the image is the wcs of the original image @@ -768,13 +776,13 @@ def test_draw(): # We can center the star at an arbitrary location on the image. # The default is equivalent to center=(x,y). So check that this is equivalent. # Also, 48 is the default stamp size, so that can be omitted here. - im2 = psf.draw(x, y, chipnum, center=(x,y), apodize=None) + im2 = psf.draw(x, y, chipnum, center=(x,y)) assert im2.bounds == im1.bounds np.testing.assert_allclose(im2.array, im1.array, rtol=1.e-14, atol=1.e-14) # Moving by an integer number of pixels should be very close to the same image # over a different slice of the array. - im3 = psf.draw(x, y, chipnum, center=(x+1, y+3), apodize=None) + im3 = psf.draw(x, y, chipnum, center=(x+1, y+3)) assert im3.bounds == im1.bounds # (Remember -- numpy indexing is y,x!) # Also, the FFTs will be different in detail, so only match to 1.e-6. @@ -788,14 +796,14 @@ def test_draw(): # Can center at other locations, and the hsm centroids should come out centered pretty # close to that location. # (Of course the array will be different here, so can't test that.) - im4 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), apodize=None) + im4 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8)) assert im4.bounds == im1.bounds hsm = im4.FindAdaptiveMom() np.testing.assert_allclose(hsm.moments_centroid.x, x+1.3, atol=0.01) np.testing.assert_allclose(hsm.moments_centroid.y, y-0.8, atol=0.01) # Also allowed is center=True to place in the center of the image. - im5 = psf.draw(x, y, chipnum, center=True, apodize=None) + im5 = psf.draw(x, y, chipnum, center=True) assert im5.bounds == im1.bounds assert im5.array.shape == (48,48) np.testing.assert_allclose(im5.bounds.true_center.x, x, atol=0.5) @@ -814,7 +822,7 @@ def test_draw(): # then center=True works fine to draw in the center of that image. im6 = im5.copy() im6.setCenter(0,0) - psf.draw(x, y, chipnum, center=True, image=im6, apodize=None) + psf.draw(x, y, chipnum, center=True, image=im6) assert im6.bounds.center == galsim.PositionI(0,0) np.testing.assert_allclose(im6.array.sum(), 1., rtol=1.e-3) hsm = im6.FindAdaptiveMom() @@ -824,7 +832,7 @@ def test_draw(): np.testing.assert_allclose(im6.array, im5.array, rtol=1.e-14, atol=1.e-14) # Check non-even stamp size. Also, not unit flux while we're at it. - im7 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), stamp_size=43, flux=23.7, apodize=None) + im7 = psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), stamp_size=43, flux=23.7) assert im7.array.shape == (43,43) np.testing.assert_allclose(im7.bounds.true_center.x, x, atol=0.5) np.testing.assert_allclose(im7.bounds.true_center.y, y, atol=0.5) @@ -836,7 +844,7 @@ def test_draw(): # Can't do mixed even/odd shape with stamp_size, but it will respect a provided image. im8 = galsim.Image(43,44) im8.setCenter(x,y) # It will respect the given bounds, so put it near the right place. - psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), image=im8, flux=23.7, apodize=None) + psf.draw(x, y, chipnum, center=(x+1.3,y-0.8), image=im8, flux=23.7) assert im8.array.shape == (44,43) np.testing.assert_allclose(im8.array.sum(), 23.7, rtol=1.e-3) hsm = im8.FindAdaptiveMom() @@ -845,7 +853,7 @@ def test_draw(): # The offset parameter can add an additional to whatever center is used. # Here center=None, so this is equivalent to im4 above. - im9 = psf.draw(x, y, chipnum, offset=(1.3,-0.8), apodize=None) + im9 = psf.draw(x, y, chipnum, offset=(1.3,-0.8)) assert im9.bounds == im1.bounds hsm = im9.FindAdaptiveMom() np.testing.assert_allclose(im9.array, im4.array, rtol=1.e-14, atol=1.e-14) @@ -854,7 +862,7 @@ def test_draw(): # use for this, but it's allowed. (The above with default center is used in unit # tests a number of times, so that version at least is useful if only for us. # I'm hard pressed to imaging end users wanting to specify things this way though.) - im10 = psf.draw(x, y, chipnum, center=(x+0.8, y-0.3), offset=(0.5,-0.5), apodize=None) + im10 = psf.draw(x, y, chipnum, center=(x+0.8, y-0.3), offset=(0.5,-0.5)) assert im10.bounds == im1.bounds np.testing.assert_allclose(im10.array, im4.array, rtol=1.e-14, atol=1.e-14) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index eaed388e..d5fd7d9a 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -253,6 +253,11 @@ def test_wrongwcs(): def test_single(): """Same as test_focal, but using the SingleCCD PSF type, which does a separate fit on each CCD. """ + + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + wcs1 = galsim.TanWCS( galsim.AffineTransform(0.26, 0.05, -0.08, -0.24, galsim.PositionD(1024,1024)), galsim.CelestialCoord(-5 * galsim.arcmin, -25 * galsim.degrees) @@ -366,7 +371,7 @@ def test_single(): # This is the more user-friendly way to do this. # Equivalent to ~machine precision. - im = psf.draw(x, y, chipnum=chipnum, apodize=None) + im = psf.draw(x, y, chipnum=chipnum) print('im = ',im) print('star im = ',star.data.image) print('max diff = ',np.max(np.abs(im.array - star.data.image.array))) @@ -601,6 +606,10 @@ def test_olddes(): else: logger = piff.config.setup_logger(log_file='output/test_olddes.log') + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + fname = os.path.join('input', 'D00240560_r_c01_r2362p01_piff.fits') psf = piff.PSF.read(fname, logger=logger) @@ -612,7 +621,7 @@ def test_olddes(): print('area at 0,0 = ',psf.wcs[0].pixelArea(galsim.PositionD(0,0)),' = %f**2'%( psf.wcs[0].pixelArea(galsim.PositionD(0,0))**0.5)) assert np.isclose(psf.wcs[0].pixelArea(galsim.PositionD(0,0)), 0.2628**2, rtol=1.e-3) - image = psf.draw(x=103.3, y=592.0, logger=logger, apodize=None) + image = psf.draw(x=103.3, y=592.0, logger=logger) print('image shape = ',image.array.shape) print('image near center = ',image.array[23:26,23:26]) print('image sum = ',image.array.sum()) @@ -628,7 +637,7 @@ def test_olddes(): # Also check that it is picklable. psf2 = copy.deepcopy(psf) - image2 = psf2.draw(x=103.3, y=592.0, apodize=None) + image2 = psf2.draw(x=103.3, y=592.0) np.testing.assert_equal(image2.array, image.array) @timer @@ -651,6 +660,10 @@ def test_newdes(): else: logger = piff.config.setup_logger(log_file='output/test_newdes.log') + # turn off apodization + import piff.pixelgrid + piff.pixelgrid.APODIZE_PARAMS = None + fname = os.path.join('input', 'D00232418_i_c19_r5006p01_piff-model.fits') with warnings.catch_warnings(): # This file was written with GalSim 2.1, and now raises a deprecation warning for 2.2. @@ -668,7 +681,7 @@ def test_newdes(): print('area at 0,0 = ',psf.wcs[0].pixelArea(galsim.PositionD(0,0)),' = %f**2'%( psf.wcs[0].pixelArea(galsim.PositionD(0,0))**0.5)) assert np.isclose(psf.wcs[0].pixelArea(galsim.PositionD(0,0)), 0.263021**2, rtol=1.e-3) - image = psf.draw(x=103.3, y=592.0, logger=logger, apodize=None) + image = psf.draw(x=103.3, y=592.0, logger=logger) print('image shape = ',image.array.shape) print('image near center = ',image.array[23:26,23:26]) print('image sum = ',image.array.sum()) @@ -684,7 +697,7 @@ def test_newdes(): # Also check that it is picklable. psf2 = copy.deepcopy(psf) - image2 = psf2.draw(x=103.3, y=592.0, apodize=None) + image2 = psf2.draw(x=103.3, y=592.0) np.testing.assert_equal(image2.array, image.array) @timer From 122ef5ce0ddaaf3b66654196e9701b612e0a18ac Mon Sep 17 00:00:00 2001 From: beckermr Date: Tue, 16 Jul 2024 11:47:03 -0500 Subject: [PATCH 23/27] refactor: make code clear it is not a copy --- piff/pixelgrid.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/piff/pixelgrid.py b/piff/pixelgrid.py index 213c174a..0c10bb8c 100644 --- a/piff/pixelgrid.py +++ b/piff/pixelgrid.py @@ -474,10 +474,10 @@ def getProfile(self, params): r2 = dx**2 + dy**2 apwidth, aprad = APODIZE_PARAMS # in arcsec - apwidth /= self.scale # convert to pixels - aprad /= self.scale # convert to pixels + _apwidth = apwidth / self.scale # convert to pixels + _aprad = aprad / self.scale # convert to pixels - apim = im._array * _ap_kern_kern(aprad, np.sqrt(r2), apwidth) + apim = im._array * _ap_kern_kern(_aprad, np.sqrt(r2), _apwidth) im._array = apim / np.sum(apim) * np.sum(im._array) return galsim.InterpolatedImage(im, x_interpolant=self.interp, From 44e8e2f3f436709fe8520f895c72a1b53b3e1384 Mon Sep 17 00:00:00 2001 From: "Matthew R. Becker" Date: Tue, 16 Jul 2024 11:47:33 -0500 Subject: [PATCH 24/27] Update piff/psf.py --- piff/psf.py | 1 - 1 file changed, 1 deletion(-) diff --git a/piff/psf.py b/piff/psf.py index b8f877f1..f1b526d5 100644 --- a/piff/psf.py +++ b/piff/psf.py @@ -24,7 +24,6 @@ from .star import Star, StarData from .util import write_kwargs, read_kwargs - class PSF(object): """The base class for describing a PSF model across a field of view. From a7165d4112647d896c4280809dceb8272a621fb3 Mon Sep 17 00:00:00 2001 From: beckermr Date: Tue, 16 Jul 2024 11:57:52 -0500 Subject: [PATCH 25/27] test: try this --- tests/conftest.py | 4 ++++ tests/test_pixel.py | 4 ---- tests/test_simple.py | 8 -------- tests/test_wcs.py | 13 ------------- 4 files changed, 4 insertions(+), 25 deletions(-) create mode 100644 tests/conftest.py diff --git a/tests/conftest.py b/tests/conftest.py new file mode 100644 index 00000000..079d67ae --- /dev/null +++ b/tests/conftest.py @@ -0,0 +1,4 @@ +# turn off apodization +import piff.pixelgrid + +piff.pixelgrid.APODIZE_PARAMS = None diff --git a/tests/test_pixel.py b/tests/test_pixel.py index 43799979..2cd33e73 100644 --- a/tests/test_pixel.py +++ b/tests/test_pixel.py @@ -1436,10 +1436,6 @@ def test_var(): def test_color(): """Test having a color-dependent size. """ - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - pixel_scale = 0.3 image = galsim.ImageF(1024,1024, scale=pixel_scale) diff --git a/tests/test_simple.py b/tests/test_simple.py index 9170d818..57e23910 100644 --- a/tests/test_simple.py +++ b/tests/test_simple.py @@ -139,10 +139,6 @@ def test_single_image(): else: logger = piff.config.setup_logger(log_file='output/test_single_image.log') - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - # Make the image image = galsim.Image(2048, 2048, scale=0.26) @@ -712,10 +708,6 @@ def test_draw(): else: logger = piff.config.setup_logger(log_file='output/test_draw.log') - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - # Use an existing Piff solution to match as closely as possible how users would actually # use this function. psf = piff.read('input/test_single_py27.piff', logger=logger) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index d5fd7d9a..2a373671 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -253,11 +253,6 @@ def test_wrongwcs(): def test_single(): """Same as test_focal, but using the SingleCCD PSF type, which does a separate fit on each CCD. """ - - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - wcs1 = galsim.TanWCS( galsim.AffineTransform(0.26, 0.05, -0.08, -0.24, galsim.PositionD(1024,1024)), galsim.CelestialCoord(-5 * galsim.arcmin, -25 * galsim.degrees) @@ -606,10 +601,6 @@ def test_olddes(): else: logger = piff.config.setup_logger(log_file='output/test_olddes.log') - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - fname = os.path.join('input', 'D00240560_r_c01_r2362p01_piff.fits') psf = piff.PSF.read(fname, logger=logger) @@ -660,10 +651,6 @@ def test_newdes(): else: logger = piff.config.setup_logger(log_file='output/test_newdes.log') - # turn off apodization - import piff.pixelgrid - piff.pixelgrid.APODIZE_PARAMS = None - fname = os.path.join('input', 'D00232418_i_c19_r5006p01_piff-model.fits') with warnings.catch_warnings(): # This file was written with GalSim 2.1, and now raises a deprecation warning for 2.2. From 0271d4fc7e0da09b04d979e1207f1d866cdc9de3 Mon Sep 17 00:00:00 2001 From: beckermr Date: Tue, 16 Jul 2024 12:55:55 -0500 Subject: [PATCH 26/27] test: add test --- piff/pixelgrid.py | 5 +++++ tests/conftest.py | 2 +- tests/test_wcs.py | 42 ++++++++++++++++++++++++++++++++++++++++++ 3 files changed, 48 insertions(+), 1 deletion(-) diff --git a/piff/pixelgrid.py b/piff/pixelgrid.py index 0c10bb8c..a7765afd 100644 --- a/piff/pixelgrid.py +++ b/piff/pixelgrid.py @@ -28,6 +28,11 @@ APODIZE_PARAMS = (1.0 * 0.263, 4.25 * 0.263) +def set_apodize_params(pars): + global APODIZE_PARAMS + APODIZE_PARAMS = pars + + def _ap_kern_kern(x, m, h): # cumulative triweight kernel y = (x - m) / h + 3 diff --git a/tests/conftest.py b/tests/conftest.py index 079d67ae..cf6d85c5 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,4 +1,4 @@ # turn off apodization import piff.pixelgrid -piff.pixelgrid.APODIZE_PARAMS = None +piff.pixelgrid.set_apodize_params(None) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index 2a373671..9307b09b 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -747,6 +747,48 @@ def test_des_wcs(): np.testing.assert_allclose(wcs3.toWorld(im.center).y, wcs1.toWorld(im.center).y, rtol=0.04) +@timer +def test_newdes_apodize(): + # This is a DES Y6 PSF file made by Robert Gruendl using python 2, so + # check that this also works correctly. + try: + import pixmappy + except ImportError: + print('pixmappy not installed. Skipping test_newdes()') + return + # Also make sure pixmappy is recent enough to work. + if 'exposure_file' not in pixmappy.GalSimWCS._opt_params: + print('pixmappy not recent enough version. Skipping test_newdes()') + return + + import piff + import piff.pixelgrid + + if __name__ == '__main__': + logger = piff.config.setup_logger(verbose=2) + else: + logger = piff.config.setup_logger(log_file='output/test_newdes.log') + + fname = os.path.join('input', 'D00232418_i_c19_r5006p01_piff-model.fits') + with warnings.catch_warnings(): + # This file was written with GalSim 2.1, and now raises a deprecation warning for 2.2. + warnings.simplefilter("ignore", galsim.GalSimDeprecationWarning) + warnings.simplefilter("ignore", DeprecationWarning) + psf = piff.PSF.read(fname, logger=logger) + + ims = [] + for appars in [None, (1.0 * 0.263, 4.25 * 0.263)]: + piff.pixelgrid.set_apodize_params(appars) + ims.append(psf.draw(x=103.3, y=592.0, logger=logger)) + + assert not np.allclose(ims[0].array, ims[1].array) + assert np.allclose(ims[1].array[0, :], 0, rtol=1.e-2) + assert np.allclose(ims[1].array[-1, :], 0, rtol=1.e-2) + assert np.allclose(ims[1].array[:, 0], 0, rtol=1.e-2) + assert np.allclose(ims[1].array[:, -1], 0, rtol=1.e-2) + assert ims[1].array.sum() > 0 + print('sum = ',ims[1].array.sum()) + if __name__ == '__main__': #import cProfile, pstats #pr = cProfile.Profile() From c1c8eba46706f6ed6a2d217403a5de5f377da820 Mon Sep 17 00:00:00 2001 From: beckermr Date: Tue, 16 Jul 2024 15:41:27 -0500 Subject: [PATCH 27/27] test: add one more test --- tests/test_wcs.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/tests/test_wcs.py b/tests/test_wcs.py index 9307b09b..68f54ab3 100644 --- a/tests/test_wcs.py +++ b/tests/test_wcs.py @@ -781,13 +781,18 @@ def test_newdes_apodize(): piff.pixelgrid.set_apodize_params(appars) ims.append(psf.draw(x=103.3, y=592.0, logger=logger)) + print('sum = ',ims[1].array.sum()) assert not np.allclose(ims[0].array, ims[1].array) assert np.allclose(ims[1].array[0, :], 0, rtol=1.e-2) assert np.allclose(ims[1].array[-1, :], 0, rtol=1.e-2) assert np.allclose(ims[1].array[:, 0], 0, rtol=1.e-2) assert np.allclose(ims[1].array[:, -1], 0, rtol=1.e-2) assert ims[1].array.sum() > 0 - print('sum = ',ims[1].array.sum()) + np.testing.assert_allclose( + ims[0].array[23:26,22:25] / ims[0].array[23:26,22:25].sum(), + ims[1].array[23:26,22:25] / ims[1].array[23:26,22:25].sum(), + rtol=1.e-5, + ) if __name__ == '__main__': #import cProfile, pstats