lsst.cp.pipe  21.0.0-5-gb7080ec+f575c71e2b
astierCovPtcFit.py
Go to the documentation of this file.
1 # This file is part of cp_pipe.
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 import numpy as np
23 import copy
24 import itertools
25 from scipy.signal import fftconvolve
26 from scipy.optimize import leastsq
27 from .astierCovFitParameters import FitParameters
28 
29 import lsst.log as lsstLog
30 
31 __all__ = ["CovFit"]
32 
33 
34 def computeApproximateAcoeffs(covModel, muEl, gain):
35  """Compute the "a" coefficients of the Antilogus+14 (1402.0725) model as in
36  Guyonnet+15 (1501.01577, eq. 16, the slope of cov/var at a given flux mu in electrons).
37 
38  Eq. 16 of 1501.01577 is an approximation to the more complete model in Astier+19 (1905.08677).
39 
40  Parameters
41  ---------
42  covModel : `list`
43  Covariance model from Eq. 20 in Astier+19.
44 
45  muEl : `np.array`
46  Mean signal in electrons
47 
48  gain : `float`
49  Gain in e-/ADU.
50 
51  Returns
52  -------
53  aCoeffsOld: `numpy.array`
54  Slope of cov/var at a given flux mu in electrons.
55 
56  Notes
57  -----
58  Returns the "a" array, computed this way, to be compared to the actual a_array from the full model
59  (fit.geA()).
60  """
61  covModel = np.array(covModel)
62  var = covModel[0, 0, 0] # ADU^2
63  # For a result in electrons^-1, we have to use mu in electrons.
64  return covModel[0, :, :]/(var*muEl)
65 
66 
67 def makeCovArray(inputTuple, maxRangeFromTuple=8):
68  """Make covariances array from tuple.
69 
70  Parameters
71  ----------
72  inputTuple: `numpy.recarray`
73  Recarray with rows with at least (mu, cov, var, i, j, npix), where:
74  mu : 0.5*(m1 + m2), where:
75  mu1: mean value of flat1
76  mu2: mean value of flat2
77  cov: covariance value at lag(i, j)
78  var: variance(covariance value at lag(0, 0))
79  i: lag dimension
80  j: lag dimension
81  npix: number of pixels used for covariance calculation.
82 
83  maxRangeFromTuple: `int`
84  Maximum range to select from tuple.
85 
86  Returns
87  -------
88  cov: `numpy.array`
89  Covariance arrays, indexed by mean signal mu.
90 
91  vCov: `numpy.array`
92  Variance arrays, indexed by mean signal mu.
93 
94  muVals: `numpy.array`
95  List of mean signal values.
96 
97  Notes
98  -----
99 
100  The input tuple should contain the following rows:
101  (mu, cov, var, i, j, npix), with one entry per lag, and image pair.
102  Different lags(i.e. different i and j) from the same
103  image pair have the same values of mu1 and mu2. When i==j==0, cov
104  = var.
105 
106  If the input tuple contains several video channels, one should
107  select the data of a given channel *before* entering this
108  routine, as well as apply(e.g.) saturation cuts.
109 
110  The routine returns cov[k_mu, j, i], vcov[(same indices)], and mu[k]
111  where the first index of cov matches the one in mu.
112 
113  This routine implements the loss of variance due to
114  clipping cuts when measuring variances and covariance, but this should happen inside
115  the measurement code, where the cuts are readily available.
116 
117  """
118  if maxRangeFromTuple is not None:
119  cut = (inputTuple['i'] < maxRangeFromTuple) & (inputTuple['j'] < maxRangeFromTuple)
120  cutTuple = inputTuple[cut]
121  else:
122  cutTuple = inputTuple
123  # increasing mu order, so that we can group measurements with the same mu
124  muTemp = cutTuple['mu']
125  ind = np.argsort(muTemp)
126 
127  cutTuple = cutTuple[ind]
128  # should group measurements on the same image pairs(same average)
129  mu = cutTuple['mu']
130  xx = np.hstack(([mu[0]], mu))
131  delta = xx[1:] - xx[:-1]
132  steps, = np.where(delta > 0)
133  ind = np.zeros_like(mu, dtype=int)
134  ind[steps] = 1
135  ind = np.cumsum(ind) # this acts as an image pair index.
136  # now fill the 3-d cov array(and variance)
137  muVals = np.array(np.unique(mu))
138  i = cutTuple['i'].astype(int)
139  j = cutTuple['j'].astype(int)
140  c = 0.5*cutTuple['cov']
141  n = cutTuple['npix']
142  v = 0.5*cutTuple['var']
143  # book and fill
144  cov = np.ndarray((len(muVals), np.max(i)+1, np.max(j)+1))
145  var = np.zeros_like(cov)
146  cov[ind, i, j] = c
147  var[ind, i, j] = v**2/n
148  var[:, 0, 0] *= 2 # var(v) = 2*v**2/N
149  # compensate for loss of variance and covariance due to outlier elimination(sigma clipping)
150  # when computing variances(cut to 4 sigma): 1 per mill for variances and twice as
151  # much for covariances:
152  fact = 1.0 # 1.00107
153  cov *= fact*fact
154  cov[:, 0, 0] /= fact
155 
156  return cov, var, muVals
157 
158 
159 def symmetrize(inputArray):
160  """ Copy array over 4 quadrants prior to convolution.
161 
162  Parameters
163  ----------
164  inputarray: `numpy.array`
165  Input array to symmetrize.
166 
167  Returns
168  -------
169  aSym: `numpy.array`
170  Symmetrized array.
171  """
172 
173  targetShape = list(inputArray.shape)
174  r1, r2 = inputArray.shape[-1], inputArray.shape[-2]
175  targetShape[-1] = 2*r1-1
176  targetShape[-2] = 2*r2-1
177  aSym = np.ndarray(tuple(targetShape))
178  aSym[..., r2-1:, r1-1:] = inputArray
179  aSym[..., r2-1:, r1-1::-1] = inputArray
180  aSym[..., r2-1::-1, r1-1::-1] = inputArray
181  aSym[..., r2-1::-1, r1-1:] = inputArray
182 
183  return aSym
184 
185 
186 class Pol2d:
187  """A class to calculate 2D polynomials"""
188 
189  def __init__(self, x, y, z, order, w=None):
190  self.orderx = min(order, x.shape[0]-1)
191  self.ordery = min(order, x.shape[1]-1)
192  G = self.monomials(x.ravel(), y.ravel())
193  if w is None:
194  self.coeff, _, rank, _ = np.linalg.lstsq(G, z.ravel())
195  else:
196  self.coeff, _, rank, _ = np.linalg.lstsq((w.ravel()*G.T).T, z.ravel()*w.ravel())
197 
198  def monomials(self, x, y):
199  ncols = (self.orderx+1)*(self.ordery+1)
200  G = np.zeros(x.shape + (ncols,))
201  ij = itertools.product(range(self.orderx+1), range(self.ordery+1))
202  for k, (i, j) in enumerate(ij):
203  G[..., k] = x**i * y**j
204 
205  return G
206 
207  def eval(self, x, y):
208  G = self.monomials(x, y)
209 
210  return np.dot(G, self.coeff)
211 
212 
213 class CovFit:
214  """A class to fit the models in Astier+19 to flat covariances.
215 
216  This code implements the model(and the fit thereof) described in
217  Astier+19: https://arxiv.org/pdf/1905.08677.pdf
218  For the time being it uses as input a numpy recarray (tuple with named tags) which
219  contains one row per covariance and per pair: see the routine makeCovArray.
220 
221  Parameters
222  ----------
223  inputTuple: `numpy.recarray`
224  Tuple with at least (mu, cov, var, i, j, npix), where:
225  mu : 0.5*(m1 + m2), where:
226  mu1: mean value of flat1
227  mu2: mean value of flat2
228  cov: covariance value at lag(i, j)
229  var: variance(covariance value at lag(0, 0))
230  i: lag dimension
231  j: lag dimension
232  npix: number of pixels used for covariance calculation.
233 
234  maxRangeFromTuple: `int`, optional
235  Maximum range to select from tuple.
236 
237  meanSignalMask: `list`[`bool`], optional
238  Mask of mean signal 1D array. Use all entries if empty.
239  """
240 
241  def __init__(self, inputTuple, maxRangeFromTuple=8, meanSignalMask=[]):
242  self.cov, self.vcov, self.mu = makeCovArray(inputTuple, maxRangeFromTuple)
243  # make it nan safe, replacing nan's with 0 in weights
244  self.sqrtW = np.nan_to_num(1./np.sqrt(self.vcov))
245  self.r = self.cov.shape[1]
246  self.logger = lsstLog.Log.getDefaultLogger()
247  if len(meanSignalMask):
248  self.maskMu = meanSignalMask
249  else:
250  self.maskMu = np.repeat(True, len(self.mu))
251 
252  def subtractDistantOffset(self, maxLag=8, startLag=5, polDegree=1):
253  """Subtract a background/offset to the measured covariances.
254 
255  Parameters
256  ---------
257  maxLag: `int`
258  Maximum lag considered
259 
260  startLag: `int`
261  First lag from where to start the offset subtraction.
262 
263  polDegree: `int`
264  Degree of 2D polynomial to fit to covariance to define offse to be subtracted.
265  """
266  assert(startLag < self.r)
267  for k in range(len(self.mu)):
268  # Make a copy because it is going to be altered
269  w = self.sqrtW[k, ...] + 0.
270  sh = w.shape
271  i, j = np.meshgrid(range(sh[0]), range(sh[1]), indexing='ij')
272  # kill the core for the fit
273  w[:startLag, :startLag] = 0
274  poly = Pol2d(i, j, self.cov[k, ...], polDegree+1, w=w)
275  back = poly.eval(i, j)
276  self.cov[k, ...] -= back
277  self.r = maxLag
278  self.cov = self.cov[:, :maxLag, :maxLag]
279  self.vcov = self.vcov[:, :maxLag, :maxLag]
280  self.sqrtW = self.sqrtW[:, :maxLag, :maxLag]
281 
282  return
283 
284  def copy(self):
285  """Make a copy of params"""
286  cop = copy.deepcopy(self)
287  # deepcopy does not work for FitParameters.
288  if hasattr(self, 'params'):
289  cop.params = self.params.copy()
290  return cop
291 
292  def initFit(self):
293  """ Performs a crude parabolic fit of the data in order to start
294  the full fit close to the solution.
295  """
296  # number of parameters for 'a' array.
297  lenA = self.r*self.r
298  # define parameters: c corresponds to a*b in Astier+19 (Eq. 20).
299  self.params = FitParameters([('a', lenA), ('c', lenA), ('noise', lenA), ('gain', 1)])
300  self.params['gain'] = 1.
301  # c=0 in a first go.
302  self.params['c'].fix(val=0.)
303  # plumbing: extract stuff from the parameter structure
304  a = self.params['a'].full.reshape(self.r, self.r)
305  noise = self.params['noise'].full.reshape(self.r, self.r)
306  gain = self.params['gain'].full[0]
307 
308  # iterate the fit to account for higher orders
309  # the chi2 does not necessarily go down, so one could
310  # stop when it increases
311  oldChi2 = 1e30
312  for _ in range(5):
313  model = self.evalCovModel() # this computes the full model.
314  # loop on lags
315  for i in range(self.r):
316  for j in range(self.r):
317  # fit a parabola for a given lag
318  parsFit = np.polyfit(self.mu, self.cov[:, i, j] - model[:, i, j],
319  2, w=self.sqrtW[:, i, j])
320  # model equation(Eq. 20) in Astier+19:
321  a[i, j] += parsFit[0]
322  noise[i, j] += parsFit[2]*gain*gain
323  if(i + j == 0):
324  gain = 1./(1/gain+parsFit[1])
325  self.params['gain'].full[0] = gain
326  chi2 = self.chi2()
327  if chi2 > oldChi2:
328  break
329  oldChi2 = chi2
330 
331  return
332 
333  def getParamValues(self):
334  """Return an array of free parameter values (it is a copy)."""
335  return self.params.free + 0.
336 
337  def setParamValues(self, p):
338  """Set parameter values."""
339  self.params.free = p
340  return
341 
342  def evalCovModel(self, mu=None):
343  """Computes full covariances model (Eq. 20 of Astier+19).
344 
345  Parameters
346  ----------
347  mu: `numpy.array`, optional
348  List of mean signals.
349 
350  Returns
351  -------
352  covModel: `numpy.array`
353  Covariances model.
354 
355  Notes
356  -----
357  By default, computes the covModel for the mu's stored(self.mu).
358 
359  Returns cov[Nmu, self.r, self.r]. The variance for the PTC is cov[:, 0, 0].
360  mu and cov are in ADUs and ADUs squared. To use electrons for both,
361  the gain should be set to 1. This routine implements the model in Astier+19 (1905.08677).
362 
363  The parameters of the full model for C_ij(mu) ("C_ij" and "mu" in ADU^2 and ADU, respectively)
364  in Astier+19 (Eq. 20) are:
365 
366  "a" coefficients (r by r matrix), units: 1/e
367  "b" coefficients (r by r matrix), units: 1/e
368  noise matrix (r by r matrix), units: e^2
369  gain, units: e/ADU
370 
371  "b" appears in Eq. 20 only through the "ab" combination, which is defined in this code as "c=ab".
372  """
373  sa = (self.r, self.r)
374  a = self.params['a'].full.reshape(sa)
375  c = self.params['c'].full.reshape(sa)
376  gain = self.params['gain'].full[0]
377  noise = self.params['noise'].full.reshape(sa)
378  # pad a with zeros and symmetrize
379  aEnlarged = np.zeros((int(sa[0]*1.5)+1, int(sa[1]*1.5)+1))
380  aEnlarged[0:sa[0], 0:sa[1]] = a
381  aSym = symmetrize(aEnlarged)
382  # pad c with zeros and symmetrize
383  cEnlarged = np.zeros((int(sa[0]*1.5)+1, int(sa[1]*1.5)+1))
384  cEnlarged[0:sa[0], 0:sa[1]] = c
385  cSym = symmetrize(cEnlarged)
386  a2 = fftconvolve(aSym, aSym, mode='same')
387  a3 = fftconvolve(a2, aSym, mode='same')
388  ac = fftconvolve(aSym, cSym, mode='same')
389  (xc, yc) = np.unravel_index(np.abs(aSym).argmax(), a2.shape)
390  range = self.r
391  a1 = a[np.newaxis, :, :]
392  a2 = a2[np.newaxis, xc:xc + range, yc:yc + range]
393  a3 = a3[np.newaxis, xc:xc + range, yc:yc + range]
394  ac = ac[np.newaxis, xc:xc + range, yc:yc + range]
395  c1 = c[np.newaxis, ::]
396  if mu is None:
397  mu = self.mu
398  # assumes that mu is 1d
399  bigMu = mu[:, np.newaxis, np.newaxis]*gain
400  # c(=a*b in Astier+19) also has a contribution to the last term, that is absent for now.
401  covModel = (bigMu/(gain*gain)*(a1*bigMu+2./3.*(bigMu*bigMu)*(a2 + c1) +
402  (1./3.*a3 + 5./6.*ac)*(bigMu*bigMu*bigMu)) + noise[np.newaxis, :, :]/gain**2)
403  # add the Poisson term, and the read out noise (variance)
404  covModel[:, 0, 0] += mu/gain
405 
406  return covModel
407 
408  def getA(self):
409  """'a' matrix from Astier+19(e.g., Eq. 20)"""
410  return self.params['a'].full.reshape(self.r, self.r)
411 
412  def getB(self):
413  """'b' matrix from Astier+19(e.g., Eq. 20)"""
414  return self.params['c'].full.reshape(self.r, self.r)/self.getA()
415 
416  def getC(self):
417  """'c'='ab' matrix from Astier+19(e.g., Eq. 20)"""
418  return np.array(self.params['c'].full.reshape(self.r, self.r))
419 
420  def _getCovParams(self, what):
421  """Get covariance matrix of parameters from fit"""
422  indices = self.params[what].indexof()
423  i1 = indices[:, np.newaxis]
424  i2 = indices[np.newaxis, :]
425  if self.covParams is not None:
426  covp = self.covParams[i1, i2]
427  else:
428  covp = None
429  return covp
430 
431  def getACov(self):
432  """Get covariance matrix of "a" coefficients from fit"""
433  if self._getCovParams('a') is not None:
434  cova = self._getCovParams('a').reshape((self.r, self.r, self.r, self.r))
435  else:
436  cova = None
437  return cova
438 
439  def getASig(self):
440  """Square root of diagonal of the parameter covariance of the fitted "a" matrix"""
441  if self._getCovParams('a') is not None:
442  sigA = np.sqrt(self._getCovParams('a').diagonal()).reshape((self.r, self.r))
443  else:
444  sigA = None
445  return sigA
446 
447  def getBCov(self):
448  """Get covariance matrix of "a" coefficients from fit
449  b = c /a
450  """
451  covb = self._getCovParams('c')
452  aval = self.getA().flatten()
453  factor = np.outer(aval, aval)
454  covb /= factor
455  return covb.reshape((self.r, self.r, self.r, self.r))
456 
457  def getCCov(self):
458  """Get covariance matrix of "c" coefficients from fit"""
459  cova = self._getCovParams('c')
460  return cova.reshape((self.r, self.r, self.r, self.r))
461 
462  def getGainErr(self):
463  """Get error on fitted gain parameter"""
464  if self._getCovParams('gain') is not None:
465  gainErr = np.sqrt(self._getCovParams('gain')[0][0])
466  else:
467  gainErr = 0.0
468  return gainErr
469 
470  def getNoiseCov(self):
471  """Get covariances of noise matrix from fit"""
472  covNoise = self._getCovParams('noise')
473  return covNoise.reshape((self.r, self.r, self.r, self.r))
474 
475  def getNoiseSig(self):
476  """Square root of diagonal of the parameter covariance of the fitted "noise" matrix"""
477  if self._getCovParams('noise') is not None:
478  covNoise = self._getCovParams('noise')
479  noise = np.sqrt(covNoise.diagonal()).reshape((self.r, self.r))
480  else:
481  noise = None
482  return noise
483 
484  def getGain(self):
485  """Get gain (e/ADU)"""
486  return self.params['gain'].full[0]
487 
488  def getRon(self):
489  """Get readout noise (e^2)"""
490  return self.params['noise'].full[0]
491 
492  def getRonErr(self):
493  """Get error on readout noise parameter"""
494  ronSqrt = np.sqrt(np.fabs(self.getRon()))
495  if self.getNoiseSig() is not None:
496  noiseSigma = self.getNoiseSig()[0][0]
497  ronErr = 0.5*(noiseSigma/np.fabs(self.getRon()))*ronSqrt
498  else:
499  ronErr = np.nan
500  return ronErr
501 
502  def getNoise(self):
503  """Get noise matrix"""
504  return self.params['noise'].full.reshape(self.r, self.r)
505 
506  def getMaskCov(self, i, j):
507  """Get mask of Cov[i,j]"""
508  weights = self.sqrtW[:, i, j]
509  mask = weights != 0
510  return mask
511 
512  def setAandB(self, a, b):
513  """Set "a" and "b" coeffcients forfull Astier+19 model (Eq. 20). "c=a*b"."""
514  self.params['a'].full = a.flatten()
515  self.params['c'].full = a.flatten()*b.flatten()
516  return
517 
518  def chi2(self):
519  """Calculate weighted chi2 of full-model fit."""
520  return(self.weightedRes()**2).sum()
521 
522  def wres(self, params=None):
523  """To be used in weightedRes"""
524  if params is not None:
525  self.setParamValues(params)
526  covModel = self.evalCovModel()
527  weightedRes = (covModel-self.cov)*self.sqrtW
528  maskedWeightedRes = weightedRes[self.maskMu]
529 
530  return maskedWeightedRes
531 
532  def weightedRes(self, params=None):
533  """Weighted residuas.
534 
535  Notes
536  -----
537  To be used via:
538  c = CovFit(nt)
539  c.initFit()
540  coeffs, cov, _, mesg, ierr = leastsq(c.weightedRes, c.getParamValues(), full_output=True)
541  """
542  return self.wres(params).flatten()
543 
544  def fitFullModel(self, pInit=None):
545  """Fit measured covariances to full model in Astier+19 (Eq. 20)
546 
547  Parameters
548  ----------
549  pInit : `list`
550  Initial parameters of the fit.
551  len(pInit) = #entries(a) + #entries(c) + #entries(noise) + 1
552  len(pInit) = r^2 + r^2 + r^2 + 1, where "r" is the maximum lag considered for the
553  covariances calculation, and the extra "1" is the gain.
554  If "b" is 0, then "c" is 0, and len(pInit) will have r^2 fewer entries.
555 
556  Returns
557  -------
558  params : `np.array`
559  Fit parameters (see "Notes" below).
560 
561  Notes
562  -----
563  The parameters of the full model for C_ij(mu) ("C_ij" and "mu" in ADU^2 and ADU, respectively)
564  in Astier+19 (Eq. 20) are:
565 
566  "a" coefficients (r by r matrix), units: 1/e
567  "b" coefficients (r by r matrix), units: 1/e
568  noise matrix (r by r matrix), units: e^2
569  gain, units: e/ADU
570 
571  "b" appears in Eq. 20 only through the "ab" combination, which is defined in this code as "c=ab".
572  """
573 
574  if pInit is None:
575  pInit = self.getParamValues()
576 
577  params, paramsCov, _, mesg, ierr = leastsq(self.weightedRes, pInit, full_output=True)
578  self.covParams = paramsCov
579 
580  return params
581 
582  def ndof(self):
583  """Number of degrees of freedom
584 
585  Returns
586  -------
587  mask.sum() - len(self.params.free): `int`
588  Number of usable pixels - number of parameters of fit.
589  """
590  mask = self.sqrtW != 0
591 
592  return mask.sum() - len(self.params.free)
593 
594  def getFitData(self, i, j, divideByMu=False, unitsElectrons=False, returnMasked=False):
595  """Get measured signal and covariance, cov model, weigths, and mask at covariance lag (i, j).
596 
597  Parameters
598  ---------
599  i: `int`
600  Lag for covariance matrix.
601 
602  j: `int`
603  Lag for covariance matrix.
604 
605  divideByMu: `bool`, optional
606  Divide covariance, model, and weights by signal mu?
607 
608  unitsElectrons : `bool`, optional
609  mu, covariance, and model are in ADU (or powers of ADU) If tthis parameter is true, these are
610  multiplied by the adequte factors of the gain to return quantities in electrons
611  (or powers of electrons).
612 
613  returnMasked : `bool`, optional
614  Use mask (based on weights) in returned arrays (mu, covariance, and model)?
615 
616  Returns
617  -------
618  mu: `numpy.array`
619  list of signal values (mu).
620 
621  covariance: `numpy.array`
622  Covariance arrays, indexed by mean signal mu (self.cov[:, i, j]).
623 
624  covarianceModel: `numpy.array`
625  Covariance model (model).
626 
627  weights: `numpy.array`
628  Weights (self.sqrtW)
629 
630  mask : `numpy.array`, optional
631  Boolean mask of the covariance at (i,j).
632 
633  Notes
634  -----
635  Using a CovFit object, selects from (i, j) and returns
636  mu*gain, self.cov[:, i, j]*gain**2 model*gain**2, and self.sqrtW/gain**2
637  in electrons or ADU if unitsElectrons=False.
638  """
639  if unitsElectrons:
640  gain = self.getGain()
641  else:
642  gain = 1.0
643 
644  mu = self.mu*gain
645  covariance = self.cov[:, i, j]*(gain**2)
646  covarianceModel = self.evalCovModel()[:, i, j]*(gain**2)
647  weights = self.sqrtW[:, i, j]/(gain**2)
648 
649  # select data used for the fit:
650  mask = self.getMaskCov(i, j)
651  if returnMasked:
652  weights = weights[mask]
653  covarianceModel = covarianceModel[mask]
654  mu = mu[mask]
655  covariance = covariance[mask]
656 
657  if divideByMu:
658  covariance /= mu
659  covarianceModel /= mu
660  weights *= mu
661 
662  return mu, covariance, covarianceModel, weights, mask
663 
664  def __call__(self, params):
665  self.setParamValues(params)
666  chi2 = self.chi2()
667 
668  return chi2
lsst.cp.pipe.astierCovPtcFit.CovFit.mu
mu
Definition: astierCovPtcFit.py:242
lsst.cp.pipe.astierCovPtcFit.CovFit.r
r
Definition: astierCovPtcFit.py:245
lsst.cp.pipe.astierCovPtcFit.CovFit.fitFullModel
def fitFullModel(self, pInit=None)
Definition: astierCovPtcFit.py:544
lsst.cp.pipe.astierCovPtcFit.CovFit.getNoiseSig
def getNoiseSig(self)
Definition: astierCovPtcFit.py:475
lsst.cp.pipe.astierCovPtcFit.CovFit.copy
def copy(self)
Definition: astierCovPtcFit.py:284
lsst.cp.pipe.astierCovPtcFit.makeCovArray
def makeCovArray(inputTuple, maxRangeFromTuple=8)
Definition: astierCovPtcFit.py:67
lsst.cp.pipe.astierCovPtcFit.CovFit.getBCov
def getBCov(self)
Definition: astierCovPtcFit.py:447
lsst.cp.pipe.astierCovPtcFit.Pol2d.__init__
def __init__(self, x, y, z, order, w=None)
Definition: astierCovPtcFit.py:189
lsst.cp.pipe.astierCovPtcFit.CovFit.getGain
def getGain(self)
Definition: astierCovPtcFit.py:484
lsst.cp.pipe.astierCovPtcFit.CovFit.getNoise
def getNoise(self)
Definition: astierCovPtcFit.py:502
lsst.cp.pipe.astierCovPtcFit.computeApproximateAcoeffs
def computeApproximateAcoeffs(covModel, muEl, gain)
Definition: astierCovPtcFit.py:34
lsst.cp.pipe.astierCovPtcFit.CovFit.params
params
Definition: astierCovPtcFit.py:299
lsst.cp.pipe.astierCovPtcFit.CovFit.logger
logger
Definition: astierCovPtcFit.py:246
lsst.cp.pipe.astierCovPtcFit.CovFit.getRon
def getRon(self)
Definition: astierCovPtcFit.py:488
lsst.cp.pipe.astierCovPtcFit.Pol2d
Definition: astierCovPtcFit.py:186
lsst.cp.pipe.astierCovPtcFit.CovFit.cov
cov
Definition: astierCovPtcFit.py:278
lsst.cp.pipe.astierCovPtcFit.CovFit.ndof
def ndof(self)
Definition: astierCovPtcFit.py:582
lsst.cp.pipe.astierCovPtcFit.CovFit.evalCovModel
def evalCovModel(self, mu=None)
Definition: astierCovPtcFit.py:342
lsst.cp.pipe.astierCovPtcFit.Pol2d.orderx
orderx
Definition: astierCovPtcFit.py:190
lsst.cp.pipe.astierCovPtcFit.CovFit.setAandB
def setAandB(self, a, b)
Definition: astierCovPtcFit.py:512
lsst.cp.pipe.astierCovPtcFit.CovFit.getASig
def getASig(self)
Definition: astierCovPtcFit.py:439
lsst.cp.pipe.astierCovPtcFit.CovFit.getParamValues
def getParamValues(self)
Definition: astierCovPtcFit.py:333
lsst.cp.pipe.astierCovPtcFit.CovFit.__init__
def __init__(self, inputTuple, maxRangeFromTuple=8, meanSignalMask=[])
Definition: astierCovPtcFit.py:241
lsst.cp.pipe.astierCovPtcFit.CovFit.subtractDistantOffset
def subtractDistantOffset(self, maxLag=8, startLag=5, polDegree=1)
Definition: astierCovPtcFit.py:252
lsst.cp.pipe.astierCovPtcFit.CovFit
Definition: astierCovPtcFit.py:213
lsst.cp.pipe.astierCovPtcFit.CovFit.getC
def getC(self)
Definition: astierCovPtcFit.py:416
lsst.cp.pipe.astierCovPtcFit.CovFit.setParamValues
def setParamValues(self, p)
Definition: astierCovPtcFit.py:337
lsst.cp.pipe.astierCovPtcFit.CovFit.getGainErr
def getGainErr(self)
Definition: astierCovPtcFit.py:462
lsst.cp.pipe.astierCovPtcFit.CovFit.__call__
def __call__(self, params)
Definition: astierCovPtcFit.py:664
lsst.cp.pipe.astierCovPtcFit.CovFit.chi2
def chi2(self)
Definition: astierCovPtcFit.py:518
lsst.cp.pipe.astierCovPtcFit.CovFit.weightedRes
def weightedRes(self, params=None)
Definition: astierCovPtcFit.py:532
lsst.cp.pipe.astierCovPtcFit.CovFit.getA
def getA(self)
Definition: astierCovPtcFit.py:408
lsst.cp.pipe.astierCovPtcFit.CovFit.getMaskCov
def getMaskCov(self, i, j)
Definition: astierCovPtcFit.py:506
lsst.cp.pipe.astierCovPtcFit.CovFit.getFitData
def getFitData(self, i, j, divideByMu=False, unitsElectrons=False, returnMasked=False)
Definition: astierCovPtcFit.py:594
lsst.cp.pipe.astierCovPtcFit.CovFit.getRonErr
def getRonErr(self)
Definition: astierCovPtcFit.py:492
lsst.cp.pipe.astierCovPtcFit.CovFit.getCCov
def getCCov(self)
Definition: astierCovPtcFit.py:457
lsst.cp.pipe.astierCovPtcFit.CovFit.wres
def wres(self, params=None)
Definition: astierCovPtcFit.py:522
lsst.cp.pipe.astierCovPtcFit.CovFit.maskMu
maskMu
Definition: astierCovPtcFit.py:248
lsst.cp.pipe.astierCovPtcFit.CovFit.getNoiseCov
def getNoiseCov(self)
Definition: astierCovPtcFit.py:470
lsst.cp.pipe.astierCovPtcFit.CovFit.getB
def getB(self)
Definition: astierCovPtcFit.py:412
lsst.cp.pipe.astierCovPtcFit.CovFit.covParams
covParams
Definition: astierCovPtcFit.py:578
lsst.cp.pipe.astierCovPtcFit.CovFit.sqrtW
sqrtW
Definition: astierCovPtcFit.py:244
lsst.cp.pipe.astierCovPtcFit.CovFit.initFit
def initFit(self)
Definition: astierCovPtcFit.py:292
lsst.cp.pipe.astierCovPtcFit.CovFit.vcov
vcov
Definition: astierCovPtcFit.py:279
lsst.cp.pipe.astierCovPtcFit.Pol2d.monomials
def monomials(self, x, y)
Definition: astierCovPtcFit.py:198
lsst.cp.pipe.astierCovFitParameters.FitParameters
Definition: astierCovFitParameters.py:164
lsst::log
lsst.cp.pipe.astierCovPtcFit.CovFit.getACov
def getACov(self)
Definition: astierCovPtcFit.py:431
lsst.cp.pipe.astierCovPtcFit.Pol2d.eval
def eval(self, x, y)
Definition: astierCovPtcFit.py:207
lsst.cp.pipe.astierCovPtcFit.symmetrize
def symmetrize(inputArray)
Definition: astierCovPtcFit.py:159
lsst.cp.pipe.astierCovPtcFit.Pol2d.ordery
ordery
Definition: astierCovPtcFit.py:191
lsst.cp.pipe.astierCovPtcFit.CovFit._getCovParams
def _getCovParams(self, what)
Definition: astierCovPtcFit.py:420