LSST Applications  21.0.0-172-gfb10e10a+18fedfabac,22.0.0+297cba6710,22.0.0+80564b0ff1,22.0.0+8d77f4f51a,22.0.0+a28f4c53b1,22.0.0+dcf3732eb2,22.0.1-1-g7d6de66+2a20fdde0d,22.0.1-1-g8e32f31+297cba6710,22.0.1-1-geca5380+7fa3b7d9b6,22.0.1-12-g44dc1dc+2a20fdde0d,22.0.1-15-g6a90155+515f58c32b,22.0.1-16-g9282f48+790f5f2caa,22.0.1-2-g92698f7+dcf3732eb2,22.0.1-2-ga9b0f51+7fa3b7d9b6,22.0.1-2-gd1925c9+bf4f0e694f,22.0.1-24-g1ad7a390+a9625a72a8,22.0.1-25-g5bf6245+3ad8ecd50b,22.0.1-25-gb120d7b+8b5510f75f,22.0.1-27-g97737f7+2a20fdde0d,22.0.1-32-gf62ce7b1+aa4237961e,22.0.1-4-g0b3f228+2a20fdde0d,22.0.1-4-g243d05b+871c1b8305,22.0.1-4-g3a563be+32dcf1063f,22.0.1-4-g44f2e3d+9e4ab0f4fa,22.0.1-42-gca6935d93+ba5e5ca3eb,22.0.1-5-g15c806e+85460ae5f3,22.0.1-5-g58711c4+611d128589,22.0.1-5-g75bb458+99c117b92f,22.0.1-6-g1c63a23+7fa3b7d9b6,22.0.1-6-g50866e6+84ff5a128b,22.0.1-6-g8d3140d+720564cf76,22.0.1-6-gd805d02+cc5644f571,22.0.1-8-ge5750ce+85460ae5f3,master-g6e05de7fdc+babf819c66,master-g99da0e417a+8d77f4f51a,w.2021.48
LSST Data Management Base Package
Public Member Functions | Public Attributes | List of all members
lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm Class Reference
Inheritance diagram for lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm:

Public Member Functions

def __init__ (self, diffim, posImage=None, negImage=None)
 
def fitDipoleImpl (self, source, tol=1e-7, rel_weight=0.5, fitBackground=1, bgGradientOrder=1, maxSepInSigma=5., separateNegParams=True, verbose=False)
 
def fitDipole (self, source, tol=1e-7, rel_weight=0.1, fitBackground=1, maxSepInSigma=5., separateNegParams=True, bgGradientOrder=1, verbose=False, display=False)
 
def displayFitResults (self, footprint, result)
 

Public Attributes

 diffim
 
 posImage
 
 negImage
 
 psfSigma
 
 log
 
 debug
 

Detailed Description

Fit a dipole model using an image difference.

See also:
`DMTN-007: Dipole characterization for image differencing  <https://dmtn-007.lsst.io>`_.

Definition at line 484 of file dipoleFitTask.py.

Constructor & Destructor Documentation

◆ __init__()

def lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.__init__ (   self,
  diffim,
  posImage = None,
  negImage = None 
)
Algorithm to run dipole measurement on a diaSource

Parameters
----------
diffim : `lsst.afw.image.Exposure`
    Exposure on which the diaSources were detected
posImage : `lsst.afw.image.Exposure`
    "Positive" exposure from which the template was subtracted
negImage : `lsst.afw.image.Exposure`
    "Negative" exposure which was subtracted from the posImage

Definition at line 510 of file dipoleFitTask.py.

510  def __init__(self, diffim, posImage=None, negImage=None):
511  """Algorithm to run dipole measurement on a diaSource
512 
513  Parameters
514  ----------
515  diffim : `lsst.afw.image.Exposure`
516  Exposure on which the diaSources were detected
517  posImage : `lsst.afw.image.Exposure`
518  "Positive" exposure from which the template was subtracted
519  negImage : `lsst.afw.image.Exposure`
520  "Negative" exposure which was subtracted from the posImage
521  """
522 
523  self.diffim = diffim
524  self.posImage = posImage
525  self.negImage = negImage
526  self.psfSigma = None
527  if diffim is not None:
528  self.psfSigma = diffim.getPsf().computeShape().getDeterminantRadius()
529 
530  self.log = Log.getLogger(__name__)
531 
532  import lsstDebug
533  self.debug = lsstDebug.Info(__name__).debug
534 

