Logo ROOT   6.08/07
Reference Guide
List of all members | Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Static Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
TUnfold Class Reference

TUnfold is used to decompose a measurement y into several sources x given the measurement uncertainties and a matrix of migrations A.

For most applications, it is better to use TUnfoldDensity instead of using TUnfoldSys or TUnfold**

If you use this software, please consider the following citation S.Schmitt, JINST 7 (2012) T10003 [arXiv:1205.6201]

More documentation and updates are available on http://www.desy.de/~sschmitt

A short summary of the algorithm is given in the following: the "best" x matching the measurement y within errors is determined by minimizing the function L1+L2+L3

where L1 = (y-Ax)# Vyy^-1 (y-Ax) L2 = tau^2 (L(x-x0))# L(x-x0) L3 = lambda sum_i(y_i -(Ax)_i)

[the notation # means that the matrix is transposed ]

The term L1 is familiar from a least-square minimisation The term L2 defines the regularisation (smootheness condition on x), where the parameter tau^2 gives the strength of teh regularisation The term L3 is an optional area constraint with Lagrangian parameter lambda, ensuring that that the normalisation of x is consistent with the normalisation of y

The method can be applied to a very large number of problems, where the measured distribution y is a linear superposition of several Monte Carlo shapes

Input from measurement:

y: vector of measured quantities (dimension ny) Vyy: covariance matrix for y (dimension ny x ny) in many cases V is diagonal and calculated from the errors of y

From simulation:

A: migration matrix (dimension ny x nx)

Result

x: unknown underlying distribution (dimension nx) The error matrix of x, V_xx, is also determined

Regularisation

tau: parameter, defining the regularisation strength L: matrix of regularisation conditions (dimension nl x nx) depends on the structure of the input data x0: bias distribution, from simulation

Preservation of the area

lambda: lagrangian multiplier y_i: one component of the vector y (Ax)_i: one component of the vector Ax

Determination of the unfolding result x:

The constraint can be useful to reduce biases on the result x in cases where the vector y follows non-Gaussian probability densities (example: Poisson statistics at counting experiments in particle physics)

Some random examples:

  1. measure a cross-section as a function of, say, E_T(detector) and unfold it to obtain the underlying distribution E_T(generator)
  2. measure a lifetime distribution and unfold the contributions from different flavours
  3. measure the transverse mass and decay angle and unfold for the true mass distribution plus background

Documentation

Some technical documentation is available here: http://www.desy.de/~sschmitt

Note:

For most applications it is better to use the derived class TUnfoldSys or even better TUnfoldDensity

TUnfoldSys extends the functionality of TUnfold such that systematic errors are propagated to teh result and that the unfolding can be done with proper background subtraction

TUnfoldDensity extends further the functionality of TUnfoldSys complex binning schemes are supported The binning of input histograms is handeld consistently: (1) the regularisation may be done by density, i.e respecting the bin widths (2) methods are provided which preserve the proper binning of the result histograms

Implementation

The result of the unfolding is calculated as follows:

Lsquared = L::L regularisation conditions squared

epsilon_j = sum_i A_ij vector of efficiencies

E^-1 = ((A# Vyy^-1 A)+tau^2 Lsquared)

x = E (A# Vyy^-1 y + tau^2 Lsquared x0 +lambda/2 * epsilon) is the result

The derivatives dx_k/dy_i dx_k/dA_ij dx_k/d(tau^2) are calculated for further usage.

The covariance matrix V_xx is calculated as: Vxx_ij = sum_kl dx_i/dy_k Vyy_kl dx_j/dy_l

Warning:

The algorithm is based on "standard" matrix inversion, with the known limitations in numerical accuracy and computing cost for matrices with large dimensions.

Thus the algorithm should not used for large dimensions of x and y nx should not be much larger than 200 ny should not be much larger than 1000

Proper choice of tau

One of the difficult questions is about the choice of tau. The method implemented in TUnfold is the L-curve method: a two-dimensional curve is plotted x-axis: log10(chisquare) y-axis: log10(regularisation condition) In many cases this curve has an L-shape. The best choice of tau is in the kink of the L

Within TUnfold a simple version of the L-curve analysis is available. It tests a given number of points in a predefined tau-range and searches for the maximum of the curvature in the L-curve (kink position). if no tau range is given, the range of the scan is determined automatically

A nice overview of the L-curve method is given in: The L-curve and Its Use in the Numerical Treatment of Inverse Problems (2000) by P. C. Hansen, in Computational Inverse Problems in Electrocardiology, ed. P. Johnston, Advances in Computational Bioengineering http://www.imm.dtu.dk/~pch/TR/Lcurve.ps

Alternative Regularisation conditions

Regularisation is needed for most unfolding problems, in order to avoid large oscillations and large correlations on the output bins. It means that some extra conditions are applied on the output bins

Within TUnfold these conditions are posed on the difference (x-x0), where x: unfolding output x0: the bias distribution, by default calculated from the input matrix A. There is a method SetBias() to change the bias distribution. The 3rd argument to DoUnfold() is a scale factor applied to the bias bias_default[j] = sum_i A[i][j] x0[j] = scaleBias*bias[j] The scale factor can be used to (a) completely suppress the bias by setting it to zero (b) compensate differences in the normalisation between data and Monte Carlo

If the regularisation is strong, i.e. large parameter tau, then the distribution x or its derivatives will look like the bias distribution. If the parameter tau is small, the distribution x is independent of the bias.

Three basic types of regularisation are implemented in TUnfold

condition regularisation
kRegModeNone none
kRegModeSize minimize the size of (x-x0)
kRegModeDerivative minimize the 1st derivative of (x-x0)
kRegModeCurvature minimize the 2nd derivative of (x-x0)

kRegModeSize is the regularisation scheme which often is found in literature. The second derivative is often named curvature. Sometimes the bias is not discussed, equivalent to a bias scale factor of zero.

The regularisation schemes kRegModeDerivative and kRegModeCurvature have the nice feature that they create correlations between x-bins, whereas the non-regularized unfolding tends to create negative correlations between bins. For these regularisation schemes the parameter tau could be tuned such that the correlations are smallest, as an alternative to the L-curve method.

If kRegModeSize is chosen or if x is a smooth function through all bins, the regularisation condition can be set on all bins together by giving the appropriate argument in the constructor (see examples above).

If x is composed of independent groups of bins (for example, signal and background binning in two variables), it may be necessary to set regularisation conditions for the individual groups of bins. In this case, give kRegModeNone in the constructor and specify the bin grouping with calls to RegularizeBins() specify a 1-dimensional group of bins RegularizeBins2D() specify a 2-dimensional group of bins

Note, the class TUnfoldDensity provides an automatic setup of complex regularisation schemes

For ultimate flexibility, the regularisation condition can be set on each bin individually -> give kRegModeNone in the constructor and use RegularizeSize() regularize one bin RegularizeDerivative() regularize the slope given by two bins RegularizeCurvature() regularize the curvature given by three bins AddRegularisationCondition() define an arbitrary regulatisation condition

Definition at line 99 of file TUnfold.h.

Public Types

enum  EConstraint { kEConstraintNone =0, kEConstraintArea =1 }
 
enum  EHistMap { kHistMapOutputHoriz = 0, kHistMapOutputVert = 1 }
 
enum  ERegMode {
  kRegModeNone = 0, kRegModeSize = 1, kRegModeDerivative = 2, kRegModeCurvature = 3,
  kRegModeMixed = 4
}
 
- Public Types inherited from TObject
enum  { kIsOnHeap = 0x01000000, kNotDeleted = 0x02000000, kZombie = 0x04000000, kBitMask = 0x00ffffff }
 
enum  { kSingleKey = BIT(0), kOverwrite = BIT(1), kWriteDelete = BIT(2) }
 
enum  EStatusBits {
  kCanDelete = BIT(0), kMustCleanup = BIT(3), kObjInCanvas = BIT(3), kIsReferenced = BIT(4),
  kHasUUID = BIT(5), kCannotPick = BIT(6), kNoContextMenu = BIT(8), kInvalidObject = BIT(13)
}
 

Public Member Functions

 TUnfold (const TH2 *hist_A, EHistMap histmap, ERegMode regmode=kRegModeSize, EConstraint constraint=kEConstraintArea)
 
virtual ~TUnfold (void)
 
Double_t DoUnfold (Double_t tau, const TH1 *hist_y, Double_t scaleBias=0.0)
 
virtual Double_t DoUnfold (Double_t tau)
 
void GetBias (TH1 *bias, const Int_t *binMap=0) const
 
Double_t GetChi2A (void) const
 
Double_t GetChi2L (void) const
 
void GetEmatrix (TH2 *ematrix, const Int_t *binMap=0) const
 
Double_t GetEpsMatrix (void) const
 
void GetFoldedOutput (TH1 *folded, const Int_t *binMap=0) const
 
void GetInput (TH1 *inputData, const Int_t *binMap=0) const
 
void GetInputInverseEmatrix (TH2 *ematrix)
 
void GetL (TH2 *l) const
 
virtual Double_t GetLcurveX (void) const
 
virtual Double_t GetLcurveY (void) const
 
void GetLsquared (TH2 *lsquared) const
 
Int_t GetNdf (void) const
 
void GetNormalisationVector (TH1 *s, const Int_t *binMap=0) const
 
Int_t GetNpar (void) const
 
Int_t GetNr (void) const
 
void GetOutput (TH1 *output, const Int_t *binMap=0) const
 
void GetProbabilityMatrix (TH2 *A, EHistMap histmap) const
 
Double_t GetRhoAvg (void) const
 
Double_t GetRhoI (TH1 *rhoi, const Int_t *binMap=0, TH2 *invEmat=0) const
 
void GetRhoIJ (TH2 *rhoij, const Int_t *binMap=0) const
 
Double_t GetRhoMax (void) const
 
Double_t GetTau (void) const
 
Int_t RegularizeBins (int start, int step, int nbin, ERegMode regmode)
 
Int_t RegularizeBins2D (int start_bin, int step1, int nbin1, int step2, int nbin2, ERegMode regmode)
 
Int_t RegularizeCurvature (int left_bin, int center_bin, int right_bin, Double_t scale_left=1.0, Double_t scale_right=1.0)
 
Int_t RegularizeDerivative (int left_bin, int right_bin, Double_t scale=1.0)
 
Int_t RegularizeSize (int bin, Double_t scale=1.0)
 
virtual Int_t ScanLcurve (Int_t nPoint, Double_t tauMin, Double_t tauMax, TGraph **lCurve, TSpline **logTauX=0, TSpline **logTauY=0)
 
void SetBias (const TH1 *bias)
 
void SetConstraint (EConstraint constraint)
 
void SetEpsMatrix (Double_t eps)
 
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)
 
- Public Member Functions inherited from TObject
 TObject ()
 TObject constructor. More...
 
 TObject (const TObject &object)
 TObject copy ctor. More...
 
virtual ~TObject ()
 TObject destructor. More...
 
void AbstractMethod (const char *method) const
 Use this method to implement an "abstract" method that you don't want to leave purely abstract. More...
 
virtual void AppendPad (Option_t *option="")
 Append graphics object to current pad. More...
 
virtual void Browse (TBrowser *b)
 Browse object. May be overridden for another default action. More...
 
virtual const char * ClassName () const
 Returns name of class to which the object belongs. More...
 
virtual void Clear (Option_t *="")
 
virtual TObjectClone (const char *newname="") const
 Make a clone of an object using the Streamer facility. More...
 
virtual Int_t Compare (const TObject *obj) const
 Compare abstract method. More...
 
virtual void Copy (TObject &object) const
 Copy this to obj. More...
 
virtual void Delete (Option_t *option="")
 Delete this object. More...
 
virtual Int_t DistancetoPrimitive (Int_t px, Int_t py)
 Computes distance from point (px,py) to the object. More...
 
virtual void Draw (Option_t *option="")
 Default Draw method for all objects. More...
 
virtual void DrawClass () const
 Draw class inheritance tree of the class to which this object belongs. More...
 
virtual TObjectDrawClone (Option_t *option="") const
 Draw a clone of this object in the current pad. More...
 
virtual void Dump () const
 Dump contents of object on stdout. More...
 
virtual void Error (const char *method, const char *msgfmt,...) const
 Issue error message. More...
 
virtual void Execute (const char *method, const char *params, Int_t *error=0)
 Execute method on this object with the given parameter string, e.g. More...
 
virtual void Execute (TMethod *method, TObjArray *params, Int_t *error=0)
 Execute method on this object with parameters stored in the TObjArray. More...
 
virtual void ExecuteEvent (Int_t event, Int_t px, Int_t py)
 Execute action corresponding to an event at (px,py). More...
 
virtual void Fatal (const char *method, const char *msgfmt,...) const
 Issue fatal error message. More...
 
virtual TObjectFindObject (const char *name) const
 Must be redefined in derived classes. More...
 
virtual TObjectFindObject (const TObject *obj) const
 Must be redefined in derived classes. More...
 
virtual Option_tGetDrawOption () const
 Get option used by the graphics system to draw this object. More...
 
virtual const char * GetIconName () const
 Returns mime type name of object. More...
 
virtual const char * GetName () const
 Returns name of object. More...
 
virtual char * GetObjectInfo (Int_t px, Int_t py) const
 Returns string containing info about the object at position (px,py). More...
 
virtual Option_tGetOption () const
 
virtual const char * GetTitle () const
 Returns title of object. More...
 
virtual UInt_t GetUniqueID () const
 Return the unique object id. More...
 
virtual Bool_t HandleTimer (TTimer *timer)
 Execute action in response of a timer timing out. More...
 
virtual ULong_t Hash () const
 Return hash value for this object. More...
 
virtual void Info (const char *method, const char *msgfmt,...) const
 Issue info message. More...
 
virtual Bool_t InheritsFrom (const char *classname) const
 Returns kTRUE if object inherits from class "classname". More...
 
virtual Bool_t InheritsFrom (const TClass *cl) const
 Returns kTRUE if object inherits from TClass cl. More...
 
virtual void Inspect () const
 Dump contents of this object in a graphics canvas. More...
 
void InvertBit (UInt_t f)
 
virtual Bool_t IsEqual (const TObject *obj) const
 Default equal comparison (objects are equal if they have the same address in memory). More...
 
virtual Bool_t IsFolder () const
 Returns kTRUE in case object contains browsable objects (like containers or lists of other objects). More...
 
Bool_t IsOnHeap () const
 
virtual Bool_t IsSortable () const
 
Bool_t IsZombie () const
 
virtual void ls (Option_t *option="") const
 The ls function lists the contents of a class on stdout. More...
 
void MayNotUse (const char *method) const
 Use this method to signal that a method (defined in a base class) may not be called in a derived class (in principle against good design since a child class should not provide less functionality than its parent, however, sometimes it is necessary). More...
 
virtual Bool_t Notify ()
 This method must be overridden to handle object notification. More...
 
void Obsolete (const char *method, const char *asOfVers, const char *removedFromVers) const
 Use this method to declare a method obsolete. More...
 
void operator delete (void *ptr)
 Operator delete. More...
 
void operator delete[] (void *ptr)
 Operator delete []. More...
 
voidoperator new (size_t sz)
 
voidoperator new (size_t sz, void *vp)
 
voidoperator new[] (size_t sz)
 
voidoperator new[] (size_t sz, void *vp)
 
TObjectoperator= (const TObject &rhs)
 TObject assignment operator. More...
 
virtual void Paint (Option_t *option="")
 This method must be overridden if a class wants to paint itself. More...
 
virtual void Pop ()
 Pop on object drawn in a pad to the top of the display list. More...
 
virtual void Print (Option_t *option="") const
 This method must be overridden when a class wants to print itself. More...
 
virtual Int_t Read (const char *name)
 Read contents of object with specified name from the current directory. More...
 
virtual void RecursiveRemove (TObject *obj)
 Recursively remove this object from a list. More...
 
void ResetBit (UInt_t f)
 
virtual void SaveAs (const char *filename="", Option_t *option="") const
 Save this object in the file specified by filename. More...
 
virtual void SavePrimitive (std::ostream &out, Option_t *option="")
 Save a primitive as a C++ statement(s) on output stream "out". More...
 
void SetBit (UInt_t f, Bool_t set)
 Set or unset the user status bits as specified in f. More...
 
void SetBit (UInt_t f)
 
virtual void SetDrawOption (Option_t *option="")
 Set drawing option for object. More...
 
virtual void SetUniqueID (UInt_t uid)
 Set the unique object id. More...
 
virtual void SysError (const char *method, const char *msgfmt,...) const
 Issue system error message. More...
 
Bool_t TestBit (UInt_t f) const
 
Int_t TestBits (UInt_t f) const
 
virtual void UseCurrentStyle ()
 Set current style settings in this object This function is called when either TCanvas::UseCurrentStyle or TROOT::ForceStyle have been invoked. More...
 
virtual void Warning (const char *method, const char *msgfmt,...) const
 Issue warning message. More...
 
virtual Int_t Write (const char *name=0, Int_t option=0, Int_t bufsize=0)
 Write this object to the current directory. More...
 
virtual Int_t Write (const char *name=0, Int_t option=0, Int_t bufsize=0) const
 Write this object to the current directory. More...
 

Static Public Member Functions

static const char * GetTUnfoldVersion (void)
 
- Static Public Member Functions inherited from TObject
static Long_t GetDtorOnly ()
 Return destructor only flag. More...
 
static Bool_t GetObjectStat ()
 Get status of object stat flag. More...
 
static void SetDtorOnly (void *obj)
 Set destructor only flag. More...
 
static void SetObjectStat (Bool_t stat)
 Turn on/off tracking of objects in the TObjectTable. More...
 

Protected Member Functions

 TUnfold (void)
 
void AddMSparse (TMatrixDSparse *dest, Double_t f, const TMatrixDSparse *src) const
 
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.)
 
Bool_t AddRegularisationCondition (Int_t nEle, const Int_t *indices, const Double_t *rowData)
 
void ClearHistogram (TH1 *h, Double_t x=0.) const
 
virtual void ClearResults (void)
 
TMatrixDSparseCreateSparseMatrix (Int_t nrow, Int_t ncol, Int_t nele, Int_t *row, Int_t *col, Double_t *data) const
 
virtual Double_t DoUnfold (void)
 
void ErrorMatrixToHist (TH2 *ematrix, const TMatrixDSparse *emat, const Int_t *binMap, Bool_t doClear) const
 
const TMatrixDSparseGetAx (void) const
 
Int_t GetBinFromRow (int ix) const
 
const TMatrixDSparseGetDXDAM (int i) const
 
const TMatrixDSparseGetDXDAZ (int i) const
 
const TMatrixDSparseGetDXDtauSquared (void) const
 
const TMatrixDSparseGetDXDY (void) const
 
const TMatrixDSparseGetE (void) const
 
const TMatrixDSparseGetEinv (void) const
 
Int_t GetNx (void) const
 
Int_t GetNy (void) const
 
virtual TString GetOutputBinName (Int_t iBinX) const
 
Double_t GetRhoIFromMatrix (TH1 *rhoi, const TMatrixDSparse *eOrig, const Int_t *binMap, TH2 *invEmat) const
 
Int_t GetRowFromBin (int ix) const
 
const TMatrixDSparseGetVxx (void) const
 
const TMatrixDSparseGetVxxInv (void) const
 
const TMatrixDSparseGetVyyInv (void) const
 
const TMatrixDGetX (void) const
 
TMatrixDSparseInvertMSparseSymmPos (const TMatrixDSparse *A, Int_t *rank) const
 
TMatrixDSparseMultiplyMSparseM (const TMatrixDSparse *a, const TMatrixD *b) const
 
TMatrixDSparseMultiplyMSparseMSparse (const TMatrixDSparse *a, const TMatrixDSparse *b) const
 
TMatrixDSparseMultiplyMSparseMSparseTranspVector (const TMatrixDSparse *m1, const TMatrixDSparse *m2, const TMatrixTBase< Double_t > *v) const
 
TMatrixDSparseMultiplyMSparseTranspMSparse (const TMatrixDSparse *a, const TMatrixDSparse *b) const
 
- Protected Member Functions inherited from TObject
virtual void DoError (int level, const char *location, const char *fmt, va_list va) const
 Interface to ErrorHandler (protected). More...
 
void MakeZombie ()
 

Static Protected Member Functions

static void DeleteMatrix (TMatrixD **m)
 
static void DeleteMatrix (TMatrixDSparse **m)
 

Protected Attributes

TMatrixDSparsefA
 
Double_t fBiasScale
 
EConstraint fConstraint
 
TArrayI fHistToX
 
TMatrixDSparsefL
 
ERegMode fRegMode
 
TArrayD fSumOverY
 
Double_t fTauSquared
 
TMatrixDSparsefVyy
 
TMatrixDfX0
 
TArrayI fXToHist
 
TMatrixDfY
 

Private Member Functions

void InitTUnfold (void)
 

Private Attributes

TMatrixDSparsefAx
 
Double_t fChi2A
 
TMatrixDSparsefDXDAM [2]
 
TMatrixDSparsefDXDAZ [2]
 
TMatrixDSparsefDXDtauSquared
 
TMatrixDSparsefDXDY
 
TMatrixDSparsefE
 
TMatrixDSparsefEinv
 
Double_t fEpsMatrix
 
Int_t fIgnoredBins
 
Double_t fLXsquared
 
Int_t fNdf
 
Double_t fRhoAvg
 
Double_t fRhoMax
 
TMatrixDSparsefVxx
 
TMatrixDSparsefVxxInv
 
TMatrixDSparsefVyyInv
 
TMatrixDfX
 

#include <TUnfold.h>

Inheritance diagram for TUnfold:
[legend]

Member Enumeration Documentation

◆ EConstraint

Enumerator
kEConstraintNone 
kEConstraintArea 

Definition at line 103 of file TUnfold.h.

◆ EHistMap

Enumerator
kHistMapOutputHoriz 
kHistMapOutputVert 

Definition at line 188 of file TUnfold.h.

◆ ERegMode

Enumerator
kRegModeNone 
kRegModeSize 
kRegModeDerivative 
kRegModeCurvature 
kRegModeMixed 

Definition at line 107 of file TUnfold.h.

Constructor & Destructor Documentation

◆ TUnfold() [1/2]

TUnfold::TUnfold ( void  )
protected

Definition at line 372 of file TUnfold.cxx.

◆ TUnfold() [2/2]

TUnfold::TUnfold ( const TH2 hist_A,
EHistMap  histmap,
ERegMode  regmode = kRegModeSize,
EConstraint  constraint = kEConstraintArea 
)

Definition at line 1745 of file TUnfold.cxx.

◆ ~TUnfold()

TUnfold::~TUnfold ( void  )
virtual

Definition at line 1948 of file TUnfold.cxx.

Member Function Documentation

◆ AddMSparse()

void TUnfold::AddMSparse ( TMatrixDSparse dest,
Double_t  f,
const TMatrixDSparse src 
) const
protected

Definition at line 1001 of file TUnfold.cxx.

◆ AddRegularisationCondition() [1/2]

Bool_t TUnfold::AddRegularisationCondition ( Int_t  i0,
Double_t  f0,
Int_t  i1 = -1,
Double_t  f1 = 0.,
Int_t  i2 = -1,
Double_t  f2 = 0. 
)
protected

Definition at line 1973 of file TUnfold.cxx.

◆ AddRegularisationCondition() [2/2]

Bool_t TUnfold::AddRegularisationCondition ( Int_t  nEle,
const Int_t indices,
const Double_t rowData 
)
protected

Definition at line 2007 of file TUnfold.cxx.

◆ ClearHistogram()

void TUnfold::ClearHistogram ( TH1 h,
Double_t  x = 0. 
) const
protected

Definition at line 3394 of file TUnfold.cxx.

◆ ClearResults()

void TUnfold::ClearResults ( void  )
protectedvirtual

Reimplemented in TUnfoldSys.

Definition at line 345 of file TUnfold.cxx.

◆ CreateSparseMatrix()

TMatrixDSparse * TUnfold::CreateSparseMatrix ( Int_t  nrow,
Int_t  ncol,
Int_t  nele,
Int_t row,
Int_t col,
Double_t data 
) const
protected

Definition at line 698 of file TUnfold.cxx.

◆ DeleteMatrix() [1/2]

void TUnfold::DeleteMatrix ( TMatrixD **  m)
staticprotected

Definition at line 333 of file TUnfold.cxx.

◆ DeleteMatrix() [2/2]

void TUnfold::DeleteMatrix ( TMatrixDSparse **  m)
staticprotected

Definition at line 339 of file TUnfold.cxx.

◆ DoUnfold() [1/3]

Double_t TUnfold::DoUnfold ( void  )
protectedvirtual

Definition at line 378 of file TUnfold.cxx.

◆ DoUnfold() [2/3]

Double_t TUnfold::DoUnfold ( Double_t  tau,
const TH1 hist_y,
Double_t  scaleBias = 0.0 
)

Definition at line 2186 of file TUnfold.cxx.

◆ DoUnfold() [3/3]

Double_t TUnfold::DoUnfold ( Double_t  tau)
virtual

Definition at line 2407 of file TUnfold.cxx.

◆ ErrorMatrixToHist()

void TUnfold::ErrorMatrixToHist ( TH2 ematrix,
const TMatrixDSparse emat,
const Int_t binMap,
Bool_t  doClear 
) const
protected

Definition at line 3109 of file TUnfold.cxx.

◆ GetAx()

const TMatrixDSparse* TUnfold::GetAx ( void  ) const
inlineprotected

Definition at line 174 of file TUnfold.h.

◆ GetBias()

void TUnfold::GetBias ( TH1 bias,
const Int_t binMap = 0 
) const

Definition at line 2794 of file TUnfold.cxx.

◆ GetBinFromRow()

Int_t TUnfold::GetBinFromRow ( int  ix) const
inlineprotected

Definition at line 182 of file TUnfold.h.

◆ GetChi2A()

Double_t TUnfold::GetChi2A ( void  ) const
inline

Definition at line 239 of file TUnfold.h.

◆ GetChi2L()

Double_t TUnfold::GetChi2L ( void  ) const

Definition at line 3009 of file TUnfold.cxx.

◆ GetDXDAM()

const TMatrixDSparse* TUnfold::GetDXDAM ( int  i) const
inlineprotected

Definition at line 171 of file TUnfold.h.

◆ GetDXDAZ()

const TMatrixDSparse* TUnfold::GetDXDAZ ( int  i) const
inlineprotected

Definition at line 172 of file TUnfold.h.

◆ GetDXDtauSquared()

const TMatrixDSparse* TUnfold::GetDXDtauSquared ( void  ) const
inlineprotected

Definition at line 173 of file TUnfold.h.

◆ GetDXDY()

const TMatrixDSparse* TUnfold::GetDXDY ( void  ) const
inlineprotected

Definition at line 170 of file TUnfold.h.

◆ GetE()

const TMatrixDSparse* TUnfold::GetE ( void  ) const
inlineprotected

Definition at line 176 of file TUnfold.h.

◆ GetEinv()

const TMatrixDSparse* TUnfold::GetEinv ( void  ) const
inlineprotected

Definition at line 175 of file TUnfold.h.

◆ GetEmatrix()

void TUnfold::GetEmatrix ( TH2 ematrix,
const Int_t binMap = 0 
) const

Definition at line 3175 of file TUnfold.cxx.

◆ GetEpsMatrix()

Double_t TUnfold::GetEpsMatrix ( void  ) const
inline

Definition at line 246 of file TUnfold.h.

◆ GetFoldedOutput()

void TUnfold::GetFoldedOutput ( TH1 folded,
const Int_t binMap = 0 
) const

Definition at line 2815 of file TUnfold.cxx.

◆ GetInput()

void TUnfold::GetInput ( TH1 inputData,
const Int_t binMap = 0 
) const

Definition at line 2884 of file TUnfold.cxx.

◆ GetInputInverseEmatrix()

void TUnfold::GetInputInverseEmatrix ( TH2 ematrix)

Definition at line 2917 of file TUnfold.cxx.

◆ GetL()

void TUnfold::GetL ( TH2 l) const

Definition at line 2977 of file TUnfold.cxx.

◆ GetLcurveX()

Double_t TUnfold::GetLcurveX ( void  ) const
virtual

Definition at line 3021 of file TUnfold.cxx.

◆ GetLcurveY()

Double_t TUnfold::GetLcurveY ( void  ) const
virtual

Definition at line 3027 of file TUnfold.cxx.

◆ GetLsquared()

void TUnfold::GetLsquared ( TH2 lsquared) const

Definition at line 2958 of file TUnfold.cxx.

◆ GetNdf()

Int_t TUnfold::GetNdf ( void  ) const
inline

Definition at line 243 of file TUnfold.h.

◆ GetNormalisationVector()

void TUnfold::GetNormalisationVector ( TH1 s,
const Int_t binMap = 0 
) const

Definition at line 2774 of file TUnfold.cxx.

◆ GetNpar()

Int_t TUnfold::GetNpar ( void  ) const

Definition at line 3015 of file TUnfold.cxx.

◆ GetNr()

Int_t TUnfold::GetNr ( void  ) const
inline

Definition at line 231 of file TUnfold.h.

◆ GetNx()

Int_t TUnfold::GetNx ( void  ) const
inlineprotected

Definition at line 162 of file TUnfold.h.

◆ GetNy()

Int_t TUnfold::GetNy ( void  ) const
inlineprotected

Definition at line 165 of file TUnfold.h.

◆ GetOutput()

void TUnfold::GetOutput ( TH1 output,
const Int_t binMap = 0 
) const

Definition at line 3033 of file TUnfold.cxx.

◆ GetOutputBinName()

TString TUnfold::GetOutputBinName ( Int_t  iBinX) const
protectedvirtual

Reimplemented in TUnfoldDensity.

Definition at line 1737 of file TUnfold.cxx.

◆ GetProbabilityMatrix()

void TUnfold::GetProbabilityMatrix ( TH2 A,
EHistMap  histmap 
) const

Definition at line 2863 of file TUnfold.cxx.

◆ GetRhoAvg()

Double_t TUnfold::GetRhoAvg ( void  ) const
inline

Definition at line 238 of file TUnfold.h.

◆ GetRhoI()

Double_t TUnfold::GetRhoI ( TH1 rhoi,
const Int_t binMap = 0,
TH2 invEmat = 0 
) const

Definition at line 3216 of file TUnfold.cxx.

◆ GetRhoIFromMatrix()

Double_t TUnfold::GetRhoIFromMatrix ( TH1 rhoi,
const TMatrixDSparse eOrig,
const Int_t binMap,
TH2 invEmat 
) const
protected

Definition at line 3275 of file TUnfold.cxx.

◆ GetRhoIJ()

void TUnfold::GetRhoIJ ( TH2 rhoij,
const Int_t binMap = 0 
) const

Definition at line 3188 of file TUnfold.cxx.

◆ GetRhoMax()

Double_t TUnfold::GetRhoMax ( void  ) const
inline

Definition at line 237 of file TUnfold.h.

◆ GetRowFromBin()

Int_t TUnfold::GetRowFromBin ( int  ix) const
inlineprotected

Definition at line 181 of file TUnfold.h.

◆ GetTau()

Double_t TUnfold::GetTau ( void  ) const

Definition at line 3003 of file TUnfold.cxx.

◆ GetTUnfoldVersion()

const char * TUnfold::GetTUnfoldVersion ( void  )
static

Definition at line 289 of file TUnfold.cxx.

◆ GetVxx()

const TMatrixDSparse* TUnfold::GetVxx ( void  ) const
inlineprotected

Definition at line 177 of file TUnfold.h.

◆ GetVxxInv()

const TMatrixDSparse* TUnfold::GetVxxInv ( void  ) const
inlineprotected

Definition at line 178 of file TUnfold.h.

◆ GetVyyInv()

const TMatrixDSparse* TUnfold::GetVyyInv ( void  ) const
inlineprotected

Definition at line 179 of file TUnfold.h.

◆ GetX()

const TMatrixD* TUnfold::GetX ( void  ) const
inlineprotected

Definition at line 180 of file TUnfold.h.

◆ InitTUnfold()

void TUnfold::InitTUnfold ( void  )
private

Definition at line 295 of file TUnfold.cxx.

◆ InvertMSparseSymmPos()

TMatrixDSparse * TUnfold::InvertMSparseSymmPos ( const TMatrixDSparse A,
Int_t rank 
) const
protected

Definition at line 1064 of file TUnfold.cxx.

◆ MultiplyMSparseM()

TMatrixDSparse * TUnfold::MultiplyMSparseM ( const TMatrixDSparse a,
const TMatrixD b 
) const
protected

Definition at line 858 of file TUnfold.cxx.

◆ MultiplyMSparseMSparse()

TMatrixDSparse * TUnfold::MultiplyMSparseMSparse ( const TMatrixDSparse a,
const TMatrixDSparse b 
) const
protected

Definition at line 711 of file TUnfold.cxx.

◆ MultiplyMSparseMSparseTranspVector()

TMatrixDSparse * TUnfold::MultiplyMSparseMSparseTranspVector ( const TMatrixDSparse m1,
const TMatrixDSparse m2,
const TMatrixTBase< Double_t > *  v 
) const
protected

Definition at line 911 of file TUnfold.cxx.

◆ MultiplyMSparseTranspMSparse()

TMatrixDSparse * TUnfold::MultiplyMSparseTranspMSparse ( const TMatrixDSparse a,
const TMatrixDSparse b 
) const
protected

Definition at line 780 of file TUnfold.cxx.

◆ RegularizeBins()

Int_t TUnfold::RegularizeBins ( int  start,
int  step,
int  nbin,
ERegMode  regmode 
)

Definition at line 2123 of file TUnfold.cxx.

◆ RegularizeBins2D()

Int_t TUnfold::RegularizeBins2D ( int  start_bin,
int  step1,
int  nbin1,
int  step2,
int  nbin2,
ERegMode  regmode 
)

Definition at line 2163 of file TUnfold.cxx.

◆ RegularizeCurvature()

Int_t TUnfold::RegularizeCurvature ( int  left_bin,
int  center_bin,
int  right_bin,
Double_t  scale_left = 1.0,
Double_t  scale_right = 1.0 
)

Definition at line 2099 of file TUnfold.cxx.

◆ RegularizeDerivative()

Int_t TUnfold::RegularizeDerivative ( int  left_bin,
int  right_bin,
Double_t  scale = 1.0 
)

Definition at line 2083 of file TUnfold.cxx.

◆ RegularizeSize()

Int_t TUnfold::RegularizeSize ( int  bin,
Double_t  scale = 1.0 
)

Definition at line 2069 of file TUnfold.cxx.

◆ ScanLcurve()

Int_t TUnfold::ScanLcurve ( Int_t  nPoint,
Double_t  tauMin,
Double_t  tauMax,
TGraph **  lCurve,
TSpline **  logTauX = 0,
TSpline **  logTauY = 0 
)
virtual

Definition at line 2424 of file TUnfold.cxx.

◆ SetBias()

void TUnfold::SetBias ( const TH1 bias)

Definition at line 1961 of file TUnfold.cxx.

◆ SetConstraint()

void TUnfold::SetConstraint ( EConstraint  constraint)

Definition at line 2995 of file TUnfold.cxx.

◆ SetEpsMatrix()

void TUnfold::SetEpsMatrix ( Double_t  eps)

Definition at line 3418 of file TUnfold.cxx.

◆ SetInput()

Int_t TUnfold::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 
)
virtual

Reimplemented in TUnfoldSys.

Definition at line 2208 of file TUnfold.cxx.

Member Data Documentation

◆ fA

TMatrixDSparse* TUnfold::fA
protected

Definition at line 115 of file TUnfold.h.

◆ fAx

TMatrixDSparse* TUnfold::fAx
private

Definition at line 134 of file TUnfold.h.

◆ fBiasScale

Double_t TUnfold::fBiasScale
protected

Definition at line 121 of file TUnfold.h.

◆ fChi2A

Double_t TUnfold::fChi2A
private

Definition at line 135 of file TUnfold.h.

◆ fConstraint

EConstraint TUnfold::fConstraint
protected

Definition at line 125 of file TUnfold.h.

◆ fDXDAM

TMatrixDSparse* TUnfold::fDXDAM[2]
private

Definition at line 140 of file TUnfold.h.

◆ fDXDAZ

TMatrixDSparse* TUnfold::fDXDAZ[2]
private

Definition at line 141 of file TUnfold.h.

◆ fDXDtauSquared

TMatrixDSparse* TUnfold::fDXDtauSquared
private

Definition at line 142 of file TUnfold.h.

◆ fDXDY

TMatrixDSparse* TUnfold::fDXDY
private

Definition at line 143 of file TUnfold.h.

◆ fE

TMatrixDSparse* TUnfold::fE
private

Definition at line 145 of file TUnfold.h.

◆ fEinv

TMatrixDSparse* TUnfold::fEinv
private

Definition at line 144 of file TUnfold.h.

◆ fEpsMatrix

Double_t TUnfold::fEpsMatrix
private

Definition at line 129 of file TUnfold.h.

◆ fHistToX

TArrayI TUnfold::fHistToX
protected

Definition at line 123 of file TUnfold.h.

◆ fIgnoredBins

Int_t TUnfold::fIgnoredBins
private

Definition at line 128 of file TUnfold.h.

◆ fL

TMatrixDSparse* TUnfold::fL
protected

Definition at line 116 of file TUnfold.h.

◆ fLXsquared

Double_t TUnfold::fLXsquared
private

Definition at line 136 of file TUnfold.h.

◆ fNdf

Int_t TUnfold::fNdf
private

Definition at line 139 of file TUnfold.h.

◆ fRegMode

ERegMode TUnfold::fRegMode
protected

Definition at line 126 of file TUnfold.h.

◆ fRhoAvg

Double_t TUnfold::fRhoAvg
private

Definition at line 138 of file TUnfold.h.

◆ fRhoMax

Double_t TUnfold::fRhoMax
private

Definition at line 137 of file TUnfold.h.

◆ fSumOverY

TArrayD TUnfold::fSumOverY
protected

Definition at line 124 of file TUnfold.h.

◆ fTauSquared

Double_t TUnfold::fTauSquared
protected

Definition at line 120 of file TUnfold.h.

◆ fVxx

TMatrixDSparse* TUnfold::fVxx
private

Definition at line 132 of file TUnfold.h.

◆ fVxxInv

TMatrixDSparse* TUnfold::fVxxInv
private

Definition at line 133 of file TUnfold.h.

◆ fVyy

TMatrixDSparse* TUnfold::fVyy
protected

Definition at line 117 of file TUnfold.h.

◆ fVyyInv

TMatrixDSparse* TUnfold::fVyyInv
private

Definition at line 131 of file TUnfold.h.

◆ fX

TMatrixD* TUnfold::fX
private

Definition at line 130 of file TUnfold.h.

◆ fX0

TMatrixD* TUnfold::fX0
protected

Definition at line 119 of file TUnfold.h.

◆ fXToHist

TArrayI TUnfold::fXToHist
protected

Definition at line 122 of file TUnfold.h.

◆ fY

TMatrixD* TUnfold::fY
protected

Definition at line 118 of file TUnfold.h.


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