RecoLocalTracker/SiStripRecHitConverter/test/Resolutionwc.C - master - githubmirror/cmssw

RecoLocalTracker/SiStripRecHitConverter/test/Resolutionwc.C
160 строк6 KB
Giulio Eulisse
Snapshot of CMSSW_6_2_0_pre8.
26 июн 2013, 18:12
26 июн 2013, 18:12214ab73
// Double_t fit_func(Double_t *x, Double_t *par){
//   Float_t xx =x[0];
//   if(xx<par[0])f = par[1]*x+par[2];
//   else f=par[3]*x+par[0]*(par[1]-par[3])+par[2];
//   return f;
// }
{
  int NBINS=32;
    TFile *resfile=new TFile("resolution_2.5.root");
 // TFile *resfile=new TFile("resolution_good.root");
  float tibrms[NBINS],tobrms[NBINS],tidrms[NBINS],tecrms[NBINS];
  float tibsigma[NBINS],tobsigma[NBINS],tidsigma[NBINS],tecsigma[NBINS];
  float tibrmserror[NBINS],tobrmserror[NBINS],tidrmserror[NBINS],tecrmserror[NBINS];
  float tibsigmaerror[NBINS],tobsigmaerror[NBINS],tidsigmaerror[NBINS],tecsigmaerror[NBINS];
  TH1F *tibres[NBINS],*tobres[NBINS],*tidres[NBINS],*tecres[NBINS];
  TH1F *allres[NBINS];
  float allrms[NBINS],allrmserror[NBINS];
  float allsigma[NBINS],allsigmaerror[NBINS];
  for(int i=0;i<NBINS;i++){
    TH1F *tibres[i]=(TH1F*)resfile->Get(Form("TIBreswc_%d_strip",i));
    TH1F *tobres[i]=(TH1F*)resfile->Get(Form("TOBreswc_%d_strip",i));
    TH1F *tidres[i]=(TH1F*)resfile->Get(Form("TIDreswc_%d_strip",i));
    TH1F *tecres[i]=(TH1F*)resfile->Get(Form("TECreswc_%d_strip",i));
    allres[i]=(TH1F*)tibres[i]->Clone(Form("Allres_%d_strip",i));
    allres[i]->Add(tobres[i]);
    allres[i]->Add(tidres[i]);
    allres[i]->Add(tecres[i]);
//     tibres[i]->Fit("gaus");
//     tobres[i]->Fit("gaus");
//     tidres[i]->Fit("gaus");
//     tecres[i]->Fit("gaus");
//     allres[i]->Fit("gaus");
//     tibrms[i]=tibres[i]->GetFunction("gaus")->GetParameter(2);
//     tibrmserror[i]=tibres[i]->GetFunction("gaus")->GetParError(2);
//     tobrms[i]=tobres[i]->GetFunction("gaus")->GetParameter(2);
//     tobrmserror[i]=tobres[i]->GetFunction("gaus")->GetParError(2);
//     tidrms[i]=tidres[i]->GetFunction("gaus")->GetParameter(2);
//     tidrmserror[i]=tidres[i]->GetFunction("gaus")->GetParError(2);
//     tecrms[i]=tecres[i]->GetFunction("gaus")->GetParameter(2);
//     tecrmserror[i]=tecres[i]->GetFunction("gaus")->GetParError(2);
//     allrms[i]=allres[i]->GetFunction("gaus")->GetParameter(2);
//     allrmserror[i]=allres[i]->GetFunction("gaus")->GetParError(2);
       TF1 *newparam_old=new TF1("error_param_old","[0]*x*exp(-x*[1])+[2]",0.,8);
       //         TF1 *newparam_old=new TF1("error_param_old","1/sqrt(12)",0.,8);
    newparam_old->SetParameter(0,-3.18e-01);
    newparam_old->SetParameter(1,6.67e-01);
    newparam_old->SetParameter(2,3.19e-01);
    //    float range=newparam_old->Eval(float(i)/4+0.125)*sqrt(30);
    float range=5;
    tibres[i]->GetXaxis()->SetRangeUser(-range,range);
    tidres[i]->GetXaxis()->SetRangeUser(-range,range);
    tobres[i]->GetXaxis()->SetRangeUser(-range,range);
    tecres[i]->GetXaxis()->SetRangeUser(-range,range);
    allres[i]->GetXaxis()->SetRangeUser(-range,range);
    tibrms[i]=tibres[i]->GetRMS();
    tibrmserror[i]=tibres[i]->GetRMSError();
    tobrms[i]=tobres[i]->GetRMS();
    tobrmserror[i]=tobres[i]->GetRMSError();
    tidrms[i]=tidres[i]->GetRMS();
    tidrmserror[i]=tidres[i]->GetRMSError();
    tecrms[i]=tecres[i]->GetRMS();
    tecrmserror[i]=tecres[i]->GetRMSError();
    allrms[i]=allres[i]->GetRMS();
    allrmserror[i]=allres[i]->GetRMSError();
  }
  float proj[NBINS],proje[NBINS];
  for (int i=0;i<NBINS;i++){proj[i]=i; proje[i]=0.;};
  TGraphErrors *TIB=new TGraphErrors(NBINS,proj,tibrms,proje,tibrmserror);
  TGraphErrors *TOB=new TGraphErrors(NBINS,proj,tobrms,proje,tobrmserror);
  TGraphErrors *TID=new TGraphErrors(NBINS,proj,tidrms,proje,tidrmserror);
  TGraphErrors *TEC=new TGraphErrors(NBINS,proj,tecrms,proje,tecrmserror);
  TGraphErrors *ALL=new TGraphErrors(NBINS,proj,allrms,proje,allrmserror);

  TF1 *errparamtib=new TF1("error_paramtib","((x-[0])*(x-[0])*([1]-[2])/([0]*[0])+[2])",0.01,4);
  errparamtib->SetParameter(0,38.07*0.032);
  errparamtib->SetParameter(1,0.3184);
  errparamtib->SetParameter(2,0.09828);
  TF1 *errparamtob=new TF1("error_paramtob","((x-[0])*(x-[0])*([1]-[2])/([0]*[0])+[2])",0.01,4);
  errparamtob->SetParameter(0,38.07*0.05);
  errparamtob->SetParameter(1,0.3184);
  errparamtob->SetParameter(2,0.09828);

  //  TF1 *newparam=new TF1("error_param","(((x-[0])*(x-[0])*([1]-[2])/([0]*[0]))*(([3]-x)/[3])+[2])*(((x-[0])*(x-[0])*([1]-[2])/([0]*[0]))*(([3]-x)/[3])+[2])",0,8); 
  TF1 *newparam=new TF1("error_param","[0]*x*exp(-x*[1])+[2]",0.,4);
  //TF1 *newparam=new TF1("error_param","[0]*exp(-0.5*((x-[3])/[1])**2)+[2]*x",0.0001,8);
  // TF1 *newparam=new TF1("error_param",fit_func,0.0001,8,4);
  newparam->SetParameter(0,-0.411);
  newparam->SetParameter(1,0.740);
  newparam->SetParameter(2,0.325);
  // newparam->SetParameter(3,0.25);
 
 TCanvas * plot=new TCanvas("resolution","resolution");
  plot->Divide(2,2);
  plot->cd(1);
  TIB->SetMaximum(4);
  TIB->SetMinimum(0.);
  TIB->GetHistogram()->GetXaxis()->SetTitle("Track projection");
  TIB->GetHistogram()->SetTitle("TIB");
  TIB->Draw("ap");
  //  TIB->Fit(newparam);
  TLine *l=new TLine(4,4/sqrt(12),32,32./sqrt(12));
  l->SetLineColor(3);
  l->SetLineWidth(2);
  l->Draw("same");
  errparamtib->SetLineColor(4);
  //errparamtib->Draw("same");
  // newparam->Draw("same");

  plot->cd(2);
  TOB->SetMaximum(4);
  TOB->SetMinimum(0.);
  TOB->GetHistogram()->GetXaxis()->SetTitle("Track projection");
  TOB->GetHistogram()->SetTitle("TOB");
  TOB->Draw("ap");
  //  TOB->Fit(newparam);
  TLine *l=new TLine(4,4./sqrt(12),32,32./sqrt(12));
  l->SetLineColor(3);
  l->SetLineWidth(2);
  l->Draw("same");
  errparamtob->SetLineColor(4);
  // errparamtob->Draw("same");
  // newparam->Draw("same");

  plot->cd(3);
  TID->SetMaximum(4);
  TID->SetMinimum(0.);
  TID->GetHistogram()->GetXaxis()->SetTitle("Track projection");
  TID->GetHistogram()->SetTitle("TID");
  TID->Draw("ap");
  //  TID->Fit(newparam);
  TLine *l=new TLine(4,4./sqrt(12),32,32./sqrt(12));
  l->SetLineColor(3);
  l->Draw("same");
  errparamtib->SetLineColor(4);
  //errparamtib->Draw("same");
  // newparam->Draw("same");

  plot->cd(4);
  TEC->SetMaximum(4);
  TEC->SetMinimum(0.);
  TEC->GetHistogram()->GetXaxis()->SetTitle("Track projection");
  TEC->GetHistogram()->SetTitle("TEC");
  TEC->Draw("ap");
  //  TEC->Fit(newparam);
  TLine *l=new TLine(4,4./sqrt(12),32,32./sqrt(12));
  l->SetLineColor(3);
  l->SetLineWidth(2);
  l->Draw("same");
  // newparam->Draw("same");

 TCanvas * plot2=new TCanvas("resolution2","resolution2");
 ALL->SetMaximum(4);
 ALL->SetMinimum(0.);
 ALL->Draw("ap");
 // ALL->Fit(newparam);
 TLine *l=new TLine(4,2./sqrt(12),32,16./sqrt(12));
 l->SetLineColor(3);
 l->Draw("same");
 // newparam_old->Draw("same");
}