Member Function Documentation

◆ displayFitResults()

def lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.displayFitResults (   self,
  footprint,
  result 
)
Display data, model fits and residuals (currently uses matplotlib display functions).

Parameters
----------
footprint : TODO: DM-17458
    Footprint containing the dipole that was fit
result : `lmfit.MinimizerResult`
    `lmfit.MinimizerResult` object returned by `lmfit` optimizer

Returns
-------
fig : `matplotlib.pyplot.plot`

Definition at line 885 of file dipoleFitTask.py.

885  def displayFitResults(self, footprint, result):
886  """Display data, model fits and residuals (currently uses matplotlib display functions).
887 
888  Parameters
889  ----------
890  footprint : TODO: DM-17458
891  Footprint containing the dipole that was fit
892  result : `lmfit.MinimizerResult`
893  `lmfit.MinimizerResult` object returned by `lmfit` optimizer
894 
895  Returns
896  -------
897  fig : `matplotlib.pyplot.plot`
898  """
899  try:
900  import matplotlib.pyplot as plt
901  except ImportError as err:
902  self.log.warning('Unable to import matplotlib: %s', err)
903  raise err
904 
905  def display2dArray(arr, title='Data', extent=None):
906  """Use `matplotlib.pyplot.imshow` to display a 2-D array with a given coordinate range.
907  """
908  fig = plt.imshow(arr, origin='lower', interpolation='none', cmap='gray', extent=extent)
909  plt.title(title)
910  plt.colorbar(fig, cmap='gray')
911  return fig
912 
913  z = result.data
914  fit = result.best_fit
915  bbox = footprint.getBBox()
916  extent = (bbox.getBeginX(), bbox.getEndX(), bbox.getBeginY(), bbox.getEndY())
917  if z.shape[0] == 3:
918  fig = plt.figure(figsize=(8, 8))
919  for i in range(3):
920  plt.subplot(3, 3, i*3+1)
921  display2dArray(z[i, :], 'Data', extent=extent)
922  plt.subplot(3, 3, i*3+2)
923  display2dArray(fit[i, :], 'Model', extent=extent)
924  plt.subplot(3, 3, i*3+3)
925  display2dArray(z[i, :] - fit[i, :], 'Residual', extent=extent)
926  return fig
927  else:
928  fig = plt.figure(figsize=(8, 2.5))
929  plt.subplot(1, 3, 1)
930  display2dArray(z, 'Data', extent=extent)
931  plt.subplot(1, 3, 2)
932  display2dArray(fit, 'Model', extent=extent)
933  plt.subplot(1, 3, 3)
934  display2dArray(z - fit, 'Residual', extent=extent)
935  return fig
936 
937  plt.show()
938 
939 
940 @measBase.register("ip_diffim_DipoleFit")

◆ fitDipole()

def lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.fitDipole (   self,
  source,
  tol = 1e-7,
  rel_weight = 0.1,
  fitBackground = 1,
  maxSepInSigma = 5.,
  separateNegParams = True,
  bgGradientOrder = 1,
  verbose = False,
  display = False 
)
Fit a dipole model to an input ``diaSource`` (wraps `fitDipoleImpl`).

Actually, fits the subimage bounded by the input source's
footprint) and optionally constrain the fit using the
pre-subtraction images self.posImage (science) and
self.negImage (template). Wraps the output into a
`pipeBase.Struct` named tuple after computing additional
statistics such as orientation and SNR.

Parameters
----------
source : `lsst.afw.table.SourceRecord`
    Record containing the (merged) dipole source footprint detected on the diffim
tol : `float`, optional
    Tolerance parameter for scipy.leastsq() optimization
rel_weight : `float`, optional
    Weighting of posImage/negImage relative to the diffim in the fit
fitBackground : `int`, {0, 1, 2}, optional
    How to fit linear background gradient in posImage/negImage

        - 0: do not fit background at all
        - 1 (default): pre-fit the background using linear least squares and then do not fit it
          as part of the dipole fitting optimization
        - 2: pre-fit the background using linear least squares (as in 1), and use the parameter
          estimates from that fit as starting parameters for an integrated "re-fit" of the
          background as part of the overall dipole fitting optimization.
