using namespace std;
{
do {
do {
} while(t>=1.0);
} while(t<=0.0);
return t;
} else {
do {
do {
} while(t>=1.0);
} while(t>=0.0);
return t;
}
}
return rnd->
Gaus(mTrue+smallBias,smallSigma);
} else {
return rnd->
Gaus(mTrue+wideBias,wideSigma);
}
}
int testUnfold2()
{
TH1D *histMgenMC=
new TH1D(
"MgenMC",
";mass(gen)",nGen,xminGen,xmaxGen);
TH1D *histMdetMC=
new TH1D(
"MdetMC",
";mass(det)",nDet,xminDet,xmaxDet);
TH2D *histMdetGenMC=
new TH2D(
"MdetgenMC",
";mass(det);mass(gen)",nDet,xminDet,xmaxDet,
nGen,xminGen,xmaxGen);
for(
Int_t i=0;i<neventMC;i++) {
4.0,
0.2);
histMgenMC->
Fill(mGen,luminosityData/luminosityMC);
histMdetMC->
Fill(mDet,luminosityData/luminosityMC);
histMdetGenMC->
Fill(mDet,mGen,luminosityData/luminosityMC);
}
TH1D *histMgenData=
new TH1D(
"MgenData",
";mass(gen)",nGen,xminGen,xmaxGen);
TH1D *histMdetData=
new TH1D(
"MdetData",
";mass(det)",nDet,xminDet,xmaxDet);
for(
Int_t i=0;i<neventData;i++) {
3.8,
0.15);
histMgenData->
Fill(mGen);
histMdetData->
Fill(mDet);
}
Int_t iPeek=(
Int_t)(nGen*(estimatedPeakPosition-xminGen)/(xmaxGen-xminGen)
+1.5);
unfold.RegularizeBins(1,1,iPeek-nPeek,regMode);
unfold.RegularizeBins(iPeek+nPeek,1,nGen-(iPeek+nPeek),regMode);
if(unfold.SetInput(histMdetData,0.0)>=10000) {
std::cout<<"Unfolding result may be wrong\n";
}
iBest=unfold.ScanLcurve(nScan,tauMin,tauMax,&lCurve,&logTauX,&logTauY);
std::cout<<"tau="<<unfold.GetTau()<<"\n";
std::cout<<"chi**2="<<unfold.GetChi2A()<<"+"<<unfold.GetChi2L()
<<" / "<<unfold.GetNdf()<<"\n";
for(
Int_t i=1;i<=nGen;i++) binMap[i]=i;
binMap[0]=-1;
binMap[nGen+1]=-1;
TH1D *histMunfold=
new TH1D(
"Unfolded",
";mass(gen)",nGen,xminGen,xmaxGen);
unfold.GetOutput(histMunfold,binMap);
TH1D *histMdetFold=
new TH1D(
"FoldedBack",
"mass(det)",nDet,xminDet,xmaxDet);
unfold.GetFoldedOutput(histMdetFold);
TH1D *histRhoi=
new TH1D(
"rho_I",
"mass",nGen,xminGen,xmaxGen);
unfold.GetRhoI(histRhoi,binMap);
delete[] binMap;
binMap=0;
histMdetGenMC->
Draw(
"BOX");
histMgenData->
Draw(
"SAME");
histMgenMC->
Draw(
"SAME HIST");
histMdetData->
Draw(
"SAME");
histMdetMC->
Draw(
"SAME HIST");
output.SaveAs(
"testUnfold2.ps");
return 0;
}
virtual void SetLineColor(Color_t lcolor)
Set the line color.
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
A TGraph is an object made of two arrays X and Y with npoints each.
virtual void Draw(Option_t *chopt="")
Draw this graph with its current attributes.
1-D histogram with a double per channel (see TH1 documentation)}
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
static void SetDefaultSumw2(Bool_t sumw2=kTRUE)
When this static function is called with sumw2=kTRUE, all new histograms will automatically activate ...
virtual void Draw(Option_t *option="")
Draw this histogram with options.
2-D histogram with a double per channel (see TH1 documentation)}
Int_t Fill(Double_t)
Invalid Fill method.
Random number generator class based on M.
This is the base class for the ROOT Random number generators.
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
virtual Int_t Poisson(Double_t mean)
Generates a random integer N according to a Poisson law.
virtual Double_t Rndm()
Machine independent random number generator.
Base class for spline implementation containing the Draw/Paint methods.
virtual void Draw(Option_t *option="")
Draw this function with its current attributes.
virtual void GetKnot(Int_t i, Double_t &x, Double_t &y) const =0
An algorithm to unfold distributions from detector to truth level.
ERegMode
choice of regularisation scheme
@ kRegModeNone
no regularisation, or defined later by RegularizeXXX() methods
@ kRegModeCurvature
regularize the 2nd derivative of the output distribution
@ kHistMapOutputVert
truth level on y-axis of the response matrix
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
static void output(int code)