lsst.ip.diffim  18.1.0-6-gdda7f3e+7
diffimTools.py
Go to the documentation of this file.
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/>.
21 
22 __all__ = ["backgroundSubtract", "writeKernelCellSet", "sourceToFootprintList", "NbasisEvaluator"]
23 
24 # python
25 import time
26 import os
27 from collections import Counter
28 import numpy as np
29 
30 # all the c++ level classes and routines
31 from . import diffimLib
32 
33 # all the other LSST packages
34 import lsst.afw.geom as afwGeom
35 import lsst.afw.image as afwImage
36 import lsst.afw.table as afwTable
37 import lsst.afw.detection as afwDetect
38 import lsst.afw.math.mathLib as afwMath
39 import lsst.geom as geom
40 from lsst.log import Log
41 import lsst.pex.config as pexConfig
42 from .makeKernelBasisList import makeKernelBasisList
43 
44 # Helper functions for ipDiffim; mostly viewing of results and writing
45 # debugging info to disk.
46 
47 
50 
51 
52 def makeFlatNoiseImage(mi, seedStat=afwMath.MAX):
53  img = mi.getImage()
54  seed = int(10.*afwMath.makeStatistics(mi.getImage(), seedStat).getValue() + 1)
55  rdm = afwMath.Random(afwMath.Random.MT19937, seed)
56  rdmImage = img.Factory(img.getDimensions())
57  afwMath.randomGaussianImage(rdmImage, rdm)
58  return rdmImage
59 
60 
62  """Return a Poisson noise image based on im
63 
64  Parameters
65  ----------
66  im : `lsst.afw.image.Image`
67  image; the output image has the same dtype, dimensions, and shape
68  and its expectation value is the value of ``im`` at each pixel
69 
70  Returns
71  -------
72  noiseIm : `lsst.afw.image.Image`
73  Newly constructed image instance, same type as ``im``.
74 
75  Notes
76  -----
77  - Warning: This uses an undocumented numpy API (the documented API
78  uses a single float expectation value instead of an array).
79 
80  - Uses numpy.random; you may wish to call numpy.random.seed first.
81  """
82  import numpy.random as rand
83  imArr = im.getArray()
84  noiseIm = im.Factory(im.getBBox())
85  noiseArr = noiseIm.getArray()
86 
87  with np.errstate(invalid='ignore'):
88  intNoiseArr = rand.poisson(imArr)
89 
90  noiseArr[:, :] = intNoiseArr.astype(noiseArr.dtype)
91  return noiseIm
92 
93 
98 
99 
101  kCoeffs = ((1.0, 0.0, 0.0),
102  (0.005, -0.000001, 0.000001),
103  (0.005, 0.000004, 0.000004),
104  (-0.001, -0.000030, 0.000030),
105  (-0.001, 0.000015, 0.000015),
106  (-0.005, -0.000050, 0.000050))
107  return kCoeffs
108 
109 
110 def makeFakeKernelSet(sizeCell=128, nCell=3,
111  deltaFunctionCounts=1.e4, tGaussianWidth=1.0,
112  addNoise=True, bgValue=100., display=False):
113  """Generate test template and science images with sources.
114 
115  Parameters
116  ----------
117  sizeCell : `int`, optional
118  Size of the square spatial cells in pixels.
119  nCell : `int`, optional
120  Number of adjacent spatial cells in both direction in both images.
121  deltaFunctionCounts : `float`, optional
122  Flux value for the template image sources.
123  tGaussianWidth : `float`, optional
124  Sigma of the generated Gaussian PSF sources in the template image.
125  addNoise : `bool`, optional
126  If `True`, Poisson noise is added to both the generated template
127  and science images.
128  bgValue : `float`, optional
129  Background level to be added to the generated science image.
130  display : `bool`, optional
131  If `True` displays the generated template and science images by
132  `lsst.afw.display.Display`.
133 
134  Notes
135  -----
136  - The generated images consist of adjacent ``nCell x nCell`` cells, each
137  of pixel size ``sizeCell x sizeCell``.
138  - The sources in the science image are generated by convolving the
139  template by ``sKernel``. ``sKernel`` is a spatial `LinearCombinationKernel`
140  of hard wired kernel bases functions. The linear combination has first
141  order polynomial spatial dependence with polynomial parameters from ``fakeCoeffs()``.
142  - The template image sources are generated in the center of each spatial
143  cell from one pixel, set to `deltaFunctionCounts` counts, then convolved
144  by a 2D Gaussian with sigma of `tGaussianWidth` along each axis.
145  - The sources are also returned in ``kernelCellSet`` each source is "detected"
146  exactly at the center of a cell.
147 
148  Returns
149  -------
150  tMi : `lsst.afw.image.MaskedImage`
151  Generated template image.
152  sMi : `lsst.afw.image.MaskedImage`
153  Generated science image.
154  sKernel : `lsst.afw.math.LinearCombinationKernel`
155  The spatial kernel used to generate the sources in the science image.
156  kernelCellSet : `lsst.afw.math.SpatialCellSet`
157  Cell grid of `lsst.afw.math.SpatialCell` instances, containing
158  `lsst.ip.diffim.KernelCandidate` instances around all the generated sources
159  in the science image.
160  configFake : `lsst.ip.diffim.ImagePsfMatchConfig`
161  Config instance used in the image generation.
162  """
163  from . import imagePsfMatch
164  configFake = imagePsfMatch.ImagePsfMatchConfig()
165  configFake.kernel.name = "AL"
166  subconfigFake = configFake.kernel.active
167  subconfigFake.alardNGauss = 1
168  subconfigFake.alardSigGauss = [2.5, ]
169  subconfigFake.alardDegGauss = [2, ]
170  subconfigFake.sizeCellX = sizeCell
171  subconfigFake.sizeCellY = sizeCell
172  subconfigFake.spatialKernelOrder = 1
173  subconfigFake.spatialModelType = "polynomial"
174  subconfigFake.singleKernelClipping = False # variance is a hack
175  subconfigFake.spatialKernelClipping = False # variance is a hack
176  if bgValue > 0.0:
177  subconfigFake.fitForBackground = True
178 
179  policyFake = pexConfig.makePolicy(subconfigFake)
180 
181  basisList = makeKernelBasisList(subconfigFake)
182  kSize = subconfigFake.kernelSize
183 
184  # This sets the final extent of each convolved delta function
185  gaussKernelWidth = sizeCell//2
186 
187  # This sets the scale over which pixels are correlated in the
188  # spatial convolution; should be at least as big as the kernel you
189  # are trying to fit for
190  spatialKernelWidth = kSize
191 
192  # Number of bad pixels due to convolutions
193  border = (gaussKernelWidth + spatialKernelWidth)//2
194 
195  # Make a fake image with a matrix of delta functions
196  totalSize = nCell*sizeCell + 2*border
197  tim = afwImage.ImageF(geom.Extent2I(totalSize, totalSize))
198  for x in range(nCell):
199  for y in range(nCell):
200  tim[x*sizeCell + sizeCell//2 + border - 1,
201  y*sizeCell + sizeCell//2 + border - 1,
202  afwImage.LOCAL] = deltaFunctionCounts
203 
204  # Turn this into stars with a narrow width; conserve counts
205  gaussFunction = afwMath.GaussianFunction2D(tGaussianWidth, tGaussianWidth)
206  gaussKernel = afwMath.AnalyticKernel(gaussKernelWidth, gaussKernelWidth, gaussFunction)
207  cim = afwImage.ImageF(tim.getDimensions())
208  afwMath.convolve(cim, tim, gaussKernel, True)
209  tim = cim
210 
211  # Trim off border pixels
212  bbox = gaussKernel.shrinkBBox(tim.getBBox(afwImage.LOCAL))
213  tim = afwImage.ImageF(tim, bbox, afwImage.LOCAL)
214 
215  # Now make a science image which is this convolved with some
216  # spatial function. Use input basis list.
217  polyFunc = afwMath.PolynomialFunction2D(1)
218  kCoeffs = fakeCoeffs()
219  nToUse = min(len(kCoeffs), len(basisList))
220 
221  # Make the full convolved science image
222  sKernel = afwMath.LinearCombinationKernel(basisList[:nToUse], polyFunc)
223  sKernel.setSpatialParameters(kCoeffs[:nToUse])
224  sim = afwImage.ImageF(tim.getDimensions())
225  afwMath.convolve(sim, tim, sKernel, True)
226 
227  # Get the good subregion
228  bbox = sKernel.shrinkBBox(sim.getBBox(afwImage.LOCAL))
229 
230  # Add background
231  sim += bgValue
232 
233  # Watch out for negative values
234  tim += 2*np.abs(np.min(tim.getArray()))
235 
236  # Add noise?
237  if addNoise:
238  sim = makePoissonNoiseImage(sim)
239  tim = makePoissonNoiseImage(tim)
240 
241  # And turn into MaskedImages
242  sim = afwImage.ImageF(sim, bbox, afwImage.LOCAL)
243  svar = afwImage.ImageF(sim, True)
244  smask = afwImage.Mask(sim.getDimensions())
245  smask.set(0x0)
246  sMi = afwImage.MaskedImageF(sim, smask, svar)
247 
248  tim = afwImage.ImageF(tim, bbox, afwImage.LOCAL)
249  tvar = afwImage.ImageF(tim, True)
250  tmask = afwImage.Mask(tim.getDimensions())
251  tmask.set(0x0)
252  tMi = afwImage.MaskedImageF(tim, tmask, tvar)
253 
254  if display:
255  import lsst.afw.display as afwDisplay
256  afwDisplay.Display(frame=1).mtv(tMi)
257  afwDisplay.Display(frame=2).mtv(sMi)
258 
259  # Finally, make a kernelSet from these 2 images
260  kernelCellSet = afwMath.SpatialCellSet(geom.Box2I(geom.Point2I(0, 0),
261  geom.Extent2I(sizeCell*nCell,
262  sizeCell*nCell)),
263  sizeCell,
264  sizeCell)
265  stampHalfWidth = 2*kSize
266  for x in range(nCell):
267  for y in range(nCell):
268  xCoord = x*sizeCell + sizeCell//2
269  yCoord = y*sizeCell + sizeCell//2
270  p0 = geom.Point2I(xCoord - stampHalfWidth,
271  yCoord - stampHalfWidth)
272  p1 = geom.Point2I(xCoord + stampHalfWidth,
273  yCoord + stampHalfWidth)
274  bbox = geom.Box2I(p0, p1)
275  tsi = afwImage.MaskedImageF(tMi, bbox, origin=afwImage.LOCAL)
276  ssi = afwImage.MaskedImageF(sMi, bbox, origin=afwImage.LOCAL)
277 
278  kc = diffimLib.makeKernelCandidate(xCoord, yCoord, tsi, ssi, policyFake)
279  kernelCellSet.insertCandidate(kc)
280 
281  tMi.setXY0(0, 0)
282  sMi.setXY0(0, 0)
283  return tMi, sMi, sKernel, kernelCellSet, configFake
284 
285 
286 
289 
290 def backgroundSubtract(config, maskedImages):
291  """Subtract the background from masked images.
292 
293  Parameters
294  ----------
295  config : TODO: DM-17458
296  TODO: DM-17458
297  maskedImages : `list` of `lsst.afw.image.MaskedImage`
298  TODO: DM-17458
299 
300  Returns
301  -------
302  TODO: DM-17458
303  TODO: DM-17458
304  """
305  backgrounds = []
306  t0 = time.time()
307  algorithm = config.algorithm
308  binsize = config.binSize
309  undersample = config.undersampleStyle
310  bctrl = afwMath.BackgroundControl(algorithm)
311  bctrl.setUndersampleStyle(undersample)
312  for maskedImage in maskedImages:
313  bctrl.setNxSample(maskedImage.getWidth()//binsize + 1)
314  bctrl.setNySample(maskedImage.getHeight()//binsize + 1)
315  image = maskedImage.getImage()
316  backobj = afwMath.makeBackground(image, bctrl)
317 
318  image -= backobj.getImageF()
319  backgrounds.append(backobj.getImageF())
320  del backobj
321 
322  t1 = time.time()
323  logger = Log.getLogger("ip.diffim.backgroundSubtract")
324  logger.debug("Total time for background subtraction : %.2f s", (t1 - t0))
325  return backgrounds
326 
327 
330 
331 
332 def writeKernelCellSet(kernelCellSet, psfMatchingKernel, backgroundModel, outdir):
333  """TODO: DM-17458
334 
335  Parameters
336  ----------
337  kernelCellSet : TODO: DM-17458
338  TODO: DM-17458
339  psfMatchingKernel : TODO: DM-17458
340  TODO: DM-17458
341  backgroundModel : TODO: DM-17458
342  TODO: DM-17458
343  outdir : TODO: DM-17458
344  TODO: DM-17458
345  """
346  if not os.path.isdir(outdir):
347  os.makedirs(outdir)
348 
349  for cell in kernelCellSet.getCellList():
350  for cand in cell.begin(False): # False = include bad candidates
351  if cand.getStatus() == afwMath.SpatialCellCandidate.GOOD:
352  xCand = int(cand.getXCenter())
353  yCand = int(cand.getYCenter())
354  idCand = cand.getId()
355  diffIm = cand.getDifferenceImage(diffimLib.KernelCandidateF.ORIG)
356  kernel = cand.getKernelImage(diffimLib.KernelCandidateF.ORIG)
357  diffIm.writeFits(os.path.join(outdir, 'diffim_c%d_x%d_y%d.fits' % (idCand, xCand, yCand)))
358  kernel.writeFits(os.path.join(outdir, 'kernel_c%d_x%d_y%d.fits' % (idCand, xCand, yCand)))
359 
360  # Diffim from spatial model
361  ski = afwImage.ImageD(kernel.getDimensions())
362  psfMatchingKernel.computeImage(ski, False, xCand, yCand)
363  sk = afwMath.FixedKernel(ski)
364  sbg = backgroundModel(xCand, yCand)
365  sdmi = cand.getDifferenceImage(sk, sbg)
366  sdmi.writeFits(os.path.join(outdir, 'sdiffim_c%d_x%d_y%d.fits' % (idCand, xCand, yCand)))
367 
368 
371 
372 
373 def sourceToFootprintList(candidateInList, templateExposure, scienceExposure, kernelSize, config, log):
374  """Convert a list of sources for the PSF-matching Kernel to Footprints.
375 
376  Parameters
377  ----------
378  candidateInList : TODO: DM-17458
379  Input list of Sources
380  templateExposure : TODO: DM-17458
381  Template image, to be checked for Mask bits in Source Footprint
382  scienceExposure : TODO: DM-17458
383  Science image, to be checked for Mask bits in Source Footprint
384  kernelSize : TODO: DM-17458
385  TODO: DM-17458
386  config : TODO: DM-17458
387  Config that defines the Mask planes that indicate an invalid Source and Bbox grow radius
388  log : TODO: DM-17458
389  Log for output
390 
391  Returns
392  -------
393  candidateOutList : `list`
394  a list of dicts having a "source" and "footprint" field, to be used for Psf-matching
395 
396  Raises
397  ------
398  RuntimeError
399  TODO: DM-17458
400 
401  Notes
402  -----
403  Takes an input list of Sources that were selected to constrain
404  the Psf-matching Kernel and turns them into a List of Footprints,
405  which are used to seed a set of KernelCandidates. The function
406  checks both the template and science image for masked pixels,
407  rejecting the Source if certain Mask bits (defined in config) are
408  set within the Footprint.
409  """
410 
411  candidateOutList = []
412  fsb = diffimLib.FindSetBitsU()
413  badBitMask = 0
414  for mp in config.badMaskPlanes:
415  badBitMask |= afwImage.Mask.getPlaneBitMask(mp)
416  bbox = scienceExposure.getBBox()
417 
418  # Size to grow Sources
419  if config.scaleByFwhm:
420  fpGrowPix = int(config.fpGrowKernelScaling*kernelSize + 0.5)
421  else:
422  fpGrowPix = config.fpGrowPix
423  log.info("Growing %d kernel candidate stars by %d pixels", len(candidateInList), fpGrowPix)
424 
425  for kernelCandidate in candidateInList:
426  if not type(kernelCandidate) == afwTable.SourceRecord:
427  raise RuntimeError("Candiate not of type afwTable.SourceRecord")
428  bm1 = 0
429  bm2 = 0
430  center = geom.Point2I(scienceExposure.getWcs().skyToPixel(kernelCandidate.getCoord()))
431  if center[0] < bbox.getMinX() or center[0] > bbox.getMaxX():
432  continue
433  if center[1] < bbox.getMinY() or center[1] > bbox.getMaxY():
434  continue
435 
436  xmin = center[0] - fpGrowPix
437  xmax = center[0] + fpGrowPix
438  ymin = center[1] - fpGrowPix
439  ymax = center[1] + fpGrowPix
440 
441  # Keep object centered
442  if (xmin - bbox.getMinX()) < 0:
443  xmax += (xmin - bbox.getMinX())
444  xmin -= (xmin - bbox.getMinX())
445  if (ymin - bbox.getMinY()) < 0:
446  ymax += (ymin - bbox.getMinY())
447  ymin -= (ymin - bbox.getMinY())
448  if (bbox.getMaxX() - xmax) < 0:
449  xmin -= (bbox.getMaxX() - xmax)
450  xmax += (bbox.getMaxX() - xmax)
451  if (bbox.getMaxY() - ymax) < 0:
452  ymin -= (bbox.getMaxY() - ymax)
453  ymax += (bbox.getMaxY() - ymax)
454  if xmin > xmax or ymin > ymax:
455  continue
456 
457  kbbox = geom.Box2I(geom.Point2I(xmin, ymin), geom.Point2I(xmax, ymax))
458  try:
459  fsb.apply(afwImage.MaskedImageF(templateExposure.getMaskedImage(), kbbox, deep=False).getMask())
460  bm1 = fsb.getBits()
461  fsb.apply(afwImage.MaskedImageF(scienceExposure.getMaskedImage(), kbbox, deep=False).getMask())
462  bm2 = fsb.getBits()
463  except Exception:
464  pass
465  else:
466  if not((bm1 & badBitMask) or (bm2 & badBitMask)):
467  candidateOutList.append({'source': kernelCandidate,
468  'footprint': afwDetect.Footprint(afwGeom.SpanSet(kbbox))})
469  log.info("Selected %d / %d sources for KernelCandidacy", len(candidateOutList), len(candidateInList))
470  return candidateOutList
471 
472 
473 def sourceTableToCandidateList(sourceTable, templateExposure, scienceExposure, kConfig, dConfig, log,
474  basisList, doBuild=False):
475  """Convert a list of Sources into KernelCandidates.
476 
477  The KernelCandidates are used for fitting the Psf-matching kernel.
478 
479  Parameters
480  ----------
481  sourceTable : TODO: DM-17458
482  TODO: DM-17458
483  templateExposure : TODO: DM-17458
484  TODO: DM-17458
485  scienceExposure : TODO: DM-17458
486  TODO: DM-17458
487  kConfig : TODO: DM-17458
488  TODO: DM-17458
489  dConfig : TODO: DM-17458
490  TODO: DM-17458
491  log : TODO: DM-17458
492  TODO: DM-17458
493  basisList : TODO: DM-17458
494  TODO: DM-17458
495  doBuild : `bool`, optional
496  TODO: DM-17458
497 
498  Returns
499  -------
500  TODO: DM-17458
501  TODO: DM-17458
502  """
503  kernelSize = basisList[0].getWidth()
504  footprintList = sourceToFootprintList(list(sourceTable), templateExposure, scienceExposure,
505  kernelSize, dConfig, log)
506  candList = []
507 
508  if doBuild and not basisList:
509  doBuild = False
510  else:
511  policy = pexConfig.makePolicy(kConfig)
512  visitor = diffimLib.BuildSingleKernelVisitorF(basisList, policy)
513 
514  policy = pexConfig.makePolicy(kConfig)
515  for cand in footprintList:
516  bbox = cand['footprint'].getBBox()
517  tmi = afwImage.MaskedImageF(templateExposure.getMaskedImage(), bbox)
518  smi = afwImage.MaskedImageF(scienceExposure.getMaskedImage(), bbox)
519  kCand = diffimLib.makeKernelCandidate(cand['source'], tmi, smi, policy)
520  if doBuild:
521  visitor.processCandidate(kCand)
522  kCand.setStatus(afwMath.SpatialCellCandidate.UNKNOWN)
523  candList.append(kCand)
524  return candList
525 
526 
527 
530 
531 
532 class NbasisEvaluator(object):
533  """A functor to evaluate the Bayesian Information Criterion for the number of basis sets
534  going into the kernel fitting"""
535 
536  def __init__(self, psfMatchConfig, psfFwhmPixTc, psfFwhmPixTnc):
537  self.psfMatchConfig = psfMatchConfig
538  self.psfFwhmPixTc = psfFwhmPixTc
539  self.psfFwhmPixTnc = psfFwhmPixTnc
540  if not self.psfMatchConfig.kernelBasisSet == "alard-lupton":
541  raise RuntimeError("BIC only implemnted for AL (alard lupton) basis")
542 
543  def __call__(self, kernelCellSet, log):
544  d1, d2, d3 = self.psfMatchConfig.alardDegGauss
545  bicArray = {}
546  for d1i in range(1, d1 + 1):
547  for d2i in range(1, d2 + 1):
548  for d3i in range(1, d3 + 1):
549  dList = [d1i, d2i, d3i]
550  bicConfig = type(self.psfMatchConfig)(self.psfMatchConfig, alardDegGauss=dList)
551  kList = makeKernelBasisList(bicConfig, self.psfFwhmPixTc, self.psfFwhmPixTnc)
552  k = len(kList)
553  visitor = diffimLib.BuildSingleKernelVisitorF(kList, pexConfig.makePolicy(bicConfig))
554  visitor.setSkipBuilt(False)
555  kernelCellSet.visitCandidates(visitor, bicConfig.nStarPerCell)
556 
557  for cell in kernelCellSet.getCellList():
558  for cand in cell.begin(False): # False = include bad candidates
559  if cand.getStatus() != afwMath.SpatialCellCandidate.GOOD:
560  continue
561  diffIm = cand.getDifferenceImage(diffimLib.KernelCandidateF.RECENT)
562  bbox = cand.getKernel(diffimLib.KernelCandidateF.RECENT).shrinkBBox(
563  diffIm.getBBox(afwImage.LOCAL))
564  diffIm = type(diffIm)(diffIm, bbox, True)
565  chi2 = diffIm.getImage().getArray()**2/diffIm.getVariance().getArray()
566  n = chi2.shape[0]*chi2.shape[1]
567  bic = np.sum(chi2) + k*np.log(n)
568  if cand.getId() not in bicArray:
569  bicArray[cand.getId()] = {}
570  bicArray[cand.getId()][(d1i, d2i, d3i)] = bic
571 
572  bestConfigs = []
573  for candId in bicArray:
574  cconfig, cvals = list(bicArray[candId].keys()), list(bicArray[candId].values())
575  idx = np.argsort(cvals)
576  bestConfig = cconfig[idx[0]]
577  bestConfigs.append(bestConfig)
578 
579  counter = Counter(bestConfigs).most_common(3)
580  log.info("B.I.C. prefers basis complexity %s %d times; %s %d times; %s %d times",
581  counter[0][0], counter[0][1],
582  counter[1][0], counter[1][1],
583  counter[2][0], counter[2][1])
584  return counter[0][0], counter[1][0], counter[2][0]
def makeFlatNoiseImage(mi, seedStat=afwMath.MAX)
Add noise.
Definition: diffimTools.py:52
def makeKernelBasisList(config, targetFwhmPix=None, referenceFwhmPix=None, basisDegGauss=None, metadata=None)
def fakeCoeffs()
Make fake images for testing; one is a delta function (or narrow gaussian) and the other is a convolu...
Definition: diffimTools.py:100
std::shared_ptr< Background > makeBackground(ImageT const &img, BackgroundControl const &bgCtrl)
def backgroundSubtract(config, maskedImages)
Background subtraction for ip_diffim.
Definition: diffimTools.py:290
def makeFakeKernelSet(sizeCell=128, nCell=3, deltaFunctionCounts=1.e4, tGaussianWidth=1.0, addNoise=True, bgValue=100., display=False)
Definition: diffimTools.py:112
Statistics makeStatistics(lsst::afw::math::MaskedVector< EntryT > const &mv, std::vector< WeightPixel > const &vweights, int const flags, StatisticsControl const &sctrl=StatisticsControl())
def writeKernelCellSet(kernelCellSet, psfMatchingKernel, backgroundModel, outdir)
More coarse debugging.
Definition: diffimTools.py:332
void randomGaussianImage(ImageT *image, Random &rand)
def sourceTableToCandidateList(sourceTable, templateExposure, scienceExposure, kConfig, dConfig, log, basisList, doBuild=False)
Definition: diffimTools.py:474
def __init__(self, psfMatchConfig, psfFwhmPixTc, psfFwhmPixTnc)
Definition: diffimTools.py:536
void convolve(OutImageT &convolvedImage, InImageT const &inImage, KernelT const &kernel, bool doNormalize, bool doCopyEdge=false)
def __call__(self, kernelCellSet, log)
Definition: diffimTools.py:543
def sourceToFootprintList(candidateInList, templateExposure, scienceExposure, kernelSize, config, log)
Converting types.
Definition: diffimTools.py:373