maxSepInSigma : `float`, optional
    Allowed window of centroid parameters relative to peak in input source footprint
separateNegParams : `bool`, optional
    Fit separate parameters to the flux and background gradient in
bgGradientOrder : `int`, {0, 1, 2}, optional
    Desired polynomial order of background gradient
verbose: `bool`, optional
    Be verbose
display
    Display input data, best fit model(s) and residuals in a matplotlib window.

Returns
-------
result : `struct`
    `pipeBase.Struct` object containing the fit parameters and other information.

result : `callable`
    `lmfit.MinimizerResult` object for debugging and error estimation, etc.

Notes
-----
Parameter `fitBackground` has three options, thus it is an integer:

Definition at line 772 of file dipoleFitTask.py.

774  bgGradientOrder=1, verbose=False, display=False):
775  """Fit a dipole model to an input ``diaSource`` (wraps `fitDipoleImpl`).
776 
777  Actually, fits the subimage bounded by the input source's
778  footprint) and optionally constrain the fit using the
779  pre-subtraction images self.posImage (science) and
780  self.negImage (template). Wraps the output into a
781  `pipeBase.Struct` named tuple after computing additional
782  statistics such as orientation and SNR.
783 
784  Parameters
785  ----------
786  source : `lsst.afw.table.SourceRecord`
787  Record containing the (merged) dipole source footprint detected on the diffim
788  tol : `float`, optional
789  Tolerance parameter for scipy.leastsq() optimization
790  rel_weight : `float`, optional
791  Weighting of posImage/negImage relative to the diffim in the fit
792  fitBackground : `int`, {0, 1, 2}, optional
793  How to fit linear background gradient in posImage/negImage
794 
795  - 0: do not fit background at all
796  - 1 (default): pre-fit the background using linear least squares and then do not fit it
797  as part of the dipole fitting optimization
798  - 2: pre-fit the background using linear least squares (as in 1), and use the parameter
799  estimates from that fit as starting parameters for an integrated "re-fit" of the
800  background as part of the overall dipole fitting optimization.
801  maxSepInSigma : `float`, optional
802  Allowed window of centroid parameters relative to peak in input source footprint
803  separateNegParams : `bool`, optional
804  Fit separate parameters to the flux and background gradient in
805  bgGradientOrder : `int`, {0, 1, 2}, optional
806  Desired polynomial order of background gradient
807  verbose: `bool`, optional
808  Be verbose
809  display
810  Display input data, best fit model(s) and residuals in a matplotlib window.
811 
812  Returns
813  -------
814  result : `struct`
815  `pipeBase.Struct` object containing the fit parameters and other information.
816 
817  result : `callable`
818  `lmfit.MinimizerResult` object for debugging and error estimation, etc.
819 
820  Notes
821  -----
822  Parameter `fitBackground` has three options, thus it is an integer:
823 
824  """
825 
826  fitResult = self.fitDipoleImpl(
827  source, tol=tol, rel_weight=rel_weight, fitBackground=fitBackground,
828  maxSepInSigma=maxSepInSigma, separateNegParams=separateNegParams,
829  bgGradientOrder=bgGradientOrder, verbose=verbose)
830 
831  # Display images, model fits and residuals (currently uses matplotlib display functions)
832  if display:
833  fp = source.getFootprint()
834  self.displayFitResults(fp, fitResult)
835 
836  fitParams = fitResult.best_values
837  if fitParams['flux'] <= 1.: # usually around 0.1 -- the minimum flux allowed -- i.e. bad fit.
838  out = Struct(posCentroidX=np.nan, posCentroidY=np.nan,
839  negCentroidX=np.nan, negCentroidY=np.nan,
840  posFlux=np.nan, negFlux=np.nan, posFluxErr=np.nan, negFluxErr=np.nan,
841  centroidX=np.nan, centroidY=np.nan, orientation=np.nan,
842  signalToNoise=np.nan, chi2=np.nan, redChi2=np.nan)
843  return out, fitResult
844 
845  centroid = ((fitParams['xcenPos'] + fitParams['xcenNeg']) / 2.,
846  (fitParams['ycenPos'] + fitParams['ycenNeg']) / 2.)
847  dx, dy = fitParams['xcenPos'] - fitParams['xcenNeg'], fitParams['ycenPos'] - fitParams['ycenNeg']
848  angle = np.arctan2(dy, dx) / np.pi * 180. # convert to degrees (should keep as rad?)
849 
850  # Exctract flux value, compute signalToNoise from flux/variance_within_footprint
851  # Also extract the stderr of flux estimate.
852  def computeSumVariance(exposure, footprint):
853  box = footprint.getBBox()
854  subim = afwImage.MaskedImageF(exposure.getMaskedImage(), box, origin=afwImage.PARENT)
855  return np.sqrt(np.nansum(subim.getArrays()[1][:, :]))
856 
857  fluxVal = fluxVar = fitParams['flux']
858  fluxErr = fluxErrNeg = fitResult.params['flux'].stderr
859  if self.posImage is not None:
860  fluxVar = computeSumVariance(self.posImage, source.getFootprint())
861  else:
862  fluxVar = computeSumVariance(self.diffim, source.getFootprint())
863 
864  fluxValNeg, fluxVarNeg = fluxVal, fluxVar
865  if separateNegParams:
866  fluxValNeg = fitParams['fluxNeg']
867  fluxErrNeg = fitResult.params['fluxNeg'].stderr
868  if self.negImage is not None:
869  fluxVarNeg = computeSumVariance(self.negImage, source.getFootprint())
870 
871  try:
872  signalToNoise = np.sqrt((fluxVal/fluxVar)**2 + (fluxValNeg/fluxVarNeg)**2)
873  except ZeroDivisionError: # catch divide by zero - should never happen.
874  signalToNoise = np.nan
875 
876  out = Struct(posCentroidX=fitParams['xcenPos'], posCentroidY=fitParams['ycenPos'],
877  negCentroidX=fitParams['xcenNeg'], negCentroidY=fitParams['ycenNeg'],
878  posFlux=fluxVal, negFlux=-fluxValNeg, posFluxErr=fluxErr, negFluxErr=fluxErrNeg,
879  centroidX=centroid[0], centroidY=centroid[1], orientation=angle,
880  signalToNoise=signalToNoise, chi2=fitResult.chisqr, redChi2=fitResult.redchi)
881 
882  # fitResult may be returned for debugging
883  return out, fitResult
884 

