DoubleCrystalBall.cc File Reference

#include <cstdlib>
#include "Core/JEC/interface/resolution.h"
#include "Core/JEC/interface/DoubleCrystalBall.h"
#include <TROOT.h>
#include <TCanvas.h>
#include <TLine.h>
#include <TLegend.h>
#include <TF1.h>

Include dependency graph for DoubleCrystalBall.cc:

Functions
int	main (int argc, char *argv[])

Function Documentation

main()

int main	(	int	argc,
		char *	argv[]
	)

{
    gROOT->SetBatch();
    auto c = new TCanvas("DoubleCrystalBall", "DoubleCrystalBall", 600, 800);
 
    double x1 = -0.8, 
           x2 = 0.8,
           y1 = 1e-6,
           y2 = 100;
    c->DrawFrame(x1, y1, x2, y2, ";#Delta;");
    c->SetTicks(1,1);
    c->SetLogy();
 
    double N = 1,
           mu = 0.,
           sigma = 0.1,
           kL = -0.05,
           nL = 100,
           kR = 0.3,
           nR = 3;
 
    double p[] = { N, mu, sigma, kL, nL, kR, nR };
 
    TF1 * f = new TF1("DoubleCrystalBall", DAS::DoubleCrystalBall::Distribution(), x1, x2, 7);
    f->SetParameters(p);
    f->Draw("same");
    cout << "Integral of CB = " << f->Integral(x1, x2) << endl;
 
    TF1 * g = new TF1("gaus", "gaus(0)", x1, x2);
    double q[] = {N/(sigma*sqrt(2*M_PI)), mu, sigma};
    g->SetParameters(q);
    g->SetLineColor(kBlue);
    g->SetLineStyle(2);
    g->SetLineWidth(1);
    g->Draw("same");
    cout << "Integral of Gaussian = " << g->Integral(x1, x2) << endl;
 
    TLegend * legend = new TLegend(0.6, 0.7, 0.89, 0.85);
    legend->SetBorderSize(0);
    legend->AddEntry(f, "Double Crystal-Ball", "l");
    legend->AddEntry(g, "Gaussian with same core parameters", "l");
 
    auto line = [&](TF1 * f, Style_t s, double x, TString title = "") {
        TLine * line = new TLine;
        line->SetLineStyle(s);
        line->DrawLine(x, y1, x, f->Eval(x));
        if (title != "") legend->AddEntry(line, title, "l");
    };
    line(f, 1, mu, Form("mean (#mu = %.2f)", mu));
    line(f, 2, mu+sigma, Form("mean #pm width (#sigma = %.2f)", sigma));
    line(f, 2, mu-sigma);
    line(f, 3, kL, Form("k_{L} = %.2f, n_{L} = %.0f", kL, nL));
    line(f, 3, kR, Form("k_{R} = %.2f, n_{R} = %.0f", kR, nR));
 
    legend->Draw();
    c->RedrawAxis();
    c->Print("DoubleCrystalBall.pdf(");
 
 
    auto c2 = new TCanvas("IntDoubleCrystalBall", "IntDoubleCrystalBall", 600, 800);
 
    y1 = -0.5;//1e-2; //1e-10;
    y2 = 1.5;//10*N/sigma; //20;
    c2->DrawFrame(x1, y1, x2, y2, ";#Delta;");
    c2->SetTicks(1,1);
    //c2->SetLogy();
 
    TF1 * f2 = new TF1("IntDoubleCrystalBall", DAS::DoubleCrystalBall::Integral(), x1, x2, 7);
    f2->SetParameters(p);
    f2->Draw("same");
 
    auto myErf = [](double * x, double *p) {
        double N = p[0], mu = p[1], sigma = p[2];
        double z = (*x-mu)/sigma;
        return N*(1+erf(z/sqrt(2)))/2;
    };
 
    TF1 * g2 = new TF1("myErf", myErf, x1, x2, 3);
    g2->SetParameters(p);
    g2->SetLineColor(kBlue);
    g2->SetLineStyle(2);
    g2->SetLineWidth(1);
    g2->Draw("same");
 
    legend = new TLegend(0.2, 0.7, 0.49, 0.85);
    legend->SetBorderSize(0);
    legend->AddEntry(f, "Integrated DCB", "l");
    legend->AddEntry(g, "Error function", "l");
 
    line(f2, 1, mu, Form("mean (#mu = %.2f)", mu));
    line(f2, 2, mu+sigma, Form("mean #pm width (#sigma = %.2f)", sigma));
    line(f2, 2, mu-sigma);
    line(f2, 3, kL, Form("k_{L} = %.2f, n_{L} = %.0f", kL, nL));
    line(f2, 3, kR, Form("k_{R} = %.2f, n_{R} = %.0f", kR, nR));
 
    legend->Draw();
    c2->RedrawAxis();
    c2->Print("DoubleCrystalBall.pdf");
 
 
    auto c3 = new TCanvas("DerivLogDoubleCrystalBall", "DerivLogDoubleCrystalBall", 600, 800);
 
    y1 = -50;
    y2 = 50;
    c3->DrawFrame(x1, y1, x2, y2, ";#Delta;");
    c3->SetTicks(1,1);
 
    TF1 * f3 = new TF1("DerivLogDoubleCrystalBall", DAS::DoubleCrystalBall::DLog(), x1, x2, 7);
    f3->SetParameters(p);
    f3->Draw("same");
 
    TF1 * g3 = new TF1("derivLogGauss", "-(x-[0])/([1]*[1])", x1, x2);
    g3->SetParameters(&p[1]);
    g3->SetLineColor(kBlue);
    g3->SetLineStyle(2);
    g3->SetLineWidth(1);
    g3->Draw("same");
 
    legend = new TLegend(0.6, 0.7, 0.89, 0.85);
    legend->SetBorderSize(0);
    legend->AddEntry(f, "Derivative of Log of DCB", "l");
    legend->AddEntry(g, "Same for pure Gaussian", "l");
 
    line(f3, 1, mu, Form("mean (#mu = %.2f)", mu));
    line(f3, 2, mu+sigma, Form("mean #pm width (#sigma = %.2f)", sigma));
    line(f3, 2, mu-sigma);
    line(f3, 3, kL, Form("k_{L} = %.2f, n_{L} = %.0f", kL, nL));
    line(f3, 3, kR, Form("k_{R} = %.2f, n_{R} = %.0f", kR, nR));
 
    legend->Draw();
 
    c3->RedrawAxis();
    c3->Print("DoubleCrystalBall.pdf)");
 
    return EXIT_SUCCESS;
}