test/macro/FitEta.C

0001 #include "CrystalBall.C"
0002 #include "effSigma.C"
0003
0004 #include "TChain.h"
0005
0006 TChain* chain = new TChain("energyScale");
0007 void load_noET() {
0008   chain->Add("../SuperClusters.root");
0009   //Need to be added more files if there are multiple output root files
0010 }
0011
0012 double fit0, fit1, fit2, fit3, fit4;
0013
0014 void plotEnergy(TString var, TString Eta, int etaBin, bool newE) {
0015
0016   float mean_value = 0;
0017   float mean_error = 0;
0018   float sigma_value = 0;
0019   float ChiSquare = 0;
0020
0021   int nBin = 1000;
0022   float xMin = 0.6;
0023   float xMax = 1.1;
0024
0025   gStyle->SetOptFit(1);
0026   gStyle->SetStatX(0.48);
0027   gStyle->SetStatY(0.82);
0028   gStyle->SetStatW(0.2);
0029   gStyle->SetStatH(0.15);
0030
0031   TString Cut = Eta + "&&em_scType==1";
0032   Cut = Cut + "&&em_et>10&&em_et<300&&!em_isInCrack";
0033   TH1F *h = new TH1F("h"," Et^{RECO}/Et^{GEN} GeV",nBin, xMin, xMax);
0034   chain->Draw(var + "/mc_et>>h",Cut);
0035
0036   if ( newE ) {
0037     fit0 = 1800;
0038     fit1 = 0.9;
0039     fit2 = 5;
0040     fit3 = 1.000;
0041     fit4 = 0.006;
0042   }
0043
0044   //TF1* g1 = new TF1(cb,"","",Emin, Emax);
0045   TF1* g1 = new TF1("g1",CrystalBall,0.85, 1.2, 5);
0046   g1->SetParameters(fit0, fit1, fit2, fit3, fit4);
0047   g1->SetLineColor(4);
0048   g1->SetLineWidth(3);
0049   float max = 1.02;
0050   if ( etaBin < 3 )
0051     max = 1.02;
0052   else if ( etaBin < 5 )
0053     max = 1.03;
0054   else
0055     max = 1.05;
0056
0057   h->Fit(g1,"","",0.6,max);
0058   cout << endl;
0059   mean_value  = g1->GetParameter(3);
0060   mean_error  = g1->GetParError(3);
0061   sigma_value = g1->GetParameter(4);
0062   ChiSquare = g1->GetChisquare();
0063
0064   fit0 = g1->GetParameter(0);
0065   fit1 = g1->GetParameter(1);
0066   fit2 = g1->GetParameter(2);
0067   fit3 = g1->GetParameter(3);
0068   fit4 = g1->GetParameter(4);
0069
0070   cout << "Result:\t" << var << "\t" << Eta << "\t" << mean_value << "\t" << mean_error << "\t" << effSigma(h) << endl;
0071
0072   delete g1;
0073   delete h;
0074 }
0075
0076 void derivationEnergyCorrection_phiWidth() {
0077
0078   // MODULE 1: 0.018 < ETA < 0.423: 0.018, 0.119, 0.220, 0.320, 0.423
0079   // MODULE 2: 0.461 < ETA < 0.770: 0.461, 0.538, 0.616, 0.693, 0.770
0080   // MODULE 3: 0.806 < ETA < 1.127: 0.806, 0.886, 0.967, 1.047, 1.127
0081   // MODULE 4: 1.163 < ETA < 1.460: 1.163, 1.237, 1.312, 1.386, 1.460
0082
0083   TString leftEta[24] = {"0.018", "0.0855", "0.153", "0.2205", "0.288", "0.3555",
0084              "0.461", "0.5125", "0.564", "0.6155", "0.667", "0.7185",
0085              "0.806", "0.8595", "0.913", "0.9665", "1.02",  "1.0735",
0086              "1.163", "1.2125", "1.262", "1.3115", "1.361", "1.4105"};
0087
0088   TString rightEta[24] = {"0.0855", "0.153", "0.2205", "0.288", "0.3555", "0.423",
0089               "0.5125", "0.564", "0.6155", "0.667", "0.7185", "0.770",
0090               "0.8595", "0.913",  "0.9665", "1.02", "1.0735", "1.127",
0091               "1.2125", "1.262",  "1.3115","1.361","1.4105",  "1.46"};
0092
0093
0094   TString vars[2] = {"emCorr_et", "em_et"};
0095
0096   for (int i=0; i < 2; ++i) {
0097     bool newE = true;
0098     for (int j=0; j < 24; ++j) {
0099       TString etaCut = "abs(em_eta)>" + leftEta[j] + "&&abs(em_eta)<" + rightEta[j];
0100       plotEnergy(vars[i], etaCut, j, newE);
0101       newE = false;
0102     }
0103   }
0104 }
0105
0106 void FitEta(){
0107   load_noET();
0108   derivationEnergyCorrection_phiWidth();
0109   gROOT->ProcessLine(".q");
0110 }
0111
0112 TString IntToString(int number){
0113   ostringstream oss;
0114   oss << number;
0115   return oss.str();
0116 }
0117
0118 TString FloatToString(float number){
0119   ostringstream oss;
0120   oss << number;
0121   return oss.str();
0122 }