Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:33:34

0001 #include "BSvsPVPlots.h"
0002 #include <vector>
0003 #include <algorithm>
0004 #include <string>
0005 #include <map>
0006 #include <iostream>
0007 #include "TPad.h"
0008 #include "TH1D.h"
0009 #include "TFile.h"
0010 #include "TH2F.h"
0011 #include "TH1F.h"
0012 #include "TF1.h"
0013 #include "TProfile.h"
0014 #include "TProfile2D.h"
0015 #include "DPGAnalysis/SiStripTools/interface/CommonAnalyzer.h"
0016 #include "TCanvas.h"
0017 #include "TSystem.h"
0018 #include "TStyle.h"
0019 #include "TText.h"
0020 #include "TLegend.h"
0021 
0022 void BSvsPVPlots(const char* fullname,
0023                  const char* module,
0024                  const char* label,
0025                  const char* postfix,
0026                  const char* shortname,
0027                  const char* outtrunk) {
0028   char modfull[300];
0029   sprintf(modfull, "%s%s", module, postfix);
0030   char labfull[300];
0031   sprintf(labfull, "%s%s", label, postfix);
0032 
0033   char dirname[400];
0034   sprintf(dirname, "%s", shortname);
0035 
0036   //  char fullname[300];
0037   //  if(strlen(family)==0) {  sprintf(fullname,"rootfiles/Tracking_PFG_%s.root",filename);}
0038   //  else {  sprintf(fullname,"rootfiles/%s.root",dirname); }
0039 
0040   std::string workdir = outtrunk;
0041   workdir += dirname;
0042   gSystem->cd(workdir.c_str());
0043   gSystem->MakeDirectory(labfull);
0044   //  gSystem->cd("/afs/cern.ch/cms/tracking/output");
0045 
0046   TFile ff(fullname);
0047 
0048   gStyle->SetOptStat(111111);
0049   // Colliding events
0050 
0051   CommonAnalyzer castat(&ff, "", modfull);
0052   {
0053     TH1F* deltax = (TH1F*)castat.getObject("deltax");
0054     if (deltax && deltax->GetEntries() > 0) {
0055       deltax->Draw();
0056       gPad->SetLogy(1);
0057       std::string plotfilename;
0058       plotfilename += outtrunk;
0059       plotfilename += dirname;
0060       plotfilename += "/";
0061       plotfilename += labfull;
0062       plotfilename += "/deltax_";
0063       plotfilename += labfull;
0064       plotfilename += "_";
0065       plotfilename += dirname;
0066       plotfilename += ".gif";
0067       gPad->Print(plotfilename.c_str());
0068       gPad->SetLogy(0);
0069       delete deltax;
0070     }
0071     TH1F* deltay = (TH1F*)castat.getObject("deltay");
0072     if (deltay && deltay->GetEntries() > 0) {
0073       deltay->Draw();
0074       gPad->SetLogy(1);
0075       std::string plotfilename;
0076       plotfilename += outtrunk;
0077       plotfilename += dirname;
0078       plotfilename += "/";
0079       plotfilename += labfull;
0080       plotfilename += "/deltay_";
0081       plotfilename += labfull;
0082       plotfilename += "_";
0083       plotfilename += dirname;
0084       plotfilename += ".gif";
0085       gPad->Print(plotfilename.c_str());
0086       gPad->SetLogy(0);
0087       delete deltay;
0088     }
0089     TH1F* deltaz = (TH1F*)castat.getObject("deltaz");
0090     if (deltaz && deltaz->GetEntries() > 0) {
0091       deltaz->Draw();
0092       gPad->SetLogy(1);
0093       std::string plotfilename;
0094       plotfilename += outtrunk;
0095       plotfilename += dirname;
0096       plotfilename += "/";
0097       plotfilename += labfull;
0098       plotfilename += "/deltaz_";
0099       plotfilename += labfull;
0100       plotfilename += "_";
0101       plotfilename += dirname;
0102       plotfilename += ".gif";
0103       gPad->Print(plotfilename.c_str());
0104       gPad->SetLogy(0);
0105       delete deltaz;
0106     }
0107     gStyle->SetOptStat(111);
0108     gStyle->SetOptFit(111);
0109     TProfile* deltaxvsz = (TProfile*)castat.getObject("deltaxvsz");
0110     if (deltaxvsz && deltaxvsz->GetEntries() > 0) {
0111       //    deltaxvsz->Draw();
0112       deltaxvsz->Fit("pol1", "", "", -10., 10.);
0113       if (deltaxvsz->GetFunction("pol1")) {
0114         deltaxvsz->GetFunction("pol1")->SetLineColor(kRed);
0115         deltaxvsz->GetFunction("pol1")->SetLineWidth(1);
0116       }
0117       deltaxvsz->GetYaxis()->SetRangeUser(-0.001, 0.001);
0118       std::string plotfilename;
0119       plotfilename += outtrunk;
0120       plotfilename += dirname;
0121       plotfilename += "/";
0122       plotfilename += labfull;
0123       plotfilename += "/deltaxvsz_";
0124       plotfilename += labfull;
0125       plotfilename += "_";
0126       plotfilename += dirname;
0127       plotfilename += ".gif";
0128       gPad->Print(plotfilename.c_str());
0129       delete deltaxvsz;
0130     }
0131     TProfile* deltayvsz = (TProfile*)castat.getObject("deltayvsz");
0132     if (deltayvsz && deltayvsz->GetEntries() > 0) {
0133       //    deltayvsz->Draw();
0134       deltayvsz->Fit("pol1", "", "", -10., 10.);
0135       if (deltayvsz->GetFunction("pol1")) {
0136         deltayvsz->GetFunction("pol1")->SetLineColor(kRed);
0137         deltayvsz->GetFunction("pol1")->SetLineWidth(1);
0138       }
0139       deltayvsz->GetYaxis()->SetRangeUser(-0.001, 0.001);
0140       std::string plotfilename;
0141       plotfilename += outtrunk;
0142       plotfilename += dirname;
0143       plotfilename += "/";
0144       plotfilename += labfull;
0145       plotfilename += "/deltayvsz_";
0146       plotfilename += labfull;
0147       plotfilename += "_";
0148       plotfilename += dirname;
0149       plotfilename += ".gif";
0150       gPad->Print(plotfilename.c_str());
0151       delete deltayvsz;
0152     }
0153   }
0154   gStyle->SetOptStat(111111);
0155   gStyle->SetOptFit(1111);
0156 
0157   // Define fitting functions
0158   TF1* fdoubleg = new TF1("doubleg",
0159                           "[1]*exp(-0.5*((x-[0])/[2])**2)+[3]*exp(-0.5*((x-[0])/[4])**2)+[5]*exp(sqrt((x-[0])**2)*[6])",
0160                           -.2,
0161                           .2);
0162   fdoubleg->SetLineColor(kRed);
0163   fdoubleg->SetLineWidth(1);
0164   /*
0165   fdoubleg->SetParLimits(1,0.,1e9);
0166   fdoubleg->SetParLimits(3,0.,1e9);
0167   fdoubleg->SetParLimits(5,0.,1e9);
0168   */
0169 
0170   gStyle->SetOptFit(1111);
0171 
0172   // Summary histograms
0173   TH1D* deltaxsum = new TH1D("deltaxsum", "(PV-BS) Fitted X position vs run", 10, 0., 10.);
0174   deltaxsum->SetCanExtend(TH1::kAllAxes);
0175   TH1D* deltaysum = new TH1D("deltaysum", "(PV-BS) Fitted Y position vs run", 10, 0., 10.);
0176   deltaysum->SetCanExtend(TH1::kAllAxes);
0177 
0178   TH1D* deltaxmeansum = new TH1D("deltaxmeansum", "(PV-BS) Mean X position vs run", 10, 0., 10.);
0179   deltaxmeansum->SetCanExtend(TH1::kAllAxes);
0180   TH1D* deltaymeansum = new TH1D("deltaymeansum", "(PV-BS) Mean Y position vs run", 10, 0., 10.);
0181   deltaymeansum->SetCanExtend(TH1::kAllAxes);
0182 
0183   TH1D* deltazmeansum = new TH1D("deltazmeansum", "(PV-BS) Mean Z position vs run", 10, 0., 10.);
0184   deltazmeansum->SetCanExtend(TH1::kAllAxes);
0185 
0186   std::vector<unsigned int> runs = castat.getRunList();
0187   std::sort(runs.begin(), runs.end());
0188 
0189   {
0190     std::cout << "Found " << runs.size() << " runs" << std::endl;
0191 
0192     for (unsigned int i = 0; i < runs.size(); ++i) {
0193       char runlabel[100];
0194       sprintf(runlabel, "%d", runs[i]);
0195       char runpath[100];
0196       sprintf(runpath, "run_%d", runs[i]);
0197       castat.setPath(runpath);
0198       std::cout << runpath << std::endl;
0199 
0200       TH1F* deltax = (TH1F*)castat.getObject("deltaxrun");
0201       if (deltax && deltax->GetEntries() > 0) {
0202         //  deltax->Draw();
0203         fdoubleg->SetParameter(0, deltax->GetMean());
0204         fdoubleg->SetParameter(2, deltax->GetRMS());
0205         fdoubleg->SetParameter(4, deltax->GetRMS());
0206         fdoubleg->SetParameter(1, deltax->GetMaximum());
0207         fdoubleg->SetParameter(3, 0.1 * deltax->GetMaximum());
0208         fdoubleg->SetParameter(5, 0.1 * deltax->GetMaximum());
0209         const int result = deltax->Fit(fdoubleg, "b", "", -.05, .05);
0210         gPad->SetLogy(1);
0211         char tresult[100];
0212         sprintf(tresult, "%d", result);
0213         TText res;
0214         res.SetTextColor(kRed);
0215         if (result != 0)
0216           res.DrawTextNDC(.2, .8, tresult);
0217 
0218         int bin = deltaxsum->Fill(runlabel, fdoubleg->GetParameter(0));
0219         deltaxsum->SetBinError(bin, fdoubleg->GetParError(0));
0220 
0221         bin = deltaxmeansum->Fill(runlabel, deltax->GetMean());
0222         deltaxmeansum->SetBinError(bin, deltax->GetMeanError());
0223 
0224         std::string plotfilename;
0225         plotfilename += outtrunk;
0226         plotfilename += dirname;
0227         plotfilename += "/";
0228         plotfilename += labfull;
0229         plotfilename += "/deltaxrun_";
0230         plotfilename += labfull;
0231         plotfilename += "_";
0232         plotfilename += dirname;
0233         plotfilename += "_";
0234         plotfilename += runpath;
0235         plotfilename += ".gif";
0236         gPad->Print(plotfilename.c_str());
0237         gPad->SetLogy(0);
0238         delete deltax;
0239       }
0240       TH1F* deltay = (TH1F*)castat.getObject("deltayrun");
0241       if (deltay && deltay->GetEntries() > 0) {
0242         deltay->Draw();
0243         fdoubleg->SetParameter(0, deltay->GetMean());
0244         fdoubleg->SetParameter(2, deltay->GetRMS());
0245         fdoubleg->SetParameter(4, deltay->GetRMS());
0246         fdoubleg->SetParameter(1, deltay->GetMaximum());
0247         fdoubleg->SetParameter(3, 0.1 * deltay->GetMaximum());
0248         fdoubleg->SetParameter(5, 0.1 * deltay->GetMaximum());
0249         const int result = deltay->Fit(fdoubleg, "b", "", -.05, .05);
0250         gPad->SetLogy(1);
0251         char tresult[100];
0252         sprintf(tresult, "%d", result);
0253         TText res;
0254         res.SetTextColor(kRed);
0255         if (result != 0)
0256           res.DrawTextNDC(.2, .8, tresult);
0257 
0258         int bin = deltaysum->Fill(runlabel, fdoubleg->GetParameter(0));
0259         deltaysum->SetBinError(bin, fdoubleg->GetParError(0));
0260 
0261         bin = deltaymeansum->Fill(runlabel, deltay->GetMean());
0262         deltaymeansum->SetBinError(bin, deltay->GetMeanError());
0263 
0264         std::string plotfilename;
0265         plotfilename += outtrunk;
0266         plotfilename += dirname;
0267         plotfilename += "/";
0268         plotfilename += labfull;
0269         plotfilename += "/deltayrun_";
0270         plotfilename += labfull;
0271         plotfilename += "_";
0272         plotfilename += dirname;
0273         plotfilename += "_";
0274         plotfilename += runpath;
0275         plotfilename += ".gif";
0276         gPad->Print(plotfilename.c_str());
0277         gPad->SetLogy(0);
0278         delete deltay;
0279       }
0280       TH1F* deltaz = (TH1F*)castat.getObject("deltazrun");
0281       if (deltaz && deltaz->GetEntries() > 0) {
0282         deltaz->Draw();
0283         gPad->SetLogy(1);
0284 
0285         int bin = deltazmeansum->Fill(runlabel, deltaz->GetMean());
0286         deltazmeansum->SetBinError(bin, deltaz->GetMeanError());
0287 
0288         std::string plotfilename;
0289         plotfilename += outtrunk;
0290         plotfilename += dirname;
0291         plotfilename += "/";
0292         plotfilename += labfull;
0293         plotfilename += "/deltazrun_";
0294         plotfilename += labfull;
0295         plotfilename += "_";
0296         plotfilename += dirname;
0297         plotfilename += "_";
0298         plotfilename += runpath;
0299         plotfilename += ".gif";
0300         gPad->Print(plotfilename.c_str());
0301         gPad->SetLogy(0);
0302         delete deltaz;
0303       }
0304 
0305       TProfile* deltaxvsz = (TProfile*)castat.getObject("deltaxvszrun");
0306       if (deltaxvsz && deltaxvsz->GetEntries() > 0) {
0307         deltaxvsz->Draw();
0308         std::string plotfilename;
0309         plotfilename += outtrunk;
0310         plotfilename += dirname;
0311         plotfilename += "/";
0312         plotfilename += labfull;
0313         plotfilename += "/deltaxvszrun_";
0314         plotfilename += labfull;
0315         plotfilename += "_";
0316         plotfilename += dirname;
0317         plotfilename += "_";
0318         plotfilename += runpath;
0319         plotfilename += ".gif";
0320         gPad->Print(plotfilename.c_str());
0321         delete deltaxvsz;
0322       }
0323 
0324       TProfile* deltayvsz = (TProfile*)castat.getObject("deltayvszrun");
0325       if (deltayvsz && deltayvsz->GetEntries() > 0) {
0326         deltayvsz->Draw();
0327         std::string plotfilename;
0328         plotfilename += outtrunk;
0329         plotfilename += dirname;
0330         plotfilename += "/";
0331         plotfilename += labfull;
0332         plotfilename += "/deltayvszrun_";
0333         plotfilename += labfull;
0334         plotfilename += "_";
0335         plotfilename += dirname;
0336         plotfilename += "_";
0337         plotfilename += runpath;
0338         plotfilename += ".gif";
0339         gPad->Print(plotfilename.c_str());
0340         delete deltayvsz;
0341       }
0342 
0343       TH1F* deltaxvsorb = (TH1F*)castat.getObject("deltaxvsorbrun");
0344       if (deltaxvsorb && deltaxvsorb->GetEntries() > 0) {
0345         deltaxvsorb->Draw();
0346         std::string plotfilename;
0347         plotfilename += outtrunk;
0348         plotfilename += dirname;
0349         plotfilename += "/";
0350         plotfilename += labfull;
0351         plotfilename += "/deltaxvsorb_";
0352         plotfilename += labfull;
0353         plotfilename += "_";
0354         plotfilename += dirname;
0355         plotfilename += "_";
0356         plotfilename += runpath;
0357         plotfilename += ".gif";
0358         gPad->Print(plotfilename.c_str());
0359         delete deltaxvsorb;
0360       }
0361       TH1F* deltayvsorb = (TH1F*)castat.getObject("deltayvsorbrun");
0362       if (deltayvsorb && deltayvsorb->GetEntries() > 0) {
0363         deltayvsorb->Draw();
0364         std::string plotfilename;
0365         plotfilename += outtrunk;
0366         plotfilename += dirname;
0367         plotfilename += "/";
0368         plotfilename += labfull;
0369         plotfilename += "/deltayvsorb_";
0370         plotfilename += labfull;
0371         plotfilename += "_";
0372         plotfilename += dirname;
0373         plotfilename += "_";
0374         plotfilename += runpath;
0375         plotfilename += ".gif";
0376         gPad->Print(plotfilename.c_str());
0377         delete deltayvsorb;
0378       }
0379       TH1F* deltazvsorb = (TH1F*)castat.getObject("deltazvsorbrun");
0380       if (deltazvsorb && deltazvsorb->GetEntries() > 0) {
0381         deltazvsorb->Draw();
0382         std::string plotfilename;
0383         plotfilename += outtrunk;
0384         plotfilename += dirname;
0385         plotfilename += "/";
0386         plotfilename += labfull;
0387         plotfilename += "/deltazvsorb_";
0388         plotfilename += labfull;
0389         plotfilename += "_";
0390         plotfilename += dirname;
0391         plotfilename += "_";
0392         plotfilename += runpath;
0393         plotfilename += ".gif";
0394         gPad->Print(plotfilename.c_str());
0395         delete deltazvsorb;
0396       }
0397     }
0398   }
0399 
0400   gStyle->SetOptStat(1111);
0401   gStyle->SetOptFit(0);
0402 
0403   if (!runs.empty()) {
0404     std::string plotfilename;
0405 
0406     plotfilename = outtrunk;
0407     plotfilename += dirname;
0408     plotfilename += "/";
0409     plotfilename += labfull;
0410     plotfilename += "/deltaxsum_";
0411     plotfilename += labfull;
0412     plotfilename += "_";
0413     plotfilename += dirname;
0414     plotfilename += ".gif";
0415 
0416     TCanvas* cwidedeltax = new TCanvas(plotfilename.c_str(), plotfilename.c_str(), 1500, 500);
0417 
0418     deltaxsum->SetLineColor(kRed);
0419     deltaxsum->SetMarkerColor(kRed);
0420     deltaxsum->GetYaxis()->SetRangeUser(-.002, .002);
0421     deltaxsum->GetYaxis()->SetTitle("#Delta x (cm)");
0422     deltaxsum->Draw();
0423     deltaxmeansum->Draw("esame");
0424     TLegend deltaxleg(.7, .8, .85, .9, "#Delta(x)");
0425     deltaxleg.AddEntry(deltaxsum, "fitted mean", "l");
0426     deltaxleg.AddEntry(deltaxmeansum, "aritm. mean", "l");
0427     deltaxleg.Draw();
0428     gPad->Print(plotfilename.c_str());
0429     delete cwidedeltax;
0430 
0431     plotfilename = outtrunk;
0432     plotfilename += dirname;
0433     plotfilename += "/";
0434     plotfilename += labfull;
0435     plotfilename += "/deltaysum_";
0436     plotfilename += labfull;
0437     plotfilename += "_";
0438     plotfilename += dirname;
0439     plotfilename += ".gif";
0440 
0441     TCanvas* cwidedeltay = new TCanvas(plotfilename.c_str(), plotfilename.c_str(), 1500, 500);
0442 
0443     deltaysum->SetLineColor(kRed);
0444     deltaysum->SetMarkerColor(kRed);
0445     deltaysum->GetYaxis()->SetRangeUser(-.002, .002);
0446     deltaysum->GetYaxis()->SetTitle("#Delta y (cm)");
0447     deltaysum->Draw();
0448     deltaymeansum->Draw("esame");
0449     TLegend deltayleg(.7, .8, .85, .9, "#Delta(y)");
0450     deltayleg.AddEntry(deltaysum, "fitted mean", "l");
0451     deltayleg.AddEntry(deltaymeansum, "aritm. mean", "l");
0452     deltayleg.Draw();
0453     gPad->Print(plotfilename.c_str());
0454     delete cwidedeltay;
0455 
0456     plotfilename = outtrunk;
0457     plotfilename += dirname;
0458     plotfilename += "/";
0459     plotfilename += labfull;
0460     plotfilename += "/deltazsum_";
0461     plotfilename += labfull;
0462     plotfilename += "_";
0463     plotfilename += dirname;
0464     plotfilename += ".gif";
0465 
0466     TCanvas* cwidedeltaz = new TCanvas(plotfilename.c_str(), plotfilename.c_str(), 1500, 500);
0467 
0468     deltazmeansum->GetYaxis()->SetRangeUser(-2., 2.);
0469     deltazmeansum->GetYaxis()->SetTitle("#Delta z (cm)");
0470     deltazmeansum->Draw();
0471     gPad->Print(plotfilename.c_str());
0472     delete cwidedeltaz;
0473   }
0474   delete deltaxsum;
0475   delete deltaysum;
0476   delete deltaxmeansum;
0477   delete deltaymeansum;
0478   delete deltazmeansum;
0479 
0480   ff.Close();
0481   delete fdoubleg;
0482 }