136 if(rnd->
Rndm()>frac) {
137 return rnd->
Gaus(mTrue+smallBias,smallSigma);
139 return rnd->
Gaus(mTrue+wideBias,wideSigma);
152 Double_t const luminosityData=100000;
153 Double_t const luminosityMC=1000000;
156 Int_t const nDet=250;
157 Int_t const nGen=100;
166 TH1D *histMgenMC=
new TH1D(
"MgenMC",
";mass(gen)",nGen,xminGen,xmaxGen);
167 TH1D *histMdetMC=
new TH1D(
"MdetMC",
";mass(det)",nDet,xminDet,xmaxDet);
168 TH2D *histMdetGenMC=
new TH2D(
"MdetgenMC",
";mass(det);mass(gen)",nDet,xminDet,xmaxDet,
169 nGen,xminGen,xmaxGen);
171 for(
Int_t i=0;i<neventMC;i++) {
184 histMgenMC->
Fill(mGen,luminosityData/luminosityMC);
186 histMdetMC->
Fill(mDet,luminosityData/luminosityMC);
201 histMdetGenMC->
Fill(mDet,mGen,luminosityData/luminosityMC);
207 TH1D *histMgenData=
new TH1D(
"MgenData",
";mass(gen)",nGen,xminGen,xmaxGen);
208 TH1D *histMdetData=
new TH1D(
"MdetData",
";mass(det)",nDet,xminDet,xmaxDet);
209 Int_t neventData=rnd->
Poisson(luminosityData*crossSection);
210 for(
Int_t i=0;i<neventData;i++) {
217 histMgenData->
Fill(mGen);
220 histMdetData->
Fill(mDet);
245 Int_t iPeek=(
Int_t)(nGen*(estimatedPeakPosition-xminGen)/(xmaxGen-xminGen)
251 unfold.RegularizeBins(1,1,iPeek-nPeek,regMode);
253 unfold.RegularizeBins(iPeek+nPeek,1,nGen-(iPeek+nPeek),regMode);
259 if(unfold.SetInput(histMdetData,0.0)>=10000) {
260 std::cout<<
"Unfolding result may be wrong\n";
272 iBest=unfold.ScanLcurve(nScan,tauMin,tauMax,&lCurve,&logTauX,&logTauY);
273 std::cout<<
"tau="<<unfold.GetTau()<<
"\n";
274 std::cout<<
"chi**2="<<unfold.GetChi2A()<<
"+"<<unfold.GetChi2L()
275 <<
" / "<<unfold.GetNdf()<<
"\n";
291 for(
Int_t i=1;i<=nGen;i++) binMap[i]=i;
295 TH1D *histMunfold=
new TH1D(
"Unfolded",
";mass(gen)",nGen,xminGen,xmaxGen);
296 unfold.GetOutput(histMunfold,binMap);
297 TH1D *histMdetFold=
new TH1D(
"FoldedBack",
"mass(det)",nDet,xminDet,xmaxDet);
298 unfold.GetFoldedOutput(histMdetFold);
301 TH1D *histRhoi=
new TH1D(
"rho_I",
"mass",nGen,xminGen,xmaxGen);
302 unfold.GetRhoI(histRhoi,binMap);
321 histMdetGenMC->
Draw(
"BOX");
331 histMgenData->
Draw(
"SAME");
332 histMgenMC->
Draw(
"SAME HIST");
340 histMdetFold->
Draw();
342 histMdetData->
Draw(
"SAME");
343 histMdetMC->
Draw(
"SAME HIST");
357 bestLogTauX->
Draw(
"*");
362 bestLcurve->
Draw(
"*");
364 output.SaveAs(
"testUnfold2.ps");
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)