From dbd9e0762f9e394902e2d0beb96915ac1e30d210 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 14 Mar 2024 21:21:36 -0400 Subject: [PATCH 01/11] working on centers --- imsim/stamp.py | 58 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 58 insertions(+) diff --git a/imsim/stamp.py b/imsim/stamp.py index fcb945ff..8943d07c 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -783,8 +783,66 @@ def draw(self, prof, image, method, offset, config, base, logger): poisson_flux=False) base['realized_flux'] = image.added_flux + xcen, ycen = _get_centroid_from_photons( + gal=gal, + rng=self.rng, + bp_for_drawImage=bp_for_drawImage, + image=image, + sensor=sensor, + photon_ops=photon_ops, + offset=offset, + ) + + # these can be saved in the config using @xcentroid, @ycentroid + base['xcentroid'] = xcen + base['ycentroid'] = ycen + return image + +def _get_centroid_from_photons( + gal, + rng, + bp_for_drawImage, + image, + + sensor, + photon_ops, + offset, +): + """ + draw another image with fixed number of photons + and use to get centroid. We can't do this above + because with maxN the photons cannot be saved + + 100_000 photons should give centroid to a few miliarcsec + """ + + timage = image.copy() + gal.drawImage( + bp_for_drawImage, + image=timage, + + method='phot', + n_photons=1_000_000, + sensor=sensor, + photon_ops=photon_ops, + poisson_flux=False, + save_photons=True, + + offset=offset, + rng=rng, + ) + + # print('imcen:', imcen, 'bounds:', image.bounds, 'off:', xoff, yoff) + + # location of true center, actually in big image + imcen = image.true_center + xcen = imcen.x + timage.photons.x.mean() + ycen = imcen.y + timage.photons.y.mean() + return xcen, ycen + + # Pick the right function to be _fix_seds. if galsim.__version_info__ < (2,5): LSST_SiliconBuilder._fix_seds = LSST_SiliconBuilder._fix_seds_24 From e7286628251c4062f65a0bc09c51eb3eb3a89dd5 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 14 Mar 2024 22:25:50 -0400 Subject: [PATCH 02/11] 100000 is enough --- imsim/stamp.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 8943d07c..07c3850e 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -824,7 +824,7 @@ def _get_centroid_from_photons( image=timage, method='phot', - n_photons=1_000_000, + n_photons=100_000, sensor=sensor, photon_ops=photon_ops, poisson_flux=False, From 04fdf2ba1ae6bf2eb1a73dce112dcf2d874a8e16 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Sat, 23 Mar 2024 22:35:15 -0400 Subject: [PATCH 03/11] refactor cen --- imsim/stamp.py | 56 ++++++++++++++++++++++++++++++++++++++++++-------- 1 file changed, 48 insertions(+), 8 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 07c3850e..7f1a682f 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -758,9 +758,10 @@ def draw(self, prof, image, method, offset, config, base, logger): else: bp_for_drawImage = bandpass - # Put the psfs at the start of the photon_ops. - # Probably a little better to put them a bit later than the start in some cases - # (e.g. after TimeSampler, PupilAnnulusSampler), but leave that as a todo for now. + # Put the psfs at the start of the photon_ops. Probably a little + # better to put them a bit later than the start in some cases (e.g. + # after TimeSampler, PupilAnnulusSampler), but leave that as a todo + # for now. photon_ops = psfs + photon_ops if faint: @@ -783,7 +784,7 @@ def draw(self, prof, image, method, offset, config, base, logger): poisson_flux=False) base['realized_flux'] = image.added_flux - xcen, ycen = _get_centroid_from_photons( + xvals, yvals = _get_photon_positions( gal=gal, rng=self.rng, bp_for_drawImage=bp_for_drawImage, @@ -792,6 +793,20 @@ def draw(self, prof, image, method, offset, config, base, logger): photon_ops=photon_ops, offset=offset, ) + xcen, ycen = _get_robust_centroids( + xvals, yvals, stat_type='sigma_clip', + ) + + print( + 'shapeyx:', image.array.shape[0], image.array.shape[1], + 'true_cen:', image.true_center.y, image.true_center.x, + 'offset:', offset.y, offset.x, + 'cen:', ycen, xcen, + 'cenoff:', ycen-image.true_center.y, xcen-image.true_center.x, + 'maxoff:', + np.abs(yvals-image.true_center.y).max(), + np.abs(xcen-image.true_center.x).max(), + ) # these can be saved in the config using @xcentroid, @ycentroid base['xcentroid'] = xcen @@ -800,7 +815,7 @@ def draw(self, prof, image, method, offset, config, base, logger): return image -def _get_centroid_from_photons( +def _get_photon_positions( gal, rng, bp_for_drawImage, @@ -836,15 +851,40 @@ def _get_centroid_from_photons( # print('imcen:', imcen, 'bounds:', image.bounds, 'off:', xoff, yoff) + photx = timage.photons.x + photy = timage.photons.y + + logic = np.isnan(photx) | np.isnan(photy) + wnan, = np.where(logic) + wgood, = np.where(~logic) + if wnan.size > 0: + print(f'found {wnan.size} nan in photon positions') + # location of true center, actually in big image imcen = image.true_center - xcen = imcen.x + timage.photons.x.mean() - ycen = imcen.y + timage.photons.y.mean() + # xcen = imcen.x + timage.photons.x.mean() + # ycen = imcen.y + timage.photons.y.mean() + xvals = imcen.x + photx[wgood] + yvals = imcen.y + photy[wgood] + return xvals, yvals + + +def _get_robust_centroids(xvals, yvals, stat_type='sigma_clip'): + if stat_type == 'median': + xcen = np.median(xvals) + ycen = np.median(yvals) + elif stat_type == 'sigma_clip': + from esutil.stat import sigma_clip + xcen, xsig = sigma_clip(xvals) + ycen, ysig = sigma_clip(yvals) + else: + raise ValueError(f'Bad stat_type {stat_type}') + return xcen, ycen # Pick the right function to be _fix_seds. -if galsim.__version_info__ < (2,5): +if galsim.__version_info__ < (2, 5): LSST_SiliconBuilder._fix_seds = LSST_SiliconBuilder._fix_seds_24 else: LSST_SiliconBuilder._fix_seds = LSST_SiliconBuilder._fix_seds_25 From a1bf7e239b91c42e1d234c0b480dae4de8f8adb8 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Sun, 24 Mar 2024 12:09:49 -0400 Subject: [PATCH 04/11] quite --- imsim/stamp.py | 22 ++++++++++------------ 1 file changed, 10 insertions(+), 12 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 7f1a682f..c4ab479e 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -797,16 +797,16 @@ def draw(self, prof, image, method, offset, config, base, logger): xvals, yvals, stat_type='sigma_clip', ) - print( - 'shapeyx:', image.array.shape[0], image.array.shape[1], - 'true_cen:', image.true_center.y, image.true_center.x, - 'offset:', offset.y, offset.x, - 'cen:', ycen, xcen, - 'cenoff:', ycen-image.true_center.y, xcen-image.true_center.x, - 'maxoff:', - np.abs(yvals-image.true_center.y).max(), - np.abs(xcen-image.true_center.x).max(), - ) + # print( + # 'shapeyx:', image.array.shape[0], image.array.shape[1], + # 'true_cen:', image.true_center.y, image.true_center.x, + # 'offset:', offset.y, offset.x, + # 'cen:', ycen, xcen, + # 'cenoff:', ycen-image.true_center.y, xcen-image.true_center.x, + # 'maxoff:', + # np.abs(yvals-image.true_center.y).max(), + # np.abs(xcen-image.true_center.x).max(), + # ) # these can be saved in the config using @xcentroid, @ycentroid base['xcentroid'] = xcen @@ -849,8 +849,6 @@ def _get_photon_positions( rng=rng, ) - # print('imcen:', imcen, 'bounds:', image.bounds, 'off:', xoff, yoff) - photx = timage.photons.x photy = timage.photons.y From 036c4a66a2db96b9eabb832c169ef773600023a2 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Mon, 25 Mar 2024 11:30:47 -0400 Subject: [PATCH 05/11] more testing --- imsim/stamp.py | 20 +++++++++++++++++--- 1 file changed, 17 insertions(+), 3 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 7f1a682f..77915f26 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -805,7 +805,7 @@ def draw(self, prof, image, method, offset, config, base, logger): 'cenoff:', ycen-image.true_center.y, xcen-image.true_center.x, 'maxoff:', np.abs(yvals-image.true_center.y).max(), - np.abs(xcen-image.true_center.x).max(), + np.abs(xvals-image.true_center.x).max(), ) # these can be saved in the config using @xcentroid, @ycentroid @@ -857,8 +857,6 @@ def _get_photon_positions( logic = np.isnan(photx) | np.isnan(photy) wnan, = np.where(logic) wgood, = np.where(~logic) - if wnan.size > 0: - print(f'found {wnan.size} nan in photon positions') # location of true center, actually in big image imcen = image.true_center @@ -866,6 +864,22 @@ def _get_photon_positions( # ycen = imcen.y + timage.photons.y.mean() xvals = imcen.x + photx[wgood] yvals = imcen.y + photy[wgood] + + ymaxoff = np.abs(photy-image.true_center.y).max() + xmaxoff = np.abs(photx-image.true_center.x).max() + + if wnan.size > 0 or xmaxoff > 3_000_000 or ymaxoff > 3_000_000: + if wnan.size > 0: + print(f'found {wnan.size} nan in photon positions') + else: + print('maxoff:', ymaxoff, xmaxoff) + + print('gal') + print(repr(gal)) + # print('photon_ops') + # print(repr(photon_ops)) + # import IPython; IPython.embed() + return xvals, yvals From af8fa4ede0a9142b7fac864c3a53e705c4f94a21 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Tue, 26 Mar 2024 11:40:03 -0400 Subject: [PATCH 06/11] more checks on photons --- imsim/stamp.py | 62 ++++++++++++++++++++++++++++---------------------- 1 file changed, 35 insertions(+), 27 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 77915f26..f27ae18f 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -797,16 +797,17 @@ def draw(self, prof, image, method, offset, config, base, logger): xvals, yvals, stat_type='sigma_clip', ) - print( - 'shapeyx:', image.array.shape[0], image.array.shape[1], - 'true_cen:', image.true_center.y, image.true_center.x, - 'offset:', offset.y, offset.x, - 'cen:', ycen, xcen, - 'cenoff:', ycen-image.true_center.y, xcen-image.true_center.x, - 'maxoff:', - np.abs(yvals-image.true_center.y).max(), - np.abs(xvals-image.true_center.x).max(), - ) + # print( + # 'shapeyx:', image.array.shape[0], image.array.shape[1], + # 'true_cen:', image.true_center.y, image.true_center.x, + # 'offset:', offset.y, offset.x, + # 'cen:', ycen, xcen, + # 'cenoff:', + # ycen-image.true_center.y, xcen-image.true_center.x, + # 'maxoff:', + # np.abs(yvals-image.true_center.y).max(), + # np.abs(xvals-image.true_center.x).max(), + # ) # these can be saved in the config using @xcentroid, @ycentroid base['xcentroid'] = xcen @@ -854,31 +855,38 @@ def _get_photon_positions( photx = timage.photons.x photy = timage.photons.y + flux = timage.photons.flux + wvig, = np.where(flux == 0) + logic = np.isnan(photx) | np.isnan(photy) wnan, = np.where(logic) - wgood, = np.where(~logic) + + # flux == 0 means it was vignetted + wgood, = np.where( + np.isfinite(photx) + & np.isfinite(photy) + & (timage.photons.flux > 0) + ) # location of true center, actually in big image imcen = image.true_center - # xcen = imcen.x + timage.photons.x.mean() - # ycen = imcen.y + timage.photons.y.mean() xvals = imcen.x + photx[wgood] yvals = imcen.y + photy[wgood] - ymaxoff = np.abs(photy-image.true_center.y).max() - xmaxoff = np.abs(photx-image.true_center.x).max() - - if wnan.size > 0 or xmaxoff > 3_000_000 or ymaxoff > 3_000_000: - if wnan.size > 0: - print(f'found {wnan.size} nan in photon positions') - else: - print('maxoff:', ymaxoff, xmaxoff) - - print('gal') - print(repr(gal)) - # print('photon_ops') - # print(repr(photon_ops)) - # import IPython; IPython.embed() + # ymaxoff = np.abs(photy-image.true_center.y).max() + # xmaxoff = np.abs(photx-image.true_center.x).max() + + # if wnan.size > 0 or xmaxoff > 3_000_000 or ymaxoff > 3_000_000: + # if wnan.size > 0: + # print(f'found {wnan.size} nan in photon positions') + # else: + # print('maxoff:', ymaxoff, xmaxoff) + # + # # print('gal') + # # print(repr(gal)) + # # print('photon_ops') + # # print(repr(photon_ops)) + # # import IPython; IPython.embed() return xvals, yvals From 57472bb1adcfaf70142072937df431c244120f8a Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 11 Apr 2024 12:32:26 -0400 Subject: [PATCH 07/11] put robust stat code here was using esutil --- imsim/stamp.py | 90 ++++++++++++++++++++++++++++++++++++++++---------- 1 file changed, 73 insertions(+), 17 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index f46909ad..12991a56 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -793,13 +793,13 @@ def draw(self, prof, image, method, offset, config, base, logger): photon_ops=photon_ops, offset=offset, ) - xcen, ycen = _get_robust_centroids( - xvals, yvals, stat_type='sigma_clip', - ) + cenres = _get_robust_centroids(xvals, yvals) # these can be saved in the config using @xcentroid, @ycentroid - base['xcentroid'] = xcen - base['ycentroid'] = ycen + base['xcentroid'] = cenres['x'] + base['xcentroid_err'] = cenres['xerr'] + base['ycentroid'] = cenres['y'] + base['ycentroid_err'] = cenres['yerr'] return image @@ -813,6 +813,7 @@ def _get_photon_positions( sensor, photon_ops, offset, + n_photons=100_000, ): """ draw another image with fixed number of photons @@ -828,7 +829,7 @@ def _get_photon_positions( image=timage, method='phot', - n_photons=100_000, + n_photons=n_photons, sensor=sensor, photon_ops=photon_ops, poisson_flux=False, @@ -862,18 +863,73 @@ def _get_photon_positions( return xvals, yvals -def _get_robust_centroids(xvals, yvals, stat_type='sigma_clip'): - if stat_type == 'median': - xcen = np.median(xvals) - ycen = np.median(yvals) - elif stat_type == 'sigma_clip': - from esutil.stat import sigma_clip - xcen, xsig = sigma_clip(xvals) - ycen, ysig = sigma_clip(yvals) - else: - raise ValueError(f'Bad stat_type {stat_type}') +def _get_robust_centroids(xvals, yvals): + xcen, _, xerr = sigma_clip(xvals) + ycen, _, yerr = sigma_clip(yvals) + + return { + 'x': xcen, + 'xerr': xerr, + 'y': ycen, + 'yerr': yerr, + } + + +def sigma_clip(arrin, niter=4, nsig=4): + """ + Calculate the mean, sigma, error of an array with sigma clipping. + + parameters + ---------- + arr: array or sequence + A numpy array or sequence + niter: int, optional + number of iterations, defaults to 4 + nsig: float, optional + number of sigma, defaults to 4 + + returns + ------- + mean, stdev, err + """ + + arr = np.array(arrin, ndmin=1, copy=False) + + if len(arr.shape) > 1: + raise ValueError( + 'only 1-dimensional arrays suppored, got {arr.shape}' + ) + + indices = np.arange(arr.size) + nold = arr.size + + mn, sig, err = _get_sigma_clip_stats(arr) + + for i in range(1, niter + 1): + + w, = np.where((np.abs(arr[indices] - mn)) < nsig * sig) + + if w.size == 0: + # everything clipped, nothing to do but report latest + # statistics + break + + if w.size == nold: + break + + indices = indices[w] + nold = w.size + + mn, sig, err = _get_sigma_clip_stats(arr[indices]) + + return mn, sig, err + - return xcen, ycen +def _get_sigma_clip_stats(arr): + mn = arr.mean() + sig = arr.std() + err = sig / np.sqrt(arr.size) + return mn, sig, err # Pick the right function to be _fix_seds. From 4c1259a97d7f2d471d25b8cf165123c69ff53a07 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 11 Apr 2024 13:42:28 -0400 Subject: [PATCH 08/11] adding test but this and other tests fail because nothing is actually drawn --- tests/test_stamp.py | 47 ++++++++++++++++++++++++++++++++++++++++++--- 1 file changed, 44 insertions(+), 3 deletions(-) diff --git a/tests/test_stamp.py b/tests/test_stamp.py index f8c4d0ba..2bfe31e9 100644 --- a/tests/test_stamp.py +++ b/tests/test_stamp.py @@ -571,7 +571,48 @@ def new_toplevel_only(self, object_types): np.testing.assert_allclose(realized_X, predict, rtol=rtol) +def test_lsst_silicon_builder_calculates_centroid() -> None: + """ + LSST_SiliconBuilder.draw passes the correct list of photon_ops to + prof.drawImage. + + This fails because no flux is actually put in the image, but I'm + not sure why + """ + lsst_silicon = create_test_lsst_silicon(faint=False) + image = galsim.Image(ncol=256, nrow=256) + offset = galsim.PositionD(0,0) + config = create_test_config() + + prof = galsim.Gaussian(sigma=2) * galsim.SED('vega.txt', 'nm', 'flambda') + method = 'phot' + logger = galsim.config.LoggerWrapper(None) + + if galsim.__version_info__ < (2,5): + phot_type = 'galsim.ChromaticTransformation' + else: + # In GalSim 2.5+, this is now a SimpleChromaticTransformation + phot_type = 'galsim.SimpleChromaticTransformation' + + with mock.patch(phot_type+'.drawImage', return_value=image) as mock_drawImage: + image.added_flux = lsst_silicon.phot_flux + lsst_silicon.draw( + prof, + image, + method, + offset, + config=config["stamp"], + base=config, + logger=logger, + ) + for c in ('x', 'y'): + for n in ('', '_err'): + name = f'{c}centroid{n}' + assert name in base + + if __name__ == "__main__": - testfns = [v for k, v in vars().items() if k[:5] == 'test_' and callable(v)] - for testfn in testfns: - testfn() + test_lsst_silicon_builder_calculates_centroid() + # testfns = [v for k, v in vars().items() if k[:5] == 'test_' and callable(v)] + # for testfn in testfns: + # testfn() From 700a8d5d286011a647c72729b6ac438783efbc14 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 11 Apr 2024 13:56:11 -0400 Subject: [PATCH 09/11] make centroid calculation optional The only way is through the stamp sub-config --- imsim/stamp.py | 39 ++++++++++++++++++++++----------------- tests/test_stamp.py | 11 +++++++---- 2 files changed, 29 insertions(+), 21 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 12991a56..373f13f8 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -784,22 +784,24 @@ def draw(self, prof, image, method, offset, config, base, logger): poisson_flux=False) base['realized_flux'] = image.added_flux - xvals, yvals = _get_photon_positions( - gal=gal, - rng=self.rng, - bp_for_drawImage=bp_for_drawImage, - image=image, - sensor=sensor, - photon_ops=photon_ops, - offset=offset, - ) - cenres = _get_robust_centroids(xvals, yvals) + if 'centroid' in config: + xvals, yvals = _get_photon_positions( + gal=gal, + rng=self.rng, + bp_for_drawImage=bp_for_drawImage, + image=image, + sensor=sensor, + photon_ops=photon_ops, + offset=offset, + config=config['centroid'], + ) + cenres = _get_robust_centroids(xvals, yvals) - # these can be saved in the config using @xcentroid, @ycentroid - base['xcentroid'] = cenres['x'] - base['xcentroid_err'] = cenres['xerr'] - base['ycentroid'] = cenres['y'] - base['ycentroid_err'] = cenres['yerr'] + # these can be saved in the config using @xcentroid, @ycentroid + base['xcentroid'] = cenres['x'] + base['xcentroid_err'] = cenres['xerr'] + base['ycentroid'] = cenres['y'] + base['ycentroid_err'] = cenres['yerr'] return image @@ -813,7 +815,7 @@ def _get_photon_positions( sensor, photon_ops, offset, - n_photons=100_000, + config, ): """ draw another image with fixed number of photons @@ -823,13 +825,16 @@ def _get_photon_positions( 100_000 photons should give centroid to a few miliarcsec """ + if 'n_photons' not in config: + raise ValueError('n_photons must be in centroid config') + timage = image.copy() gal.drawImage( bp_for_drawImage, image=timage, method='phot', - n_photons=n_photons, + n_photons=config['n_photons'], sensor=sensor, photon_ops=photon_ops, poisson_flux=False, diff --git a/tests/test_stamp.py b/tests/test_stamp.py index 2bfe31e9..5be8a706 100644 --- a/tests/test_stamp.py +++ b/tests/test_stamp.py @@ -583,6 +583,9 @@ def test_lsst_silicon_builder_calculates_centroid() -> None: image = galsim.Image(ncol=256, nrow=256) offset = galsim.PositionD(0,0) config = create_test_config() + config['stamp']['centroid'] = { + 'n_photons': 1000, + } prof = galsim.Gaussian(sigma=2) * galsim.SED('vega.txt', 'nm', 'flambda') method = 'phot' @@ -612,7 +615,7 @@ def test_lsst_silicon_builder_calculates_centroid() -> None: if __name__ == "__main__": - test_lsst_silicon_builder_calculates_centroid() - # testfns = [v for k, v in vars().items() if k[:5] == 'test_' and callable(v)] - # for testfn in testfns: - # testfn() + # test_lsst_silicon_builder_calculates_centroid() + testfns = [v for k, v in vars().items() if k[:5] == 'test_' and callable(v)] + for testfn in testfns: + testfn() From ea1922de4bf8e02c4ccb3a14f3b209c8bc584d15 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Thu, 11 Apr 2024 14:07:50 -0400 Subject: [PATCH 10/11] move centroid config check to LSST_SiliconBuilder --- imsim/stamp.py | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/imsim/stamp.py b/imsim/stamp.py index 373f13f8..3df3fe79 100644 --- a/imsim/stamp.py +++ b/imsim/stamp.py @@ -92,7 +92,8 @@ def setup(self, config, base, xsize, ysize, ignore, logger): req = {} opt = {'camera': str, 'diffraction_fft': dict, - 'airmass': float, 'rawSeeing': float, 'band': str} + 'airmass': float, 'rawSeeing': float, 'band': str, + 'centroid': dict} if self.vignetting: req['det_name'] = str else: @@ -247,6 +248,9 @@ def setup(self, config, base, xsize, ysize, ignore, logger): else: world_pos = None + if 'centroid' in params and 'n_photons' not in params['centroid']: + raise ValueError('n_photons not found in centroid config') + return xsize, ysize, image_pos, world_pos def _getGoodPhotImageSize(self, obj_list, keep_sb_level, pixel_scale): @@ -825,9 +829,6 @@ def _get_photon_positions( 100_000 photons should give centroid to a few miliarcsec """ - if 'n_photons' not in config: - raise ValueError('n_photons must be in centroid config') - timage = image.copy() gal.drawImage( bp_for_drawImage, From f1c69356c3b06c78eebb1d97e364c52926880846 Mon Sep 17 00:00:00 2001 From: Erin Sheldon Date: Fri, 12 Apr 2024 11:43:23 -0400 Subject: [PATCH 11/11] fix bug where pupil or sampler is present There was no time sampler or pupil. --- tests/test_stamp.py | 56 ++++++++++++++++++++++----------------------- 1 file changed, 27 insertions(+), 29 deletions(-) diff --git a/tests/test_stamp.py b/tests/test_stamp.py index 5be8a706..cfffac53 100644 --- a/tests/test_stamp.py +++ b/tests/test_stamp.py @@ -571,51 +571,49 @@ def new_toplevel_only(self, object_types): np.testing.assert_allclose(realized_X, predict, rtol=rtol) -def test_lsst_silicon_builder_calculates_centroid() -> None: +def test_centroid_smoke() -> None: """ - LSST_SiliconBuilder.draw passes the correct list of photon_ops to - prof.drawImage. + Test that centroids are calculated. - This fails because no flux is actually put in the image, but I'm - not sure why + We are not checking accuracy here """ lsst_silicon = create_test_lsst_silicon(faint=False) image = galsim.Image(ncol=256, nrow=256) - offset = galsim.PositionD(0,0) + offset = galsim.PositionD(0, 0) config = create_test_config() config['stamp']['centroid'] = { 'n_photons': 1000, } + # order matters here + config['stamp']['photon_ops'] = [ + {'type': 'TimeSampler', 't0': 0.0, 'exptime': 30.0}, + {'type': 'PupilAnnulusSampler', 'R_outer': 4.18, 'R_inner': 2.55}, + ] + config['stamp']['photon_ops'] + prof = galsim.Gaussian(sigma=2) * galsim.SED('vega.txt', 'nm', 'flambda') method = 'phot' logger = galsim.config.LoggerWrapper(None) - if galsim.__version_info__ < (2,5): - phot_type = 'galsim.ChromaticTransformation' - else: - # In GalSim 2.5+, this is now a SimpleChromaticTransformation - phot_type = 'galsim.SimpleChromaticTransformation' - - with mock.patch(phot_type+'.drawImage', return_value=image) as mock_drawImage: - image.added_flux = lsst_silicon.phot_flux - lsst_silicon.draw( - prof, - image, - method, - offset, - config=config["stamp"], - base=config, - logger=logger, - ) - for c in ('x', 'y'): - for n in ('', '_err'): - name = f'{c}centroid{n}' - assert name in base + image.added_flux = lsst_silicon.phot_flux + lsst_silicon.draw( + prof, + image, + method, + offset, + config=config["stamp"], + base=config, + logger=logger, + ) + for c in ('x', 'y'): + for n in ('', '_err'): + name = f'{c}centroid{n}' + assert name in config if __name__ == "__main__": - # test_lsst_silicon_builder_calculates_centroid() - testfns = [v for k, v in vars().items() if k[:5] == 'test_' and callable(v)] + testfns = [ + v for k, v in vars().items() if k[:5] == 'test_' and callable(v) + ] for testfn in testfns: testfn()