LSSTApplications  19.0.0-10-g920eed2,19.0.0-11-g48a0200+2,19.0.0-18-gfc4e62b+13,19.0.0-2-g3b2f90d+2,19.0.0-2-gd671419+5,19.0.0-20-g5a5a17ab+11,19.0.0-21-g2644856+13,19.0.0-23-g84eeccb+1,19.0.0-24-g878c510+1,19.0.0-25-g6c8df7140,19.0.0-25-gb330496+1,19.0.0-3-g2b32d65+5,19.0.0-3-g8227491+12,19.0.0-3-g9c54d0d+12,19.0.0-3-gca68e65+8,19.0.0-3-gcfc5f51+5,19.0.0-3-ge110943+11,19.0.0-3-ge74d124,19.0.0-3-gfe04aa6+13,19.0.0-30-g9c3fd16+1,19.0.0-4-g06f5963+5,19.0.0-4-g3d16501+13,19.0.0-4-g4a9c019+5,19.0.0-4-g5a8b323,19.0.0-4-g66397f0+1,19.0.0-4-g8278b9b+1,19.0.0-4-g8557e14,19.0.0-4-g8964aba+13,19.0.0-4-ge404a01+12,19.0.0-5-g40f3a5a,19.0.0-5-g4db63b3,19.0.0-5-gfb03ce7+13,19.0.0-6-gbaebbfb+12,19.0.0-61-gec4c6e08+1,19.0.0-7-g039c0b5+11,19.0.0-7-gbea9075+4,19.0.0-7-gc567de5+13,19.0.0-71-g41c0270,19.0.0-9-g2f02add+1,19.0.0-9-g463f923+12,w.2020.22
LSSTDataManagementBasePackage
utils.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 """Support utilities for Measuring sources"""
23 
24 # Export DipoleTestImage to expose fake image generating funcs
25 __all__ = ["DipoleTestImage"]
26 
27 import numpy as np
28 
29 import lsst.geom as geom
30 import lsst.afw.detection as afwDet
31 import lsst.afw.display as afwDisplay
32 import lsst.afw.geom as afwGeom
33 import lsst.afw.image as afwImage
34 import lsst.afw.math as afwMath
35 import lsst.afw.table as afwTable
36 from lsst.log import Log
37 import lsst.meas.algorithms as measAlg
38 import lsst.meas.base as measBase
39 from .dipoleFitTask import DipoleFitAlgorithm
40 from . import diffimLib
41 from . import diffimTools
42 
43 afwDisplay.setDefaultMaskTransparency(75)
44 keptPlots = False # Have we arranged to keep spatial plots open?
45 
46 
47 def showSourceSet(sSet, xy0=(0, 0), frame=0, ctype=afwDisplay.GREEN, symb="+", size=2):
48  """Draw the (XAstrom, YAstrom) positions of a set of Sources.
49 
50  Image has the given XY0.
51  """
52  disp = afwDisplay.afwDisplay(frame=frame)
53  with disp.Buffering():
54  for s in sSet:
55  xc, yc = s.getXAstrom() - xy0[0], s.getYAstrom() - xy0[1]
56 
57  if symb == "id":
58  disp.dot(str(s.getId()), xc, yc, ctype=ctype, size=size)
59  else:
60  disp.dot(symb, xc, yc, ctype=ctype, size=size)
61 
62 
63 # Kernel display utilities
64 #
65 
66 
67 def showKernelSpatialCells(maskedIm, kernelCellSet, showChi2=False, symb="o",
68  ctype=None, ctypeUnused=None, ctypeBad=None, size=3,
69  frame=None, title="Spatial Cells"):
70  """Show the SpatialCells.
71 
72  If symb is something that display.dot understands (e.g. "o"), the top
73  nMaxPerCell candidates will be indicated with that symbol, using ctype
74  and size.
75  """
76  disp = afwDisplay.Display(frame=frame)
77  disp.mtv(maskedIm, title=title)
78  with disp.Buffering():
79  origin = [-maskedIm.getX0(), -maskedIm.getY0()]
80  for cell in kernelCellSet.getCellList():
81  afwDisplay.utils.drawBBox(cell.getBBox(), origin=origin, display=disp)
82 
83  goodies = ctypeBad is None
84  for cand in cell.begin(goodies):
85  xc, yc = cand.getXCenter() + origin[0], cand.getYCenter() + origin[1]
86  if cand.getStatus() == afwMath.SpatialCellCandidate.BAD:
87  color = ctypeBad
88  elif cand.getStatus() == afwMath.SpatialCellCandidate.GOOD:
89  color = ctype
90  elif cand.getStatus() == afwMath.SpatialCellCandidate.UNKNOWN:
91  color = ctypeUnused
92  else:
93  continue
94 
95  if color:
96  disp.dot(symb, xc, yc, ctype=color, size=size)
97 
98  if showChi2:
99  rchi2 = cand.getChi2()
100  if rchi2 > 1e100:
101  rchi2 = np.nan
102  disp.dot("%d %.1f" % (cand.getId(), rchi2),
103  xc - size, yc - size - 4, ctype=color, size=size)
104 
105 
106 def showDiaSources(sources, exposure, isFlagged, isDipole, frame=None):
107  """Display Dia Sources.
108  """
109  #
110  # Show us the ccandidates
111  #
112  # Too many mask planes in diffims
113  disp = afwDisplay.Display(frame=frame)
114  for plane in ("BAD", "CR", "EDGE", "INTERPOlATED", "INTRP", "SAT", "SATURATED"):
115  disp.setMaskPlaneColor(plane, color="ignore")
116 
117  mos = afwDisplay.utils.Mosaic()
118  for i in range(len(sources)):
119  source = sources[i]
120  badFlag = isFlagged[i]
121  dipoleFlag = isDipole[i]
122  bbox = source.getFootprint().getBBox()
123  stamp = exposure.Factory(exposure, bbox, True)
124  im = afwDisplay.utils.Mosaic(gutter=1, background=0, mode="x")
125  im.append(stamp.getMaskedImage())
126  lab = "%.1f,%.1f:" % (source.getX(), source.getY())
127  if badFlag:
128  ctype = afwDisplay.RED
129  lab += "BAD"
130  if dipoleFlag:
131  ctype = afwDisplay.YELLOW
132  lab += "DIPOLE"
133  if not badFlag and not dipoleFlag:
134  ctype = afwDisplay.GREEN
135  lab += "OK"
136  mos.append(im.makeMosaic(), lab, ctype)
137  title = "Dia Sources"
138  mosaicImage = mos.makeMosaic(display=disp, title=title)
139  return mosaicImage
140 
141 
142 def showKernelCandidates(kernelCellSet, kernel, background, frame=None, showBadCandidates=True,
143  resids=False, kernels=False):
144  """Display the Kernel candidates.
145 
146  If kernel is provided include spatial model and residuals;
147  If chi is True, generate a plot of residuals/sqrt(variance), i.e. chi.
148  """
149  #
150  # Show us the ccandidates
151  #
152  if kernels:
153  mos = afwDisplay.utils.Mosaic(gutter=5, background=0)
154  else:
155  mos = afwDisplay.utils.Mosaic(gutter=5, background=-1)
156  #
157  candidateCenters = []
158  candidateCentersBad = []
159  candidateIndex = 0
160  for cell in kernelCellSet.getCellList():
161  for cand in cell.begin(False): # include bad candidates
162  # Original difference image; if does not exist, skip candidate
163  try:
164  resid = cand.getDifferenceImage(diffimLib.KernelCandidateF.ORIG)
165  except Exception:
166  continue
167 
168  rchi2 = cand.getChi2()
169  if rchi2 > 1e100:
170  rchi2 = np.nan
171 
172  if not showBadCandidates and cand.isBad():
173  continue
174 
175  im_resid = afwDisplay.utils.Mosaic(gutter=1, background=-0.5, mode="x")
176 
177  try:
178  im = cand.getScienceMaskedImage()
179  im = im.Factory(im, True)
180  im.setXY0(cand.getScienceMaskedImage().getXY0())
181  except Exception:
182  continue
183  if (not resids and not kernels):
184  im_resid.append(im.Factory(im, True))
185  try:
186  im = cand.getTemplateMaskedImage()
187  im = im.Factory(im, True)
188  im.setXY0(cand.getTemplateMaskedImage().getXY0())
189  except Exception:
190  continue
191  if (not resids and not kernels):
192  im_resid.append(im.Factory(im, True))
193 
194  # Difference image with original basis
195  if resids:
196  var = resid.getVariance()
197  var = var.Factory(var, True)
198  np.sqrt(var.getArray(), var.getArray()) # inplace sqrt
199  resid = resid.getImage()
200  resid /= var
201  bbox = kernel.shrinkBBox(resid.getBBox())
202  resid = resid.Factory(resid, bbox, deep=True)
203  elif kernels:
204  kim = cand.getKernelImage(diffimLib.KernelCandidateF.ORIG).convertF()
205  resid = kim.Factory(kim, True)
206  im_resid.append(resid)
207 
208  # residuals using spatial model
209  ski = afwImage.ImageD(kernel.getDimensions())
210  kernel.computeImage(ski, False, int(cand.getXCenter()), int(cand.getYCenter()))
211  sk = afwMath.FixedKernel(ski)
212  sbg = 0.0
213  if background:
214  sbg = background(int(cand.getXCenter()), int(cand.getYCenter()))
215  sresid = cand.getDifferenceImage(sk, sbg)
216  resid = sresid
217  if resids:
218  resid = sresid.getImage()
219  resid /= var
220  bbox = kernel.shrinkBBox(resid.getBBox())
221  resid = resid.Factory(resid, bbox, deep=True)
222  elif kernels:
223  kim = ski.convertF()
224  resid = kim.Factory(kim, True)
225  im_resid.append(resid)
226 
227  im = im_resid.makeMosaic()
228 
229  lab = "%d chi^2 %.1f" % (cand.getId(), rchi2)
230  ctype = afwDisplay.RED if cand.isBad() else afwDisplay.GREEN
231 
232  mos.append(im, lab, ctype)
233 
234  if False and np.isnan(rchi2):
235  disp = afwDisplay.Display(frame=1)
236  disp.mtv(cand.getScienceMaskedImage.getImage(), title="candidate")
237  print("rating", cand.getCandidateRating())
238 
239  im = cand.getScienceMaskedImage()
240  center = (candidateIndex, cand.getXCenter() - im.getX0(), cand.getYCenter() - im.getY0())
241  candidateIndex += 1
242  if cand.isBad():
243  candidateCentersBad.append(center)
244  else:
245  candidateCenters.append(center)
246 
247  if resids:
248  title = "chi Diffim"
249  elif kernels:
250  title = "Kernels"
251  else:
252  title = "Candidates & residuals"
253 
254  disp = afwDisplay.Display(frame=frame)
255  mosaicImage = mos.makeMosaic(display=disp, title=title)
256 
257  return mosaicImage
258 
259 
260 def showKernelBasis(kernel, frame=None):
261  """Display a Kernel's basis images.
262  """
263  mos = afwDisplay.utils.Mosaic()
264 
265  for k in kernel.getKernelList():
266  im = afwImage.ImageD(k.getDimensions())
267  k.computeImage(im, False)
268  mos.append(im)
269 
270  disp = afwDisplay.Display(frame=frame)
271  mos.makeMosaic(display=disp, title="Kernel Basis Images")
272 
273  return mos
274 
275 
276 
277 
278 def plotKernelSpatialModel(kernel, kernelCellSet, showBadCandidates=True,
279  numSample=128, keepPlots=True, maxCoeff=10):
280  """Plot the Kernel spatial model.
281  """
282  try:
283  import matplotlib.pyplot as plt
284  import matplotlib.colors
285  except ImportError as e:
286  print("Unable to import numpy and matplotlib: %s" % e)
287  return
288 
289  x0 = kernelCellSet.getBBox().getBeginX()
290  y0 = kernelCellSet.getBBox().getBeginY()
291 
292  candPos = list()
293  candFits = list()
294  badPos = list()
295  badFits = list()
296  candAmps = list()
297  badAmps = list()
298  for cell in kernelCellSet.getCellList():
299  for cand in cell.begin(False):
300  if not showBadCandidates and cand.isBad():
301  continue
302  candCenter = geom.PointD(cand.getXCenter(), cand.getYCenter())
303  try:
304  im = cand.getTemplateMaskedImage()
305  except Exception:
306  continue
307 
308  targetFits = badFits if cand.isBad() else candFits
309  targetPos = badPos if cand.isBad() else candPos
310  targetAmps = badAmps if cand.isBad() else candAmps
311 
312  # compare original and spatial kernel coefficients
313  kp0 = np.array(cand.getKernel(diffimLib.KernelCandidateF.ORIG).getKernelParameters())
314  amp = cand.getCandidateRating()
315 
316  targetFits = badFits if cand.isBad() else candFits
317  targetPos = badPos if cand.isBad() else candPos
318  targetAmps = badAmps if cand.isBad() else candAmps
319 
320  targetFits.append(kp0)
321  targetPos.append(candCenter)
322  targetAmps.append(amp)
323 
324  xGood = np.array([pos.getX() for pos in candPos]) - x0
325  yGood = np.array([pos.getY() for pos in candPos]) - y0
326  zGood = np.array(candFits)
327 
328  xBad = np.array([pos.getX() for pos in badPos]) - x0
329  yBad = np.array([pos.getY() for pos in badPos]) - y0
330  zBad = np.array(badFits)
331  numBad = len(badPos)
332 
333  xRange = np.linspace(0, kernelCellSet.getBBox().getWidth(), num=numSample)
334  yRange = np.linspace(0, kernelCellSet.getBBox().getHeight(), num=numSample)
335 
336  if maxCoeff:
337  maxCoeff = min(maxCoeff, kernel.getNKernelParameters())
338  else:
339  maxCoeff = kernel.getNKernelParameters()
340 
341  for k in range(maxCoeff):
342  func = kernel.getSpatialFunction(k)
343  dfGood = zGood[:, k] - np.array([func(pos.getX(), pos.getY()) for pos in candPos])
344  yMin = dfGood.min()
345  yMax = dfGood.max()
346  if numBad > 0:
347  dfBad = zBad[:, k] - np.array([func(pos.getX(), pos.getY()) for pos in badPos])
348  # Can really screw up the range...
349  yMin = min([yMin, dfBad.min()])
350  yMax = max([yMax, dfBad.max()])
351  yMin -= 0.05*(yMax - yMin)
352  yMax += 0.05*(yMax - yMin)
353 
354  fRange = np.ndarray((len(xRange), len(yRange)))
355  for j, yVal in enumerate(yRange):
356  for i, xVal in enumerate(xRange):
357  fRange[j][i] = func(xVal, yVal)
358 
359  fig = plt.figure(k)
360 
361  fig.clf()
362  try:
363  fig.canvas._tkcanvas._root().lift() # == Tk's raise, but raise is a python reserved word
364  except Exception: # protect against API changes
365  pass
366 
367  fig.suptitle('Kernel component %d' % k)
368 
369  # LL
370  ax = fig.add_axes((0.1, 0.05, 0.35, 0.35))
371  vmin = fRange.min() # - 0.05*np.fabs(fRange.min())
372  vmax = fRange.max() # + 0.05*np.fabs(fRange.max())
373  norm = matplotlib.colors.Normalize(vmin=vmin, vmax=vmax)
374  im = ax.imshow(fRange, aspect='auto', norm=norm,
375  extent=[0, kernelCellSet.getBBox().getWidth() - 1,
376  0, kernelCellSet.getBBox().getHeight() - 1])
377  ax.set_title('Spatial polynomial')
378  plt.colorbar(im, orientation='horizontal', ticks=[vmin, vmax])
379 
380  # UL
381  ax = fig.add_axes((0.1, 0.55, 0.35, 0.35))
382  ax.plot(-2.5*np.log10(candAmps), zGood[:, k], 'b+')
383  if numBad > 0:
384  ax.plot(-2.5*np.log10(badAmps), zBad[:, k], 'r+')
385  ax.set_title("Basis Coefficients")
386  ax.set_xlabel("Instr mag")
387  ax.set_ylabel("Coeff")
388 
389  # LR
390  ax = fig.add_axes((0.55, 0.05, 0.35, 0.35))
391  ax.set_autoscale_on(False)
392  ax.set_xbound(lower=0, upper=kernelCellSet.getBBox().getHeight())
393  ax.set_ybound(lower=yMin, upper=yMax)
394  ax.plot(yGood, dfGood, 'b+')
395  if numBad > 0:
396  ax.plot(yBad, dfBad, 'r+')
397  ax.axhline(0.0)
398  ax.set_title('dCoeff (indiv-spatial) vs. y')
399 
400  # UR
401  ax = fig.add_axes((0.55, 0.55, 0.35, 0.35))
402  ax.set_autoscale_on(False)
403  ax.set_xbound(lower=0, upper=kernelCellSet.getBBox().getWidth())
404  ax.set_ybound(lower=yMin, upper=yMax)
405  ax.plot(xGood, dfGood, 'b+')
406  if numBad > 0:
407  ax.plot(xBad, dfBad, 'r+')
408  ax.axhline(0.0)
409  ax.set_title('dCoeff (indiv-spatial) vs. x')
410 
411  fig.show()
412 
413  global keptPlots
414  if keepPlots and not keptPlots:
415  # Keep plots open when done
416  def show():
417  print("%s: Please close plots when done." % __name__)
418  try:
419  plt.show()
420  except Exception:
421  pass
422  print("Plots closed, exiting...")
423  import atexit
424  atexit.register(show)
425  keptPlots = True
426 
427 
428 def plotKernelCoefficients(spatialKernel, kernelCellSet, showBadCandidates=False, keepPlots=True):
429  """Plot the individual kernel candidate and the spatial kernel solution coefficients.
430 
431  Parameters
432  ----------
433 
434  spatialKernel : `lsst.afw.math.LinearCombinationKernel`
435  The spatial spatialKernel solution model which is a spatially varying linear combination
436  of the spatialKernel basis functions.
437  Typically returned by `lsst.ip.diffim.SpatialKernelSolution.getSolutionPair()`.
438 
439  kernelCellSet : `lsst.afw.math.SpatialCellSet`
440  The spatial cells that was used for solution for the spatialKernel. They contain the
441  local solutions of the AL kernel for the selected sources.
442 
443  showBadCandidates : `bool`, optional
444  If True, plot the coefficient values for kernel candidates where the solution was marked
445  bad by the numerical algorithm. Defaults to False.
446 
447  keepPlots: `bool`, optional
448  If True, sets ``plt.show()`` to be called before the task terminates, so that the plots
449  can be explored interactively. Defaults to True.
450 
451  Notes
452  -----
453  This function produces 3 figures per image subtraction operation.
454  * A grid plot of the local solutions. Each grid cell corresponds to a proportional area in
455  the image. In each cell, local kernel solution coefficients are plotted of kernel candidates (color)
456  that fall into this area as a function of the kernel basis function number.
457  * A grid plot of the spatial solution. Each grid cell corresponds to a proportional area in
458  the image. In each cell, the spatial solution coefficients are evaluated for the center of the cell.
459  * Histogram of the local solution coefficients. Red line marks the spatial solution value at
460  center of the image.
461 
462  This function is called if ``lsst.ip.diffim.psfMatch.plotKernelCoefficients==True`` in lsstDebug. This
463  function was implemented as part of DM-17825.
464  """
465  try:
466  import matplotlib.pyplot as plt
467  except ImportError as e:
468  print("Unable to import matplotlib: %s" % e)
469  return
470 
471  # Image dimensions
472  imgBBox = kernelCellSet.getBBox()
473  x0 = imgBBox.getBeginX()
474  y0 = imgBBox.getBeginY()
475  wImage = imgBBox.getWidth()
476  hImage = imgBBox.getHeight()
477  imgCenterX = imgBBox.getCenterX()
478  imgCenterY = imgBBox.getCenterY()
479 
480  # Plot the local solutions
481  # ----
482 
483  # Grid size
484  nX = 8
485  nY = 8
486  wCell = wImage / nX
487  hCell = hImage / nY
488 
489  fig = plt.figure()
490  fig.suptitle("Kernel candidate parameters on an image grid")
491  arrAx = fig.subplots(nrows=nY, ncols=nX, sharex=True, sharey=True, gridspec_kw=dict(
492  wspace=0, hspace=0))
493 
494  # Bottom left panel is for bottom left part of the image
495  arrAx = arrAx[::-1, :]
496 
497  allParams = []
498  for cell in kernelCellSet.getCellList():
499  cellBBox = afwGeom.Box2D(cell.getBBox())
500  # Determine which panel this spatial cell belongs to
501  iX = int((cellBBox.getCenterX() - x0)//wCell)
502  iY = int((cellBBox.getCenterY() - y0)//hCell)
503 
504  for cand in cell.begin(False):
505  try:
506  kernel = cand.getKernel(cand.ORIG)
507  except Exception:
508  continue
509 
510  if not showBadCandidates and cand.isBad():
511  continue
512 
513  nKernelParams = kernel.getNKernelParameters()
514  kernelParams = np.array(kernel.getKernelParameters())
515  allParams.append(kernelParams)
516 
517  if cand.isBad():
518  color = 'red'
519  else:
520  color = None
521  arrAx[iY, iX].plot(np.arange(nKernelParams), kernelParams, '.-',
522  color=color, drawstyle='steps-mid', linewidth=0.1)
523  for ax in arrAx.ravel():
524  ax.grid(True, axis='y')
525 
526  # Plot histogram of the local parameters and the global solution at the image center
527  # ----
528 
529  spatialFuncs = spatialKernel.getSpatialFunctionList()
530  nKernelParams = spatialKernel.getNKernelParameters()
531  nX = 8
532  fig = plt.figure()
533  fig.suptitle("Hist. of parameters marked with spatial solution at img center")
534  arrAx = fig.subplots(nrows=int(nKernelParams//nX)+1, ncols=nX)
535  arrAx = arrAx[::-1, :]
536  allParams = np.array(allParams)
537  for k in range(nKernelParams):
538  ax = arrAx.ravel()[k]
539  ax.hist(allParams[:, k], bins=20, edgecolor='black')
540  ax.set_xlabel('P{}'.format(k))
541  valueParam = spatialFuncs[k](imgCenterX, imgCenterY)
542  ax.axvline(x=valueParam, color='red')
543  ax.text(0.1, 0.9, '{:.1f}'.format(valueParam),
544  transform=ax.transAxes, backgroundcolor='lightsteelblue')
545 
546  # Plot grid of the spatial solution
547  # ----
548 
549  nX = 8
550  nY = 8
551  wCell = wImage / nX
552  hCell = hImage / nY
553  x0 += wCell / 2
554  y0 += hCell / 2
555 
556  fig = plt.figure()
557  fig.suptitle("Spatial solution of kernel parameters on an image grid")
558  arrAx = fig.subplots(nrows=nY, ncols=nX, sharex=True, sharey=True, gridspec_kw=dict(
559  wspace=0, hspace=0))
560  arrAx = arrAx[::-1, :]
561  kernelParams = np.zeros(nKernelParams, dtype=float)
562 
563  for iX in range(nX):
564  for iY in range(nY):
565  x = x0 + iX * wCell
566  y = y0 + iY * hCell
567  # Evaluate the spatial solution functions for this x,y location
568  kernelParams = [f(x, y) for f in spatialFuncs]
569  arrAx[iY, iX].plot(np.arange(nKernelParams), kernelParams, '.-', drawstyle='steps-mid')
570  arrAx[iY, iX].grid(True, axis='y')
571 
572  global keptPlots
573  if keepPlots and not keptPlots:
574  # Keep plots open when done
575  def show():
576  print("%s: Please close plots when done." % __name__)
577  try:
578  plt.show()
579  except Exception:
580  pass
581  print("Plots closed, exiting...")
582  import atexit
583  atexit.register(show)
584  keptPlots = True
585 
586 
587 def showKernelMosaic(bbox, kernel, nx=7, ny=None, frame=None, title=None,
588  showCenter=True, showEllipticity=True):
589  """Show a mosaic of Kernel images.
590  """
591  mos = afwDisplay.utils.Mosaic()
592 
593  x0 = bbox.getBeginX()
594  y0 = bbox.getBeginY()
595  width = bbox.getWidth()
596  height = bbox.getHeight()
597 
598  if not ny:
599  ny = int(nx*float(height)/width + 0.5)
600  if not ny:
601  ny = 1
602 
603  schema = afwTable.SourceTable.makeMinimalSchema()
604  centroidName = "base_SdssCentroid"
605  shapeName = "base_SdssShape"
606  control = measBase.SdssCentroidControl()
607  schema.getAliasMap().set("slot_Centroid", centroidName)
608  schema.getAliasMap().set("slot_Centroid_flag", centroidName + "_flag")
609  centroider = measBase.SdssCentroidAlgorithm(control, centroidName, schema)
610  sdssShape = measBase.SdssShapeControl()
611  shaper = measBase.SdssShapeAlgorithm(sdssShape, shapeName, schema)
612  table = afwTable.SourceTable.make(schema)
613  table.defineCentroid(centroidName)
614  table.defineShape(shapeName)
615 
616  centers = []
617  shapes = []
618  for iy in range(ny):
619  for ix in range(nx):
620  x = int(ix*(width - 1)/(nx - 1)) + x0
621  y = int(iy*(height - 1)/(ny - 1)) + y0
622 
623  im = afwImage.ImageD(kernel.getDimensions())
624  ksum = kernel.computeImage(im, False, x, y)
625  lab = "Kernel(%d,%d)=%.2f" % (x, y, ksum) if False else ""
626  mos.append(im, lab)
627 
628  # SdssCentroidAlgorithm.measure requires an exposure of floats
629  exp = afwImage.makeExposure(afwImage.makeMaskedImage(im.convertF()))
630 
631  w, h = im.getWidth(), im.getHeight()
632  centerX = im.getX0() + w//2
633  centerY = im.getY0() + h//2
634  src = table.makeRecord()
635  spans = afwGeom.SpanSet(exp.getBBox())
636  foot = afwDet.Footprint(spans)
637  foot.addPeak(centerX, centerY, 1)
638  src.setFootprint(foot)
639 
640  try: # The centroider requires a psf, so this will fail if none is attached to exp
641  centroider.measure(src, exp)
642  centers.append((src.getX(), src.getY()))
643 
644  shaper.measure(src, exp)
645  shapes.append((src.getIxx(), src.getIxy(), src.getIyy()))
646  except Exception:
647  pass
648 
649  disp = afwDisplay.Display(frame=frame)
650  mos.makeMosaic(display=disp, title=title if title else "Model Kernel", mode=nx)
651 
652  if centers and frame is not None:
653  disp = afwDisplay.Display(frame=frame)
654  i = 0
655  with disp.Buffering():
656  for cen, shape in zip(centers, shapes):
657  bbox = mos.getBBox(i)
658  i += 1
659  xc, yc = cen[0] + bbox.getMinX(), cen[1] + bbox.getMinY()
660  if showCenter:
661  disp.dot("+", xc, yc, ctype=afwDisplay.BLUE)
662 
663  if showEllipticity:
664  ixx, ixy, iyy = shape
665  disp.dot("@:%g,%g,%g" % (ixx, ixy, iyy), xc, yc, ctype=afwDisplay.RED)
666 
667  return mos
668 
669 
670 def plotPixelResiduals(exposure, warpedTemplateExposure, diffExposure, kernelCellSet,
671  kernel, background, testSources, config,
672  origVariance=False, nptsFull=1e6, keepPlots=True, titleFs=14):
673  """Plot diffim residuals for LOCAL and SPATIAL models.
674  """
675  candidateResids = []
676  spatialResids = []
677  nonfitResids = []
678 
679  for cell in kernelCellSet.getCellList():
680  for cand in cell.begin(True): # only look at good ones
681  # Be sure
682  if not (cand.getStatus() == afwMath.SpatialCellCandidate.GOOD):
683  continue
684 
685  diffim = cand.getDifferenceImage(diffimLib.KernelCandidateF.ORIG)
686  orig = cand.getScienceMaskedImage()
687 
688  ski = afwImage.ImageD(kernel.getDimensions())
689  kernel.computeImage(ski, False, int(cand.getXCenter()), int(cand.getYCenter()))
690  sk = afwMath.FixedKernel(ski)
691  sbg = background(int(cand.getXCenter()), int(cand.getYCenter()))
692  sdiffim = cand.getDifferenceImage(sk, sbg)
693 
694  # trim edgs due to convolution
695  bbox = kernel.shrinkBBox(diffim.getBBox())
696  tdiffim = diffim.Factory(diffim, bbox)
697  torig = orig.Factory(orig, bbox)
698  tsdiffim = sdiffim.Factory(sdiffim, bbox)
699 
700  if origVariance:
701  candidateResids.append(np.ravel(tdiffim.getImage().getArray()
702  / np.sqrt(torig.getVariance().getArray())))
703  spatialResids.append(np.ravel(tsdiffim.getImage().getArray()
704  / np.sqrt(torig.getVariance().getArray())))
705  else:
706  candidateResids.append(np.ravel(tdiffim.getImage().getArray()
707  / np.sqrt(tdiffim.getVariance().getArray())))
708  spatialResids.append(np.ravel(tsdiffim.getImage().getArray()
709  / np.sqrt(tsdiffim.getVariance().getArray())))
710 
711  fullIm = diffExposure.getMaskedImage().getImage().getArray()
712  fullMask = diffExposure.getMaskedImage().getMask().getArray()
713  if origVariance:
714  fullVar = exposure.getMaskedImage().getVariance().getArray()
715  else:
716  fullVar = diffExposure.getMaskedImage().getVariance().getArray()
717 
718  bitmaskBad = 0
719  bitmaskBad |= afwImage.Mask.getPlaneBitMask('NO_DATA')
720  bitmaskBad |= afwImage.Mask.getPlaneBitMask('SAT')
721  idx = np.where((fullMask & bitmaskBad) == 0)
722  stride = int(len(idx[0])//nptsFull)
723  sidx = idx[0][::stride], idx[1][::stride]
724  allResids = fullIm[sidx]/np.sqrt(fullVar[sidx])
725 
726  testFootprints = diffimTools.sourceToFootprintList(testSources, warpedTemplateExposure,
727  exposure, config, Log.getDefaultLogger())
728  for fp in testFootprints:
729  subexp = diffExposure.Factory(diffExposure, fp["footprint"].getBBox())
730  subim = subexp.getMaskedImage().getImage()
731  if origVariance:
732  subvar = afwImage.ExposureF(exposure, fp["footprint"].getBBox()).getMaskedImage().getVariance()
733  else:
734  subvar = subexp.getMaskedImage().getVariance()
735  nonfitResids.append(np.ravel(subim.getArray()/np.sqrt(subvar.getArray())))
736 
737  candidateResids = np.ravel(np.array(candidateResids))
738  spatialResids = np.ravel(np.array(spatialResids))
739  nonfitResids = np.ravel(np.array(nonfitResids))
740 
741  try:
742  import pylab
743  from matplotlib.font_manager import FontProperties
744  except ImportError as e:
745  print("Unable to import pylab: %s" % e)
746  return
747 
748  fig = pylab.figure()
749  fig.clf()
750  try:
751  fig.canvas._tkcanvas._root().lift() # == Tk's raise, but raise is a python reserved word
752  except Exception: # protect against API changes
753  pass
754  if origVariance:
755  fig.suptitle("Diffim residuals: Normalized by sqrt(input variance)", fontsize=titleFs)
756  else:
757  fig.suptitle("Diffim residuals: Normalized by sqrt(diffim variance)", fontsize=titleFs)
758 
759  sp1 = pylab.subplot(221)
760  sp2 = pylab.subplot(222, sharex=sp1, sharey=sp1)
761  sp3 = pylab.subplot(223, sharex=sp1, sharey=sp1)
762  sp4 = pylab.subplot(224, sharex=sp1, sharey=sp1)
763  xs = np.arange(-5, 5.05, 0.1)
764  ys = 1./np.sqrt(2*np.pi)*np.exp(-0.5*xs**2)
765 
766  sp1.hist(candidateResids, bins=xs, normed=True, alpha=0.5, label="N(%.2f, %.2f)"
767  % (np.mean(candidateResids), np.var(candidateResids)))
768  sp1.plot(xs, ys, "r-", lw=2, label="N(0,1)")
769  sp1.set_title("Candidates: basis fit", fontsize=titleFs - 2)
770  sp1.legend(loc=1, fancybox=True, shadow=True, prop=FontProperties(size=titleFs - 6))
771 
772  sp2.hist(spatialResids, bins=xs, normed=True, alpha=0.5, label="N(%.2f, %.2f)"
773  % (np.mean(spatialResids), np.var(spatialResids)))
774  sp2.plot(xs, ys, "r-", lw=2, label="N(0,1)")
775  sp2.set_title("Candidates: spatial fit", fontsize=titleFs - 2)
776  sp2.legend(loc=1, fancybox=True, shadow=True, prop=FontProperties(size=titleFs - 6))
777 
778  sp3.hist(nonfitResids, bins=xs, normed=True, alpha=0.5, label="N(%.2f, %.2f)"
779  % (np.mean(nonfitResids), np.var(nonfitResids)))
780  sp3.plot(xs, ys, "r-", lw=2, label="N(0,1)")
781  sp3.set_title("Control sample: spatial fit", fontsize=titleFs - 2)
782  sp3.legend(loc=1, fancybox=True, shadow=True, prop=FontProperties(size=titleFs - 6))
783 
784  sp4.hist(allResids, bins=xs, normed=True, alpha=0.5, label="N(%.2f, %.2f)"
785  % (np.mean(allResids), np.var(allResids)))
786  sp4.plot(xs, ys, "r-", lw=2, label="N(0,1)")
787  sp4.set_title("Full image (subsampled)", fontsize=titleFs - 2)
788  sp4.legend(loc=1, fancybox=True, shadow=True, prop=FontProperties(size=titleFs - 6))
789 
790  pylab.setp(sp1.get_xticklabels() + sp1.get_yticklabels(), fontsize=titleFs - 4)
791  pylab.setp(sp2.get_xticklabels() + sp2.get_yticklabels(), fontsize=titleFs - 4)
792  pylab.setp(sp3.get_xticklabels() + sp3.get_yticklabels(), fontsize=titleFs - 4)
793  pylab.setp(sp4.get_xticklabels() + sp4.get_yticklabels(), fontsize=titleFs - 4)
794 
795  sp1.set_xlim(-5, 5)
796  sp1.set_ylim(0, 0.5)
797  fig.show()
798 
799  global keptPlots
800  if keepPlots and not keptPlots:
801  # Keep plots open when done
802  def show():
803  print("%s: Please close plots when done." % __name__)
804  try:
805  pylab.show()
806  except Exception:
807  pass
808  print("Plots closed, exiting...")
809  import atexit
810  atexit.register(show)
811  keptPlots = True
812 
813 
814 def calcCentroid(arr):
815  """Calculate first moment of a (kernel) image.
816  """
817  y, x = arr.shape
818  sarr = arr*arr
819  xarr = np.asarray([[el for el in range(x)] for el2 in range(y)])
820  yarr = np.asarray([[el2 for el in range(x)] for el2 in range(y)])
821  narr = xarr*sarr
822  sarrSum = sarr.sum()
823  centx = narr.sum()/sarrSum
824  narr = yarr*sarr
825  centy = narr.sum()/sarrSum
826  return centx, centy
827 
828 
829 def calcWidth(arr, centx, centy):
830  """Calculate second moment of a (kernel) image.
831  """
832  y, x = arr.shape
833  # Square the flux so we don't have to deal with negatives
834  sarr = arr*arr
835  xarr = np.asarray([[el for el in range(x)] for el2 in range(y)])
836  yarr = np.asarray([[el2 for el in range(x)] for el2 in range(y)])
837  narr = sarr*np.power((xarr - centx), 2.)
838  sarrSum = sarr.sum()
839  xstd = np.sqrt(narr.sum()/sarrSum)
840  narr = sarr*np.power((yarr - centy), 2.)
841  ystd = np.sqrt(narr.sum()/sarrSum)
842  return xstd, ystd
843 
844 
845 def printSkyDiffs(sources, wcs):
846  """Print differences in sky coordinates.
847 
848  The difference is that between the source Position and its Centroid mapped
849  through Wcs.
850  """
851  for s in sources:
852  sCentroid = s.getCentroid()
853  sPosition = s.getCoord().getPosition(geom.degrees)
854  dra = 3600*(sPosition.getX() - wcs.pixelToSky(sCentroid).getPosition(geom.degrees).getX())/0.2
855  ddec = 3600*(sPosition.getY() - wcs.pixelToSky(sCentroid).getPosition(geom.degrees).getY())/0.2
856  if np.isfinite(dra) and np.isfinite(ddec):
857  print(dra, ddec)
858 
859 
860 def makeRegions(sources, outfilename, wcs=None):
861  """Create regions file for display from input source list.
862  """
863  fh = open(outfilename, "w")
864  fh.write("global color=red font=\"helvetica 10 normal\" "
865  "select=1 highlite=1 edit=1 move=1 delete=1 include=1 fixed=0 source\nfk5\n")
866  for s in sources:
867  if wcs:
868  (ra, dec) = wcs.pixelToSky(s.getCentroid()).getPosition(geom.degrees)
869  else:
870  (ra, dec) = s.getCoord().getPosition(geom.degrees)
871  if np.isfinite(ra) and np.isfinite(dec):
872  fh.write("circle(%f,%f,2\")\n"%(ra, dec))
873  fh.flush()
874  fh.close()
875 
876 
877 def showSourceSetSky(sSet, wcs, xy0, frame=0, ctype=afwDisplay.GREEN, symb="+", size=2):
878  """Draw the (RA, Dec) positions of a set of Sources. Image has the XY0.
879  """
880  disp = afwDisplay.Display(frame=frame)
881  with disp.Buffering():
882  for s in sSet:
883  (xc, yc) = wcs.skyToPixel(s.getCoord().getRa(), s.getCoord().getDec())
884  xc -= xy0[0]
885  yc -= xy0[1]
886  disp.dot(symb, xc, yc, ctype=ctype, size=size)
887 
888 
889 def plotWhisker(results, newWcs):
890  """Plot whisker diagram of astromeric offsets between results.matches.
891  """
892  refCoordKey = results.matches[0].first.getTable().getCoordKey()
893  inCentroidKey = results.matches[0].second.getTable().getCentroidKey()
894  positions = [m.first.get(refCoordKey) for m in results.matches]
895  residuals = [m.first.get(refCoordKey).getOffsetFrom(
896  newWcs.pixelToSky(m.second.get(inCentroidKey))) for
897  m in results.matches]
898  import matplotlib.pyplot as plt
899  fig = plt.figure()
900  sp = fig.add_subplot(1, 1, 0)
901  xpos = [x[0].asDegrees() for x in positions]
902  ypos = [x[1].asDegrees() for x in positions]
903  xpos.append(0.02*(max(xpos) - min(xpos)) + min(xpos))
904  ypos.append(0.98*(max(ypos) - min(ypos)) + min(ypos))
905  xidxs = np.isfinite(xpos)
906  yidxs = np.isfinite(ypos)
907  X = np.asarray(xpos)[xidxs]
908  Y = np.asarray(ypos)[yidxs]
909  distance = [x[1].asArcseconds() for x in residuals]
910  distance.append(0.2)
911  distance = np.asarray(distance)[xidxs]
912  # NOTE: This assumes that the bearing is measured positive from +RA through North.
913  # From the documentation this is not clear.
914  bearing = [x[0].asRadians() for x in residuals]
915  bearing.append(0)
916  bearing = np.asarray(bearing)[xidxs]
917  U = (distance*np.cos(bearing))
918  V = (distance*np.sin(bearing))
919  sp.quiver(X, Y, U, V)
920  sp.set_title("WCS Residual")
921  plt.show()
922 
923 
925  """Utility class for dipole measurement testing.
926 
927  Generate an image with simulated dipoles and noise; store the original
928  "pre-subtraction" images and catalogs as well.
929  Used to generate test data for DMTN-007 (http://dmtn-007.lsst.io).
930  """
931 
932  def __init__(self, w=101, h=101, xcenPos=[27.], ycenPos=[25.], xcenNeg=[23.], ycenNeg=[25.],
933  psfSigma=2., flux=[30000.], fluxNeg=None, noise=10., gradientParams=None):
934  self.w = w
935  self.h = h
936  self.xcenPos = xcenPos
937  self.ycenPos = ycenPos
938  self.xcenNeg = xcenNeg
939  self.ycenNeg = ycenNeg
940  self.psfSigma = psfSigma
941  self.flux = flux
942  self.fluxNeg = fluxNeg
943  if fluxNeg is None:
944  self.fluxNeg = self.flux
945  self.noise = noise
946  self.gradientParams = gradientParams
947  self._makeDipoleImage()
948 
949  def _makeDipoleImage(self):
950  """Generate an exposure and catalog with the given dipole source(s).
951  """
952  # Must seed the pos/neg images with different values to ensure they get different noise realizations
953  posImage, posCatalog = self._makeStarImage(
954  xc=self.xcenPos, yc=self.ycenPos, flux=self.flux, randomSeed=111)
955 
956  negImage, negCatalog = self._makeStarImage(
957  xc=self.xcenNeg, yc=self.ycenNeg, flux=self.fluxNeg, randomSeed=222)
958 
959  dipole = posImage.clone()
960  di = dipole.getMaskedImage()
961  di -= negImage.getMaskedImage()
962 
963  # Carry through pos/neg detection masks to new planes in diffim
964  dm = di.getMask()
965  posDetectedBits = posImage.getMaskedImage().getMask().getArray() == dm.getPlaneBitMask("DETECTED")
966  negDetectedBits = negImage.getMaskedImage().getMask().getArray() == dm.getPlaneBitMask("DETECTED")
967  pos_det = dm.addMaskPlane("DETECTED_POS") # new mask plane -- different from "DETECTED"
968  neg_det = dm.addMaskPlane("DETECTED_NEG") # new mask plane -- different from "DETECTED_NEGATIVE"
969  dma = dm.getArray()
970  # set the two custom mask planes to these new masks
971  dma[:, :] = posDetectedBits*pos_det + negDetectedBits*neg_det
972  self.diffim, self.posImage, self.posCatalog, self.negImage, self.negCatalog \
973  = dipole, posImage, posCatalog, negImage, negCatalog
974 
975  def _makeStarImage(self, xc=[15.3], yc=[18.6], flux=[2500], schema=None, randomSeed=None):
976  """Generate an exposure and catalog with the given stellar source(s).
977  """
978  from lsst.meas.base.tests import TestDataset
979  bbox = geom.Box2I(geom.Point2I(0, 0), geom.Point2I(self.w - 1, self.h - 1))
980  dataset = TestDataset(bbox, psfSigma=self.psfSigma, threshold=1.)
981 
982  for i in range(len(xc)):
983  dataset.addSource(instFlux=flux[i], centroid=geom.Point2D(xc[i], yc[i]))
984 
985  if schema is None:
986  schema = TestDataset.makeMinimalSchema()
987  exposure, catalog = dataset.realize(noise=self.noise, schema=schema, randomSeed=randomSeed)
988 
989  if self.gradientParams is not None:
990  y, x = np.mgrid[:self.w, :self.h]
991  gp = self.gradientParams
992  gradient = gp[0] + gp[1]*x + gp[2]*y
993  if len(self.gradientParams) > 3: # it includes a set of 2nd-order polynomial params
994  gradient += gp[3]*x*y + gp[4]*x*x + gp[5]*y*y
995  imgArr = exposure.getMaskedImage().getArrays()[0]
996  imgArr += gradient
997 
998  return exposure, catalog
999 
1000  def fitDipoleSource(self, source, **kwds):
1001  alg = DipoleFitAlgorithm(self.diffim, self.posImage, self.negImage)
1002  fitResult = alg.fitDipole(source, **kwds)
1003  return fitResult
1004 
1005  def detectDipoleSources(self, doMerge=True, diffim=None, detectSigma=5.5, grow=3, minBinSize=32):
1006  """Utility function for detecting dipoles.
1007 
1008  Detect pos/neg sources in the diffim, then merge them. A
1009  bigger "grow" parameter leads to a larger footprint which
1010  helps with dipole measurement for faint dipoles.
1011 
1012  Parameters
1013  ----------
1014  doMerge : `bool`
1015  Whether to merge the positive and negagive detections into a single
1016  source table.
1017  diffim : `lsst.afw.image.exposure.exposure.ExposureF`
1018  Difference image on which to perform detection.
1019  detectSigma : `float`
1020  Threshold for object detection.
1021  grow : `int`
1022  Number of pixels to grow the footprints before merging.
1023  minBinSize : `int`
1024  Minimum bin size for the background (re)estimation (only applies if
1025  the default leads to min(nBinX, nBinY) < fit order so the default
1026  config parameter needs to be decreased, but not to a value smaller
1027  than ``minBinSize``, in which case the fitting algorithm will take
1028  over and decrease the fit order appropriately.)
1029 
1030  Returns
1031  -------
1032  sources : `lsst.afw.table.SourceCatalog`
1033  If doMerge=True, the merged source catalog is returned OR
1034  detectTask : `lsst.meas.algorithms.SourceDetectionTask`
1035  schema : `lsst.afw.table.Schema`
1036  If doMerge=False, the source detection task and its schema are
1037  returned.
1038  """
1039  if diffim is None:
1040  diffim = self.diffim
1041 
1042  # Start with a minimal schema - only the fields all SourceCatalogs need
1043  schema = afwTable.SourceTable.makeMinimalSchema()
1044 
1045  # Customize the detection task a bit (optional)
1046  detectConfig = measAlg.SourceDetectionConfig()
1047  detectConfig.returnOriginalFootprints = False # should be the default
1048 
1049  psfSigma = diffim.getPsf().computeShape().getDeterminantRadius()
1050 
1051  # code from imageDifference.py:
1052  detectConfig.thresholdPolarity = "both"
1053  detectConfig.thresholdValue = detectSigma
1054  # detectConfig.nSigmaToGrow = psfSigma
1055  detectConfig.reEstimateBackground = True # if False, will fail often for faint sources on gradients?
1056  detectConfig.thresholdType = "pixel_stdev"
1057  # Test images are often quite small, so may need to adjust background binSize
1058  while ((min(diffim.getWidth(), diffim.getHeight()))//detectConfig.background.binSize
1059  < detectConfig.background.approxOrderX and detectConfig.background.binSize > minBinSize):
1060  detectConfig.background.binSize = max(minBinSize, detectConfig.background.binSize//2)
1061 
1062  # Create the detection task. We pass the schema so the task can declare a few flag fields
1063  detectTask = measAlg.SourceDetectionTask(schema, config=detectConfig)
1064 
1065  table = afwTable.SourceTable.make(schema)
1066  catalog = detectTask.run(table, diffim, sigma=psfSigma)
1067 
1068  # Now do the merge.
1069  if doMerge:
1070  fpSet = catalog.fpSets.positive
1071  fpSet.merge(catalog.fpSets.negative, grow, grow, False)
1072  sources = afwTable.SourceCatalog(table)
1073  fpSet.makeSources(sources)
1074 
1075  return sources
1076 
1077  else:
1078  return detectTask, schema
lsst::afw::image
Backwards-compatibility support for depersisting the old Calib (FluxMag0/FluxMag0Err) objects.
Definition: imageAlgorithm.dox:1
lsst::ip::diffim.utils.DipoleTestImage.h
h
Definition: utils.py:934
lsst::ip::diffim.utils.DipoleTestImage.xcenPos
xcenPos
Definition: utils.py:935
lsst::ip::diffim.utils.plotWhisker
def plotWhisker(results, newWcs)
Definition: utils.py:889
lsst::ip::diffim.utils.DipoleTestImage.gradientParams
gradientParams
Definition: utils.py:945
lsst::ip::diffim.utils.showKernelMosaic
def showKernelMosaic(bbox, kernel, nx=7, ny=None, frame=None, title=None, showCenter=True, showEllipticity=True)
Definition: utils.py:587
lsst::ip::diffim.utils.DipoleTestImage.ycenNeg
ycenNeg
Definition: utils.py:938
lsst::afw::image::makeExposure
std::shared_ptr< Exposure< ImagePixelT, MaskPixelT, VariancePixelT > > makeExposure(MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > &mimage, std::shared_ptr< geom::SkyWcs const > wcs=std::shared_ptr< geom::SkyWcs const >())
A function to return an Exposure of the correct type (cf.
Definition: Exposure.h:442
lsst::ip::diffim.utils.showSourceSetSky
def showSourceSetSky(sSet, wcs, xy0, frame=0, ctype=afwDisplay.GREEN, symb="+", size=2)
Definition: utils.py:877
lsst::ip::diffim.utils.calcCentroid
def calcCentroid(arr)
Definition: utils.py:814
pex.config.history.format
def format(config, name=None, writeSourceLine=True, prefix="", verbose=False)
Definition: history.py:174
lsst::afw::math::FixedKernel
A kernel created from an Image.
Definition: Kernel.h:518
lsst::afw.display
Definition: __init__.py:1
lsst::ip::diffim.utils.DipoleTestImage.__init__
def __init__(self, w=101, h=101, xcenPos=[27.], ycenPos=[25.], xcenNeg=[23.], ycenNeg=[25.], psfSigma=2., flux=[30000.], fluxNeg=None, noise=10., gradientParams=None)
Definition: utils.py:932
lsst::meas::algorithms.objectSizeStarSelector.plot
def plot(mag, width, centers, clusterId, marker="o", markersize=2, markeredgewidth=0, ltype='-', magType="model", clear=True)
Definition: objectSizeStarSelector.py:264
lsst::ip::diffim.utils.DipoleTestImage.fluxNeg
fluxNeg
Definition: utils.py:941
lsst::ip::diffim.utils.DipoleTestImage.noise
noise
Definition: utils.py:944
lsst::ip::diffim.utils.showKernelSpatialCells
def showKernelSpatialCells(maskedIm, kernelCellSet, showChi2=False, symb="o", ctype=None, ctypeUnused=None, ctypeBad=None, size=3, frame=None, title="Spatial Cells")
Definition: utils.py:67
lsst::meas::base.tests
Definition: tests.py:1
lsst::afw::geom::SpanSet
A compact representation of a collection of pixels.
Definition: SpanSet.h:78
lsst::ip::diffim.utils.DipoleTestImage.w
w
Definition: utils.py:933
lsst::ip::diffim.utils.makeRegions
def makeRegions(sources, outfilename, wcs=None)
Definition: utils.py:860
lsst.pipe.drivers.visualizeVisit.background
background
Definition: visualizeVisit.py:37
lsst::meas::base
Definition: Algorithm.h:37
lsst::afw::table._source.SourceCatalog
Definition: _source.py:33
lsst::ip::diffim.utils.plotPixelResiduals
def plotPixelResiduals(exposure, warpedTemplateExposure, diffExposure, kernelCellSet, kernel, background, testSources, config, origVariance=False, nptsFull=1e6, keepPlots=True, titleFs=14)
Definition: utils.py:670
lsst::log
Definition: Log.h:706
lsst::ip::diffim.utils.DipoleTestImage.ycenPos
ycenPos
Definition: utils.py:936
max
int max
Definition: BoundedField.cc:104
lsst::afw::table
Definition: table.dox:3
lsst::ip::diffim.utils.showKernelCandidates
def showKernelCandidates(kernelCellSet, kernel, background, frame=None, showBadCandidates=True, resids=False, kernels=False)
Definition: utils.py:142
lsst::ip::diffim.dipoleFitTask.DipoleFitAlgorithm
Definition: dipoleFitTask.py:483
lsst::ip::diffim.utils.DipoleTestImage._makeStarImage
def _makeStarImage(self, xc=[15.3], yc=[18.6], flux=[2500], schema=None, randomSeed=None)
Definition: utils.py:975
object
lsst::ip::diffim.utils.plotKernelCoefficients
def plotKernelCoefficients(spatialKernel, kernelCellSet, showBadCandidates=False, keepPlots=True)
Definition: utils.py:428
lsst::afw::detection
Definition: Footprint.h:50
lsst::geom
Definition: geomOperators.dox:4
lsst::ip::diffim.utils.DipoleTestImage.xcenNeg
xcenNeg
Definition: utils.py:937
lsst::ip::diffim.utils.printSkyDiffs
def printSkyDiffs(sources, wcs)
Definition: utils.py:845
list
daf::base::PropertyList * list
Definition: fits.cc:913
lsst::ip::diffim.utils.DipoleTestImage
Definition: utils.py:924
min
int min
Definition: BoundedField.cc:103
lsst::afw::math
Definition: statistics.dox:6
lsst::geom::Point< double, 2 >
lsst::ip::diffim.utils.plotKernelSpatialModel
def plotKernelSpatialModel(kernel, kernelCellSet, showBadCandidates=True, numSample=128, keepPlots=True, maxCoeff=10)
Definition: utils.py:278
lsst::geom::Box2I
An integer coordinate rectangle.
Definition: Box.h:55
lsst::ip::diffim.utils.DipoleTestImage.psfSigma
psfSigma
Definition: utils.py:939
lsst::ip::diffim.utils.DipoleTestImage.flux
flux
Definition: utils.py:940
lsst::ip::diffim.utils.DipoleTestImage.detectDipoleSources
def detectDipoleSources(self, doMerge=True, diffim=None, detectSigma=5.5, grow=3, minBinSize=32)
Definition: utils.py:1005
lsst::afw::image::makeMaskedImage
MaskedImage< ImagePixelT, MaskPixelT, VariancePixelT > * makeMaskedImage(typename std::shared_ptr< Image< ImagePixelT >> image, typename std::shared_ptr< Mask< MaskPixelT >> mask=Mask< MaskPixelT >(), typename std::shared_ptr< Image< VariancePixelT >> variance=Image< VariancePixelT >())
A function to return a MaskedImage of the correct type (cf.
Definition: MaskedImage.h:1279
lsst::afw.display.ds9.show
def show(frame=None)
Definition: ds9.py:89
lsst::geom::Box2D
A floating-point coordinate rectangle geometry.
Definition: Box.h:413
lsst::afw::detection::Footprint
Class to describe the properties of a detected object from an image.
Definition: Footprint.h:63
lsst::ip::diffim.utils.showKernelBasis
def showKernelBasis(kernel, frame=None)
Definition: utils.py:260
lsst::meas::algorithms
Fit spatial kernel using approximate fluxes for candidates, and solving a linear system of equations.
Definition: CoaddBoundedField.h:34
lsst::ip::diffim.utils.calcWidth
def calcWidth(arr, centx, centy)
Definition: utils.py:829
lsst::ip::diffim.utils.showSourceSet
def showSourceSet(sSet, xy0=(0, 0), frame=0, ctype=afwDisplay.GREEN, symb="+", size=2)
Definition: utils.py:47
lsst::ip::diffim.utils.showDiaSources
def showDiaSources(sources, exposure, isFlagged, isDipole, frame=None)
Definition: utils.py:106
set
daf::base::PropertySet * set
Definition: fits.cc:912
lsst::ip::diffim.utils.DipoleTestImage.fitDipoleSource
def fitDipoleSource(self, source, **kwds)
Definition: utils.py:1000
lsst::ip::diffim.utils.DipoleTestImage._makeDipoleImage
def _makeDipoleImage(self)
Definition: utils.py:949
lsst::afw::geom
Definition: frameSetUtils.h:40