132 for (
Int_t i=0;i<nPoints;i++) {
139 for (
Int_t i=0;i<nPoints;i++) {
140 for (
Int_t j=0;j<nPoints;j++) {
150 return par[0]/(dm*dm+par[2]*par[2]);
163 if(rnd->
Rndm()>bgr) {
203 if(rnd->
Rndm()>frac) {
204 return rnd->
Gaus(mTrue+smallBias,smallSigma);
206 return rnd->
Gaus(mTrue+wideBias,wideSigma);
222 Double_t const luminosityData=100000;
223 Double_t const luminosityMC=1000000;
226 Int_t const nDet=250;
227 Int_t const nGen=100;
236 TH1D *histMgenMC=
new TH1D(
"MgenMC",
";mass(gen)",nGen,xminGen,xmaxGen);
237 TH1D *histMdetMC=
new TH1D(
"MdetMC",
";mass(det)",nDet,xminDet,xmaxDet);
238 TH2D *histMdetGenMC=
new TH2D(
"MdetgenMC",
";mass(det);mass(gen)",
239 nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
241 for(
Int_t i=0;i<neventMC;i++) {
254 histMgenMC->
Fill(mGen,luminosityData/luminosityMC);
256 histMdetMC->
Fill(mDet,luminosityData/luminosityMC);
271 histMdetGenMC->
Fill(mDet,mGen,luminosityData/luminosityMC);
278 TH2D *histMdetGenSysMC=
new TH2D(
"MdetgenSysMC",
";mass(det);mass(gen)",
279 nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
280 neventMC=rnd->
Poisson(luminosityMC*crossSection);
281 for(
Int_t i=0;i<neventMC;i++) {
287 histMdetGenSysMC->
Fill(mDet,mGen,luminosityData/luminosityMC);
293 TH1D *histMgenData=
new TH1D(
"MgenData",
";mass(gen)",nGen,xminGen,xmaxGen);
294 TH1D *histMdetData=
new TH1D(
"MdetData",
";mass(det)",nDet,xminDet,xmaxDet);
295 Int_t neventData=rnd->
Poisson(luminosityData*crossSection);
296 for(
Int_t i=0;i<neventData;i++) {
303 histMgenData->
Fill(mGen);
306 histMdetData->
Fill(mDet);
311 TH1D *histDensityGenData=
new TH1D(
"DensityGenData",
";mass(gen)",
312 nGen,xminGen,xmaxGen);
313 TH1D *histDensityGenMC=
new TH1D(
"DensityGenMC",
";mass(gen)",
314 nGen,xminGen,xmaxGen);
315 for(
Int_t i=1;i<=nGen;i++) {
333 if(unfold.SetInput(histMdetData)>=10000) {
334 std::cout<<
"Unfolding result may be wrong\n";
350#ifdef VERBOSE_LCURVE_SCAN
356 iBest=unfold.ScanLcurve(nScan,tauMin,tauMax,&lCurve,&logTauX,&logTauY);
359#ifdef VERBOSE_LCURVE_SCAN
369 TH2D *histMdetGenSys1=
new TH2D(
"Mdetgensys1",
";mass(det);mass(gen)",
370 nDet,xminDet,xmaxDet,nGen,xminGen,xmaxGen);
371 for(
Int_t i=0;i<=nDet+1;i++) {
373 for(
Int_t j=0;j<=nGen+1;j++) {
386 std::cout<<
"tau="<<unfold.GetTau()<<
"\n";
387 std::cout<<
"chi**2="<<unfold.GetChi2A()<<
"+"<<unfold.GetChi2L()
388 <<
" / "<<unfold.GetNdf()<<
"\n";
389 std::cout<<
"chi**2(sys)="<<unfold.GetChi2Sys()<<
"\n";
405 TH1 *histMunfold=unfold.GetOutput(
"Unfolded");
408 TH1 *histMdetFold=unfold.GetFoldedOutput(
"FoldedBack");
416 TH2 *histEmatTotal=unfold.GetEmatrixTotal(
"EmatTotal");
419 TH1D *histTotalError=
420 new TH1D(
"TotalError",
";mass(gen)",nGen,xminGen,xmaxGen);
421 for(
Int_t bin=1;bin<=nGen;bin++) {
432 TH1 *histRhoi=unfold.GetRhoItotal(
"rho_I",
446 gFitter->
SetFCN(chisquare_corr);
448 TF1 *bw=
new TF1(
"bw",bw_func,xminGen,xmaxGen,3);
455 histMunfold->
Fit(bw,
"UE");
469 histMdetGenMC->
Draw(
"BOX");
477 histTotalError->
Draw(
"E");
479 histMunfold->
Draw(
"SAME E1");
481 histDensityGenData->
Draw(
"SAME");
482 histDensityGenMC->
Draw(
"SAME HIST");
490 histMdetFold->
Draw();
491 histMdetMC->
Draw(
"SAME HIST");
493 TH1 *histInput=unfold.GetInput(
"Minput",
";mass(det)");
496 histInput->
Draw(
"SAME");
506 bestLogTauLogChi2->
Draw(
"*");
512 bestLcurve->
Draw(
"*");
514 output.SaveAs(
"testUnfold1.ps");
R__EXTERN Int_t gErrorIgnoreLevel
R__EXTERN TStyle * gStyle
virtual void SetLineColor(Color_t lcolor)
Set the line color.
virtual void SetMarkerColor(Color_t mcolor=1)
Set the marker color.
static void RejectPoint(Bool_t reject=kTRUE)
Static function to set the global flag to reject points the fgRejectPoint global flag is tested by al...
virtual Int_t GetNpar() const
virtual void SetNumberFitPoints(Int_t npfits)
virtual Double_t EvalPar(const Double_t *x, const Double_t *params=0)
Evaluate function with given coordinates and parameters.
virtual void InitArgs(const Double_t *x, const Double_t *params)
Initialize parameters addresses.
static Bool_t RejectedPoint()
See TF1::RejectPoint above.
virtual void SetParameter(Int_t param, Double_t value)
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 Double_t GetBinCenter(Int_t bin) const
Return bin center for 1D histogram.
virtual TFitResultPtr Fit(const char *formula, Option_t *option="", Option_t *goption="", Double_t xmin=0, Double_t xmax=0)
Fit histogram with function fname.
virtual Int_t GetNbinsX() const
virtual void SetBinError(Int_t bin, Double_t error)
Set the bin Error Note that this resets the bin eror option to be of Normal Type and for the non-empt...
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 SetBinContent(Int_t bin, Double_t content)
Set bin content see convention for numbering bins in TH1::GetBin In case the bin number is greater th...
virtual void Draw(Option_t *option="")
Draw this histogram with options.
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
virtual Double_t GetBinWidth(Int_t bin) const
Return bin width for 1D histogram.
2-D histogram with a double per channel (see TH1 documentation)}
Service class for 2-Dim histogram classes.
Int_t Fill(Double_t)
Invalid Fill method.
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
virtual void SetBinContent(Int_t bin, Double_t content)
Set bin content.
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
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
An algorithm to unfold distributions from detector to truth level.
@ kSysErrModeRelative
matrix gives the relative shifts
@ kSysErrModeMatrix
matrix is an alternative to the default matrix, the errors are the difference to the original matrix
@ kHistMapOutputVert
truth level on y-axis of the response matrix
Abstract Base Class for Fitting.
virtual TObject * GetObjectFit() const
virtual void SetFCN(void(*fcn)(Int_t &, Double_t *, Double_t &f, Double_t *, Int_t))
To set the address of the minimization objective function called by the native compiler (see function...
virtual TObject * GetUserFunc() const
static TVirtualFitter * Fitter(TObject *obj, Int_t maxpar=25)
Static function returning a pointer to the current fitter.
Double_t Sqrt(Double_t x)
LongDouble_t Power(LongDouble_t x, LongDouble_t y)
static void output(int code)