ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/tschum/FWlite_Analysis/PlotTool.cc
Revision: 1.17
Committed: Mon Feb 8 16:51:26 2010 UTC (15 years, 2 months ago) by thomsen
Content type: text/plain
Branch: MAIN
Changes since 1.16: +72 -2 lines
Log Message:
JetContent aka Peters Plot

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