An algorithm to unfold distributions from detector to truth level. More...
#include <TUnfold.h>
Inherited by TUnfoldSysV17.
Public Types | |
enum | EConstraint { kEConstraintNone = 0, kEConstraintArea = 1 } |
type of extra constraint More... | |
enum | ERegMode { kRegModeNone = 0, kRegModeSize = 1, kRegModeDerivative = 2, kRegModeCurvature = 3, kRegModeMixed = 4 } |
choice of regularisation scheme More... | |
enum | EHistMap { kHistMapOutputHoriz = 0, kHistMapOutputVert = 1 } |
arrangement of axes for the response matrix (TH2 histogram) More... | |
Public Member Functions | |
TUnfoldV17 (const TH2 *hist_A, EHistMap histmap, ERegMode regmode=kRegModeSize, EConstraint constraint=kEConstraintArea) | |
Set up response matrix and regularisation scheme. | |
TUnfoldV17 (void) | |
only for use by root streamer or derived classes | |
virtual | ~TUnfoldV17 (void) |
virtual Int_t | SetInput (const TH1 *hist_y, Double_t scaleBias=0.0, Double_t oneOverZeroError=0.0, const TH2 *hist_vyy=0, const TH2 *hist_vyy_inv=0) |
Define input data for subsequent calls to DoUnfold(tau). | |
virtual Double_t | DoUnfold (Double_t tau) |
perform the unfolding for a given regularisation parameter tau | |
Double_t | DoUnfold (Double_t tau, const TH1 *hist_y, Double_t scaleBias=0.0) |
perform the unfolding for a given input and regularisation | |
virtual Int_t | ScanLcurve (Int_t nPoint, Double_t tauMin, Double_t tauMax, TGraph **lCurve, TSpline **logTauX=0, TSpline **logTauY=0, TSpline **logTauCurvature=0) |
scan the L curve, determine tau and unfold at the final value of tau | |
Double_t | GetTau (void) const |
return regularisation parameter | |
void | GetOutput (TH1 *output, const Int_t *binMap=0) const |
get output distribution, possibly cumulated over several bins | |
void | GetEmatrix (TH2 *ematrix, const Int_t *binMap=0) const |
get output covariance matrix, possibly cumulated over several bins | |
void | GetRhoIJ (TH2 *rhoij, const Int_t *binMap=0) const |
get correlation coefficiencts, possibly cumulated over several bins | |
Double_t | GetRhoI (TH1 *rhoi, const Int_t *binMap=0, TH2 *invEmat=0) const |
get global correlation coefficiencts, possibly cumulated over several bins | |
void | GetFoldedOutput (TH1 *folded, const Int_t *binMap=0) const |
get unfolding result on detector level | |
void | GetProbabilityMatrix (TH2 *A, EHistMap histmap) const |
get matrix of probabilities | |
void | GetNormalisationVector (TH1 *s, const Int_t *binMap=0) const |
histogram of truth bins, determined from suming over the response matrix | |
void | GetInput (TH1 *inputData, const Int_t *binMap=0) const |
Input vector of measurements. | |
void | GetInputInverseEmatrix (TH2 *ematrix) |
get inverse of the measurement's covariance matrix | |
void | GetBias (TH1 *bias, const Int_t *binMap=0) const |
get bias vector including bias scale | |
Int_t | GetNr (void) const |
get number of regularisation conditions | |
void | GetL (TH2 *l) const |
get matrix of regularisation conditions | |
void | GetLsquared (TH2 *lsquared) const |
get matrix of regularisation conditions squared | |
Double_t | GetRhoMax (void) const |
get maximum global correlation determined in recent unfolding | |
Double_t | GetRhoAvg (void) const |
get average global correlation determined in recent unfolding | |
Double_t | GetChi2A (void) const |
get 2A contribution determined in recent unfolding | |
Double_t | GetChi2L (void) const |
get 2L contribution determined in recent unfolding | |
virtual Double_t | GetLcurveX (void) const |
get value on x-axis of L-curve determined in recent unfolding | |
virtual Double_t | GetLcurveY (void) const |
get value on y-axis of L-curve determined in recent unfolding | |
Int_t | GetNdf (void) const |
get number of degrees of freedom determined in recent unfolding | |
Int_t | GetNpar (void) const |
get number of truth parameters determined in recent unfolding | |
void | SetBias (const TH1 *bias) |
set bias vector | |
void | SetConstraint (EConstraint constraint) |
set type of area constraint | |
Int_t | RegularizeSize (int bin, Double_t scale=1.0) |
add a regularisation condition on the magnitude of a truth bin | |
Int_t | RegularizeDerivative (int left_bin, int right_bin, Double_t scale=1.0) |
add a regularisation condition on the difference of two truth bin | |
Int_t | RegularizeCurvature (int left_bin, int center_bin, int right_bin, Double_t scale_left=1.0, Double_t scale_right=1.0) |
add a regularisation condition on the curvature of three truth bin | |
Int_t | RegularizeBins (int start, int step, int nbin, ERegMode regmode) |
add regularisation conditions for a group of bins | |
Int_t | RegularizeBins2D (int start_bin, int step1, int nbin1, int step2, int nbin2, ERegMode regmode) |
add regularisation conditions for 2d unfolding | |
Double_t | GetEpsMatrix (void) const |
get numerical accuracy for Eigenvalue analysis when inverting matrices with rank problems | |
void | SetEpsMatrix (Double_t eps) |
set numerical accuracy for Eigenvalue analysis when inverting matrices with rank problems | |
Static Public Member Functions | |
static const char * | GetTUnfoldVersion (void) |
return a string describing the TUnfold version | |
Protected Member Functions | |
virtual Double_t | DoUnfold (void) |
core unfolding algorithm | |
virtual void | ClearResults (void) |
reset all results | |
void | ClearHistogram (TH1 *h, Double_t x=0.) const |
Initialize bin contents and bin errors for a given histogram. | |
virtual TString | GetOutputBinName (Int_t iBinX) const |
Get bin name of an outpt bin. | |
TMatrixDSparse * | MultiplyMSparseM (const TMatrixDSparse *a, const TMatrixD *b) const |
multiply sparse matrix and a non-sparse matrix | |
TMatrixDSparse * | MultiplyMSparseMSparse (const TMatrixDSparse *a, const TMatrixDSparse *b) const |
multiply two sparse matrices | |
TMatrixDSparse * | MultiplyMSparseTranspMSparse (const TMatrixDSparse *a, const TMatrixDSparse *b) const |
multiply a transposed Sparse matrix with another Sparse matrix | |
TMatrixDSparse * | MultiplyMSparseMSparseTranspVector (const TMatrixDSparse *m1, const TMatrixDSparse *m2, const TMatrixTBase< Double_t > *v) const |
calculate a sparse matrix product M1*V*M2T where the diagonal matrix V is given by a vector | |
TMatrixDSparse * | InvertMSparseSymmPos (const TMatrixDSparse *A, Int_t *rank) const |
get the inverse or pseudo-inverse of a positive, sparse matrix | |
void | AddMSparse (TMatrixDSparse *dest, Double_t f, const TMatrixDSparse *src) const |
add a sparse matrix, scaled by a factor, to another scaled matrix | |
TMatrixDSparse * | CreateSparseMatrix (Int_t nrow, Int_t ncol, Int_t nele, Int_t *row, Int_t *col, Double_t *data) const |
create a sparse matrix, given the nonzero elements | |
Int_t | GetNx (void) const |
returns internal number of output (truth) matrix rows | |
Int_t | GetRowFromBin (int ix) const |
converts truth histogram bin number to matrix row | |
Int_t | GetBinFromRow (int ix) const |
converts matrix row to truth histogram bin number | |
Int_t | GetNy (void) const |
returns the number of measurement bins | |
const TMatrixD * | GetX (void) const |
vector of the unfolding result | |
const TMatrixDSparse * | GetVxx (void) const |
covariance matrix of the result | |
const TMatrixDSparse * | GetVxxInv (void) const |
inverse of covariance matrix of the result | |
const TMatrixDSparse * | GetAx (void) const |
vector of folded-back result | |
const TMatrixDSparse * | GetDXDY (void) const |
matrix of derivatives dx/dy | |
const TMatrixDSparse * | GetDXDAM (int i) const |
matrix contributions of the derivative dx/dA | |
const TMatrixDSparse * | GetDXDAZ (int i) const |
vector contributions of the derivative dx/dA | |
const TMatrixDSparse * | GetEinv (void) const |
matrix E-1, using internal bin counting | |
const TMatrixDSparse * | GetE (void) const |
matrix E, using internal bin counting | |
const TMatrixDSparse * | GetVyyInv (void) const |
inverse of covariance matrix of the data y | |
void | ErrorMatrixToHist (TH2 *ematrix, const TMatrixDSparse *emat, const Int_t *binMap, Bool_t doClear) const |
add up an error matrix, also respecting the bin mapping | |
Double_t | GetRhoIFromMatrix (TH1 *rhoi, const TMatrixDSparse *eOrig, const Int_t *binMap, TH2 *invEmat) const |
const TMatrixDSparse * | GetDXDtauSquared (void) const |
vector of derivative dx/dtauSquared, using internal bin counting | |
Bool_t | AddRegularisationCondition (Int_t i0, Double_t f0, Int_t i1=-1, Double_t f1=0., Int_t i2=-1, Double_t f2=0.) |
add a row of regularisation conditions to the matrix L | |
Bool_t | AddRegularisationCondition (Int_t nEle, const Int_t *indices, const Double_t *rowData) |
add a row of regularisation conditions to the matrix L | |
Static Protected Member Functions | |
static void | DeleteMatrix (TMatrixD **m) |
delete matrix and invalidate pointer | |
static void | DeleteMatrix (TMatrixDSparse **m) |
delete sparse matrix and invalidate pointer | |
Protected Attributes | |
TMatrixDSparse * | fA |
response matrix A | |
TMatrixDSparse * | fL |
regularisation conditions L | |
TMatrixD * | fY |
input (measured) data y | |
TMatrixDSparse * | fVyy |
covariance matrix Vyy corresponding to y | |
Double_t | fBiasScale |
scale factor for the bias | |
TMatrixD * | fX0 |
bias vector x0 | |
Double_t | fTauSquared |
regularisation parameter tau squared | |
TArrayI | fXToHist |
mapping of matrix indices to histogram bins | |
TArrayI | fHistToX |
mapping of histogram bins to matrix indices | |
TArrayD | fSumOverY |
truth vector calculated from the non-normalized response matrix | |
EConstraint | fConstraint |
type of constraint to use for the unfolding | |
ERegMode | fRegMode |
type of regularisation | |
Private Member Functions | |
void | InitTUnfold (void) |
initialize data menbers, for use in constructors | |
Private Attributes | |
Int_t | fIgnoredBins |
number of input bins which are dropped because they have error=0 | |
Double_t | fEpsMatrix |
machine accuracy used to determine matrix rank after eigenvalue analysis | |
TMatrixD * | fX |
unfolding result x | |
TMatrixDSparse * | fVxx |
covariance matrix Vxx | |
TMatrixDSparse * | fVxxInv |
inverse of covariance matrix Vxx-1 | |
TMatrixDSparse * | fVyyInv |
inverse of the input covariance matrix Vyy-1 | |
TMatrixDSparse * | fAx |
result x folded back A*x | |
Double_t | fChi2A |
chi**2 contribution from (y-Ax)Vyy-1(y-Ax) | |
Double_t | fLXsquared |
chi**2 contribution from (x-s*x0)TLTL(x-s*x0) | |
Double_t | fRhoMax |
maximum global correlation coefficient | |
Double_t | fRhoAvg |
average global correlation coefficient | |
Int_t | fNdf |
number of degrees of freedom | |
TMatrixDSparse * | fDXDAM [2] |
matrix contribution to the of derivative dx_k/dA_ij | |
TMatrixDSparse * | fDXDAZ [2] |
vector contribution to the of derivative dx_k/dA_ij | |
TMatrixDSparse * | fDXDtauSquared |
derivative of the result wrt tau squared | |
TMatrixDSparse * | fDXDY |
derivative of the result wrt dx/dy | |
TMatrixDSparse * | fEinv |
matrix E^(-1) | |
TMatrixDSparse * | fE |
matrix E |
An algorithm to unfold distributions from detector to truth level.
TUnfold is used to decompose a measurement y into several sources x, given the measurement uncertainties and a matrix of migrations A. The method can be applied to a large number of problems, where the measured distribution y is a linear superposition of several Monte Carlo shapes. Beyond such a simple template fit, TUnfold has an adjustable regularisation term and also supports an optional constraint on the total number of events.
For most applications, it is better to use the derived class TUnfoldDensity instead of TUnfold. TUnfoldDensity adds various features to TUnfold, such as: background subtraction, propagation of systematic uncertainties, complex multidimensional arrangements of the bins. For innocent users, the most notable improvement of TUnfoldDensity over TUnfold are the getter functions. For TUnfold, histograms have to be booked by the user and the getter functions fill the histogram bins. TUnfoldDensity simply returns a new, already filled histogram.
If you use this software, please consider the following citation
S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201]
Detailed documentation and updates are available on http://desy.de/~sschmitt
Brief recipy to use TUnfold:
Basic formulae:
2A=(Ax-y)TVyy-1(Ax-y)
2L=(x-f*x0)TLTL(x-f*x0)
2unf=2A+22L+i(Ax-y)i
x:result, A:probabilities, y:data, Vyy:data covariance, f:bias scale, x0:bias, L:regularisation conditions, :regularisation strength, :Lagrangian multiplier
Without area constraint, is set to zero, and 2unf is minimized to determine x. With area constraint, both x and are determined.
Definition at line 105 of file TUnfold.h.
enum TUnfoldV17::EHistMap |
enum TUnfoldV17::ERegMode |
choice of regularisation scheme
TUnfoldV17::TUnfoldV17 | ( | const TH2 * | hist_A, | |
EHistMap | histmap, | |||
ERegMode | regmode = kRegModeSize , |
|||
EConstraint | constraint = kEConstraintArea | |||
) |
Set up response matrix and regularisation scheme.
[in] | hist_A | matrix of MC events that describes the migrations |
[in] | histmap | mapping of the histogram axes |
[in] | regmode | (default=kRegModeSize) global regularisation mode |
[in] | constraint | (default=kEConstraintArea) type of constraint |
Treatment of overflow bins in the matrix hist_A
If unsure, do the following:
Definition at line 1696 of file TUnfoldV17.cxx.
TUnfoldV17::TUnfoldV17 | ( | void | ) |
only for use by root streamer or derived classes
Definition at line 235 of file TUnfoldV17.cxx.
virtual TUnfoldV17::~TUnfoldV17 | ( | void | ) | [virtual] |
void TUnfoldV17::AddMSparse | ( | TMatrixDSparse * | dest, | |
Double_t | f, | |||
const TMatrixDSparse * | src | |||
) | const [protected] |
add a sparse matrix, scaled by a factor, to another scaled matrix
inout] | dest destination matrix | |
[in] | f | scaling factor |
[in] | src | matrix to be added to dest |
a replacement for (*dest) += f * (*src) which suffered from a bug in old root versions
Definition at line 912 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), TUnfoldSysV17::GetSummedErrorMatrixXX(), TUnfoldSysV17::GetSummedErrorMatrixYY(), and TUnfoldSysV17::PrepareSysError().
Bool_t TUnfoldV17::AddRegularisationCondition | ( | Int_t | nEle, | |
const Int_t * | indices, | |||
const Double_t * | rowData | |||
) | [protected] |
add a row of regularisation conditions to the matrix L
[in] | nEle | number of valid entries in indeces and rowData |
[in] | indices | column numbers of L to fill |
[in] | rowData | data to fill into the new row of L |
returns true if a row was added, false otherwise
A new row k is added to the matrix L, its dimension is expanded. The new elements Lki are filled from the array rowData[] where the indices i which are taken from the array indices[].
Definition at line 1952 of file TUnfoldV17.cxx.
Bool_t TUnfoldV17::AddRegularisationCondition | ( | Int_t | i0, | |
Double_t | f0, | |||
Int_t | i1 = -1 , |
|||
Double_t | f1 = 0. , |
|||
Int_t | i2 = -1 , |
|||
Double_t | f2 = 0. | |||
) | [protected] |
add a row of regularisation conditions to the matrix L
[in] | i0 | truth histogram bin number |
[in] | f0 | entry in the matrix L, column i0 |
[in] | i1 | truth histogram bin number |
[in] | f1 | entry in the matrix L, column i1 |
[in] | i2 | truth histogram bin number |
[in] | f2 | entry in the matrix L, column i2 |
the arguments are used to form one row (k) of the matrix L, where
Lk,i0=f0 and Lk,i1=f1 and Lk,i2=f2
negative indexes i0,i1,i2 are ignored
Definition at line 1915 of file TUnfoldV17.cxx.
Referenced by RegularizeCurvature(), RegularizeDerivative(), and RegularizeSize().
void TUnfoldV17::ClearHistogram | ( | TH1 * | h, | |
Double_t | x = 0. | |||
) | const [protected] |
Initialize bin contents and bin errors for a given histogram.
[out] | h | histogram |
[in] | x | new histogram content |
all histgram errors are set to zero, all contents are set to x
Definition at line 3624 of file TUnfoldV17.cxx.
Referenced by GetBias(), GetFoldedOutput(), GetInput(), GetNormalisationVector(), GetOutput(), GetRhoI(), and TUnfoldSysV17::GetRhoItotal().
void TUnfoldV17::ClearResults | ( | void | ) | [protected, virtual] |
reset all results
Reimplemented in TUnfoldSysV17.
Definition at line 205 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), SetConstraint(), and SetInput().
TMatrixDSparse * TUnfoldV17::CreateSparseMatrix | ( | Int_t | nrow, | |
Int_t | ncol, | |||
Int_t | nel, | |||
Int_t * | row, | |||
Int_t * | col, | |||
Double_t * | data | |||
) | const [protected] |
create a sparse matrix, given the nonzero elements
[in] | nrow | number of rows |
[in] | ncol | number of columns |
[in] | nel | number of non-zero elements |
[in] | row | row indexes of non-zero elements |
[in] | col | column indexes of non-zero elements |
[in] | data | non-zero elements data |
return pointer to a new sparse matrix
shortcut to new TMatrixDSparse() followed by SetMatrixArray()
Definition at line 576 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), TUnfoldSysV17::PrepareSysError(), SetInput(), and TUnfoldV17().
void TUnfoldV17::DeleteMatrix | ( | TMatrixDSparse ** | m | ) | [static, protected] |
delete sparse matrix and invalidate pointer
inout] | m pointer to a matrix-pointer |
if the matrix pointer os non-zero, thematrix id deleted. The matrox pointer is set to zero.
Definition at line 197 of file TUnfoldV17.cxx.
void TUnfoldV17::DeleteMatrix | ( | TMatrixD ** | m | ) | [static, protected] |
delete matrix and invalidate pointer
inout] | m pointer to a matrix-pointer |
if the matrix pointer os non-zero, thematrix id deleted. The matrox pointer is set to zero.
Definition at line 185 of file TUnfoldV17.cxx.
Referenced by ClearResults(), TUnfoldSysV17::ClearResults(), TUnfoldSysV17::DoBackgroundSubtraction(), DoUnfold(), TUnfoldSysV17::GetChi2Sys(), GetFoldedOutput(), GetLsquared(), GetRhoIFromMatrix(), TUnfoldSysV17::GetRhoItotal(), TUnfoldSysV17::GetSummedErrorMatrixXX(), TUnfoldSysV17::GetSummedErrorMatrixYY(), TUnfoldSysV17::PrepareSysError(), SetBias(), SetInput(), and TUnfoldSysV17::SetTauError().
Double_t TUnfoldV17::DoUnfold | ( | Double_t | tau_reg, | |
const TH1 * | input, | |||
Double_t | scaleBias = 0.0 | |||
) |
perform the unfolding for a given input and regularisation
[in] | tau_reg | regularisation parameter |
[in] | input | input distribution with uncertainties |
[in] | scaleBias | (default=0.0) scale factor applied to the bias |
This is a shortcut for { SetInput(input,scaleBias); DoUnfold(tau); }
Definition at line 2232 of file TUnfoldV17.cxx.
Double_t TUnfoldV17::DoUnfold | ( | Double_t | tau | ) | [virtual] |
perform the unfolding for a given regularisation parameter tau
[in] | tau | regularisation parameter |
this method sets tau and then calls the core unfolding algorithm
Definition at line 2485 of file TUnfoldV17.cxx.
Double_t TUnfoldV17::DoUnfold | ( | void | ) | [protected, virtual] |
core unfolding algorithm
Definition at line 243 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), and ScanLcurve().
void TUnfoldV17::ErrorMatrixToHist | ( | TH2 * | ematrix, | |
const TMatrixDSparse * | emat, | |||
const Int_t * | binMap, | |||
Bool_t | doClear | |||
) | const [protected] |
add up an error matrix, also respecting the bin mapping
inout] | ematrix error matrix histogram | |
[in] | emat | error matrix stored with internal mapping (member fXToHist) |
[in] | binMap | mapping of histogram bins |
[in] | doClear | if true, ematrix is cleared prior to adding elements of emat to it. |
the array binMap is explained with the method GetOutput(). The matrix emat must have dimension NxN where N=fXToHist.size() The flag doClear may be used to add covariance matrices from several uncertainty sources.
Definition at line 3323 of file TUnfoldV17.cxx.
Referenced by GetEmatrix().
const TMatrixDSparse* TUnfoldV17::GetAx | ( | void | ) | const [inline, protected] |
vector of folded-back result
Definition at line 246 of file TUnfold.h.
Referenced by TUnfoldSysV17::GetChi2Sys().
void TUnfoldV17::GetBias | ( | TH1 * | out, | |
const Int_t * | binMap = 0 | |||
) | const |
get bias vector including bias scale
[out] | out | histogram to store the scaled bias vector. The bin contents are overwritten |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
This method returns the bias vector times scaling factor, f*x0
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 2902 of file TUnfoldV17.cxx.
Int_t TUnfoldV17::GetBinFromRow | ( | int | ix | ) | const [inline, protected] |
Double_t TUnfoldV17::GetChi2A | ( | void | ) | const [inline] |
Double_t TUnfoldV17::GetChi2L | ( | void | ) | const |
get 2L contribution determined in recent unfolding
Definition at line 3175 of file TUnfoldV17.cxx.
const TMatrixDSparse* TUnfoldV17::GetDXDAM | ( | int | i | ) | const [inline, protected] |
matrix contributions of the derivative dx/dA
Definition at line 250 of file TUnfold.h.
Referenced by TUnfoldSysV17::PrepareSysError().
const TMatrixDSparse* TUnfoldV17::GetDXDAZ | ( | int | i | ) | const [inline, protected] |
const TMatrixDSparse* TUnfoldV17::GetDXDtauSquared | ( | void | ) | const [inline, protected] |
vector of derivative dx/dtauSquared, using internal bin counting
Definition at line 263 of file TUnfold.h.
Referenced by TUnfoldSysV17::PrepareSysError().
const TMatrixDSparse* TUnfoldV17::GetDXDY | ( | void | ) | const [inline, protected] |
const TMatrixDSparse* TUnfoldV17::GetE | ( | void | ) | const [inline, protected] |
const TMatrixDSparse* TUnfoldV17::GetEinv | ( | void | ) | const [inline, protected] |
void TUnfoldV17::GetEmatrix | ( | TH2 * | ematrix, | |
const Int_t * | binMap = 0 | |||
) | const |
get output covariance matrix, possibly cumulated over several bins
[out] | ematrix | histogram to store the covariance. The bin contents are overwritten. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 3390 of file TUnfoldV17.cxx.
Referenced by TUnfoldSysV17::GetEmatrixTotal(), and GetRhoIJ().
Double_t TUnfoldV17::GetEpsMatrix | ( | void | ) | const [inline] |
void TUnfoldV17::GetFoldedOutput | ( | TH1 * | out, | |
const Int_t * | binMap = 0 | |||
) | const |
get unfolding result on detector level
[out] | out | histogram to store the correlation coefficiencts. The bin contents and errors are overwritten. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
This method returns the unfolding output folded by the response matrix, i.e. the vector Ax.
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 2929 of file TUnfoldV17.cxx.
void TUnfoldV17::GetInput | ( | TH1 * | out, | |
const Int_t * | binMap = 0 | |||
) | const |
Input vector of measurements.
[out] | out | histogram to store the measurements. Bin content and bin errors are overwritte. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
Bins which had an uncertainty of zero in the call to SetInput() maye acquire bin contents or bin errors different from the original settings in SetInput().
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 3013 of file TUnfoldV17.cxx.
void TUnfoldV17::GetInputInverseEmatrix | ( | TH2 * | out | ) |
get inverse of the measurement's covariance matrix
[out] | out | histogram to store the inverted covariance |
Definition at line 3042 of file TUnfoldV17.cxx.
Referenced by DoUnfold().
void TUnfoldV17::GetL | ( | TH2 * | out | ) | const |
get matrix of regularisation conditions
[out] | out | histogram to store the regularisation conditions. the bincontents are overwritten |
The histogram should have dimension nr (x-axis) times nx (y-axis). nr corresponds to the number of regularisation conditions, it can be obtained using the method GetNr(). nx corresponds to the number of histogram bins in the response matrix along the truth axis.
Definition at line 3135 of file TUnfoldV17.cxx.
Double_t TUnfoldV17::GetLcurveX | ( | void | ) | const [virtual] |
get value on x-axis of L-curve determined in recent unfolding
x=log10(GetChi2A())
Definition at line 3195 of file TUnfoldV17.cxx.
Referenced by ScanLcurve().
Double_t TUnfoldV17::GetLcurveY | ( | void | ) | const [virtual] |
get value on y-axis of L-curve determined in recent unfolding
y=log10(GetChi2L())
Definition at line 3204 of file TUnfoldV17.cxx.
Referenced by ScanLcurve().
void TUnfoldV17::GetLsquared | ( | TH2 * | out | ) | const |
get matrix of regularisation conditions squared
[out] | out | histogram to store the squared matrix of regularisation conditions. the bin contents are overwritten |
This returns the square matrix LTL as a histogram
The histogram should have dimension nx times nx, where nx corresponds to the number of histogram bins in the response matrix along the truth axis.
Definition at line 3095 of file TUnfoldV17.cxx.
Int_t TUnfoldV17::GetNdf | ( | void | ) | const [inline] |
get number of degrees of freedom determined in recent unfolding
This returns the number of valid measurements minus the number of unfolded truth bins. If the area constraint is active, one further degree of freedom is subtracted
Definition at line 325 of file TUnfold.h.
Referenced by ScanLcurve().
void TUnfoldV17::GetNormalisationVector | ( | TH1 * | out, | |
const Int_t * | binMap = 0 | |||
) | const |
histogram of truth bins, determined from suming over the response matrix
[out] | out | histogram to store the truth bins. The bin contents are overwritten |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
This vector is also used to initialize the bias x0. However, the bias vector may be changed using the SetBias() method.
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 2877 of file TUnfoldV17.cxx.
Int_t TUnfoldV17::GetNpar | ( | void | ) | const |
get number of truth parameters determined in recent unfolding
empty bins of the response matrix or bins which can not be unfolded due to rank deficits are not counted
Definition at line 3186 of file TUnfoldV17.cxx.
Referenced by GetInputInverseEmatrix().
Int_t TUnfoldV17::GetNr | ( | void | ) | const |
get number of regularisation conditions
Ths returns the number of regularisation conditions, useful for booking a histogram for a subsequent call of GetL().
Definition at line 3120 of file TUnfoldV17.cxx.
Referenced by GetL().
Int_t TUnfoldV17::GetNx | ( | void | ) | const [inline, protected] |
returns internal number of output (truth) matrix rows
Definition at line 228 of file TUnfold.h.
Referenced by DoUnfold(), GetBias(), GetLsquared(), GetNormalisationVector(), GetNpar(), GetRhoI(), GetRhoIFromMatrix(), SetBias(), and TUnfoldV17().
Int_t TUnfoldV17::GetNy | ( | void | ) | const [inline, protected] |
returns the number of measurement bins
Definition at line 236 of file TUnfold.h.
Referenced by TUnfoldSysV17::DoBackgroundSubtraction(), DoUnfold(), GetFoldedOutput(), GetInput(), GetInputInverseEmatrix(), TUnfoldSysV17::PrepareSysError(), and SetInput().
void TUnfoldV17::GetOutput | ( | TH1 * | output, | |
const Int_t * | binMap = 0 | |||
) | const |
get output distribution, possibly cumulated over several bins
[out] | output | existing output histogram. content and errors will be updated. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
If nonzero, the array binMap must have dimension n+2, where n corresponds to the number of bins on the truth axis of the response matrix (the histogram specified with the TUnfold constructor). The indexes of binMap correspond to the truth bins (including underflow and overflow) of the response matrix. The element binMap[i] specifies the histogram number in output where the corresponding truth bin will be stored. It is possible to specify the same output bin number for multiple indexes, in which case these bins are added. Set binMap[i]=-1 to ignore an unfolded truth bin. The uncertainties are calculated from the corresponding parts of the covariance matrix, properly taking care of added truth bins.
If the pointer binMap is zero, the bins are mapped one-to-one. Truth bin zero (underflow) is stored in the output underflow, truth bin 1 is stored in bin number 1, etc.
Definition at line 3233 of file TUnfoldV17.cxx.
TString TUnfoldV17::GetOutputBinName | ( | Int_t | iBinX | ) | const [protected, virtual] |
Get bin name of an outpt bin.
[in] | iBinX | bin number |
Return value: name of the bin
For TUnfold and TUnfoldSys, this function simply returns the bin number as a string. This function really only makes sense in the context of TUnfoldDensity, where binnig schemes are implemented using the class TUnfoldBinning, and non-trivial bin names are returned.
Reimplemented in TUnfoldDensityV17.
Definition at line 1664 of file TUnfoldV17.cxx.
Referenced by SetInput().
void TUnfoldV17::GetProbabilityMatrix | ( | TH2 * | A, | |
EHistMap | histmap | |||
) | const |
get matrix of probabilities
[out] | A | two-dimensional histogram to store the probabilities (normalized response matrix). The bin contents are overwritten |
[in] | histmap | specify axis along which the truth bins are oriented |
Definition at line 2977 of file TUnfoldV17.cxx.
Referenced by TUnfoldDensityV17::GetProbabilityMatrix().
Double_t TUnfoldV17::GetRhoAvg | ( | void | ) | const [inline] |
Double_t TUnfoldV17::GetRhoI | ( | TH1 * | rhoi, | |
const Int_t * | binMap = 0 , |
|||
TH2 * | invEmat = 0 | |||
) | const |
get global correlation coefficiencts, possibly cumulated over several bins
[out] | rhoi | histogram to store the global correlation coefficients. The bin contents are overwritten. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
[out] | invEmat | (default=0) histogram to store the inverted covariance matrix |
for a given bin, the global correlation coefficient is defined as
i=sqrt(1-1/(Vii*V-1ii))
such that the calculation of global correlation coefficients possibly involves the inversion of a covariance matrix.
return value: maximum global correlation coefficient
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 3448 of file TUnfoldV17.cxx.
Double_t TUnfoldV17::GetRhoIFromMatrix | ( | TH1 * | rhoi, | |
const TMatrixDSparse * | eOrig, | |||
const Int_t * | binMap, | |||
TH2 * | invEmat | |||
) | const [protected] |
Definition at line 3497 of file TUnfoldV17.cxx.
Referenced by GetRhoI(), and TUnfoldSysV17::GetRhoItotal().
void TUnfoldV17::GetRhoIJ | ( | TH2 * | rhoij, | |
const Int_t * | binMap = 0 | |||
) | const |
get correlation coefficiencts, possibly cumulated over several bins
[out] | rhoij | histogram to store the correlation coefficiencts. The bin contents are overwritten. |
[in] | binMap | (default=0) array for mapping truth bins to histogram bins |
The use of binMap is explained with the documentation of the GetOutput() method
Definition at line 3405 of file TUnfoldV17.cxx.
Double_t TUnfoldV17::GetRhoMax | ( | void | ) | const [inline] |
Int_t TUnfoldV17::GetRowFromBin | ( | int | ix | ) | const [inline, protected] |
Double_t TUnfoldV17::GetTau | ( | void | ) | const |
return regularisation parameter
Definition at line 3167 of file TUnfoldV17.cxx.
const char * TUnfoldV17::GetTUnfoldVersion | ( | void | ) | [static] |
return a string describing the TUnfold version
The version is reported in the form Vmajor.minor Changes of the minor version number typically correspond to bug-fixes. Changes of the major version may result in adding or removing data attributes, such that the streamer methods are not compatible between different major versions.
Definition at line 3661 of file TUnfoldV17.cxx.
const TMatrixDSparse* TUnfoldV17::GetVxx | ( | void | ) | const [inline, protected] |
covariance matrix of the result
Definition at line 242 of file TUnfold.h.
Referenced by TUnfoldSysV17::GetSummedErrorMatrixXX().
const TMatrixDSparse* TUnfoldV17::GetVxxInv | ( | void | ) | const [inline, protected] |
const TMatrixDSparse* TUnfoldV17::GetVyyInv | ( | void | ) | const [inline, protected] |
inverse of covariance matrix of the data y
Definition at line 258 of file TUnfold.h.
Referenced by TUnfoldSysV17::DoBackgroundSubtraction().
const TMatrixD* TUnfoldV17::GetX | ( | void | ) | const [inline, protected] |
void TUnfoldV17::InitTUnfold | ( | void | ) | [private] |
initialize data menbers, for use in constructors
Definition at line 141 of file TUnfoldV17.cxx.
Referenced by TUnfoldV17().
TMatrixDSparse * TUnfoldV17::InvertMSparseSymmPos | ( | const TMatrixDSparse * | A, | |
Int_t * | rankPtr | |||
) | const [protected] |
get the inverse or pseudo-inverse of a positive, sparse matrix
[in] | A | the sparse matrix to be inverted, has to be positive |
inout] | rankPtr if zero, suppress calculation of pseudo-inverse otherwise the rank of the matrix is returned in *rankPtr |
return value: 0 or a new sparse matrix
the matrix inversion is optimized in performance for the case where a large submatrix of A is diagonal
Definition at line 990 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), TUnfoldSysV17::GetChi2Sys(), GetInputInverseEmatrix(), and GetRhoIFromMatrix().
TMatrixDSparse * TUnfoldV17::MultiplyMSparseM | ( | const TMatrixDSparse * | a, | |
const TMatrixD * | b | |||
) | const [protected] |
multiply sparse matrix and a non-sparse matrix
[in] | a | sparse matrix |
[in] | b | matrix |
returns a new sparse matrix a*b.
A replacement for: new TMatrixDSparse(a,TMatrixDSparse::kMult,b) the root implementation had problems in older versions of root
Definition at line 757 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), and TUnfoldSysV17::GetChi2Sys().
TMatrixDSparse * TUnfoldV17::MultiplyMSparseMSparse | ( | const TMatrixDSparse * | a, | |
const TMatrixDSparse * | b | |||
) | const [protected] |
multiply two sparse matrices
[in] | a | sparse matrix |
[in] | b | sparse matrix |
returns a new sparse matrix a*b.
A replacement for: new TMatrixDSparse(a,TMatrixDSparse::kMult,b) the root implementation had problems in older versions of root
Definition at line 600 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), GetFoldedOutput(), TUnfoldSysV17::GetSummedErrorMatrixYY(), and TUnfoldSysV17::PrepareSysError().
TMatrixDSparse * TUnfoldV17::MultiplyMSparseMSparseTranspVector | ( | const TMatrixDSparse * | m1, | |
const TMatrixDSparse * | m2, | |||
const TMatrixTBase< Double_t > * | v | |||
) | const [protected] |
calculate a sparse matrix product M1*V*M2T where the diagonal matrix V is given by a vector
[in] | m1 | pointer to sparse matrix with dimension I*K |
[in] | m2 | pointer to sparse matrix with dimension J*K |
[in] | v | pointer to vector (matrix) with dimension K*1 |
returns a sparse matrix R with elements rij=kM1ikVkM2jk
Definition at line 817 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), TUnfoldSysV17::GetSummedErrorMatrixXX(), and TUnfoldSysV17::GetSummedErrorMatrixYY().
TMatrixDSparse * TUnfoldV17::MultiplyMSparseTranspMSparse | ( | const TMatrixDSparse * | a, | |
const TMatrixDSparse * | b | |||
) | const [protected] |
multiply a transposed Sparse matrix with another Sparse matrix
[in] | a | sparse matrix (to be transposed) |
[in] | b | sparse matrix |
returns a new sparse matrix aT*b
this is a replacement for the root constructors new TMatrixDSparse(TMatrixDSparse(TMatrixDSparse::kTransposed,*a),TMatrixDSparse::kMult,*b)
Definition at line 675 of file TUnfoldV17.cxx.
Referenced by DoUnfold(), GetLsquared(), and SetInput().
Int_t TUnfoldV17::RegularizeBins | ( | int | start, | |
int | step, | |||
int | nbin, | |||
ERegMode | regmode | |||
) |
add regularisation conditions for a group of bins
[in] | start | first bin number |
[in] | step | step size |
[in] | nbin | number of bins |
[in] | regmode | regularisation mode (one of: kRegModeSize, kRegModeDerivative, kRegModeCurvature) |
add regularisation conditions for a group of equidistant bins. There are nbin bins, starting with bin start and with a distance of step between bins.
Return value: number of regularisation conditions which could not be added.
Conditions which are not added typically correspond to bin numbers where the truth can not be unfolded (either response matrix is empty or the data do not constrain).
Definition at line 2140 of file TUnfoldV17.cxx.
Referenced by RegularizeBins2D(), and TUnfoldV17().
Int_t TUnfoldV17::RegularizeBins2D | ( | int | start_bin, | |
int | step1, | |||
int | nbin1, | |||
int | step2, | |||
int | nbin2, | |||
ERegMode | regmode | |||
) |
add regularisation conditions for 2d unfolding
[in] | start_bin | first bin number |
[in] | step1 | step size, 1st dimension |
[in] | nbin1 | number of bins, 1st dimension |
[in] | step2 | step size, 2nd dimension |
[in] | nbin2 | number of bins, 2nd dimension |
[in] | regmode | regularisation mode (one of: kRegModeSize, kRegModeDerivative, kRegModeCurvature) |
add regularisation conditions for a grid of bins. The start bin is start_bin. Along the first (second) dimension, there are nbin1 (nbin2) bins and adjacent bins are spaced by step1 (step2) units.
Return value: number of regularisation conditions which could not be added. Conditions which are not added typically correspond to bin numbers where the truth can not be unfolded (either response matrix is empty or the data do not constrain).
Definition at line 2201 of file TUnfoldV17.cxx.
Int_t TUnfoldV17::RegularizeCurvature | ( | int | left_bin, | |
int | center_bin, | |||
int | right_bin, | |||
Double_t | scale_left = 1.0 , |
|||
Double_t | scale_right = 1.0 | |||
) |
add a regularisation condition on the curvature of three truth bin
[in] | left_bin | bin number |
[in] | center_bin | bin number |
[in] | right_bin | bin number |
[in] | scale_left | (default=1) scale factor |
[in] | scale_right | (default=1) scale factor |
this adds one row to L, where the element left_bin takes the value -scale_left, the element right_bin takes the value -scale_right and the element center_bin takes the value scale_left+scale_right
return value: 0 if ok, 1 if the condition has not been added. Conditions which are not added typically correspond to bin numbers where the truth can not be unfolded (either response matrix is empty or the data do not constrain).
The RegularizeXXX() methods can be used to set up a custom matrix of regularisation conditions. In this case, start with an empty matrix L (argument regmode=kRegModeNone in the constructor)
Definition at line 2095 of file TUnfoldV17.cxx.
Referenced by RegularizeBins().
Int_t TUnfoldV17::RegularizeDerivative | ( | int | left_bin, | |
int | right_bin, | |||
Double_t | scale = 1.0 | |||
) |
add a regularisation condition on the difference of two truth bin
[in] | left_bin | bin number |
[in] | right_bin | bin number |
[in] | scale | (default=1) scale factor |
this adds one row to L, where the element left_bin takes the value -scale and the element right_bin takes the value +scale
return value: 0 if ok, 1 if the condition has not been added. Conditions which are not added typically correspond to bin numbers where the truth can not be unfolded (either response matrix is empty or the data do not constrain).
The RegularizeXXX() methods can be used to set up a custom matrix of regularisation conditions. In this case, start with an empty matrix L (argument regmode=kRegModeNone in the constructor)
Definition at line 2056 of file TUnfoldV17.cxx.
Referenced by RegularizeBins().
Int_t TUnfoldV17::RegularizeSize | ( | int | bin, | |
Double_t | scale = 1.0 | |||
) |
add a regularisation condition on the magnitude of a truth bin
[in] | bin | bin number |
[in] | scale | (default=1) scale factor |
this adds one row to L, where the element bin takes the value scale
return value: 0 if ok, 1 if the condition has not been added. Conditions which are not added typically correspond to bin numbers where the truth can not be unfolded (either response matrix is empty or the data do not constrain).
The RegularizeXXX() methods can be used to set up a custom matrix of regularisation conditions. In this case, start with an empty matrix L (argument regmode=kRegModeNone in the constructor)
Definition at line 2022 of file TUnfoldV17.cxx.
Referenced by RegularizeBins().
Int_t TUnfoldV17::ScanLcurve | ( | Int_t | nPoint, | |
Double_t | tauMin, | |||
Double_t | tauMax, | |||
TGraph ** | lCurve, | |||
TSpline ** | logTauX = 0 , |
|||
TSpline ** | logTauY = 0 , |
|||
TSpline ** | logTauCurvature = 0 | |||
) | [virtual] |
scan the L curve, determine tau and unfold at the final value of tau
[in] | nPoint | number of points used for the scan |
[in] | tauMin | smallest tau value to study |
[in] | tauMax | largest tau value to study. If tauMin=tauMax=0, a scan interval is determined automatically. |
[out] | lCurve | if nonzero, a new TGraph is returned, containing the L-curve |
[out] | logTauX | if nonzero, a new TSpline is returned, to parameterize the L-curve's x-coordinates as a function of log10(tau) |
[out] | logTauY | if nonzero, a new TSpline is returned, to parameterize the L-curve's y-coordinates as a function of log10(tau) |
[out] | logTauCurvature | if nonzero, a new TSpline is returned of the L-curve curvature as a function of log10(tau) |
return value: the coordinate number in the logTauX,logTauY graphs corresponding to the "final" choice of tau
Recommendation: always check logTauCurvature, it should be a peaked function (similar to a Gaussian), the maximum corresponding to the final choice of tau. Also, check the lCurve it should be approximately L-shaped. If in doubt, adjust tauMin and tauMax until the results are satisfactory.
Definition at line 2526 of file TUnfoldV17.cxx.
void TUnfoldV17::SetBias | ( | const TH1 * | bias | ) |
set bias vector
[in] | bias | histogram with new bias vector |
the initial bias vector is determined from the response matrix but may be changed by using this method
Definition at line 1892 of file TUnfoldV17.cxx.
void TUnfoldV17::SetConstraint | ( | EConstraint | constraint | ) |
set type of area constraint
results of a previous unfolding are reset
Definition at line 3155 of file TUnfoldV17.cxx.
Referenced by TUnfoldV17().
void TUnfoldV17::SetEpsMatrix | ( | Double_t | eps | ) |
set numerical accuracy for Eigenvalue analysis when inverting matrices with rank problems
Definition at line 3647 of file TUnfoldV17.cxx.
Int_t TUnfoldV17::SetInput | ( | const TH1 * | input, | |
Double_t | scaleBias = 0.0 , |
|||
Double_t | oneOverZeroError = 0.0 , |
|||
const TH2 * | hist_vyy = 0 , |
|||
const TH2 * | hist_vyy_inv = 0 | |||
) | [virtual] |
Define input data for subsequent calls to DoUnfold(tau).
[in] | input | input distribution with uncertainties |
[in] | scaleBias | (default=0) scale factor applied to the bias |
[in] | oneOverZeroError | (default=0) for bins with zero error, this number defines 1/error. |
[in] | hist_vyy | (default=0) if non-zero, this defines the data covariance matrix |
[in] | hist_vyy_inv | (default=0) if non-zero and hist_vyy is set, defines the inverse of the data covariance matrix. This feature can be useful for repeated unfoldings in cases where the inversion of the input covariance matrix is lengthy |
Return value: nError1+10000*nError2
Reimplemented in TUnfoldSysV17.
Definition at line 2271 of file TUnfoldV17.cxx.
Referenced by DoUnfold().
TMatrixDSparse* TUnfoldV17::fA [protected] |
response matrix A
Definition at line 152 of file TUnfold.h.
Referenced by DoUnfold(), GetFoldedOutput(), GetNx(), GetNy(), GetProbabilityMatrix(), TUnfoldSysV17::GetSummedErrorMatrixYY(), InitTUnfold(), TUnfoldSysV17::PrepareSysError(), SetInput(), and TUnfoldV17().
TMatrixDSparse* TUnfoldV17::fAx [private] |
result x folded back A*x
Definition at line 189 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetAx(), GetFoldedOutput(), and InitTUnfold().
Double_t TUnfoldV17::fBiasScale [protected] |
scale factor for the bias
Definition at line 160 of file TUnfold.h.
Referenced by DoUnfold(), GetBias(), InitTUnfold(), and SetInput().
Double_t TUnfoldV17::fChi2A [private] |
chi**2 contribution from (y-Ax)Vyy-1(y-Ax)
Definition at line 191 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetChi2A(), GetLcurveX(), InitTUnfold(), and ScanLcurve().
EConstraint TUnfoldV17::fConstraint [protected] |
type of constraint to use for the unfolding
Definition at line 172 of file TUnfold.h.
Referenced by DoUnfold(), InitTUnfold(), and SetConstraint().
TMatrixDSparse* TUnfoldV17::fDXDAM[2] [private] |
matrix contribution to the of derivative dx_k/dA_ij
Definition at line 201 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetDXDAM(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fDXDAZ[2] [private] |
vector contribution to the of derivative dx_k/dA_ij
Definition at line 203 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetDXDAZ(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fDXDtauSquared [private] |
derivative of the result wrt tau squared
Definition at line 205 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetDXDtauSquared(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fDXDY [private] |
derivative of the result wrt dx/dy
Definition at line 207 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetDXDY(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fE [private] |
matrix E
Definition at line 211 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetE(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fEinv [private] |
matrix E^(-1)
Definition at line 209 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetEinv(), and InitTUnfold().
Double_t TUnfoldV17::fEpsMatrix [private] |
machine accuracy used to determine matrix rank after eigenvalue analysis
Definition at line 179 of file TUnfold.h.
Referenced by GetEpsMatrix(), InitTUnfold(), and SetEpsMatrix().
TArrayI TUnfoldV17::fHistToX [protected] |
mapping of histogram bins to matrix indices
Definition at line 168 of file TUnfold.h.
Referenced by ErrorMatrixToHist(), GetOutput(), GetRhoIFromMatrix(), GetRowFromBin(), InitTUnfold(), and TUnfoldV17().
Int_t TUnfoldV17::fIgnoredBins [private] |
number of input bins which are dropped because they have error=0
Definition at line 177 of file TUnfold.h.
Referenced by GetInputInverseEmatrix(), InitTUnfold(), and SetInput().
TMatrixDSparse* TUnfoldV17::fL [protected] |
regularisation conditions L
Definition at line 154 of file TUnfold.h.
Referenced by DoUnfold(), GetL(), GetLsquared(), GetNr(), InitTUnfold(), and TUnfoldV17().
Double_t TUnfoldV17::fLXsquared [private] |
chi**2 contribution from (x-s*x0)TLTL(x-s*x0)
Definition at line 193 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetChi2L(), GetLcurveY(), and InitTUnfold().
Int_t TUnfoldV17::fNdf [private] |
number of degrees of freedom
Definition at line 199 of file TUnfold.h.
Referenced by DoUnfold(), GetInputInverseEmatrix(), GetNdf(), InitTUnfold(), and SetInput().
ERegMode TUnfoldV17::fRegMode [protected] |
type of regularisation
Definition at line 174 of file TUnfold.h.
Referenced by InitTUnfold(), RegularizeCurvature(), RegularizeDerivative(), and RegularizeSize().
Double_t TUnfoldV17::fRhoAvg [private] |
average global correlation coefficient
Definition at line 197 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetRhoAvg(), and InitTUnfold().
Double_t TUnfoldV17::fRhoMax [private] |
maximum global correlation coefficient
Definition at line 195 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetRhoMax(), and InitTUnfold().
TArrayD TUnfoldV17::fSumOverY [protected] |
truth vector calculated from the non-normalized response matrix
Definition at line 170 of file TUnfold.h.
Referenced by GetNormalisationVector(), InitTUnfold(), and TUnfoldV17().
Double_t TUnfoldV17::fTauSquared [protected] |
regularisation parameter tau squared
Definition at line 164 of file TUnfold.h.
Referenced by DoUnfold(), GetChi2L(), GetTau(), InitTUnfold(), and TUnfoldSysV17::PrepareSysError().
TMatrixDSparse* TUnfoldV17::fVxx [private] |
covariance matrix Vxx
Definition at line 183 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetEmatrix(), GetFoldedOutput(), GetOutput(), GetRhoI(), GetVxx(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fVxxInv [private] |
inverse of covariance matrix Vxx-1
Definition at line 185 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetRhoI(), GetVxxInv(), and InitTUnfold().
TMatrixDSparse* TUnfoldV17::fVyy [protected] |
covariance matrix Vyy corresponding to y
Definition at line 158 of file TUnfold.h.
Referenced by TUnfoldSysV17::DoBackgroundSubtraction(), DoUnfold(), GetInput(), GetInputInverseEmatrix(), TUnfoldSysV17::GetSummedErrorMatrixYY(), InitTUnfold(), SetInput(), and TUnfoldSysV17::SetInput().
TMatrixDSparse* TUnfoldV17::fVyyInv [private] |
inverse of the input covariance matrix Vyy-1
Definition at line 187 of file TUnfold.h.
Referenced by DoUnfold(), GetInputInverseEmatrix(), GetVyyInv(), InitTUnfold(), and SetInput().
TMatrixD* TUnfoldV17::fX [private] |
unfolding result x
Definition at line 181 of file TUnfold.h.
Referenced by ClearResults(), DoUnfold(), GetOutput(), GetX(), and InitTUnfold().
TMatrixD* TUnfoldV17::fX0 [protected] |
bias vector x0
Definition at line 162 of file TUnfold.h.
Referenced by DoUnfold(), GetBias(), InitTUnfold(), SetBias(), and TUnfoldV17().
TArrayI TUnfoldV17::fXToHist [protected] |
mapping of matrix indices to histogram bins
Definition at line 166 of file TUnfold.h.
Referenced by GetBias(), GetBinFromRow(), GetL(), GetLsquared(), GetNormalisationVector(), GetProbabilityMatrix(), GetRhoI(), GetRhoIFromMatrix(), InitTUnfold(), SetBias(), SetInput(), and TUnfoldV17().
TMatrixD* TUnfoldV17::fY [protected] |
input (measured) data y
Definition at line 156 of file TUnfold.h.
Referenced by TUnfoldSysV17::DoBackgroundSubtraction(), DoUnfold(), TUnfoldSysV17::GetChi2Sys(), GetInput(), InitTUnfold(), SetInput(), and TUnfoldSysV17::SetInput().