source: svn/trunk/routines/resolutions.C@ 25

#include "TROOT.h"
#include "TFile.h"
#include "TTree.h"
#include "TCanvas.h"
#include "TProfile.h"
#include "TF1.h"
#include "TGraph.h"

#include "interface/FuncDef.h"

void JetResol(TTree *Analyze)
{
   TCanvas *c1 = new TCanvas("c1","JET resol",100,100,600,450);
   c1->cd(); 
  
   const Int_t numBin=7;
   double bins[numBin]={0,20,60,100,220,420,860};
   TProfile *ETSoverET = new TProfile("ETSoverET","Jet resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
   Analyze->Draw("JetPTResol.SmearePT:JetPTResol.NonSmearePT>>ETSoverET","JetPTResol.NonSmearePT>1");

   double rms[numBin];
   double mean[numBin];

   TCanvas *c2 = new TCanvas("c2","JET resol",0,0,1000,650);
   c2->cd(); int frame=0;
   c2->Divide(6,2);

   float x[numBin-1];
   float y[numBin-1];
   float finval=0;//valeur a remplir puis a fitter

   for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
      {
        TAxis *xaxis = ETSoverET->GetXaxis();
        float binCenter = xaxis->GetBinCenter(i+1);
        int binMin = (int)xaxis->GetBinLowEdge(i+1);
        int binMax = (int)xaxis->GetBinUpEdge(i+1);
        cout<<"binMin "<<binMin<<" binMax "<<binMax<<" bin center "<<binCenter<<endl;
        char tempMin[500];
        sprintf(tempMin,"JetPTResol.NonSmearePT>%d",binMin);
        string mystringMin(tempMin);
        char tempMax[500];
        sprintf(tempMax,"JetPTResol.NonSmearePT<%d",binMax);
        string mystringMax(tempMax);
        char tempName[500];
        sprintf(tempName,"JetPTResol.SmearePT>>hETSoverET%d",i);
        string mystringName(tempName);

        c2->cd(++frame);
        GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
        x[i]=binCenter;
        finval=rms[i]/mean[i];
        y[i]=(finval*2);
        cout<<"finval "<<finval<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;

      }

   TCanvas *c3 = new TCanvas("c3","JET resol",100,100,600,450);
   c3->cd();

   TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,3000);

   TGraph *gr11 = new TGraph((numBin-1),x,y);
   gr11->Draw("AP");
   gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
   gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}/E_{T}^{MC})/<E_{T}^{rec}/E_{T}^{MC}>");

   //fitfun->FixParameter(0,0.0541930);
   //fitfun->FixParameter(1,0.634908);
   gr11->Fit("user","RQ");
   gr11->Fit("user","RQ");
   gr11->Fit("user","RQ");
   gr11->Fit("user","R");



 delete fitfun;
// delete Gauss;
}


void ETmisResol(TTree *Analyze)
{

   TCanvas *c4 = new TCanvas("c4","ETmis resol",100,100,600,450);
   c4->cd(); 
  
   const Int_t numBin=10;
   double bins[numBin]={0,10,20,40,80,160,320,640,1080,2160};
   TProfile *ETSoverET = new TProfile("ETSoverET","ETmis resol",(numBin-1),bins,-10,10);
   Analyze->Draw("(ETmisResol.EtSmeare-ETmisResol.Et):(ETmisResol.Et)>>ETSoverET","(ETmisResol.Et)>1");

   double rms[numBin];
   double mean[numBin];

   TCanvas *c5 = new TCanvas("c5","gerrors2",0,0,1000,650);
   c5->cd(); int frame=0;
   c5->Divide(6,2);

   float x[numBin-1];
   float y[numBin-1];

   for ( int i=0; i<(numBin-1); i++)  // premiÚre bin : i ==1 et pas i == 0
      {
        TAxis *xaxis = ETSoverET->GetXaxis();
        float binCenter = xaxis->GetBinCenter(i+1);
        int binMin = (int)xaxis->GetBinLowEdge(i+1);
        int binMax = (int)xaxis->GetBinUpEdge(i+1);
        cout<<"binMin "<<binMin<<" binMax "<<binMax<<" bin center "<<binCenter<<endl;
        char tempMin[500];
        sprintf(tempMin,"(ETmisResol.Et)>%d",binMin);
        string mystringMin(tempMin);
        char tempMax[500];
        sprintf(tempMax,"(ETmisResol.Et)<%d",binMax);
        string mystringMax(tempMax);
        char tempName[500];
        sprintf(tempName,"(ETmisResol.EtSmeare-ETmisResol.Et)>>hETSoverET%d",i);
        string mystringName(tempName);

        c5->cd(++frame);
        GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
        x[i]=binCenter;
        mean[i]=ETSoverET->GetBinContent(i+1);
        y[i]=rms[i];
        cout<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;

      }

   TCanvas *c6 = new TCanvas("c6","ETmis resolution",100,100,600,450);
   c6->cd();

   TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,3000);

   TGraph *gr11 = new TGraph((numBin-1),x,y);
   gr11->Draw("AP");
   gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
   gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}>");

   //fitfun->FixParameter(0,0.0541930);
   //fitfun->FixParameter(1,0.634908);
   gr11->Fit("user","RQ");
   gr11->Fit("user","RQ");
   gr11->Fit("user","RQ");
   gr11->Fit("user","R");



 delete fitfun;
// delete Gauss;
}

void General()
{
  setTDRStyle();
  gROOT->Reset();

  TFile *f1 = new TFile("Smearing.root","read");
  TTree *Analyze = (TTree*)f1->Get("Analysis");
  JetResol(Analyze);

}