Fork me on GitHub

Changeset 39 in svn for trunk/routines/resolutions.C


Ignore:
Timestamp:
Nov 18, 2008, 10:30:58 AM (16 years ago)
Author:
severine ovyn
Message:

ok for ETmis

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/routines/resolutions.C

    r27 r39  
    1111void JetResol()
    1212{
    13 
    14    setTDRStyle();
    15    gROOT->Reset();
    16 
    17    TFile *f1 = new TFile("JETResol.root","read");
    18    TTree *Analyze = (TTree*)f1->Get("Analysis");
    19  
    20    const Int_t numBin=7;
    21    double bins[numBin]={0,20,40,60,80,100,220};
    22    TProfile *ETSoverET = new TProfile("ETSoverET","Jet resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
    23 
    24    double rms[numBin];
    25    double mean[numBin];
    26 
    27    TCanvas *c1 = new TCanvas("c1","JET resol",0,0,1000,650);
    28    c1->cd(); int frame=0;
    29    c1->Divide(6,2);
    30 
    31    float x[numBin-1];
    32    float y[numBin-1];
    33    float ex[numBin-1];
    34    float ey[numBin-1];
    35 
    36    float finval=0;//valeur a remplir puis a fitter
    37 
    38    for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
    39       {
    40         TAxis *xaxis = ETSoverET->GetXaxis();
    41         float binCenter = xaxis->GetBinCenter(i+1);
    42         int binMin = (int)xaxis->GetBinLowEdge(i+1);
    43         int binMax = (int)xaxis->GetBinUpEdge(i+1);
    44         cout<<"binMin "<<binMin<<" binMax "<<binMax<<" bin center "<<binCenter<<endl;
    45         char tempMin[500];
    46         if(i==0)binMin=5;
    47         sprintf(tempMin,"JetPTResol.NonSmearePT > %d",binMin);
    48         string mystringMin(tempMin);
    49         char tempMax[500];
    50         sprintf(tempMax,"JetPTResol.NonSmearePT < %d",binMax);
    51         string mystringMax(tempMax);
    52         char tempName[500];
    53         sprintf(tempName,"(JetPTResol.SmearePT)>>hETSoverET%d",i);
    54         string mystringName(tempName);
    55 
    56         c1->cd(++frame);
    57         GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
    58         x[i]=binCenter;
    59         finval=rms[i]/mean[i];
    60         y[i]=(finval*100);
    61         ex[i]=0;
    62         ey[i]=0;
    63         cout<<"finval "<<finval<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;
    64       }
    65 
    66    TCanvas *c2 = new TCanvas("c2","JET resol",100,100,600,450);
    67    c2->cd();
    68 
    69    TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,170);
    70 
    71    TGraphErrors *gr11 = new TGraphErrors((numBin-1),x,y,ex,ey);
    72    gr11->Draw("AP");
    73    gr11->SetTitle("");
    74    gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
    75    gr11->GetYaxis()->SetRangeUser(0,50);
    76    gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}/E_{T}^{MC})_{fit}/<E_{T}^{rec}/E_{T}^{MC}>_{fit}");
    77 
    78    Double_t* params = fitfun->GetParameters();
    79 
    80    gr11->Fit("user","R");
    81    gr11->Fit("user","RI");
    82    gr11->Fit("user","RI");
    83    gr11->Fit("user","RI");
    84    char tempResol[500];
    85    sprintf(tempResol,"#frac{#sigma(E_{T}^{rec}/E_{T}^{MC})}{<E_{T}^{rec}/E_{T}^{MC}>} = #frac{%f}{#sqrt{E_{T}^{MC}}} #oplus #frac{%f}{E_{T}^{MC}} #oplus %f",params[1]/100,params[2]/100,params[0]/100);
    86 
    87    TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
    88    leg1->Draw();
    89 
    90    delete fitfun;
    91 }
    92 
     13  setTDRStyle();
     14  gROOT->Reset();
     15 
     16  TFile *f1 = new TFile("JET.root","read");
     17  TTree *Analyze = (TTree*)f1->Get("Analysis");
     18 
     19  const Int_t numBin=7;
     20  double bins[numBin]={0,20,40,60,80,100,220};
     21  TProfile *ETSoverET = new TProfile("ETSoverET","Jet resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
     22 
     23  double rms[numBin];
     24  double mean[numBin];
     25 
     26  TCanvas *c1 = new TCanvas("c1","JET resol",0,0,1000,650);
     27  c1->cd(); int frame=0;
     28  c1->Divide(6,2);
     29 
     30  float x[numBin-1];
     31  float y[numBin-1];
     32  float ex[numBin-1];
     33  float ey[numBin-1];
     34 
     35  float finval=0;//valeur a remplir puis a fitter
     36 
     37  for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
     38    {
     39      TAxis *xaxis = ETSoverET->GetXaxis();
     40      float binCenter = xaxis->GetBinCenter(i+1);
     41      int binMin = (int)xaxis->GetBinLowEdge(i+1);
     42      int binMax = (int)xaxis->GetBinUpEdge(i+1);
     43      char tempMin[500];
     44      if(i==0)binMin=5;
     45      sprintf(tempMin,"JetPTResol.NonSmearePT > %d",binMin);
     46      string mystringMin(tempMin);
     47      char tempMax[500];
     48      sprintf(tempMax,"JetPTResol.NonSmearePT < %d",binMax);
     49      string mystringMax(tempMax);
     50      char tempName[500];
     51      sprintf(tempName,"(JetPTResol.SmearePT)>>hETSoverET%d",i);
     52      string mystringName(tempName);
     53     
     54      c1->cd(++frame);
     55      GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
     56      x[i]=binCenter;
     57      finval=rms[i]/mean[i];
     58      y[i]=(finval*100);
     59      ex[i]=0;
     60      ey[i]=0;
     61    }
     62 
     63  TCanvas *c2 = new TCanvas("c2","JET resol",100,100,600,450);
     64  c2->cd();
     65 
     66  TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,170);
     67 
     68  TGraphErrors *gr11 = new TGraphErrors((numBin-1),x,y,ex,ey);
     69  gr11->Draw("AP");
     70  gr11->SetTitle("");
     71  gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
     72  gr11->GetYaxis()->SetRangeUser(0,50);
     73  gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}/E_{T}^{MC})_{fit}/<E_{T}^{rec}/E_{T}^{MC}>_{fit}");
     74 
     75  Double_t* params = fitfun->GetParameters();
     76 
     77  gr11->Fit("user","QR");
     78  gr11->Fit("user","QRI");
     79  gr11->Fit("user","QRI");
     80  gr11->Fit("user","QRI");
     81  char tempResol[500];
     82  if(params[2]==0.000000)sprintf(tempResol,"#frac{#sigma(E_{T}^{rec}/E_{T}^{MC})}{<E_{T}^{rec}/E_{T}^{MC}>} = #frac{%f}{#sqrt{E_{T}^{MC}}} #oplus %f",params[1]/100,params[0]/100);
     83  else sprintf(tempResol,"#frac{#sigma(E_{T}^{rec}/E_{T}^{MC})}{<E_{T}^{rec}/E_{T}^{MC}>} = #frac{%f}{#sqrt{E_{T}^{MC}}} #oplus #frac{%f}{E_{T}^{MC}} #oplus %f",params[1]/100,params[2]/100,params[0]/100);
     84 
     85  TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
     86  leg1->Draw();
     87
     88  TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
     89  Delphes->Draw();
     90 
     91 
     92  delete fitfun;
     93}
     94
     95void ElecResol()
     96{
     97 
     98  setTDRStyle();
     99  gROOT->Reset();
     100 
     101  TFile *f1 = new TFile("PTMIS.root","read");
     102  TTree *Analyze = (TTree*)f1->Get("Analysis");
     103 
     104  const Int_t numBin=9;
     105  double bins[numBin]={0,20,40,60,80,100,150,220,400};
     106  TProfile *ETSminusET = new TProfile("ETSminusET","Electron resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
     107 
     108  double rms[numBin];
     109  double mean[numBin];
     110 
     111  TCanvas *c3 = new TCanvas("c3","ELEC resol",0,0,1000,650);
     112  c3->cd(); int frame=0;
     113  c3->Divide(6,2);
     114 
     115  float x[numBin-1];
     116  float y[numBin-1];
     117  float ex[numBin-1];
     118  float ey[numBin-1];
     119 
     120  float finval=0;//valeur a remplir puis a fitter
     121 
     122  for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
     123    {
     124      TAxis *xaxis = ETSminusET->GetXaxis();
     125      float binCenter = xaxis->GetBinCenter(i+1);
     126      int binMin = (int)xaxis->GetBinLowEdge(i+1);
     127      int binMax = (int)xaxis->GetBinUpEdge(i+1);
     128      char tempMin[500];
     129      if(i==0)binMin=5;
     130      sprintf(tempMin,"ElecEResol.NonSmeareE > %d",binMin);
     131      string mystringMin(tempMin);
     132      char tempMax[500];
     133      sprintf(tempMax,"ElecEResol.NonSmeareE < %d",binMax);
     134      string mystringMax(tempMax);
     135      char tempName[500];
     136      sprintf(tempName,"(ElecEResol.NonSmeareE-ElecEResol.SmeareE)>>hETSoverET%d",i);
     137      string mystringName(tempName);
     138     
     139      c3->cd(++frame);
     140      GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
     141      x[i]=binCenter;
     142      finval=rms[i]/binCenter;
     143      y[i]=(finval*100);
     144      ex[i]=0;
     145      ey[i]=0;
     146    }
     147 
     148  TCanvas *c4 = new TCanvas("c4","ELEC resol",100,100,600,450);
     149  c4->cd();
     150 
     151  TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,400);
     152 
     153  TGraphErrors *gr11 = new TGraphErrors((numBin-1),x,y,ex,ey);
     154  gr11->Draw("AP");
     155  gr11->SetTitle("");
     156  gr11->GetXaxis()->SetTitle("E [GeV]");
     157  gr11->GetYaxis()->SetRangeUser(0,5);
     158  gr11->GetYaxis()->SetTitle("#sigma/E");
     159 
     160  Double_t* params = fitfun->GetParameters();
     161 
     162  gr11->Fit("user","QR");
     163  gr11->Fit("user","QRI");
     164  gr11->Fit("user","QRI");
     165  gr11->Fit("user","QRI");
     166  char tempResol[500];
     167  sprintf(tempResol,"#frac{#sigma}{E} = #frac{%f}{#sqrt{E}} #oplus #frac{%f}{E} #oplus %f",params[1]/100,params[2]/100,params[0]/100);
     168 
     169  TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
     170  leg1->Draw();
     171 
     172  TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
     173  Delphes->Draw();
     174 
     175  TPaveText *events = MakeTPave(0.2,0.75,0.35,0.8,"t#bar{t} #rightarrow l^{+}l^{-}#nu#bar{#nu}b#bar{b}");
     176  events->Draw();
     177 
     178  delete fitfun;
     179}
     180
     181void TauJetInfo()
     182{
     183 
     184  setTDRStyle();
     185  gROOT->Reset();
     186 
     187  TFile *f1 = new TFile("PTMIS.root","read");
     188  TTree *Analyze = (TTree*)f1->Get("Analysis");
     189 
     190  TCanvas *ct = new TCanvas("ct","Tau information",0,0,1000,425);
     191  ct->cd(); int frame=0;
     192  ct->Divide(2,1);
     193 
     194  ct->cd(++frame);
     195  TH1F *tauEnergy =MakeNormTH1F(20,0.8,1,Analyze,"TauJetPTResol.EnergieCen>>tauEnergy",1, 0, 1,2,false);
     196  tauEnergy->Draw();
     197  tauEnergy->SetTitle("");
     198  tauEnergy->GetYaxis()->SetTitle("Fraction of events");
     199  tauEnergy->GetXaxis()->SetTitle("C_{#tau}");
     200 
     201 
     202  TPaveText *Delphes1 = MakeTPave(0.3,0.85,0.45,0.9,"MG/ME + Delphes");
     203  Delphes1->Draw();
     204  TPaveText *ctau = MakeTPave(0.3,0.75,0.45,0.8,"C_{#tau} = #frac{#sum E^{towers} (#DeltaR = 0.15)}{E^{jet}}");
     205  ctau->Draw();
     206  TPaveText *events1 = MakeTPave(0.3,0.65,0.45,0.7,"Events: t#bar{t} #rightarrow l^{+}#nul^{-}#bar{#nu}b#bar{b}");
     207  events1->Draw();
     208 
     209 
     210  ct->cd(++frame);
     211  TH1F *NumTrack =MakeNormTH1F(6,0,6,Analyze,"TauJetPTResol.NumTrack>>NumTrack",1, 0, 1,2,false);
     212  NumTrack->Draw();
     213  NumTrack->SetTitle("");
     214  NumTrack->GetYaxis()->SetTitle("Fraction of events");
     215  NumTrack->GetXaxis()->SetTitle("N^{tracks}");
     216 
     217  TPaveText *Delphes = MakeTPave(0.6,0.85,0.85,0.9,"MG/ME + Delphes");
     218  Delphes->Draw();
     219  TPaveText *numtracks = MakeTPave(0.6,0.75,0.85,0.8,"#DeltaR < 0.4, p_{T}^{track} > 2 GeV");
     220  numtracks->Draw();
     221  TPaveText *events = MakeTPave(0.6,0.65,0.85,0.7,"Events: t#bar{t} #rightarrow l^{+}#nul^{-}#bar{#nu}b#bar{b}");
     222  events->Draw();
     223   
     224}
    93225
    94226void ETmisResol()
    95227{
    96 
    97    setTDRStyle();
    98    gROOT->Reset();
    99 
    100    TFile *f1 = new TFile("PTMIS.root","read");
    101    TTree *Analyze = (TTree*)f1->Get("Analysis");
    102 
    103    const Int_t numBin=6;
    104    double bins[numBin]={10,20,30,50,80,100};
    105    TProfile *ETSoverET = new TProfile("ETSoverET","ETmis resol",(numBin-1),bins,-1000,1000);
    106 
    107    double rms[numBin-1];
    108    double mean[numBin-1];
    109 
    110    TCanvas *c5 = new TCanvas("c5","gerrors2",0,0,1000,650);
    111    c5->cd(); int frame=0;
    112    c5->Divide(6,2);
    113 
    114    double x[numBin-1];
    115    double y[numBin-1];
    116 
    117    for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
    118       {
    119         TAxis *xaxis = ETSoverET->GetXaxis();
    120         float binCenter = xaxis->GetBinCenter(i+1);
    121         int binMin = (int)xaxis->GetBinLowEdge(i+1);
    122         int binMax = (int)xaxis->GetBinUpEdge(i+1);
    123         cout<<"binMin "<<binMin<<" binMax "<<binMax<<" bin center "<<binCenter<<endl;
    124         char tempMin[500];
    125         sprintf(tempMin,"ETmisResol.Et>%d",binMin);
    126         string mystringMin(tempMin);
    127         char tempMax[500];
    128         sprintf(tempMax,"ETmisResol.Et<%d",binMax);
    129         string mystringMax(tempMax);
    130 
    131         char tempName[500];
    132         sprintf(tempName,"ETmisResol.EtSmeare>>hETSoverET%d",i);
    133         string mystringName(tempName);
    134 
    135         c5->cd(++frame);
    136         GaussValuesETmis(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
    137         x[i]=binCenter;
    138         y[i]=rms[i];
    139         cout<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;
    140 
    141       }
    142 
    143    TCanvas *c6 = new TCanvas("c6","ETmis resolution",100,100,600,450);
    144    c6->cd();
    145 
    146    TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",0,200);
    147 
    148    TGraph *gr11 = new TGraph((numBin-1),x,rms);
    149    gr11->Draw("AP");
    150    gr11->GetXaxis()->SetTitle("MC^{MET} [GeV]");
    151    gr11->GetYaxis()->SetTitle("#sigma(REC^{MET}-MC^{MET})");
    152    gr11->Fit("user","RQ");
    153    gr11->Fit("user","RQ");
    154    gr11->Fit("user","RQ");
    155    gr11->Fit("user","R");
    156 
    157    TCanvas *c7 = new TCanvas("c7","ETmis resolution",100,100,600,450);
    158    c7->cd();
    159 
    160    TGraph *gr12 = new TGraph((numBin-1),x,mean);
    161    gr12->Draw("AP");
    162    gr12->GetXaxis()->SetTitle("MC^{MET} [GeV]");
    163    gr12->GetYaxis()->SetTitle("<(REC^{MET}-MC^{MET})>");
    164 
    165 
    166  delete fitfun;
     228 
     229  setTDRStyle();
     230  gROOT->Reset();
     231 
     232  TFile *f1 = new TFile("test2.root","read");
     233  TTree *Analyze = (TTree*)f1->Get("Analysis");
     234 
     235  TF1 *fitfun = new TF1("user","[0]*sqrt(x)",0,600);
     236  const Int_t numBin=7;
     237//  double bins[numBin]={10,100,150,200,250,300,350,400,450,500};
     238  double bins[numBin]={200,250,300,350,400,450,500};
     239//  double bins[numBin]={0,10,20,30,40,50,60,70};
     240  TProfile *ETSoverET = new TProfile("ETSoverET","ETmis resol",(numBin-1),bins,-1000,1000);
     241  Analyze->Draw("ETmisResol.ExSmeare:ETmisResol.SEt>>ETSoverET");
     242  Analyze->Fit("user","RQ");
     243  Analyze->Fit("user","RQ");
     244 
     245 
     246  double rms[numBin-1];
     247 
     248  TCanvas *c5 = new TCanvas("c5","PTmis resol",0,0,1000,650);
     249  c5->cd(); int frame=0;
     250  c5->Divide(6,2);
     251 
     252  double x[numBin];
     253  double y[numBin];
     254 
     255  for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
     256    {
     257      TAxis *xaxis = ETSoverET->GetXaxis();
     258      float binCenter = xaxis->GetBinCenter(i+1);
     259      int binMin = (int)xaxis->GetBinLowEdge(i+1);
     260      int binMax = (int)xaxis->GetBinUpEdge(i+1);
     261      char tempMin[500];
     262      sprintf(tempMin,"ETmisResol.SEt>%d",binMin);
     263      string mystringMin(tempMin);
     264      char tempMax[500];
     265      sprintf(tempMax,"ETmisResol.SEt<%d",binMax);
     266      string mystringMax(tempMax);
     267     
     268      char tempName[500];
     269      sprintf(tempName,"ETmisResol.ExSmeare>>hETSoverET%d",i);
     270      string mystringName(tempName);
     271     
     272      c5->cd(++frame);
     273      GaussValuesETmis(Analyze,tempName,rms[i],mystringMin,mystringMax);
     274      x[i]=binCenter;
     275      y[i]=rms[i];
     276     
     277    }
     278 
     279  TCanvas *c6 = new TCanvas("c6","ETmis resolution",100,100,600,450);
     280  c6->cd();
     281 
     282  x[numBin]=0;
     283  y[numBin]=0;
     284  TGraph *gr11 = new TGraph((numBin),x,rms);
     285  gr11->Draw("AP");
     286  gr11->GetXaxis()->SetTitle("Offline sum of E_{T} [GeV]");
     287  gr11->GetYaxis()->SetTitle("Resolution of x-component  of MET [GeV]");
     288  gr11->Fit("user","RQ");
     289  gr11->Fit("user","RQ");
     290  gr11->Fit("user","RQ");
     291  gr11->Fit("user","RQ");
     292  gr11->GetYaxis()->SetRangeUser(0,30);
     293  gr11->GetXaxis()->SetRangeUser(1,600);
     294
     295  Double_t* params = fitfun->GetParameters();
     296 
     297  char tempResol[500];
     298  sprintf(tempResol,"%f * #sqrt{E_{T}}",params[0]);
     299
     300  TPaveText *leg1 = MakeTPave(0.4,0.6,0.8,0.65,tempResol);
     301  leg1->Draw();
     302
     303  TPaveText *leg2 = MakeTPave(0.2,0.8,0.8,0.85,"WHq'#rightarrow W#tau#tauq'#rightarrowjjl#tauq', m_{H}=150 GeV ");
     304  leg2->Draw();
     305
     306  TPaveText *Delphes = MakeTPave(0.2,0.15,0.35,0.2,"MG/ME + Delphes");
     307  Delphes->Draw();
     308 
     309 
     310  delete fitfun;
    167311}
    168312
    169313void General()
    170314{
    171   JetResol();
     315  //JetResol();
     316  //ElecResol();
    172317  ETmisResol();
    173 }
    174 
    175 
     318  //TauJetInfo();
     319 
     320}
     321
     322
Note: See TracChangeset for help on using the changeset viewer.