ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/tschum/FWlite_Analysis/PlotTool.cc
Revision: 1.16
Committed: Mon Feb 8 13:22:13 2010 UTC (15 years, 2 months ago) by thomsen
Content type: text/plain
Branch: MAIN
Changes since 1.15: +201 -90 lines
Log Message:
Track-Cluster matching

File Contents

# User Rev Content
1 tschum 1.1 #include "PlotTool.h"
2     //------------------------------------------------------------------
3     //Constructur: Set some global Variables to define canvas output
4 gebbert 1.6 PlotTool::PlotTool() {
5 tschum 1.1
6 gebbert 1.6 this->SetClass("TChain", 100);
7 tschum 1.1
8 gebbert 1.6 samePad_trees = true;
9     samePad_vars = false;
10     samePad_cuts = false;
11     sameCanv_trees = false;
12     sameCanv_vars = false;
13     sameCanv_cuts = false;
14 tschum 1.2
15 gebbert 1.6 showLegend = false;
16     logY = true;
17 thomsen 1.15 addTrackJets = true;
18     addHitDetInfo = true;
19     diTrackMass = true;
20     addTrigger = true;
21 tschum 1.11 addTrackJets = false;
22 tschum 1.13 addEventInfo = false;
23 tschum 1.14 addTower = false;
24 tschum 1.8 verbose = true;
25 tschum 1.1
26 tschum 1.12 globalCuts="";
27    
28 tschum 1.14 recreateTree = false; //force recreation of friend tree
29    
30 tschum 1.1 }
31    
32     //------------------------------------------------------------------
33     //Fill object PlotTool with Chains constructed from files from given source
34 tschum 1.7
35 gebbert 1.6 int PlotTool::init(string fileName, string dirPath, string treeName,
36 tschum 1.13 string fileLabel) {
37     this->New(this->GetEntries() );
38     int currChain = this->GetEntries() - 1;
39 gebbert 1.6
40 tschum 1.13 if (currChain < 0)
41     return currChain;
42 gebbert 1.6
43 tschum 1.13 TStopwatch timer;
44     timer.Start();
45    
46    
47     fileNames.clear();
48    
49     //check if alternative label is different from default(searchString)
50     if (fileLabel=="") {
51     fileLabel=fileName;
52     }
53     ((TChain*) this->At(currChain))->SetName(fileLabel.c_str());
54    
55     TList *files = new TList();
56     TSystemFile* sysFile = 0;
57     if(fileName.find(".") != string::npos) {
58     ifstream f(fileName.c_str());
59     if( ! f.is_open() ) return -1;
60     string line;
61    
62     while (!f.eof()) {
63     getline(f,line);
64     sysFile = new TSystemFile(line.c_str(),dirPath.c_str());
65     files->Add(sysFile);
66     }
67 tschum 1.8
68    
69 tschum 1.13 } else {
70     TSystemDirectory dir("sourceDir", dirPath.c_str());
71     files = dir.GetListOfFiles();
72     }
73     if (files->GetEntries()>0) {
74     TIter next(files);
75     TSystemFile *file;
76     TString fname;
77     string filePath;
78    
79     if(verbose && fileName.find(".") == string::npos ) cout<<"Open"<<dirPath.c_str()<<" Search for .root files that contain: "
80     <<fileName.c_str()<<endl;
81     if(verbose && fileName.find(".") != string::npos ) cout<<"Open"<<dirPath.c_str()<<" Search lines with .root in: "
82     <<fileName.c_str()<<endl;
83    
84    
85     while ((file=(TSystemFile*)next())) {
86     fname = file->GetName();
87     if (!fname.EndsWith(".root"))
88     continue;
89     if (!fname.Contains(fileName.c_str()) && fileName.find(".") == string::npos )
90     continue;
91    
92     filePath = dirPath;
93     filePath += fname.Data();
94     //fwlite::ChainEvent to lop over events, jets, etc
95     fileNames.push_back(filePath.c_str());
96    
97     if (((TChain*) this->At(currChain))->AddFile(filePath.c_str(), -1,
98     treeName.c_str()) )
99     if(verbose) cout<<"Chained "<<((TChain*) this->At(currChain))->GetNtrees()<<" file(s) with "<<((TChain*) this->At(currChain))->GetEntries()<<" events."<<endl;
100     else
101     return -1;
102 gebbert 1.6
103 tschum 1.13 }
104     } else
105     return -1;
106    
107     for (int i=0; i<((TChain*) this->At(currChain))->GetListOfBranches()->GetEntries(); ++i) {
108    
109     string s(((TChain*) this->At(currChain))->GetListOfBranches()->At(i)->GetName());
110    
111     string branch_alias = s;
112     string branch_name = s;
113     if (s.find(".", s.size()-1) != string::npos)
114     branch_name += "obj";
115    
116     size_t a = s.find("_");
117     if (a != string::npos) {
118     size_t b = s.find("_", a+1);
119     if (b != string::npos) {
120     size_t c = s.find("_", b+1);
121     if (c != string::npos) {
122     string _prod =s.substr(0,a);
123     string _label =s.substr(a+1, b-a-1);
124     string _instance =s.substr(b+1, c-b-1);
125     branch_alias = _label;
126     if(_instance.length() > 0 ) {
127     branch_alias += "_";
128     branch_alias += _instance;
129     ((TChain*) this->At(currChain))->SetAlias(_instance.c_str(),
130     branch_name.c_str());
131     }
132     string branch_alias_full = _prod + "_" + branch_alias;
133     ((TChain*) this->At(currChain))->SetAlias(branch_alias_full.c_str(),
134     branch_name.c_str());
135 gebbert 1.6 }
136 tschum 1.13 }
137     }
138    
139     ((TChain*) this->At(currChain))->SetAlias(branch_alias.c_str(),
140     branch_name.c_str());
141    
142     }
143    
144 thomsen 1.15 if( addEventInfo || addTrackJets || addTower || addTrackJets || addHitDetInfo || diTrackMass || addTrigger) {
145 tschum 1.10
146 tschum 1.14 if(verbose) cout<<"add friend tree with additional variables..."<<endl;
147 thomsen 1.16
148 tschum 1.13 //make file for tree friends (adding additional, computed branches)
149     string friendFileName("/scratch/hh/current/cms/user/");
150     friendFileName += gSystem->GetUserInfo()->fUser;
151 tschum 1.14 friendFileName += "/temp/";
152     friendFileName += fileLabel;
153 tschum 1.13 friendFileName += ".root";
154 tschum 1.14 string fileOpt = "update";
155     if( recreateTree ) fileOpt = "recreate";
156    
157     TFile *f = new TFile(friendFileName.c_str(),fileOpt.c_str());
158 thomsen 1.16
159 tschum 1.13 if( f->IsZombie() ) return -1;
160 thomsen 1.16
161 tschum 1.13 string friendTreeName("friendTree");
162     //char number = currChain;
163     //friendTreeName += number;
164     //TTree *friendTree = new TTree("friendTree","friendTree");
165    
166 tschum 1.14 if(! f->FindKey(friendTreeName.c_str())) {
167     if(verbose) cout<<"calculating additional variables..."<<endl;
168    
169     TTree *friendTree = new TTree(friendTreeName.c_str(),friendTreeName.c_str());
170     fwlite::ChainEvent ev(fileNames);
171     int _event, _run, _lumi;
172    
173     //tower data
174     const int kMAX = 10000;
175     int NobjTowCal;
176     towet = new float [ kMAX ];
177     toweta = new float [ kMAX ];
178     towphi = new float [ kMAX ];
179     towen = new float [ kMAX ];
180     towem = new float [ kMAX ];
181     towhd = new float [ kMAX ];
182     towoe = new float [ kMAX ];
183     towid_phi = new int[ kMAX ];
184     towid_eta = new int[ kMAX ];
185     towid = new int[ kMAX ];
186    
187 thomsen 1.15 int nJetsKT;
188     TrackJetKT = new float [100];
189     int nTracks;
190     TrackDetE = new float [3000];
191     TrackDetEECAL = new float [3000];
192     float DiTrackMass;
193     bool Trigger;
194     bool Trigger1;
195    
196    
197    
198 tschum 1.14 if( addEventInfo ) {
199 tschum 1.13
200    
201 tschum 1.14 friendTree->Branch("event", &_event, "event/I");
202     friendTree->Branch("run", &_run, "run/I");
203     friendTree->Branch("lumi", &_lumi, "lumi/I");
204    
205     }
206    
207     if( addTrackJets ) {
208    
209    
210     friendTree->Branch("nJetsKT", &nJetsKT, "nJetsKT/I");
211     friendTree->Branch("TrackJetKT", TrackJetKT, "TrackJetKT[nJetsKT]/F");
212     }
213    
214     if( addTower) {
215    
216     // CaloTower branches
217     friendTree->Branch( "NobjTowCal",&NobjTowCal,"NobjTowCal/I" );
218     friendTree->Branch( "TowId", towid, "TowId[NobjTowCal]/I" );
219     friendTree->Branch( "TowId_phi", towid_phi, "TowId_phi[NobjTowCal]/I" );
220     friendTree->Branch( "TowId_eta", towid_eta, "TowId_eta[NobjTowCal]/I" );
221     friendTree->Branch( "TowEt", towet, "TowEt[NobjTowCal]/F" );
222     friendTree->Branch( "TowEta", toweta, "TowEta[NobjTowCal]/F" );
223     friendTree->Branch( "TowPhi", towphi, "TowPhi[NobjTowCal]/F" );
224     friendTree->Branch( "TowE", towen, "TowE[NobjTowCal]/F" );
225     friendTree->Branch( "TowEm", towem, "TowEm[NobjTowCal]/F" );
226     friendTree->Branch( "TowHad", towhd, "TowHad[NobjTowCal]/F" );
227     friendTree ->Branch( "TowOE", towoe, "TowOE[NobjTowCal]/F" );
228    
229    
230     }
231    
232 thomsen 1.15 if( addTrackJets){
233     friendTree->Branch("nJetsKT", &nJetsKT, "nJetsKT/I");
234     friendTree->Branch("TrackJetKT", TrackJetKT, "TrackJetKT[nJetsKT]/F");
235     }
236    
237     if( addHitDetInfo){
238     friendTree->Branch("nTracks", &nTracks, "nTracks/I");
239     friendTree->Branch("TrackDetE", TrackDetE, "TrackDetE[nTracks]/F");
240     friendTree->Branch("TrackDetEECAL", TrackDetEECAL, "TrackDetEECAL[nTracks]/F");
241     }
242     if(diTrackMass) friendTree->Branch("DiTrackMass", &DiTrackMass, "DiTrackMass/F");
243     if(addTrigger) {
244     friendTree->Branch("Trigger", &Trigger, "Trigger/O");
245     friendTree->Branch("Trigger1", &Trigger1, "Trigger1/O");
246     }
247 tschum 1.14
248     int tenth = ev.size() / 10;
249     int eventNum =0;
250 thomsen 1.16
251 tschum 1.14 for (ev.toBegin(); !ev.atEnd(); ++ev, ++eventNum) {
252    
253     if( addEventInfo ) {
254 tschum 1.10
255    
256 tschum 1.14 if(eventNum==0) {
257     // fwlite::Handle<edm::TriggerResults> hltHandle;
258     // hltHandle.getByLabel(ev, "TriggerResults");
259 gebbert 1.6
260 tschum 1.14 // fwlite::Handle<edm::TriggerResults> hTriggerResults;
261 gebbert 1.6
262 tschum 1.14 // hTriggerResults.getByLabel(ev, "TriggerResults","","TEST");
263     // fwlite::TriggerNames const& triggerNames = ev.triggerNames(*hTriggerResults);
264 gebbert 1.6
265 tschum 1.14 // // std::vector<std::string> const& names = triggerNames.triggerNames();
266     // for (unsigned i = 0; i < triggerNames.size(); ++i) {
267     // std::cout << i << " " << triggerNames.triggerName(i) << std::endl;
268     // }
269     }
270     _event = ev.id().event();
271     _run = ev.id().run();
272     _lumi = ev.luminosityBlock();
273 tschum 1.7
274 tschum 1.14 }
275 gebbert 1.6
276 thomsen 1.15 if( addTrackJets)
277     {
278     fwlite::Handle<reco::CaloJetCollection> jets;
279 thomsen 1.16 jets.getByLabel(ev, "ak5CaloJets");
280 thomsen 1.15
281     fwlite::Handle<reco::TrackCollection> tracks;
282     tracks.getByLabel(ev, "generalTracks");
283    
284 thomsen 1.16
285     if (jets.isValid() && tracks.isValid()){
286     double trackSum;
287     nJetsKT = 0;
288     for (reco::CaloJetCollection::const_iterator jet = jets->begin(); jet
289     !=jets->end(); jet++) {
290     trackSum = 0;
291     for (reco::TrackCollection::const_iterator track = tracks->begin(); track
292     !=tracks->end(); track++) {
293     if (deltaR(jet->eta(), jet->phi(), track->eta(), track->phi())
294     < 0.5)
295     trackSum += track->pt();
296     }
297     TrackJetKT[nJetsKT] = trackSum;
298     nJetsKT++;
299 thomsen 1.15 }
300     }
301     }
302    
303 thomsen 1.16
304 thomsen 1.15 if(addHitDetInfo)
305     {
306     fwlite::Handle<reco::TrackCollection> tracks;
307     tracks.getByLabel(ev, "generalTracks");
308    
309     nTracks = 0;
310    
311 thomsen 1.16 if (tracks.isValid()){
312    
313     for (reco::TrackCollection::const_iterator track = tracks->begin(); track
314     !=tracks->end(); track++) {
315     fwlite::Handle<vector <reco::CaloCluster> > CaloClusters;
316     CaloClusters.getByLabel(ev, "hybridSuperClusters","hybridBarrelBasicClusters");//::hybridBarrelBasicClusters"); //add instance
317     //CaloClusters.getByLabel(ev, "multi5x5BasicClusters");
318    
319     //fwlite::Handle<vector <CaloTower> > CaloTowers;
320     //CaloTowers.getByLabel(ev, "towerMaker");
321 thomsen 1.15
322 thomsen 1.16 double DR = 100;
323     double DetE = 0;
324     double DetEECAL = 0;
325 thomsen 1.15
326 thomsen 1.16 if (CaloClusters.isValid()){
327     //if (CaloTowers.isValid()){
328     for (vector<reco::CaloCluster>::const_iterator CaloCluster = CaloClusters->begin(); CaloCluster
329     !=CaloClusters->end(); CaloCluster++) {
330     // int DetectorId = track->outerDetId(); //unfortuately only the Tracker DetId...
331     // /DetId *dId = new DetId(DetectorId);
332     //cout<<dId->det()<<endl;
333    
334    
335     if(deltaR( CaloCluster->eta(), CaloCluster->phi(), track->outerEta(), track->phi()) < DR)
336     {
337     DR = deltaR( CaloCluster->eta(), CaloCluster->phi(), track->outerEta(), track->phi());
338     DetE = CaloCluster->energy();
339     if(CaloCluster->caloID().detector(reco::CaloID::DET_ECAL_BARREL) ||CaloCluster->caloID().detector(reco::CaloID::DET_ECAL_ENDCAP))
340     DetEECAL = CaloCluster->energy();
341     }
342 thomsen 1.15 }
343 thomsen 1.16 }
344    
345     TrackDetE[nTracks] = DetE;
346     TrackDetEECAL[nTracks] = DetEECAL;
347    
348     nTracks++;
349     if(nTracks > 2990) cout<<"To many Tracks!!!!!!!!!"<<endl;
350 thomsen 1.15 }
351     }
352     }
353    
354     if(diTrackMass)
355     {
356     fwlite::Handle<reco::TrackCollection> tracks;
357     tracks.getByLabel(ev, "generalTracks");
358    
359     DiTrackMass = 0;
360    
361 thomsen 1.16 if (tracks.isValid()){
362 thomsen 1.15
363 thomsen 1.16 for (reco::TrackCollection::const_iterator track = tracks->begin(); track
364     !=tracks->end(); track++) {
365     if(!(track->d0() < 2 && track->dz() < 10 && track->numberOfValidHits() > 7)) continue;
366     for (reco::TrackCollection::const_iterator track1 = tracks->begin(); track1
367     !=tracks->end(); track1++) {
368     if(!(track1->d0() < 2 && track1->dz() < 10 && track1->numberOfValidHits() > 7) || track1 == track) continue;
369    
370     TLorentzVector tr1(track->px(), track->py(), track->pz(), track->p() );
371     TLorentzVector tr2(track1->px(), track1->py(), track1->pz(), track1->p() );
372    
373     double result = (tr1 + tr2).Mag();
374    
375     if(result > DiTrackMass) DiTrackMass = result;
376    
377     }
378 thomsen 1.15 }
379 tschum 1.14 }
380 thomsen 1.15
381     }
382    
383     if(addTrigger)
384     {
385    
386     Trigger = false;
387     Trigger1 = false;
388    
389     fwlite::Handle< L1GlobalTriggerReadoutRecord> gtReadoutRecord;
390     gtReadoutRecord.getByLabel(ev, "gtDigis");
391     //edm::Handle< L1GlobalTriggerReadoutRecord > gtReadoutRecord;
392     //event.getByLabel( edm::InputTag("gtDigis"), gtReadoutRecord);
393    
394 thomsen 1.16 if(gtReadoutRecord.isValid()){
395 thomsen 1.15
396    
397 thomsen 1.16 const TechnicalTriggerWord& technicalTriggerWordBeforeMask = gtReadoutRecord->technicalTriggerWord();
398 thomsen 1.15
399    
400 thomsen 1.16 //std::vector<bool>& technicalTriggerWordBeforeMask = gtReadoutRecord->technicalTriggerWord();
401    
402     bool techTrigger0 = technicalTriggerWordBeforeMask.at(0);
403     bool techTrigger40 = technicalTriggerWordBeforeMask.at(40);
404     bool techTrigger41 = technicalTriggerWordBeforeMask.at(41);
405     bool techTrigger36 = technicalTriggerWordBeforeMask.at(36);
406     bool techTrigger37 = technicalTriggerWordBeforeMask.at(37);
407     bool techTrigger38 = technicalTriggerWordBeforeMask.at(38);
408     bool techTrigger39 = technicalTriggerWordBeforeMask.at(39);
409    
410     if(techTrigger0 && (techTrigger40 || techTrigger41)) Trigger = true;
411     if(Trigger && !techTrigger36 && !techTrigger37 && !techTrigger38 && !techTrigger39) Trigger1 = true;
412     }
413 tschum 1.14 }
414 thomsen 1.15
415 tschum 1.14 if( addTower) {
416     fwlite::Handle<CaloTowerCollection> towers;
417     towers.getByLabel(ev, "towerMaker");
418 thomsen 1.15
419 tschum 1.14 int jtow = 0;
420     NobjTowCal=towers->size();
421     for(CaloTowerCollection::const_iterator tow = towers->begin();
422     tow != towers->end(); ++tow, ++jtow){
423     towet [jtow] = tow->et();
424     toweta[jtow] = tow->eta();
425     towphi[jtow] = tow->phi();
426     towen [jtow] = tow->energy();
427     towem [jtow] = tow->emEnergy();
428     towhd [jtow] = tow->hadEnergy();
429     towoe [jtow] = tow->outerEnergy();
430     towid_phi[jtow] = tow->id().iphi();
431     towid_eta[jtow] = tow->id().ieta();
432     towid [jtow] = tow->id().rawId();
433 tschum 1.7 }
434 gebbert 1.6 }
435 tschum 1.14 friendTree->Fill();
436 thomsen 1.15
437 tschum 1.14 if( eventNum != 0 && eventNum%tenth == 0) cout<<"Processed "<<eventNum <<" of "<<ev.size()<<" events. "<<endl;
438 thomsen 1.15
439 tschum 1.14 }
440     f->cd();
441     friendTree->Write();
442 tschum 1.13 }
443     f->Close();
444     ((TChain*) this->At(currChain))->AddFriend(friendTreeName.c_str(),
445     friendFileName.c_str());
446 thomsen 1.15
447 tschum 1.13 }
448    
449 tschum 1.8
450 tschum 1.13 timer.Stop();
451     if(verbose) timer.Print();
452 tschum 1.8
453 tschum 1.13 return this->GetEntries();
454 tschum 1.1
455     }
456     //------------------------------------------------------------------
457     //Draw given Chain with given variable and given cuts
458 gebbert 1.6 int PlotTool::plot(int chainIndex, string histName, string cutName,
459     int nEntries, string drwOpt) {
460 tschum 1.1
461 gebbert 1.6 if (chainIndex >= this->GetEntries() )
462     return -1;
463 tschum 1.1
464 tschum 1.8 TStopwatch timer;
465 tschum 1.11 if(verbose) {
466     cout<<"Plot: "<<histName<<" "<<cutName<<" from chain "<<this->At(chainIndex)->GetName()<<endl;
467     timer.Start();
468     }
469 tschum 1.8
470 gebbert 1.6 int currN = nEntries;
471     if (nEntries < 0)
472     currN = ((TChain*) this->At(chainIndex))->GetEntries(); //nEntries<0 : all entries are plotted!
473    
474     //++++ Create and name Canvases according to global variables +++++++++++++
475     ostringstream currHistName;
476     if (samePad_trees)
477     currHistName<<((TChain*) this->At(chainIndex))->GetName()<<":";
478     if (samePad_vars)
479     currHistName<<histName;
480     if (samePad_cuts)
481     currHistName<<"{"<<cutName<<"}";
482    
483     ostringstream currPadName;
484     if (!samePad_trees)
485     currPadName<<((TChain*) this->At(chainIndex))->GetName()<<":";
486     if (!samePad_vars)
487     currPadName<<histName;
488     if (!samePad_cuts)
489     currPadName<<"{"<<cutName<<"}";
490    
491     ostringstream currCanvName;
492     if (!sameCanv_trees && !samePad_trees)
493     currCanvName<<((TChain*) this->At(chainIndex))->GetName()<<":";
494     if (!sameCanv_vars && !samePad_vars)
495     currCanvName<<histName;
496     if (!sameCanv_cuts && !samePad_cuts)
497     currCanvName<<"{"<<cutName<<"}";
498    
499     if ( (sameCanv_trees || samePad_trees) && (sameCanv_vars || samePad_vars)
500     && (sameCanv_cuts || samePad_cuts))
501     currCanvName<<"All";
502    
503     string currOpt = drwOpt;
504     bool useSubPads = (sameCanv_trees && !samePad_trees) || (sameCanv_vars
505     && !samePad_vars) || (sameCanv_cuts && !samePad_cuts);
506    
507     if (useSubPads) {
508 tschum 1.7 if ( !canvases_[currCanvName.str()] || !gROOT->GetListOfCanvases()->FindObject(currCanvName.str().c_str()) ) {
509 gebbert 1.6 canvases_[currCanvName.str()] = new TCanvas(currCanvName.str().c_str(), currCanvName.str().c_str());
510     } else if (canvases_[currCanvName.str()]->GetNumber() >= 0) {
511     canvases_[currCanvName.str()]->SetNumber(canvases_[currCanvName.str()]->GetNumber()+1);
512     }
513     }
514 tschum 1.1
515 tschum 1.7 if ( !pads_[currPadName.str()] || !gROOT->GetListOfCanvases()->FindObject(currPadName.str().c_str()) ) {
516 gebbert 1.6 pads_[currPadName.str()] = new TCanvas(currPadName.str().c_str(), currPadName.str().c_str());
517 tschum 1.9 // if (logY)
518     // pads_[currPadName.str()]->SetLogy(1);
519 gebbert 1.6 } else {
520     pads_[currPadName.str()]->cd();
521     currOpt += "SAMES";
522     if (useSubPads)
523     canvases_[currCanvName.str()]->SetNumber(canvases_[currCanvName.str()]->GetNumber()-1);
524     }
525     //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
526 tschum 1.1
527    
528 gebbert 1.6 //Draw histogram:
529 tschum 1.12 string cutNameWithGlobal = cutName;
530     if(cutName!="" && globalCuts!="") cutNameWithGlobal += "&&";
531     if(globalCuts!="") cutNameWithGlobal += globalCuts;
532    
533     int draw_err = ((TChain*) this->At(chainIndex))->Draw(histName.c_str(), cutNameWithGlobal.c_str(),
534 gebbert 1.6 currOpt.c_str(), currN);
535     if (draw_err < 0)
536     return draw_err;
537    
538     //++++ Fix for histos with no entries +++++++++++++
539     if (draw_err == 0) {
540     if (currOpt.find("SAMES") == string::npos) {
541     pads_[currPadName.str()]->Close();
542     pads_.erase(currPadName.str() );
543     if (useSubPads)
544     canvases_[currCanvName.str()]->SetNumber(canvases_[currCanvName.str()]->GetNumber()-1);
545     }
546     cout<< "Warning: "<<currHistName.str().c_str()<<" in "<<currPadName.str()
547     <<" has no entries and is not drawn!"<<endl;
548     return draw_err;
549     }
550     //++++++++++++++++++++++++++++++++++++++++++++++++++
551 tschum 1.1
552 gebbert 1.6 ((TH1F*) pads_[currPadName.str()]->GetPrimitive("htemp"))->SetName(currHistName.str().c_str()); // Set Name of histogram
553    
554 tschum 1.8 // Update for 2D histograms
555     if( pads_[currPadName.str()]->GetPrimitive(currHistName.str().c_str())->InheritsFrom("TH2") ) {
556     pads_[currPadName.str()]->Draw();
557     setCanvas( pads_[currPadName.str()] );
558     }
559 gebbert 1.6
560 tschum 1.11
561 tschum 1.12 if(verbose && draw_err > 0) {
562 tschum 1.11 timer.Stop();
563 tschum 1.13 cout<<"Done: Selected "<<draw_err<<" objects in "<< currN <<" processed events."<<endl;
564 tschum 1.8 timer.Print();
565     }
566 gebbert 1.6 return draw_err;
567 tschum 1.1
568     }
569     //------------------------------------------------------------------
570     //Standard loop to draw all chains and multiple variables and cuts
571 gebbert 1.6 int PlotTool::loop(vector<string> _histName, vector<string> _cutName,
572     int nEntries, string drwOpt, bool correspond) {
573 tschum 1.1
574 gebbert 1.6 int numProcessed = 0;
575 tschum 1.1
576 gebbert 1.6 if (correspond == false) {
577 tschum 1.1
578 gebbert 1.6 for (int i=0; i<this->GetEntries(); ++i) {
579     for (vector<string>::iterator it1 =_histName.begin(); it1
580     != _histName.end(); ++it1) {
581     for (vector<string>::iterator it2 =_cutName.begin(); it2
582     != _cutName.end(); ++it2) {
583 tschum 1.1
584 gebbert 1.6 numProcessed += plot(i, *it1, *it2, nEntries, drwOpt);
585 tschum 1.8
586 tschum 1.1
587 gebbert 1.6 }
588     }
589     }
590     } else {
591 tschum 1.1
592 gebbert 1.6 if (_histName.size() != _cutName.size() )
593     return -1;
594 tschum 1.1
595 gebbert 1.6 for (int i=0; i<this->GetEntries(); ++i) {
596     for (vector<string>::iterator it1 =_histName.begin(), it2 =
597     _cutName.begin(); it1 != _cutName.end()&&it2
598     != _cutName.end(); ++it1, ++it2) {
599 tschum 1.1
600 gebbert 1.6 numProcessed += plot(i, *it1, *it2, nEntries, drwOpt);
601 tschum 1.1
602 gebbert 1.6 }
603     }
604     }
605 tschum 1.1
606 tschum 1.8 updatePads();
607 gebbert 1.6 fillCanvases();
608 tschum 1.1
609 gebbert 1.6 return numProcessed;
610 tschum 1.1 }
611    
612     //Helper methods to allow for using simple strings as input
613 gebbert 1.6 int PlotTool::loop(string histName, string cutName, int nEntries, string drwOpt) {
614 tschum 1.1
615 gebbert 1.6 vector<string> _histName;
616     _histName.push_back(histName);
617     vector<string> _cutName;
618     _cutName.push_back(cutName);
619 tschum 1.1
620 gebbert 1.6 return loop(_histName, _cutName, nEntries, drwOpt);
621 tschum 1.1
622     }
623    
624 gebbert 1.6 int PlotTool::loop(vector<string> _histName, string cutName, int nEntries,
625     string drwOpt) {
626 tschum 1.1
627 gebbert 1.6 vector<string> _cutName;
628     _cutName.push_back(cutName);
629 tschum 1.1
630 gebbert 1.6 return loop(_histName, _cutName, nEntries, drwOpt);
631 tschum 1.1
632     }
633    
634 gebbert 1.6 int PlotTool::loop(string histName, vector<string> _cutName, int nEntries,
635     string drwOpt) {
636 tschum 1.1
637 gebbert 1.6 vector<string> _histName;
638     _histName.push_back(histName);
639 tschum 1.1
640 gebbert 1.6 return loop(_histName, _cutName, nEntries, drwOpt);
641 tschum 1.1 }
642     //------------------------------------------------------------------
643     //redraw the canvas to make changes in style visible
644 gebbert 1.6 int PlotTool::updatePads() {
645 tschum 1.1
646 gebbert 1.6 for (map< string, TCanvas* >::iterator it=pads_.begin() ; it != pads_.end(); ++it) {
647 tschum 1.14 if ( ! gROOT->GetListOfCanvases()->FindObject((*it).first.c_str() ) ) {
648     pads_.erase(it);
649     continue;
650     }
651     if( (*it).second->GetListOfPrimitives()->GetEntries() == 0 ) {
652     (*it).second->Close();
653     pads_.erase(it);
654     continue;
655     }
656     (*it).second->Draw();
657     setCanvas((*it).second);
658    
659 gebbert 1.6 }
660 tschum 1.1
661 gebbert 1.6 return pads_.size();
662 tschum 1.1
663     }
664    
665     int PlotTool::fillCanvases() {
666    
667 gebbert 1.6 for (map< string, TCanvas* >::iterator it=canvases_.begin() ; it
668     != canvases_.end(); ++it) {
669 tschum 1.8 string canvName = (*it).first;
670     if (gROOT->GetListOfCanvases()->FindObject(canvName.c_str()) ) {
671 gebbert 1.6
672     int numP = (*it).second->GetNumber()+1;
673    
674     if (numP <= 0)
675     continue;
676    
677     int x = int( sqrt(numP) );
678     int y = x;
679     if (x*y < numP)
680     x += 1;
681     if (x*y < numP)
682     x += 1;
683     if (x*y < numP) {
684     x -= 1;
685     y += 1;
686     }
687    
688 tschum 1.8
689 gebbert 1.6 (*it).second->Divide(x, y);
690     int padIndex = 1;
691     for (map< string, TCanvas* >::iterator it2=pads_.begin() ; it2
692     != pads_.end(); ++it2) {
693     string padName = (*it2).first;
694 tschum 1.7 if (gROOT->GetListOfCanvases()->FindObject(padName.c_str() ) ) {
695 tschum 1.8 string n1 = canvName.substr(0,canvName.find(":"));
696     string n2("");
697     if(canvName.find("{") != string::npos ) n2 = canvName.substr(canvName.find("{"),canvName.find("}"));
698    
699     if ( (padName.find(canvName.c_str()) != string::npos ) || (padName.find(n1) != string::npos && padName.find(n2) != string::npos) || !strcmp( canvName.c_str(), "All") ) {
700 gebbert 1.6 (*it).second->cd(padIndex);
701     (*it2).second->DrawClonePad();
702     (*it2).second->Close();
703     pads_.erase(it2);
704     ++padIndex;
705     }
706     } else
707     pads_.erase(it2);
708     }
709     (*it).second->SetNumber(-1);
710 tschum 1.1
711 gebbert 1.6 } else
712     canvases_.erase(it);
713     }
714 tschum 1.1
715 gebbert 1.6 return canvases_.size();
716 tschum 1.1
717     }
718     //------------------------------------------------------------------
719     void PlotTool::setCanvas(TCanvas* thisCanvas) {
720    
721 gebbert 1.6 TH1* thisHist = 0;
722     TPaveStats* thisStatsBox = 0;
723     TPaletteAxis* palette =0;
724     TLegend* thisLeg = 0;
725     int counter =0;
726     double maxEntry=0;
727 tschum 1.1
728 gebbert 1.6 ((TFrame*) thisCanvas->GetFrame())->Delete();
729     thisCanvas->GetFrame()->SetLineWidth(gStyle->GetLineWidth() );
730 tschum 1.1
731 gebbert 1.6 thisCanvas->SetLeftMargin(0.2);
732     thisCanvas->SetRightMargin(0.06);
733     thisCanvas->SetBottomMargin(0.2);
734     thisCanvas->SetTopMargin(0.1);
735 tschum 1.1
736 gebbert 1.6 if (logY)
737     thisCanvas->SetLogy(1);
738     else
739     thisCanvas->SetLogy(0);
740 tschum 1.1
741 gebbert 1.6 if (showLegend)
742     thisLeg = new TLegend();
743 tschum 1.1
744 gebbert 1.6 for (int i = 0; i != thisCanvas->GetListOfPrimitives()->GetSize(); ++i) {
745 tschum 1.1
746 gebbert 1.6 if ( !thisCanvas->GetListOfPrimitives()->At(i)->InheritsFrom("TH1"))
747     continue;
748 tschum 1.1
749 gebbert 1.6 thisHist = ((TH1*) thisCanvas->GetListOfPrimitives()->At(i));
750     setColor(thisHist, counter);
751     setMathLabels(thisHist);
752 tschum 1.1
753 gebbert 1.6 if (thisHist->GetMaximum() > maxEntry)
754     maxEntry = thisHist->GetMaximum();
755 tschum 1.1
756 gebbert 1.6 thisStatsBox = (TPaveStats*) thisHist->GetListOfFunctions()->FindObject("stats");
757     if (thisStatsBox)
758     setStats(thisCanvas, thisStatsBox, thisHist, counter);
759 tschum 1.1
760 gebbert 1.6 palette = (TPaletteAxis*) thisHist->GetListOfFunctions()->FindObject("palette");
761     if (palette)
762     setPalette(thisCanvas, palette);
763 tschum 1.1
764 gebbert 1.6 if (thisLeg)
765     thisLeg->AddEntry(thisHist,thisHist->GetName())->SetTextSize(0.04);
766 tschum 1.1
767 gebbert 1.6 ++counter;
768 tschum 1.1
769 gebbert 1.6 }
770 tschum 1.1
771 gebbert 1.6 if (maxEntry != 0)
772     setHistMax(thisCanvas, maxEntry);
773 tschum 1.1
774 gebbert 1.6 thisCanvas->cd();
775     thisCanvas->Update();
776 tschum 1.1
777 gebbert 1.6 if (thisLeg)
778     setLegend(thisCanvas, thisLeg, counter);
779 tschum 1.1
780     }
781     //------------------------------------------------------------------
782     //private helper classes to set the canvas and hist style
783 gebbert 1.6 void PlotTool::setStats(TCanvas* thisCanvas, TPaveStats* thisStatsBox,
784     TH1* thisHist, int counter) {
785 tschum 1.1
786 gebbert 1.6 if (thisCanvas->GetRightMargin() < .2)
787     thisCanvas->SetRightMargin(.2);
788 tschum 1.1
789 gebbert 1.6 thisStatsBox->SetLineColor(thisHist->GetLineColor());
790     thisStatsBox->SetX2NDC(1.);
791     thisStatsBox->SetY2NDC(1-thisCanvas->GetTopMargin()-0.16*counter);
792     thisStatsBox->SetX1NDC(1-thisCanvas->GetRightMargin()+0.01);
793     thisStatsBox->SetY1NDC(thisStatsBox->GetY2NDC()-.15);
794 tschum 1.1
795     }
796    
797     void PlotTool::setMathLabels(TH1* thisHist) {
798    
799 gebbert 1.6 string t = thisHist ->GetTitle();
800     string x = thisHist->GetXaxis()->GetTitle();
801     string y = thisHist->GetYaxis()->GetTitle();
802    
803 tschum 1.12 // if (x.find("__") != string::npos)
804     // x = x.substr(0,x.find("__"));
805     // if (y.find("__") != string::npos)
806     // y = y.substr(0,y.find("__"));
807    
808 gebbert 1.6 if (x.find(".phi()") != string::npos)
809     x.replace(x.find(".phi()"), 6, " #phi");
810     if (y.find(".phi()") != string::npos)
811     y.replace(y.find(".phi()"), 6, " #phi");
812    
813     if (x.find(".eta()") != string::npos)
814     x.replace(x.find(".eta()"), 6, " #eta");
815     if (y.find(".eta()") != string::npos)
816     y.replace(y.find(".eta()"), 6, " #eta");
817    
818 tschum 1.12 if (x.find(".theta()") != string::npos)
819     x.replace(x.find(".theta()"), 6, " #theta");
820     if (y.find(".theta()") != string::npos)
821     y.replace(y.find(".theta()"), 6, " #theta");
822    
823     if (x.find(".chi2()") != string::npos)
824     x.replace(x.find(".chi2()"), 7, " #chi^{2}");
825     if (y.find(".chi2()") != string::npos)
826     y.replace(y.find(".chi2()"), 7, " #chi^{2}");
827    
828 gebbert 1.6 if (x.find(".pt()") != string::npos)
829 tschum 1.12 x.replace(x.find(".pt()"), 5, " p_{T} [GeV]");
830 gebbert 1.6 if (y.find(".pt()") != string::npos)
831 tschum 1.12 y.replace(y.find(".pt()"), 5, " p_{T} [GeV]");
832    
833     if (x.find(".et()") != string::npos)
834     x.replace(x.find(".et()"), 5, " E_{T} [GeV]");
835     if (y.find(".et()") != string::npos)
836     y.replace(y.find(".et()"), 5, " E_{T} [GeV]");
837    
838     //splitlines for many cuts
839     string test1= "{" + globalCuts + "}";
840     string test2= "&&" + globalCuts;
841    
842     if(t.find(test1) != string::npos) t.replace(t.find(test1),test1.length(),"");
843     if(t.find(test2) != string::npos) t.replace(t.find(test2),test2.length(),"");
844    
845     if (t.find("{") != string::npos && t.find("#splitline") == string::npos) {
846     t.replace(t.find_last_of("{"), 1, "}{");
847     string t_old = t;
848     t = "#splitline{";
849     t += t_old;
850     }
851 gebbert 1.6
852     thisHist ->SetTitle(t.c_str());
853     thisHist->GetXaxis()->SetTitle(x.c_str());
854     thisHist->GetYaxis()->SetTitle(y.c_str());
855 tschum 1.1
856     }
857    
858     void PlotTool::setColor(TH1* thisHist, int counter) {
859    
860 gebbert 1.6 if (counter == 0) {
861 tschum 1.1
862 gebbert 1.6 thisHist->SetLineColor(kRed+1);
863     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kRed+1);
864     //else
865     thisHist->SetFillStyle(0);
866     thisHist->SetMarkerColor(kRed+1);
867     } else if (counter == 1) {
868     thisHist->SetLineColor(kBlue+1);
869     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kBlue+1);
870     //else
871     thisHist->SetFillStyle(0);
872     thisHist->SetMarkerColor(kBlue+1);
873    
874     } else if (counter == 2) {
875     thisHist->SetLineColor(kGreen+2);
876     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kGreen+2);
877     //else
878     thisHist->SetFillStyle(0);
879     thisHist->SetMarkerColor(kGreen+2);
880    
881     } else if (counter == 3) {
882     thisHist->SetLineColor(kMagenta+2);
883     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kMagenta+2);
884     //else
885     thisHist->SetFillStyle(0);
886     thisHist->SetMarkerColor(kMagenta+2);
887    
888     } else if (counter == 4) {
889     thisHist->SetLineColor(kCyan+2);
890     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kCyan+2);
891     //else
892     thisHist->SetFillStyle(0);
893     thisHist->SetMarkerColor(kCyan+2);
894    
895     } else {
896     thisHist->SetLineColor(kBlack);
897     //if( gStyle->GetHistFillStyle() ) thisHist->SetFillColor (kBlack);
898     //else
899     thisHist->SetFillStyle(0);
900     thisHist->SetMarkerColor(kBlack);
901 tschum 1.1
902 gebbert 1.6 }
903 tschum 1.1
904     }
905    
906     void PlotTool::setPalette(TCanvas* thisCanvas, TPaletteAxis* palette) {
907    
908 gebbert 1.6 palette->SetLabelSize(0.045);
909     if (thisCanvas->GetRightMargin() < .15)
910     thisCanvas->SetRightMargin(.15);
911     palette->SetX1NDC(1-thisCanvas->GetRightMargin()+0.01);
912     palette->SetY1NDC(thisCanvas->GetBottomMargin());
913     palette->SetX2NDC(palette->GetX1NDC()+0.05);
914     palette->SetY2NDC(1-thisCanvas->GetTopMargin());
915 tschum 1.1
916     }
917    
918     void PlotTool::setHistMax(TCanvas* thisCanvas, double maxEntry) {
919    
920 gebbert 1.6 TH1* thisHist = 0;
921 tschum 1.1
922 gebbert 1.6 for (int i = 0; i != thisCanvas->GetListOfPrimitives()->GetSize(); ++i) {
923 tschum 1.1
924 gebbert 1.6 if ( !thisCanvas->GetListOfPrimitives()->At(i)->InheritsFrom("TH1") )
925     continue;
926     thisHist = ((TH1*) thisCanvas->GetListOfPrimitives()->At(i));
927     if (thisHist->GetMaximum() < maxEntry) {
928     double minEntry=thisHist->GetBinContent(thisHist->GetMinimumBin());
929     if (logY) {
930     int bin = (thisHist->GetMinimumBin()+1);
931     minEntry=thisHist->GetBinContent(bin);
932     }
933     thisHist->GetYaxis()->SetRangeUser(minEntry, maxEntry*1.2);
934     }
935     break;
936 tschum 1.1
937 gebbert 1.6 }
938 tschum 1.1
939     }
940    
941     void PlotTool::setLegend(TCanvas* thisCanvas, TLegend* thisLeg, int counter) {
942    
943 gebbert 1.6 thisLeg->SetFillStyle(0);
944     thisLeg->SetX1NDC(1-thisCanvas->GetRightMargin()-0.5);
945     thisLeg->SetY1NDC(1-thisCanvas->GetTopMargin()-0.05*counter);
946     thisLeg->SetX2NDC(1-thisCanvas->GetRightMargin());
947     thisLeg->SetY2NDC(1-thisCanvas->GetTopMargin());
948     thisLeg->SetEntrySeparation(0.5);
949     thisLeg->Draw("NDC");
950 tschum 1.1
951     }
952    
953 gebbert 1.6 void PlotTool::createColors() {
954     const Int_t NRGBs = 5;
955     const Int_t NCont = 255;
956    
957     Double_t stops[NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
958     Double_t red[NRGBs] = { 0.00, 0.00, 0.87, 1.00, 0.51 };
959     Double_t green[NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 };
960     Double_t blue[NRGBs] = { 0.51, 1.00, 0.12, 0.00, 0.00 };
961     TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
962     gStyle->SetNumberContours(NCont);
963 tschum 1.1
964     }
965    
966     //------------------------------------------------------------------
967    
968 tschum 1.9 int PlotTool::saveCanvases(string name, string path) {
969 tschum 1.1
970 gebbert 1.6 TSystemDirectory d("", path.c_str());
971     if (!d.GetListOfFiles())
972     return -1;
973    
974 tschum 1.9 string namePs = path +"/" + name;
975     string nameRt = path +"/" + name;
976    
977     if (name.find(".ps") == string::npos)
978     namePs += ".ps";
979    
980     if (name.find(".root") == string::npos)
981     nameRt += ".root";
982    
983     int savedCanvs =0;
984    
985 tschum 1.12
986     TIter next(gROOT->GetListOfCanvases());
987     TCanvas* canv;
988     while ((canv=(TCanvas*)next())) {
989     string modName = (*canv).GetName();
990     if(modName.find("{") != string::npos) modName = modName.substr(0,modName.find("{"));
991     (*canv).SetTitle( modName.c_str() );
992    
993     }
994    
995    
996 tschum 1.9 TPostScript ps(namePs.c_str(),112);
997     TFile rt(nameRt.c_str(),"recreate");
998 gebbert 1.6
999 tschum 1.12 next.Reset();
1000 gebbert 1.6 while ((canv=(TCanvas*)next())) {
1001 tschum 1.9 ps.NewPage();
1002 tschum 1.12
1003    
1004     (*canv).Write();
1005 tschum 1.9 (*canv).Draw();
1006     (*canv).Update();
1007 tschum 1.12
1008    
1009 tschum 1.9 ++savedCanvs;
1010 gebbert 1.6
1011 tschum 1.9 }
1012 gebbert 1.6
1013 tschum 1.9 ps.Close();
1014     rt.Close();
1015 tschum 1.1
1016 tschum 1.12
1017     if(verbose && savedCanvs) {
1018     cout<<"Saved file "<<rt.GetName()<<" with "<<savedCanvs<<" canvases."<<endl;
1019     cout<<"Saved file "<<ps.GetName()<<" with "<<savedCanvs<<" pages."<<endl;
1020     }
1021 tschum 1.9 return savedCanvs;
1022 tschum 1.1
1023     }
1024     //------------------------------------------------------------------
1025 tschum 1.11
1026     int PlotTool::clearCanvases()
1027 tschum 1.7 {
1028    
1029 tschum 1.11 while(gROOT->GetListOfCanvases()->GetEntries()) ((TCanvas*) gROOT->GetListOfCanvases()->At(0))->Close();
1030     pads_.clear();
1031     autoVars.clear();
1032     return pads_.size() + gROOT->GetListOfCanvases()->GetEntries();
1033    
1034     }
1035     //------------------------------------------------------------------
1036     int PlotTool::setVariables(string label)
1037     {
1038 tschum 1.7
1039     for (int i=0; i< this->GetEntries(); ++i) {
1040 tschum 1.11 cout<<"--------------------------------"<<endl;
1041 tschum 1.7 cout<<((TChain*) this->At(i))->GetName()<<endl;
1042 tschum 1.11 cout<<"------------"<<endl;
1043    
1044     TList* leaves = (TList*) ((TChain*) this->At(i))->GetListOfLeaves();
1045     for (int j=0; j< leaves->GetEntries(); ++j) {
1046     TString leafName ( leaves->At(j)->GetName() );
1047     if(! leafName.EndsWith(".") ) continue;
1048     if(! leafName.Contains(label.c_str() ) ) continue;
1049    
1050     TClass cl( ( (TLeafElement*) leaves->At(j) )->GetTypeName() );
1051     cout<<"++++++"<<endl;
1052     cout<< leafName.Data() <<endl;
1053     for(int k=0;k<cl.GetListOfAllPublicMethods()->GetEntries();++k) {
1054     string typeName ( ((TMethod*) cl.GetListOfAllPublicMethods()->At(k))->GetReturnTypeName() );
1055     string methName ( cl.GetListOfAllPublicMethods()->At(k)->GetName() );
1056     if( methName != "product") continue;
1057     cout<< typeName <<endl;
1058     string _type;
1059     TString testString( typeName.c_str() );
1060     if( testString.BeginsWith("vector<") ) {
1061     _type = typeName.substr(typeName.find("<")+1,typeName.find(">")-typeName.find("<")-1);
1062     string _varSize = leafName.Data();
1063     _varSize += "obj@.size()";
1064     autoVars.push_back( _varSize );
1065     }
1066     else _type = typeName.substr(0,typeName.find("*"));
1067     TClass _cl( _type.c_str() );
1068     for(int l=0;l<_cl.GetListOfAllPublicMethods()->GetEntries();++l) {
1069     string _typeName ( ((TMethod*) _cl.GetListOfAllPublicMethods()->At(l))->GetReturnTypeName() );
1070     string _methName ( _cl.GetListOfAllPublicMethods()->At(l)->GetName() );
1071     // if(_typeName.find("void") != string::npos ) continue;
1072     // cout<<" "<<_typeName<<" "<<_methName<<endl;
1073    
1074     cout<<" "<<_typeName<<" "<<_methName<<endl;
1075     if(_methName.find("operator")==string::npos&&(_typeName=="float"||_typeName=="double"||_typeName=="int"||_typeName=="unsigned int"||_typeName=="bool"||_typeName=="unsigned short"||_typeName=="unsigned long"||_typeName=="unsigned long long")) {
1076    
1077     cout<<"--> "<<_typeName<<" "<<_methName<<endl;
1078    
1079     string _varName = leafName.Data();
1080     _varName += "obj.";
1081     _varName += _methName;
1082     _varName += "()";
1083     autoVars.push_back( _varName );
1084     }
1085     }
1086    
1087    
1088     }
1089    
1090     }
1091 tschum 1.7
1092     }
1093    
1094 tschum 1.11 return autoVars.size();
1095 tschum 1.7
1096     }
1097 tschum 1.8 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1098     //Draw efficiencies
1099     int PlotTool::plotEff(int chainIndex, string histName, string cutName, int nEntries, double fitXmin, double fitXmax, string fitFormula)
1100     {
1101    
1102     if( chainIndex >= this->GetEntries() ) return -1;
1103    
1104     int currN = nEntries;
1105     if(nEntries < 0) currN = ((TChain*) this->At(chainIndex))->GetEntries(); //nEntries<0 : all entries are plotted!
1106    
1107     //++++ Create and name Canvases according to global variables +++++++++++++
1108     ostringstream currHistName;
1109     if( samePad_trees) currHistName<<((TChain*) this->At(chainIndex))->GetName()<<":";
1110     if( samePad_vars) currHistName<<histName;
1111     if( samePad_cuts) currHistName<<"{"<<cutName<<"}";
1112    
1113     ostringstream currPadName;
1114     if(! samePad_trees) currPadName<<((TChain*) this->At(chainIndex))->GetName()<<":";
1115     if(! samePad_vars) currPadName<<histName;
1116     if(! samePad_cuts) currPadName<<"{"<<cutName<<"}";
1117    
1118    
1119     //Draw total histogram:
1120     if( ! pads_["total"] || ! gROOT->FindObject("total") ) {
1121     pads_["total"] = new TCanvas("total", "total");
1122     } else {
1123     pads_["total"]->cd();
1124     }
1125     ostringstream bins_total;
1126     bins_total<<histName;
1127     bins_total<<">>total(100,0,1000)";
1128    
1129    
1130     int draw_err = ((TChain*) this->At(chainIndex))->Draw(bins_total.str().c_str(), "", "", currN);
1131     if( draw_err <= 0 ) return draw_err;
1132    
1133    
1134     TH1* total = ((TH1*) pads_["total"]->GetPrimitive("total"));
1135    
1136     ostringstream bins_bkg;
1137     bins_bkg<<histName;
1138     bins_bkg<<">>bkg(";
1139     bins_bkg<<total->GetNbinsX();
1140     bins_bkg<<",";
1141     bins_bkg<<total->GetXaxis()->GetXmin();
1142     bins_bkg<<",";
1143     bins_bkg<<total->GetXaxis()->GetXmax();
1144     bins_bkg<<")";
1145    
1146     //Draw bkg histogram:
1147     if( ! pads_["bkg"] || ! gROOT->FindObject("bkg") ) {
1148     pads_["bkg"] = new TCanvas("bkg", currPadName.str().c_str());
1149     } else {
1150     pads_["bkg"]->cd();
1151     }
1152    
1153     draw_err = ((TChain*) this->At(chainIndex))->Draw(bins_bkg.str().c_str(), cutName.c_str(),"" , currN);
1154     if( draw_err <= 0 ) return draw_err;
1155    
1156    
1157     TH1* bkg = ((TH1*) pads_["bkg"]->GetPrimitive("bkg"));
1158    
1159     //Draw pass histogram:
1160     ostringstream bins_pass;
1161     bins_pass<<histName;
1162     bins_pass<<">>pass(";
1163     bins_pass<<total->GetNbinsX();
1164     bins_pass<<",";
1165     bins_pass<<total->GetXaxis()->GetXmin();
1166     bins_pass<<",";
1167     bins_pass<<total->GetXaxis()->GetXmax();
1168     bins_pass<<")";
1169    
1170     ostringstream cut_pass;
1171     cut_pass<<"!(";
1172     cut_pass<<cutName;
1173     cut_pass<<")";
1174    
1175    
1176     if( ! pads_["pass"] || ! gROOT->FindObject("pass") ) {
1177     pads_["pass"] = new TCanvas("pass", currPadName.str().c_str());
1178     } else {
1179     pads_["pass"]->cd();
1180     }
1181    
1182     draw_err = ((TChain*) this->At(chainIndex))->Draw(bins_pass.str().c_str(),cut_pass.str().c_str() ,"" , currN);
1183     if( draw_err <= 0 ) return draw_err;
1184    
1185    
1186     TH1* pass = ((TH1*) pads_["pass"]->GetPrimitive("pass"));
1187    
1188    
1189     currPadName<<"Eff";
1190     //Draw Efficiency Graph:
1191     if( ! pads_["EffGraph"] || ! gROOT->FindObject("EffGraph") ) {
1192     pads_["EffGraph"] = new TCanvas("EffGraph", currPadName.str().c_str());
1193     } else {
1194     pads_["EffGraph"]->cd();
1195     }
1196    
1197     TGraphAsymmErrors* EffGraph = new TGraphAsymmErrors(bkg, total);
1198     EffGraph->SetName(currHistName.str().c_str());
1199    
1200     TF1* reverse = new TF2("reverse","1/y-1",total->GetXaxis()->GetXmin(),total->GetXaxis()->GetXmax());
1201     EffGraph->Apply(reverse);
1202     EffGraph->Draw("A*");
1203    
1204     TF1* fitfunc = new TF1("fitfunc",fitFormula.c_str(),fitXmin,fitXmax);
1205     EffGraph->Fit("fitfunc","R+");
1206    
1207    
1208     //Save fit function
1209    
1210     ostringstream savefuncName;
1211     savefuncName<<histName;
1212     savefuncName<<"_";
1213    
1214     if(cutName.find("<") != std::string::npos ) cutName.replace(cutName.find("<"),1,"_");
1215     if(cutName.find(">") != std::string::npos ) cutName.replace(cutName.find(">"),1,"_");
1216     savefuncName<<cutName;
1217    
1218    
1219     TFile file("ABCDFunctions.root","UPDATE");
1220     TF1* savefunc = new TF1(savefuncName.str().c_str(), fitfunc->GetExpFormula().Data(),0,10000);
1221     savefunc->SetParameters( fitfunc->GetParameters() );
1222     savefunc->SetParErrors ( fitfunc->GetParErrors() );
1223    
1224     //Show results
1225     if( ! pads_["abcd"] || ! gROOT->FindObject("abcd") ) {
1226     pads_["abcd"] = new TCanvas("abcd", "abcd");
1227     } else {
1228     pads_["abcd"]->cd();
1229     }
1230    
1231     bkg->Multiply(savefunc);
1232     bkg->Draw();
1233    
1234     // total->Add(bkg,-1);
1235     pass->SetLineColor(kRed);
1236     pass->Draw("sames,e");
1237     pads_["abcd"]->SetLogy();
1238    
1239    
1240     savefunc->Write();
1241     // file.Close();
1242    
1243    
1244     return draw_err;
1245    
1246     }
1247     //------------------------------------------------------------------
1248 thomsen 1.16
1249    
1250     //Eta-Phi to Det-komponent matching
1251    
1252     int PlotTool::getEtaBin(double eta){
1253     for(int i=-41;i<42;i++){
1254     double low = getEtaFromBin(i,true);
1255     double high = getEtaFromBin(i,false);
1256     if(eta >= low && eta <= high) return i;
1257     }
1258     cout<<"eta-bin not found"<<endl;
1259     return 0;
1260     }
1261    
1262     double PlotTool::getEtaFromBin(int etaBin, bool lowerEdge){
1263     // return eta bin - eta edge mappting
1264     switch(etaBin){
1265     case -41: return (lowerEdge ? -5.191 : -4.889); break;
1266     case -40: return (lowerEdge ? -4.889 : -4.716); break;
1267     case -39: return (lowerEdge ? -4.716 : -4.538); break;
1268     case -38: return (lowerEdge ? -4.538 : -4.363); break;
1269     case -37: return (lowerEdge ? -4.363 : -4.191); break;
1270     case -36: return (lowerEdge ? -4.191 : -4.013); break;
1271     case -35: return (lowerEdge ? -4.013 : -3.839); break;
1272     case -34: return (lowerEdge ? -3.839 : -3.664); break;
1273     case -33: return (lowerEdge ? -3.664 : -3.489); break;
1274     case -32: return (lowerEdge ? -3.489 : -3.314); break;
1275     case -31: return (lowerEdge ? -3.314 : -3.139); break;
1276     case -30: return (lowerEdge ? -3.139 : -2.964); break;
1277     case -29: return (lowerEdge ? -2.964 : -2.853); break;
1278     case -28: return (lowerEdge ? -2.853 : -2.65); break;
1279     case -27: return (lowerEdge ? -2.65 : -2.5); break;
1280     case -26: return (lowerEdge ? -2.5 : -2.322); break;
1281     case -25: return (lowerEdge ? -2.322 : -2.172); break;
1282     case -24: return (lowerEdge ? -2.172 : -2.043); break;
1283     case -23: return (lowerEdge ? -2.043 : -1.93); break;
1284     case -22: return (lowerEdge ? -1.93 : -1.83); break;
1285     case -21: return (lowerEdge ? -1.83 : -1.74); break;
1286     case -20: return (lowerEdge ? -1.74 : -1.653); break;
1287     case -19: return (lowerEdge ? -1.653 : -1.566); break;
1288     case -18: return (lowerEdge ? -1.566 : -1.479); break;
1289     case -17: return (lowerEdge ? -1.479 : -1.392); break;
1290     case -16: return (lowerEdge ? -1.392 : -1.305); break;
1291     case -15: return (lowerEdge ? -1.305 : -1.218); break;
1292     case -14: return (lowerEdge ? -1.218 : -1.131); break;
1293     case -13: return (lowerEdge ? -1.131 : -1.044); break;
1294     case -12: return (lowerEdge ? -1.044 : -0.957); break;
1295     case -11: return (lowerEdge ? -0.957 : -0.879); break;
1296     case -10: return (lowerEdge ? -0.879 : -0.783); break;
1297     case -9: return (lowerEdge ? -0.783 : -0.696); break;
1298     case -8: return (lowerEdge ? -0.696 : -0.609); break;
1299     case -7: return (lowerEdge ? -0.609 : -0.522); break;
1300     case -6: return (lowerEdge ? -0.522 : -0.435); break;
1301     case -5: return (lowerEdge ? -0.435 : -0.348); break;
1302     case -4: return (lowerEdge ? -0.348 : -0.261); break;
1303     case -3: return (lowerEdge ? -0.261 : -0.174); break;
1304     case -2: return (lowerEdge ? -0.174 : -0.087); break;
1305     case -1: return (lowerEdge ? -0.087 : 0); break;
1306     case +1: return (lowerEdge ? 0 : 0.087); break;
1307     case +2: return (lowerEdge ? 0.087 : 0.174); break;
1308     case +3: return (lowerEdge ? 0.174 : 0.261); break;
1309     case +4: return (lowerEdge ? 0.261 : 0.348); break;
1310     case +5: return (lowerEdge ? 0.348 : 0.435); break;
1311     case +6: return (lowerEdge ? 0.435 : 0.522); break;
1312     case +7: return (lowerEdge ? 0.522 : 0.609); break;
1313     case +8: return (lowerEdge ? 0.609 : 0.696); break;
1314     case +9: return (lowerEdge ? 0.696 : 0.783); break;
1315     case +10: return (lowerEdge ? 0.783 : 0.879); break;
1316     case +11: return (lowerEdge ? 0.879 : 0.957); break;
1317     case +12: return (lowerEdge ? 0.957 : 1.044); break;
1318     case +13: return (lowerEdge ? 1.044 : 1.131); break;
1319     case +14: return (lowerEdge ? 1.131 : 1.218); break;
1320     case +15: return (lowerEdge ? 1.218 : 1.305); break;
1321     case +16: return (lowerEdge ? 1.305 : 1.392); break;
1322     case +17: return (lowerEdge ? 1.392 : 1.479); break;
1323     case +18: return (lowerEdge ? 1.479 : 1.566); break;
1324     case +19: return (lowerEdge ? 1.566 : 1.653); break;
1325     case +20: return (lowerEdge ? 1.653 : 1.74); break;
1326     case +21: return (lowerEdge ? 1.74 : 1.83); break;
1327     case +22: return (lowerEdge ? 1.83 : 1.93); break;
1328     case +23: return (lowerEdge ? 1.93 : 2.043); break;
1329     case +24: return (lowerEdge ? 2.043 : 2.172); break;
1330     case +25: return (lowerEdge ? 2.172 : 2.322); break;
1331     case +26: return (lowerEdge ? 2.322 : 2.5); break;
1332     case +27: return (lowerEdge ? 2.5 : 2.65); break;
1333     case +28: return (lowerEdge ? 2.65 : 2.853); break;
1334     case +29: return (lowerEdge ? 2.853 : 2.964); break;
1335     case +30: return (lowerEdge ? 2.964 : 3.139); break;
1336     case +31: return (lowerEdge ? 3.139 : 3.314); break;
1337     case +32: return (lowerEdge ? 3.314 : 3.489); break;
1338     case +33: return (lowerEdge ? 3.489 : 3.664); break;
1339     case +34: return (lowerEdge ? 3.664 : 3.839); break;
1340     case +35: return (lowerEdge ? 3.839 : 4.013); break;
1341     case +36: return (lowerEdge ? 4.013 : 4.191); break;
1342     case +37: return (lowerEdge ? 4.191 : 4.363); break;
1343     case +38: return (lowerEdge ? 4.363 : 4.538); break;
1344     case +39: return (lowerEdge ? 4.538 : 4.716); break;
1345     case +40: return (lowerEdge ? 4.716 : 4.889); break;
1346     case +41: return (lowerEdge ? 4.889 : 5.191); break;
1347     //something went wrong;
1348     default : return -1; break;
1349     }
1350     }
1351    
1352     int PlotTool::getPhiBin(double phi){
1353    
1354     return 0;
1355     }
1356    
1357     double PlotTool::getPhiFromBin(int phiBin, bool lowerEdge){
1358    
1359     return 0;
1360     }