◆ fitDipoleImpl()

def lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.fitDipoleImpl (   self,
  source,
  tol = 1e-7,
  rel_weight = 0.5,
  fitBackground = 1,
  bgGradientOrder = 1,
  maxSepInSigma = 5.,
  separateNegParams = True,
  verbose = False 
)
Fit a dipole model to an input difference image.

Actually, fits the subimage bounded by the input source's
footprint) and optionally constrain the fit using the
pre-subtraction images posImage and negImage.

Parameters
----------
source : TODO: DM-17458
    TODO: DM-17458
tol : float, optional
    TODO: DM-17458
rel_weight : `float`, optional
    TODO: DM-17458
fitBackground : `int`, optional
    TODO: DM-17458
bgGradientOrder : `int`, optional
    TODO: DM-17458
maxSepInSigma : `float`, optional
    TODO: DM-17458
separateNegParams : `bool`, optional
    TODO: DM-17458
verbose : `bool`, optional
    TODO: DM-17458

Returns
-------
result : `lmfit.MinimizerResult`
    return `lmfit.MinimizerResult` object containing the fit
    parameters and other information.

Definition at line 535 of file dipoleFitTask.py.

537  separateNegParams=True, verbose=False):
538  """Fit a dipole model to an input difference image.
539 
540  Actually, fits the subimage bounded by the input source's
541  footprint) and optionally constrain the fit using the
542  pre-subtraction images posImage and negImage.
543 
544  Parameters
545  ----------
546  source : TODO: DM-17458
547  TODO: DM-17458
548  tol : float, optional
549  TODO: DM-17458
550  rel_weight : `float`, optional
551  TODO: DM-17458
552  fitBackground : `int`, optional
553  TODO: DM-17458
554  bgGradientOrder : `int`, optional
555  TODO: DM-17458
556  maxSepInSigma : `float`, optional
557  TODO: DM-17458
558  separateNegParams : `bool`, optional
559  TODO: DM-17458
560  verbose : `bool`, optional
561  TODO: DM-17458
562 
563  Returns
564  -------
565  result : `lmfit.MinimizerResult`
566  return `lmfit.MinimizerResult` object containing the fit
567  parameters and other information.
568  """
569 
570  # Only import lmfit if someone wants to use the new DipoleFitAlgorithm.
571  import lmfit
572 
573  fp = source.getFootprint()
574  bbox = fp.getBBox()
575  subim = afwImage.MaskedImageF(self.diffim.getMaskedImage(), bbox=bbox, origin=afwImage.PARENT)
576 
577  z = diArr = subim.getArrays()[0]
578  weights = 1. / subim.getArrays()[2] # get the weights (=1/variance)
579 
580  if rel_weight > 0. and ((self.posImage is not None) or (self.negImage is not None)):
581  if self.negImage is not None:
582  negSubim = afwImage.MaskedImageF(self.negImage.getMaskedImage(), bbox, origin=afwImage.PARENT)
583  if self.posImage is not None:
584  posSubim = afwImage.MaskedImageF(self.posImage.getMaskedImage(), bbox, origin=afwImage.PARENT)
585  if self.posImage is None: # no science image provided; generate it from diffim + negImage
586  posSubim = subim.clone()
587  posSubim += negSubim
588  if self.negImage is None: # no template provided; generate it from the posImage - diffim
589  negSubim = posSubim.clone()
590  negSubim -= subim
591 
592  z = np.append([z], [posSubim.getArrays()[0],
593  negSubim.getArrays()[0]], axis=0)
594  # Weight the pos/neg images by rel_weight relative to the diffim
595  weights = np.append([weights], [1. / posSubim.getArrays()[2] * rel_weight,
596  1. / negSubim.getArrays()[2] * rel_weight], axis=0)
597  else:
598  rel_weight = 0. # a short-cut for "don't include the pre-subtraction data"
599 
600  # It seems that `lmfit` requires a static functor as its optimized method, which eliminates
601  # the ability to pass a bound method or other class method. Here we write a wrapper which
602  # makes this possible.
603  def dipoleModelFunctor(x, flux, xcenPos, ycenPos, xcenNeg, ycenNeg, fluxNeg=None,
604  b=None, x1=None, y1=None, xy=None, x2=None, y2=None,
605  bNeg=None, x1Neg=None, y1Neg=None, xyNeg=None, x2Neg=None, y2Neg=None,
606  **kwargs):
607  """Generate dipole model with given parameters.
608 
609  It simply defers to `modelObj.makeModel()`, where `modelObj` comes
610  out of `kwargs['modelObj']`.
611  """
612  modelObj = kwargs.pop('modelObj')
613  return modelObj.makeModel(x, flux, xcenPos, ycenPos, xcenNeg, ycenNeg, fluxNeg=fluxNeg,
614  b=b, x1=x1, y1=y1, xy=xy, x2=x2, y2=y2,
615  bNeg=bNeg, x1Neg=x1Neg, y1Neg=y1Neg, xyNeg=xyNeg,
616  x2Neg=x2Neg, y2Neg=y2Neg, **kwargs)
617 
618  dipoleModel = DipoleModel()
619 
620  modelFunctor = dipoleModelFunctor # dipoleModel.makeModel does not work for now.
621  # Create the lmfit model (lmfit uses scipy 'leastsq' option by default - Levenberg-Marquardt)
622  # Note we can also tell it to drop missing values from the data.
623  gmod = lmfit.Model(modelFunctor, verbose=verbose, missing='drop')
624  # independent_vars=independent_vars) #, param_names=param_names)
625 
626  # Add the constraints for centroids, fluxes.
627  # starting constraint - near centroid of footprint
628  fpCentroid = np.array([fp.getCentroid().getX(), fp.getCentroid().getY()])
629  cenNeg = cenPos = fpCentroid
630 
631  pks = fp.getPeaks()
632 
633  if len(pks) >= 1:
634  cenPos = pks[0].getF() # if individual (merged) peaks were detected, use those
635  if len(pks) >= 2: # peaks are already sorted by centroid flux so take the most negative one
636  cenNeg = pks[-1].getF()
637 
638  # For close/faint dipoles the starting locs (min/max) might be way off, let's help them a bit.
639  # First assume dipole is not separated by more than 5*psfSigma.
640  maxSep = self.psfSigma * maxSepInSigma
641 
642  # As an initial guess -- assume the dipole is close to the center of the footprint.
643  if np.sum(np.sqrt((np.array(cenPos) - fpCentroid)**2.)) > maxSep:
644  cenPos = fpCentroid
645  if np.sum(np.sqrt((np.array(cenNeg) - fpCentroid)**2.)) > maxSep:
646  cenPos = fpCentroid
647 
648  # parameter hints/constraints: https://lmfit.github.io/lmfit-py/model.html#model-param-hints-section
649  # might make sense to not use bounds -- see http://lmfit.github.io/lmfit-py/bounds.html
650  # also see this discussion -- https://github.com/scipy/scipy/issues/3129
651  gmod.set_param_hint('xcenPos', value=cenPos[0],
652  min=cenPos[0]-maxSep, max=cenPos[0]+maxSep)
653  gmod.set_param_hint('ycenPos', value=cenPos[1],
654  min=cenPos[1]-maxSep, max=cenPos[1]+maxSep)
655  gmod.set_param_hint('xcenNeg', value=cenNeg[0],
656  min=cenNeg[0]-maxSep, max=cenNeg[0]+maxSep)
657  gmod.set_param_hint('ycenNeg', value=cenNeg[1],
658  min=cenNeg[1]-maxSep, max=cenNeg[1]+maxSep)
659 
660  # Use the (flux under the dipole)*5 for an estimate.
661  # Lots of testing showed that having startingFlux be too high was better than too low.
662  startingFlux = np.nansum(np.abs(diArr) - np.nanmedian(np.abs(diArr))) * 5.
663  posFlux = negFlux = startingFlux
664 
665  # TBD: set max. flux limit?
666  gmod.set_param_hint('flux', value=posFlux, min=0.1)
667 
668  if separateNegParams:
669  # TBD: set max negative lobe flux limit?
670  gmod.set_param_hint('fluxNeg', value=np.abs(negFlux), min=0.1)
671 
672  # Fixed parameters (don't fit for them if there are no pre-sub images or no gradient fit requested):
673  # Right now (fitBackground == 1), we fit a linear model to the background and then subtract
674  # it from the data and then don't fit the background again (this is faster).
675  # A slower alternative (fitBackground == 2) is to use the estimated background parameters as
676  # starting points in the integrated model fit. That is currently not performed by default,
677  # but might be desirable in some cases.
678  bgParsPos = bgParsNeg = (0., 0., 0.)
679  if ((rel_weight > 0.) and (fitBackground != 0) and (bgGradientOrder >= 0)):
680  pbg = 0.
681  bgFitImage = self.posImage if self.posImage is not None else self.negImage
682  # Fit the gradient to the background (linear model)
683  bgParsPos = bgParsNeg = dipoleModel.fitFootprintBackground(source, bgFitImage,
684  order=bgGradientOrder)
685 
686  # Generate the gradient and subtract it from the pre-subtraction image data
687  if fitBackground == 1:
688  in_x = dipoleModel._generateXYGrid(bbox)
689  pbg = dipoleModel.makeBackgroundModel(in_x, tuple(bgParsPos))
690  z[1, :] -= pbg
691  z[1, :] -= np.nanmedian(z[1, :])
692  posFlux = np.nansum(z[1, :])
693  gmod.set_param_hint('flux', value=posFlux*1.5, min=0.1)
694 
695  if separateNegParams and self.negImage is not None:
696  bgParsNeg = dipoleModel.fitFootprintBackground(source, self.negImage,
697  order=bgGradientOrder)
698  pbg = dipoleModel.makeBackgroundModel(in_x, tuple(bgParsNeg))
699  z[2, :] -= pbg
700  z[2, :] -= np.nanmedian(z[2, :])
701  if separateNegParams:
702  negFlux = np.nansum(z[2, :])
703  gmod.set_param_hint('fluxNeg', value=negFlux*1.5, min=0.1)
704 
705  # Do not subtract the background from the images but include the background parameters in the fit
706  if fitBackground == 2:
707  if bgGradientOrder >= 0:
708  gmod.set_param_hint('b', value=bgParsPos[0])
709  if separateNegParams:
710  gmod.set_param_hint('bNeg', value=bgParsNeg[0])
711  if bgGradientOrder >= 1:
712  gmod.set_param_hint('x1', value=bgParsPos[1])
713  gmod.set_param_hint('y1', value=bgParsPos[2])
714  if separateNegParams:
715  gmod.set_param_hint('x1Neg', value=bgParsNeg[1])
716  gmod.set_param_hint('y1Neg', value=bgParsNeg[2])
717  if bgGradientOrder >= 2:
718  gmod.set_param_hint('xy', value=bgParsPos[3])
719  gmod.set_param_hint('x2', value=bgParsPos[4])
720  gmod.set_param_hint('y2', value=bgParsPos[5])
721  if separateNegParams:
722  gmod.set_param_hint('xyNeg', value=bgParsNeg[3])
723  gmod.set_param_hint('x2Neg', value=bgParsNeg[4])
724  gmod.set_param_hint('y2Neg', value=bgParsNeg[5])
725 
726  y, x = np.mgrid[bbox.getBeginY():bbox.getEndY(), bbox.getBeginX():bbox.getEndX()]
727  in_x = np.array([x, y]).astype(np.float)
728  in_x[0, :] -= in_x[0, :].mean() # center it!
729  in_x[1, :] -= in_x[1, :].mean()
730 
731  # Instead of explicitly using a mask to ignore flagged pixels, just set the ignored pixels'
732  # weights to 0 in the fit. TBD: need to inspect mask planes to set this mask.
733  mask = np.ones_like(z, dtype=bool) # TBD: set mask values to False if the pixels are to be ignored
734 
735  # I'm not sure about the variance planes in the diffim (or convolved pre-sub. images
736  # for that matter) so for now, let's just do an un-weighted least-squares fit
737  # (override weights computed above).
738  weights = mask.astype(np.float64)
739  if self.posImage is not None and rel_weight > 0.:
740  weights = np.array([np.ones_like(diArr), np.ones_like(diArr)*rel_weight,
741  np.ones_like(diArr)*rel_weight])
742 
743  # Set the weights to zero if mask is False
744  if np.any(~mask):
745  weights[~mask] = 0.
746 
747  # Note that although we can, we're not required to set initial values for params here,
748  # since we set their param_hint's above.
749  # Can add "method" param to not use 'leastsq' (==levenberg-marquardt), e.g. "method='nelder'"
750  with warnings.catch_warnings():
751  warnings.simplefilter("ignore") # temporarily turn off silly lmfit warnings
752  result = gmod.fit(z, weights=weights, x=in_x,
753  verbose=verbose,
754  fit_kws={'ftol': tol, 'xtol': tol, 'gtol': tol,
755  'maxfev': 250}, # see scipy docs
756  psf=self.diffim.getPsf(), # hereon: kwargs that get passed to genDipoleModel()
757  rel_weight=rel_weight,
758  footprint=fp,
759  modelObj=dipoleModel)
760 
761  if verbose: # the ci_report() seems to fail if neg params are constrained -- TBD why.
762  # Never wanted in production - this takes a long time (longer than the fit!)
763  # This is how to get confidence intervals out:
764  # https://lmfit.github.io/lmfit-py/confidence.html and
765  # http://cars9.uchicago.edu/software/python/lmfit/model.html
766  print(result.fit_report(show_correl=False))
767  if separateNegParams:
768  print(result.ci_report())
769 
770  return result
771 

Member Data Documentation

◆ debug

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.debug

Definition at line 533 of file dipoleFitTask.py.

◆ diffim

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.diffim

Definition at line 523 of file dipoleFitTask.py.

◆ log

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.log

Definition at line 530 of file dipoleFitTask.py.

◆ negImage

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.negImage

Definition at line 525 of file dipoleFitTask.py.

◆ posImage

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.posImage

Definition at line 524 of file dipoleFitTask.py.

◆ psfSigma

lsst.ip.diffim.dipoleFitTask.DipoleFitAlgorithm.psfSigma

Definition at line 526 of file dipoleFitTask.py.


The documentation for this class was generated from the following file: