137 DeleteMatrix(&fVyyInv);
221 for(
Int_t i=0;i<2;i++) {
360 Fatal(
"Unfold",
"epsilon#Eepsilon has dimension %d != 1",
472 "epsilon#fDXDtauSquared has dimension %d != 1",
521 Warning(
"DoUnfold",
"rank of output covariance is %d expect %d",
607 if(
a->GetNcols()!=
b->GetNrows()) {
608 Fatal(
"MultiplyMSparseMSparse",
609 "inconsistent matrix col/ matrix row %d !=%d",
610 a->GetNcols(),
b->GetNrows());
681 if(
a->GetNrows() !=
b->GetNrows()) {
682 Fatal(
"MultiplyMSparseTranspMSparse",
683 "inconsistent matrix row numbers %d!=%d",
684 a->GetNrows(),
b->GetNrows());
698 typedef std::map<Int_t, MMatrixRow_t >
MMatrix_t;
721 n += (*irow).second.size();
764 if(
a->GetNcols()!=
b->GetNrows()) {
765 Fatal(
"MultiplyMSparseM",
"inconsistent matrix col /matrix row %d!=%d",
766 a->GetNcols(),
b->GetNrows());
825 (
v && ((m1->
GetNcols()!=
v->GetNrows())||(
v->GetNcols()!=1)))) {
827 Fatal(
"MultiplyMSparseMSparseTranspVector",
828 "matrix cols/vector rows %d!=%d!=%d or vector rows %d!=1\n",
831 Fatal(
"MultiplyMSparseMSparseTranspVector",
926 if((
dest->GetNrows()!=
src->GetNrows())||
927 (
dest->GetNcols()!=
src->GetNcols())) {
928 Fatal(
"AddMSparse",
"inconsistent matrix rows %d!=%d OR cols %d!=%d",
929 src->GetNrows(),
dest->GetNrows(),
src->GetNcols(),
dest->GetNcols());
936 for(
Int_t row=0;row<
dest->GetNrows();row++) {
957 Fatal(
"AddMSparse",
"Nan detected %d %d %d",
998 Fatal(
"InvertMSparseSymmPos",
"inconsistent matrix row/col %d!=%d",
1024 Fatal(
"InvertMSparseSymmPos",
1025 "Matrix has %d negative elements on the diagonal",
nError);
1082#ifdef CONDITION_BLOCK_PART
1085 for(
int i=
inc;i<nn;i++) {
1102 std::cout<<
" "<<i<<
" "<<swap[i]<<
" "<<
swapBack[i]<<
"\n";
1104 std::cout<<
"after sorting\n";
1106 if(i==
iDiagonal) std::cout<<
"iDiagonal="<<i<<
"\n";
1107 if(i==
iBlock) std::cout<<
"iBlock="<<i<<
"\n";
1108 std::cout<<
" "<<swap[i]<<
" "<<
aII(swap[i])<<
"\n";
1132 Fatal(
"InvertMSparseSymmPos",
"sparse matrix analysis failed %d %d %d %d %d",
1138 Info(
"InvertMSparseSymmPos",
"iDiagonal=%d iBlock=%d nRow=%d",
1209 Fatal(
"InvertMSparseSymmPos",
1210 "diagonal part 1 has rank %d != %d, matrix can not be inverted",
1238 Fatal(
"InvertMSparseSymmPos",
1239 "diagonal part 2 has rank %d != %d, matrix can not be inverted",
1292#ifndef FORCE_EIGENVALUE_DECOMPOSITION
1351 for(
Int_t k=0;k<i;k++) {
1366 std::cout<<
"dmin,dmax: "<<
dmin<<
" "<<
dmax<<
"\n";
1375 cinv(i,i)=1./
c(i,i);
1380 for(
Int_t k=i+1;k<
j;k++) {
1473 for(
Int_t iF=0;iF<
Finv->GetNrows();iF++) {
1488 Fatal(
"InvertMSparseSymmPos",
1489 "non-trivial part has rank < %d, matrix can not be inverted",
1496 Info(
"InvertMSparseSymmPos",
1497 "cholesky-decomposition failed, try eigenvalue analysis");
1499 std::cout<<
"nEV="<<
nEV<<
" iDiagonal="<<
iDiagonal<<
"\n";
1510 if((i<0)||(
j<0)||(i>=
nEV)||(
j>=
nEV)) {
1511 std::cout<<
" error "<<
nEV<<
" "<<i<<
" "<<
j<<
"\n";
1515 Fatal(
"InvertMSparseSymmPos",
1516 "non-finite number detected element %d %d\n",
1526 std::cout<<
"Eigenvalues\n";
1531 if(
Eigen.GetEigenValues()(0)<0.0) {
1533 }
else if(
Eigen.GetEigenValues()(0)>0.0) {
1541 Error(
"InvertMSparseSymmPos",
1542 "Largest Eigenvalue is negative %f",
1543 Eigen.GetEigenValues()(0));
1545 Error(
"InvertMSparseSymmPos",
1546 "Some Eigenvalues are negative (EV%d/EV0=%g epsilon=%g)",
1625 for(
Int_t i=
nullptr;i<
a.GetNrows();i++) {
1626 for(
Int_t j=
nullptr;
j<
a.GetNcols();
j++) {
1630 std::cout<<
"Ar is not symmetric Ar("<<i<<
","<<
j<<
")="<<
ar(i,
j)
1631 <<
" Ar("<<
j<<
","<<i<<
")="<<
ar(
j,i)<<
"\n";
1636 std::cout<<
"ArA is not equal A ArA("<<i<<
","<<
j<<
")="<<
ara(i,
j)
1637 <<
" A("<<i<<
","<<
j<<
")="<<
a(i,
j)<<
"\n";
1642 std::cout<<
"rAr is not equal r rAr("<<i<<
","<<
j<<
")="<<
rar(i,
j)
1643 <<
" r("<<i<<
","<<
j<<
")="<<
R(i,
j)<<
"\n";
1649 std::cout<<
"Matrix is not positive\n";
1739 for (
Int_t ix = 0; ix <
nx0; ix++) {
1744 for (
Int_t iy = 0; iy <
ny; iy++) {
1747 z =
hist_A->GetBinContent(ix, iy + 1);
1749 z =
hist_A->GetBinContent(iy + 1, ix);
1768 hist_A->GetBinContent(ix, 0) +
1772 hist_A->GetBinContent(0, ix) +
1783 for (
Int_t ix = 0; ix <
nx; ix++) {
1797 for (
Int_t ix = 0; ix <
nx0; ix++) {
1822 Info(
"TUnfold",
"underflow and overflow bin "
1823 "do not depend on the input data");
1825 Warning(
"TUnfold",
"%d output bins "
1827 static_cast<const char *
>(
binlist));
1838 for (
Int_t iy = 0; iy <
ny; iy++) {
1839 for (
Int_t ix = 0; ix <
nx; ix++) {
1856 Info(
"TUnfold",
"%d input bins and %d output bins (includes 2 underflow/overflow bins)",
ny,
nx);
1858 Info(
"TUnfold",
"%d input bins and %d output bins (includes 1 underflow bin)",
ny,
nx);
1860 Info(
"TUnfold",
"%d input bins and %d output bins (includes 1 overflow bin)",
ny,
nx);
1862 Info(
"TUnfold",
"%d input bins and %d output bins",
ny,
nx);
1866 Error(
"TUnfold",
"too few (ny=%d) input bins for nx=%d output bins",
ny,
nx);
1868 Warning(
"TUnfold",
"too few (ny=%d) input bins for nx=%d output bins",
ny,
nx);
1880 "%d regularisation conditions have been skipped",
nError);
1883 "One regularisation condition has been skipped");
2167 Error(
"RegularizeBins",
"regmode = %d is not valid",
regmode);
2346 if(
hist_vyy->GetBinContent(iy+1,
jy+1)!=0.0) {
2354 if(iy==
jy)
continue;
2367 "inverse of input covariance is taken from user input");
2389 "input covariance has elements C(X,Y)!=nullptr where V(X)==0");
2413 (*fY) (i, 0) =
input->GetBinContent(i + 1);
2420 for (
Int_t i = 0; i <
mAtV->GetNrows();i++) {
2421 if(
mAtV->GetRowIndexArray()[i]==
2422 mAtV->GetRowIndexArray()[i+1]) {
2429 Warning(
"SetInput",
"%d/%d input bins have zero error,"
2430 " 1/error set to %lf.",
2433 Warning(
"SetInput",
"One input bin has zero error,"
2439 Warning(
"SetInput",
"%d/%d input bins have zero error,"
2442 Warning(
"SetInput",
"One input bin has zero error,"
2443 " and is ignored.");
2452 for (
Int_t col = 0; col <
mAtV->GetNrows();col++) {
2453 if(
mAtV->GetRowIndexArray()[col]==
2454 mAtV->GetRowIndexArray()[col+1]) {
2471 Error(
"SetInput",
"%d/%d output bins are not constrained by any data.",
2474 Error(
"SetInput",
"One output bin [%d] is not constrained by any data.",
2524 for(
int i=0;i<
r.GetNrows();i++) {
2564 typedef std::map<Double_t,std::pair<Double_t,Double_t> >
XYtau_t;
2599 Error(
"ScanLcurve",
"too few input bins, NDF<=nullptr %d",
GetNdf());
2604 Info(
"ScanLcurve",
"logtau=-Infinity X=%lf Y=%lf",x0,
y0);
2606 Fatal(
"ScanLcurve",
"problem (too few input bins?) X=%f",x0);
2609 Fatal(
"ScanLcurve",
"problem (missing regularisation?) Y=%f",
y0);
2619 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2623 Info(
"ScanLcurve",
"logtau=%lf X=%lf Y=%lf",
2638 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2642 Info(
"ScanLcurve",
"logtau=%lf X=%lf Y=%lf",
2654 (((*
curve.
begin()).second.first-x0>0.00432)||
2656 (
curve.size()<2))) {
2660 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2664 Info(
"ScanLcurve",
"logtau=%lf X=%lf Y=%lf",
2675 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2678 Info(
"ScanLcurve",
"logtau=%lf X=%lf Y=%lf",
2685 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2688 Info(
"ScanLcurve",
"logtau=%lf X=%lf Y=%lf",
2707 const std::pair<Double_t,Double_t> &
xy0=(*i0).second;
2708 const std::pair<Double_t,Double_t> &
xy1=(*i1).second;
2714 logTau=0.5*((*i0).first+(*i1).first);
2720 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2747 lXi[
n]=(*i).second.first;
2748 lYi[
n]=(*i).second.second;
2756 for(
Int_t i=0;i<
n-1;i++) {
2854 Fatal(
"ScanLcurve",
"problem (missing regularisation?) X=%f Y=%f",
2857 Info(
"ScanLcurve",
"Result logtau=%lf X=%lf Y=%lf",
2868 if(!
curve.empty()) {
2877 x[
n]=(*i).second.first;
2878 y[
n]=(*i).second.second;
2884 (*lCurve)->SetTitle(
"L curve");
2944 out->GetBinContent(
dest));
2978 if(
destI<0)
continue;
2980 out->SetBinContent(
destI, (*
fAx) (i, 0)+ out->GetBinContent(
destI));
3080 if(
destI<0)
continue;
3082 out->SetBinContent(
destI, (*
fY) (i, 0)+out->GetBinContent(
destI));
3090 out->SetBinError(
destI,
e);
3110 Warning(
"GetInputInverseEmatrix",
"input covariance matrix has rank %d expect %d",
3114 Error(
"GetInputInverseEmatrix",
"number of parameters %d > %d (rank of input covariance). Problem can not be solved",
GetNpar(),
rank);
3115 }
else if(
fNdf==0) {
3116 Warning(
"GetInputInverseEmatrix",
"number of parameters %d = input rank %d. Problem is ill posed",
GetNpar(),
rank);
3125 for(
int i=0;i<=out->GetNbinsX()+1;i++) {
3126 for(
int j=0;
j<=out->GetNbinsY()+1;
j++) {
3127 out->SetBinContent(i,
j,0.);
3310 std::map<Int_t,Double_t>
e2;
3359 for(std::map<Int_t,Double_t>::const_iterator i=
e2.
begin();
3472 if((
e[i]>0.0)&&(
e[
j]>0.0)) {
3475 rhoij->SetBinContent(i,
j,0.0);
3614 if(
destI<0)
continue;
3624 if(
destJ<0)
continue;
3634 Warning(
"GetRhoIFormMatrix",
"Covariance matrix has rank %d expect %d",
3684 nxyz[0]=
h->GetNbinsX()+1;
3685 nxyz[1]=
h->GetNbinsY()+1;
3686 nxyz[2]=
h->GetNbinsZ()+1;
3687 for(
int i=
h->GetDimension();i<3;i++)
nxyz[i]=0;
3689 for(
int i=0;i<3;i++)
ixyz[i]=0;
3694 h->SetBinContent(
ibin,
x);
3695 h->SetBinError(
ibin,0.0);
3696 for(
Int_t i=0;i<3;i++) {
3793 std::map<double,ScanResult > scan;
3798 while((
int)scan.size()<
nPoint) {
3813 tau=1./
ev(
ev.GetNrows()-1);
3818 std::vector<double> t,s;
3821 for(std::map<double,ScanResult>::const_iterator i=scan.begin();
3822 i!=scan.end();i++) {
3823 t.push_back((*i).first);
3824 s.push_back((*i).second.SURE);
3833 double s0=0.,
s1=0.,
s2=0.;
3835 for(
size_t i=0;i<t.size()-1;i++) {
3865 for(
size_t i=2;i<t.size()-1;i++) {
3873 Info(
"ScanSURE",
"minimum near: [%f,%f,%f] -> [%f,%f,%f}",
3905 if((tau<=0.)&&(
GetNdf()<=0)) {
3906 Error(
"ScanSURE",
"too few input bins, NDF<=nullptr %d",
GetNdf());
3909 Info(
"ScanSURE",
"logtau=-Infinity Chi2A=%lf SURE=%lf DF=%lf X=%lf Y=%lf",
3910 r.chi2A,
r.SURE,
r.DF,
r.x,
r.y);
3912 Fatal(
"ScanSURE",
"problem (too few input bins?) x=%f",
r.x);
3915 Fatal(
"ScanSURE",
"problem (missing regularisation?) y=%f",
r.y);
3918 Info(
"ScanSURE",
"logtau=%lf Chi2A=%lf SURE=%lf DF=%lf X=%lf Y=%lf",
3927 for(std::map<double,ScanResult>::const_iterator i=scan.begin();
3928 i!=scan.end();i++) {
3931 if((*i).first>0.0) {
3936 double s=(*i).second.SURE;
3942 DF.push_back((*i).second.DF);
3943 chi2A.push_back((*i).second.chi2A);
3944 X.push_back((*i).second.x);
3945 Y.push_back((*i).second.y);
#define R(a, b, c, d, e, f, g, h, i)
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void input
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t dest
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t r
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t index
Option_t Option_t TPoint xy
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t src
TMatrixTSparse< Double_t > TMatrixDSparse
TMatrixT< Double_t > TMatrixD
const_iterator begin() const
const_iterator end() const
void Set(Int_t n) override
Set size of this array to n doubles.
const Double_t * GetArray() const
void Set(Int_t n) override
Set size of this array to n ints.
A TGraph is an object made of two arrays X and Y with npoints each.
TH1 is the base class of all histogram classes in ROOT.
Service class for 2-D histogram classes.
void SetBinContent(Int_t bin, Double_t content) override
Set bin content.
TMatrixTBase< Element > & SetMatrixArray(const Element *data, Option_t *="") override
Copy array data to matrix .
const Int_t * GetRowIndexArray() const override
const Int_t * GetColIndexArray() const override
const Element * GetMatrixArray() const override
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
virtual void Fatal(const char *method, const char *msgfmt,...) const
Issue fatal error message.
Class to create third splines to interpolate knots Arbitrary conditions can be introduced for first a...
Base class for spline implementation containing the Draw/Paint methods.
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString.
An algorithm to unfold distributions from detector to truth level.
TArrayI fHistToX
mapping of histogram bins to matrix indices
double GetDF(void) const
return the effecive number of degrees of freedom See e.g.
TMatrixDSparse * fE
matrix E
TMatrixDSparse * fEinv
matrix E^(-1)
virtual Double_t GetLcurveY(void) const
get value on y-axis of L-curve determined in recent unfolding
TMatrixDSparse * fAx
result x folded back A*x
TMatrixDSparse * MultiplyMSparseM(const TMatrixDSparse *a, const TMatrixD *b) const
multiply sparse matrix and a non-sparse matrix
virtual Double_t DoUnfold(void)
core unfolding algorithm
Double_t fChi2A
chi**2 contribution from (y-Ax)Vyy-1(y-Ax)
TMatrixD * fX0
bias vector x0
double GetSURE(void) const
return Stein's unbiased risk estimator See e.g.
void GetBias(TH1 *bias, const Int_t *binMap=nullptr) const
get bias vector including bias scale
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 * 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 RegularizeSize(int bin, Double_t scale=1.0)
add a regularisation condition on the magnitude of a truth bin
Double_t fEpsMatrix
machine accuracy used to determine matrix rank after eigenvalue analysis
void GetProbabilityMatrix(TH2 *A, EHistMap histmap) const
get matrix of probabilities
Double_t GetChi2L(void) const
get χ2L contribution determined in recent unfolding
TMatrixDSparse * fVxx
covariance matrix Vxx
Int_t GetNy(void) const
returns the number of measurement bins
virtual TString GetOutputBinName(Int_t iBinX) const
Get bin name of an outpt bin.
Double_t fBiasScale
scale factor for the bias
virtual Int_t ScanSURE(Int_t nPoint, Double_t tauMin, Double_t tauMax, TGraph **logTauSURE=nullptr, TGraph **df_chi2A=nullptr, TGraph **lCurve=nullptr)
minimize Stein's unbiased risk estimator "SURE" using successive calls to DoUnfold at various tau.
virtual Int_t ScanLcurve(Int_t nPoint, Double_t tauMin, Double_t tauMax, TGraph **lCurve, TSpline **logTauX=nullptr, TSpline **logTauY=nullptr, TSpline **logTauCurvature=nullptr)
scan the L curve, determine tau and unfold at the final value of tau
Double_t fRhoAvg
average global correlation coefficient
TMatrixDSparse * fDXDtauSquared
derivative of the result wrt tau squared
static void DeleteMatrix(TMatrixD **m)
delete matrix and invalidate pointer
void ClearHistogram(TH1 *h, Double_t x=0.) const
Initialize bin contents and bin errors for a given histogram.
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 GetNx(void) const
returns internal number of output (truth) matrix rows
static const char * GetTUnfoldVersion(void)
return a string describing the TUnfold version
void SetConstraint(EConstraint constraint)
set type of area constraint
void GetFoldedOutput(TH1 *folded, const Int_t *binMap=nullptr) const
get unfolding result on detector level
Int_t RegularizeBins(int start, int step, int nbin, ERegMode regmode)
add regularisation conditions for a group of bins
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
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
void SetBias(const TH1 *bias)
set bias vector
void GetL(TH2 *l) const
get matrix of regularisation conditions
ERegMode fRegMode
type of regularisation
Int_t GetNr(void) const
get number of regularisation conditions
TMatrixDSparse * fVxxInv
inverse of covariance matrix Vxx-1
TMatrixD * fX
unfolding result x
EConstraint
type of extra constraint
@ kEConstraintNone
use no extra constraint
virtual Double_t GetLcurveX(void) const
get value on x-axis of L-curve determined in recent unfolding
Double_t GetRhoI(TH1 *rhoi, const Int_t *binMap=nullptr, TH2 *invEmat=nullptr) const
get global correlation coefficiencts, possibly cumulated over several bins
TMatrixDSparse * fVyy
covariance matrix Vyy corresponding to y
Int_t fNdf
number of degrees of freedom
TArrayD fSumOverY
truth vector calculated from the non-normalized response matrix
ERegMode
choice of regularisation scheme
@ kRegModeNone
no regularisation, or defined later by RegularizeXXX() methods
@ kRegModeDerivative
regularize the 1st derivative of the output distribution
@ kRegModeSize
regularise the amplitude of the output distribution
@ kRegModeCurvature
regularize the 2nd derivative of the output distribution
@ kRegModeMixed
mixed regularisation pattern
void GetInput(TH1 *inputData, const Int_t *binMap=nullptr) const
Input vector of measurements.
void SetEpsMatrix(Double_t eps)
set numerical accuracy for Eigenvalue analysis when inverting matrices with rank problems
const TMatrixDSparse * GetE(void) const
matrix E, using internal bin counting
TVectorD GetSqrtEvEmatrix(void) const
void GetOutput(TH1 *output, const Int_t *binMap=nullptr) const
get output distribution, possibly cumulated over several bins
void GetRhoIJ(TH2 *rhoij, const Int_t *binMap=nullptr) const
get correlation coefficiencts, possibly cumulated over several bins
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
TArrayI fXToHist
mapping of matrix indices to histogram bins
TMatrixDSparse * fDXDY
derivative of the result wrt dx/dy
TMatrixD * fY
input (measured) data y
TMatrixDSparse * InvertMSparseSymmPos(const TMatrixDSparse *A, Int_t *rank) const
get the inverse or pseudo-inverse of a positive, sparse matrix
TMatrixDSparse * fVyyInv
inverse of the input covariance matrix Vyy-1
Double_t fLXsquared
chi**2 contribution from (x-s*x0)TLTL(x-s*x0)
TMatrixDSparse * fDXDAM[2]
matrix contribution to the of derivative dx_k/dA_ij
Double_t fTauSquared
regularisation parameter tau squared
Int_t GetNpar(void) const
get number of truth parameters determined in recent unfolding
virtual void ClearResults(void)
reset all results
Double_t fRhoMax
maximum global correlation coefficient
void GetEmatrix(TH2 *ematrix, const Int_t *binMap=nullptr) const
get output covariance matrix, possibly cumulated over several bins
TMatrixDSparse * MultiplyMSparseMSparse(const TMatrixDSparse *a, const TMatrixDSparse *b) const
multiply two sparse matrices
EConstraint fConstraint
type of constraint to use for the unfolding
TUnfold(void)
only for use by root streamer or derived classes
EHistMap
arrangement of axes for the response matrix (TH2 histogram)
@ kHistMapOutputHoriz
truth level on x-axis of the response matrix
void AddMSparse(TMatrixDSparse *dest, Double_t f, const TMatrixDSparse *src) const
add a sparse matrix, scaled by a factor, to another scaled matrix
void GetNormalisationVector(TH1 *s, const Int_t *binMap=nullptr) const
histogram of truth bins, determined from suming over the response matrix
TMatrixDSparse * fDXDAZ[2]
vector contribution to the of derivative dx_k/dA_ij
Double_t GetRhoIFromMatrix(TH1 *rhoi, const TMatrixDSparse *eOrig, const Int_t *binMap, TH2 *invEmat) const
void InitTUnfold(void)
initialize data menbers, for use in constructors
Double_t GetTau(void) const
return regularisation parameter
Double_t GetChi2A(void) const
get χ2A contribution determined in recent unfolding
Int_t RegularizeBins2D(int start_bin, int step1, int nbin1, int step2, int nbin2, ERegMode regmode)
add regularisation conditions for 2d unfolding
const TMatrixDSparse * GetDXDY(void) const
matrix of derivatives dx/dy
void GetLsquared(TH2 *lsquared) const
get matrix of regularisation conditions squared
void GetInputInverseEmatrix(TH2 *ematrix)
get inverse of the measurement's covariance matrix
TMatrixDSparse * fA
response matrix A
TMatrixDSparse * fL
regularisation conditions L
virtual Int_t SetInput(const TH1 *hist_y, Double_t scaleBias=0.0, Double_t oneOverZeroError=0.0, const TH2 *hist_vyy=nullptr, const TH2 *hist_vyy_inv=nullptr)
Define input data for subsequent calls to DoUnfold(tau)
Int_t fIgnoredBins
number of input bins which are dropped because they have error=nullptr
Int_t GetNdf(void) const
get number of degrees of freedom determined in recent unfolding
Short_t Max(Short_t a, Short_t b)
Returns the largest of a and b.
Int_t Finite(Double_t x)
Check if it is finite with a mask in order to be consistent in presence of fast math.
Double_t Log(Double_t x)
Returns the natural logarithm of x.
Double_t Sqrt(Double_t x)
Returns the square root of x.
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
Returns x raised to the power y.
Short_t Min(Short_t a, Short_t b)
Returns the smallest of a and b.
Double_t Log10(Double_t x)
Returns the common (base-10) logarithm of x.
Short_t Abs(Short_t d)
Returns the absolute value of parameter Short_t d.
static uint64_t sum(uint64_t i)