ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/Morgan/src/TotoAnalyzer.cc
Revision: 1.37
Committed: Mon Oct 12 08:32:36 2009 UTC (15 years, 6 months ago) by lethuill
Content type: text/plain
Branch: MAIN
Changes since 1.36: +533 -513 lines
Log Message:
Add private primary vertex reconstruction for Zee events

File Contents

# User Rev Content
1 lethuill 1.14 #include "../interface/TotoAnalyzer.h"
2 mlethuil 1.1
3     using namespace std;
4     using namespace reco;
5     using namespace edm;
6    
7    
8     TotoAnalyzer::TotoAnalyzer(const edm::ParameterSet& iConfig)
9     {
10 lethuill 1.37 myConfig_ = iConfig.getParameter<ParameterSet>("myConfig");
11     dataType_ = myConfig_.getUntrackedParameter<string>("dataType","unknown");
12     cout << "dataType: " << dataType_ << endl;
13     if( dataType_=="RECO" ) producersNames_ = iConfig.getParameter<ParameterSet>("producersNamesRECO");
14     else if( dataType_=="PAT" ) producersNames_ = iConfig.getParameter<ParameterSet>("producersNamesPAT");
15     else { cout << "TotoAnalyzer::TotoAnalyzer... dataType is unknown... exiting..." << endl; exit(1); }
16    
17     runInfos_ = 0;
18     rootEvent_ = 0;
19     hltAnalyzer_ = 0;
20     rootMCParticles_ = 0;
21     rootGenJets_ = 0;
22     rootGenMETs_ = 0;
23     rootMuMuGammaEvent_ = 0;
24     rootMCTopTop_ = 0;
25     rootMCPhotons_ = 0;
26     rootBeamSpot_ = 0;
27     rootVertices_ = 0;
28     rootZeeVertices_ = 0;
29     rootTracks_ = 0;
30     rootJets_ = 0;
31     rootMuons_ = 0;
32     rootElectrons_ = 0;
33     rootPhotons_ = 0;
34     rootBasicClusters_ = 0;
35     rootSuperClusters_ = 0;
36     rootConversionTracks_ = 0;
37     rootMETs_ = 0;
38 mlethuil 1.1 }
39    
40    
41     TotoAnalyzer::~TotoAnalyzer()
42 lethuill 1.11 {
43     }
44 mlethuil 1.1
45    
46     // ------------ method called once each job just before starting event loop ------------
47     void TotoAnalyzer::beginJob(const edm::EventSetup&)
48     {
49 lethuill 1.37
50     // Load Config parameters
51     verbosity_ = myConfig_.getUntrackedParameter<int>("verbosity", 0);
52     allowMissingCollection_ = producersNames_.getUntrackedParameter<bool>("allowMissingCollection", false);
53     rootFileName_ = myConfig_.getUntrackedParameter<string>("RootFileName","noname.root");
54     datasetXsection_ = myConfig_.getUntrackedParameter<double>("xsection");
55     datasetDesciption_ = myConfig_.getUntrackedParameter<string>("description","Pas de description");
56     doHLT_ = myConfig_.getUntrackedParameter<bool>("doHLT",false);
57     doMC_ = myConfig_.getUntrackedParameter<bool>("doMC",false);
58     doJetMC_ = myConfig_.getUntrackedParameter<bool>("doJetMC",false);
59     doMETMC_ = myConfig_.getUntrackedParameter<bool>("doMETMC",false);
60     doPDFInfo_ = myConfig_.getUntrackedParameter<bool>("doPDFInfo",false);
61     doSignalMuMuGamma_ = myConfig_.getUntrackedParameter<bool>("doSignalMuMuGamma",false);
62     doSignalTopTop_ = myConfig_.getUntrackedParameter<bool>("doSignalTopTop",false);
63     doPhotonConversionMC_ = myConfig_.getUntrackedParameter<bool>("doPhotonConversionMC",false);
64     drawMCTree_ = myConfig_.getUntrackedParameter<bool>("drawMCTree",false);
65     doBeamSpot_ = myConfig_.getUntrackedParameter<bool>("doBeamSpot",false);
66     doPrimaryVertex_ = myConfig_.getUntrackedParameter<bool>("doPrimaryVertex",false);
67     doZeePrimaryVertex_ = myConfig_.getUntrackedParameter<bool>("doZeePrimaryVertex",false);
68     doTrack_ = myConfig_.getUntrackedParameter<bool>("doTrack",false);
69     doJet_ = myConfig_.getUntrackedParameter<bool>("doJet",false);
70     doMuon_ = myConfig_.getUntrackedParameter<bool>("doMuon",false);
71     doElectron_ = myConfig_.getUntrackedParameter<bool>("doElectron",false);
72     doPhoton_ = myConfig_.getUntrackedParameter<bool>("doPhoton",false);
73     doCluster_ = myConfig_.getUntrackedParameter<bool>("doCluster",false);
74     doPhotonConversion_ = myConfig_.getUntrackedParameter<bool>("doPhotonConversion",false);
75     doPhotonIsolation_ = myConfig_.getUntrackedParameter<bool>("doPhotonIsolation",false);
76     doMET_ = myConfig_.getUntrackedParameter<bool>("doMET",false);
77    
78     nTotEvt_ = 0;
79    
80     // initialize root output file
81     if(verbosity_>0) cout << "New RootFile " << rootFileName_.c_str() << " is created" << endl;
82     rootFile_ = new TFile(rootFileName_.c_str(), "recreate");
83     rootFile_->cd();
84    
85     runInfos_ = new TRootRun();
86     runTree_ = new TTree("runTree", "Global Run Infos");
87     runTree_->Branch ("runInfos", "TRootRun", &runInfos_);
88     runInfos_->setXsection(datasetXsection_);
89     runInfos_->setDescription(datasetDesciption_);
90    
91     rootEvent_ = 0;
92     eventTree_ = new TTree("eventTree", "Event Infos");
93     eventTree_->Branch ("Event", "TRootEvent", &rootEvent_);
94    
95     if(doHLT_)
96     {
97     hltAnalyzer_ = new HLTAnalyzer(producersNames_, verbosity_);
98     }
99    
100     if(doMC_)
101     {
102     if(verbosity_>0) cout << "MC Particles info will be added to rootuple" << endl;
103     rootMCParticles_ = new TClonesArray("TRootMCParticle", 1000);
104     eventTree_->Branch ("MCParticles", "TClonesArray", &rootMCParticles_);
105     }
106    
107     if(doJetMC_)
108     {
109     if(verbosity_>0) cout << "genJets info will be added to rootuple" << endl;
110     rootGenJets_ = new TClonesArray("TRootParticle", 1000);
111     eventTree_->Branch ("genJets", "TClonesArray", &rootGenJets_);
112     }
113    
114     if(doMETMC_)
115     {
116     if(verbosity_>0) cout << "genMETs info will be added to rootuple" << endl;
117     rootGenMETs_ = new TClonesArray("TRootParticle", 1000);
118     eventTree_->Branch ("genMETs", "TClonesArray", &rootGenMETs_);
119     }
120    
121     if(doSignalMuMuGamma_)
122     {
123     if(verbosity_>0) cout << "MC info for Z -> mu mu gamma will be added to rootuple" << endl;
124     eventTree_->Branch ("MuMuGamma", "TRootSignalEvent", &rootMuMuGammaEvent_);
125     }
126    
127     if(doSignalTopTop_)
128     {
129     if(verbosity_>0) cout << "MC info for Top Top will be added to rootuple" << endl;
130     rootMCTopTop_ = new TClonesArray("TRootTopTop", 1000);
131     eventTree_->Branch ("rootMCTopTop", "TClonesArray", &rootMCTopTop_);
132     }
133    
134     if(doPhotonConversionMC_)
135     {
136     if(verbosity_>0) cout << "Converted MC Photons info will be added to rootuple" << endl;
137     rootMCPhotons_ = new TClonesArray("TRootMCPhoton", 1000);
138     eventTree_->Branch ("MCPhotons", "TClonesArray", &rootMCPhotons_);
139     }
140    
141     if(doBeamSpot_)
142     {
143     if(verbosity_>0) cout << "BeamSpot info will be added to rootuple" << endl;
144     rootBeamSpot_ = new TRootBeamSpot();
145     eventTree_->Branch ("BeamSpot", "TRootBeamSpot", &rootBeamSpot_);
146     }
147    
148     if(doPrimaryVertex_)
149     {
150     if(verbosity_>0) cout << "Vertices info will be added to rootuple" << endl;
151     rootVertices_ = new TClonesArray("TRootVertex", 1000);
152     eventTree_->Branch ("Vertices", "TClonesArray", &rootVertices_);
153     }
154    
155     if(doZeePrimaryVertex_)
156     {
157     if(verbosity_>0) cout << "Zee Primary Vertices info will be added to rootuple" << endl;
158     rootZeeVertices_ = new TClonesArray("TRootVertex", 1000);
159     eventTree_->Branch ("ZeeVertices", "TClonesArray", &rootZeeVertices_);
160     }
161    
162     if(doTrack_)
163     {
164     if(verbosity_>0) cout << "Tracks info will be added to rootuple" << endl;
165     rootTracks_ = new TClonesArray("TRootTrack", 1000);
166     eventTree_->Branch ("Tracks", "TClonesArray", &rootTracks_);
167     }
168    
169     if(doJet_)
170     {
171     if(verbosity_>0) cout << "Jets info will be added to rootuple" << endl;
172     rootJets_ = new TClonesArray("TRootJet", 1000);
173     eventTree_->Branch ("Jets", "TClonesArray", &rootJets_);
174     }
175    
176     if(doMuon_)
177     {
178     if(verbosity_>0) cout << "Muons info will be added to rootuple" << endl;
179     rootMuons_ = new TClonesArray("TRootMuon", 1000);
180     eventTree_->Branch ("Muons", "TClonesArray", &rootMuons_);
181     }
182    
183     if(doElectron_)
184     {
185     if(verbosity_>0) cout << "Electrons info will be added to rootuple" << endl;
186     rootElectrons_ = new TClonesArray("TRootElectron", 1000);
187     eventTree_->Branch ("Electrons", "TClonesArray", &rootElectrons_);
188     }
189    
190     if(doPhoton_)
191     {
192     if(verbosity_>0) cout << "Photons info will be added to rootuple" << endl;
193     rootPhotons_ = new TClonesArray("TRootPhoton", 1000);
194     eventTree_->Branch ("Photons", "TClonesArray", &rootPhotons_);
195     }
196    
197     if(doCluster_)
198     {
199     if(verbosity_>0) cout << "ECAL Clusters info will be added to rootuple" << endl;
200     rootBasicClusters_ = new TClonesArray("TRootCluster", 1000);
201     eventTree_->Branch ("BasicClusters", "TClonesArray", &rootBasicClusters_);
202     rootSuperClusters_ = new TClonesArray("TRootSuperCluster", 1000);
203     eventTree_->Branch ("SuperClusters", "TClonesArray", &rootSuperClusters_);
204     }
205    
206     if(doPhotonConversion_)
207     {
208     if(verbosity_>0) cout << "Conversion Tracks info will be added to rootuple" << endl;
209     rootConversionTracks_ = new TClonesArray("TRootTrack", 1000);
210     eventTree_->Branch ("ConversionTracks", "TClonesArray", &rootConversionTracks_);
211    
212     conversionLikelihoodCalculator_ = new ConversionLikelihoodCalculator();
213     std::string weightsString = myConfig_.getUntrackedParameter<string>("conversionLikelihoodWeightsFile","RecoEgamma/EgammaTools/data/TMVAnalysis_Likelihood.weights.txt");
214     edm::FileInPath weightsFile(weightsString.c_str() );
215     conversionLikelihoodCalculator_->setWeightsFile(weightsFile.fullPath().c_str());
216     }
217    
218     if(doMET_)
219     {
220     if(verbosity_>0) cout << "MET info will be added to rootuple" << endl;
221     rootMETs_ = new TClonesArray("TRootMET", 1000);
222     eventTree_->Branch ("MET", "TClonesArray", &rootMETs_);
223     }
224    
225 mlethuil 1.1 }
226    
227    
228     // ------------ method called once each job just after ending the event loop ------------
229     void TotoAnalyzer::endJob()
230     {
231 lethuill 1.37 // Trigger Summary Tables
232     if(doHLT_)
233     {
234     cout << "Trigger Summary Tables" << endl;
235     hltAnalyzer_ ->copySummary(runInfos_);
236     hltAnalyzer_->printStats();
237     }
238    
239     runTree_->Fill();
240    
241     std::cout << "Total number of events: " << nTotEvt_ << std::endl;
242     std::cout << "Closing rootfile " << rootFile_->GetName() << std::endl;
243     rootFile_->Write();
244     rootFile_->Close();
245    
246     if(doPhotonConversion_) delete conversionLikelihoodCalculator_;
247     if(doHLT_) delete hltAnalyzer_;
248 mlethuil 1.1 }
249    
250    
251    
252     // ------------ method called to for each event ------------
253     void TotoAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
254     {
255 lethuill 1.37
256     rootFile_->cd();
257     nTotEvt_++;
258     if( (verbosity_>1) || (verbosity_>0 && nTotEvt_%10==0 && nTotEvt_<=100) || (verbosity_>0 && nTotEvt_%100==0 && nTotEvt_>100) )
259     cout << endl << endl
260     << "####### TotoAnalyzer - Cumulated Events " << nTotEvt_
261     << " - Run " << iEvent.id().run()
262     << " - Event " << iEvent.id().event()
263     << " #######" << endl;
264    
265    
266     // Global Event Infos
267     rootEvent_ = new TRootEvent();
268     rootEvent_->setNb(nTotEvt_);
269     rootEvent_->setEventId(iEvent.id().event());
270     rootEvent_->setRunId(iEvent.id().run());
271     cout << *rootEvent_ << endl;
272    
273     // Trigger
274     rootEvent_->setGlobalHLT(true);
275     if(doHLT_)
276     {
277     if(verbosity_>1) cout << endl << "Get TriggerResults..." << endl;
278     if (nTotEvt_==1) hltAnalyzer_->init(iEvent, rootEvent_);
279     hltAnalyzer_->process(iEvent, rootEvent_);
280     }
281     //if ( ! rootEvent_->passGlobalHLT() ) return;
282    
283    
284     // MC Info
285     if(doMC_)
286     {
287     if(verbosity_>1) cout << endl << "Analysing MC info..." << endl;
288     MCAnalyzer* myMCAnalyzer = new MCAnalyzer(myConfig_, producersNames_);
289     myMCAnalyzer->setVerbosity(verbosity_);
290     if (drawMCTree_) myMCAnalyzer->drawMCTree(iEvent, iSetup, myConfig_, producersNames_);
291     if ( doPDFInfo_ ) myMCAnalyzer->pdfInfo(iEvent, rootEvent_);
292     myMCAnalyzer->processMCParticle(iEvent, rootMCParticles_);
293     if(doJetMC_) myMCAnalyzer->processGenJets(iEvent, rootGenJets_);
294     if(doMETMC_) myMCAnalyzer->processGenMETs(iEvent, rootGenMETs_);
295     if (doSignalMuMuGamma_)
296     {
297     rootMuMuGammaEvent_ = new TRootSignalEvent();
298     myMCAnalyzer->processMuMuGammaEvent(iEvent, rootMuMuGammaEvent_);
299     }
300     if (doSignalTopTop_) myMCAnalyzer->processTopTopEvent(iEvent, rootMCTopTop_);
301     if (doPhotonConversionMC_) myMCAnalyzer->processConvertedPhoton(iEvent, rootMCPhotons_);
302     delete myMCAnalyzer;
303     }
304    
305    
306     // Get Primary Vertices and Beam Spot
307     if(doPrimaryVertex_)
308     {
309     if(verbosity_>1) cout << endl << "Analysing beam spot and primary vertices collection..." << endl;
310     VertexAnalyzer* myVertexAnalyzer = new VertexAnalyzer(producersNames_, verbosity_);
311     if(doBeamSpot_) myVertexAnalyzer->getBeamSpot(iEvent, rootBeamSpot_);
312     myVertexAnalyzer->getVertices(iEvent, rootVertices_);
313     myVertexAnalyzer->selectPrimary(rootEvent_, rootVertices_);
314     delete myVertexAnalyzer;
315     }
316    
317    
318     // Reconstruct Zee Primary Vertices with and without electron pair
319     if(doZeePrimaryVertex_)
320     {
321     if(verbosity_>1) cout << endl << "Reconstruct Zee Primary Vertices with and without electron pair..." << endl;
322     ZeeVertexAnalyzer* myZeeVertexAnalyzer = new ZeeVertexAnalyzer(myConfig_, producersNames_, verbosity_);
323     myZeeVertexAnalyzer->getVertices(iEvent, iSetup, rootZeeVertices_);
324     delete myZeeVertexAnalyzer;
325     }
326    
327    
328     // Tracks
329     if(doTrack_)
330     {
331     if(verbosity_>1) cout << endl << "Analysing tracks collection..." << endl;
332     TrackAnalyzer* myTrackAnalyzer = new TrackAnalyzer(producersNames_, verbosity_);
333     myTrackAnalyzer->process(iEvent, rootTracks_);
334     delete myTrackAnalyzer;
335     }
336    
337    
338     // Calo Jet
339     if(doJet_)
340     {
341     if(verbosity_>1) cout << endl << "Analysing jets collection..." << endl;
342     JetAnalyzer* myJetAnalyzer = new JetAnalyzer(producersNames_, myConfig_, verbosity_);
343     myJetAnalyzer->process(iEvent, rootJets_);
344     delete myJetAnalyzer;
345     }
346    
347    
348     // Muons
349     if(doMuon_)
350     {
351     if(verbosity_>1) cout << endl << "Analysing muons collection..." << endl;
352     MuonAnalyzer* myMuonAnalyzer = new MuonAnalyzer(producersNames_, myConfig_, verbosity_);
353     if(doPrimaryVertex_) myMuonAnalyzer->initIPCalculator(iEvent, iSetup, rootEvent_, rootBeamSpot_);
354     myMuonAnalyzer->process(iEvent, rootBeamSpot_, rootMuons_);
355     delete myMuonAnalyzer;
356     }
357    
358    
359     // Lazy Tools to calculate Cluster shape variables
360     EcalClusterLazyTools* lazyTools = 0;
361     if( doElectron_ || doPhoton_ || doCluster_ )
362     {
363     if(verbosity_>1) cout << endl << "Loading egamma LazyTools..." << endl;
364     edm::InputTag reducedBarrelEcalRecHitCollection_ = producersNames_.getParameter<edm::InputTag>("reducedBarrelEcalRecHitCollection");
365     edm::InputTag reducedEndcapEcalRecHitCollection_ = producersNames_.getParameter<edm::InputTag>("reducedEndcapEcalRecHitCollection");
366     try
367     {
368     lazyTools = new EcalClusterLazyTools( iEvent, iSetup, reducedBarrelEcalRecHitCollection_, reducedEndcapEcalRecHitCollection_ );
369     }
370     catch (cms::Exception& exception)
371     {
372     if ( !allowMissingCollection_ )
373     {
374     cout << " ##### ERROR IN TotoAnalyzer::analyze => CaloGeometryRecord, CaloGeometryRecord or EcalRecHitCollection are missing #####"<<endl;
375     throw exception;
376     }
377     if(verbosity_>1) cout << " ===> CaloGeometryRecord, CaloGeometryRecord or EcalRecHitCollection are missing, skip calculation of cluster shape variables" << endl;
378     lazyTools = 0; // FIXME - delete authorized after an exception ?
379     }
380     }
381    
382    
383     // Electrons
384     if(doElectron_)
385     {
386     if(verbosity_>1) cout << endl << "Analysing electrons collection..." << endl;
387     ElectronAnalyzer* myElectronAnalyzer = new ElectronAnalyzer(producersNames_, myConfig_, verbosity_);
388     if(doPrimaryVertex_) myElectronAnalyzer->initIPCalculator(iEvent, iSetup, rootEvent_, rootBeamSpot_);
389     myElectronAnalyzer->process(iEvent, rootBeamSpot_, rootElectrons_, lazyTools);
390     delete myElectronAnalyzer;
391     }
392    
393    
394     // Photons
395     if(doPhoton_)
396     {
397     if(verbosity_>1) cout << endl << "Analysing photons collection..." << endl;
398     PhotonAnalyzer* myPhotonAnalyzer = new PhotonAnalyzer(producersNames_, myConfig_, verbosity_);
399     myPhotonAnalyzer->process(iEvent, iSetup, rootEvent_, rootPhotons_, rootConversionTracks_, conversionLikelihoodCalculator_, lazyTools);
400     delete myPhotonAnalyzer;
401     }
402    
403    
404     // Ecal Clusters
405     if(doCluster_)
406     {
407     if(verbosity_>1) cout << endl << "Analysing BasicClusters collection..." << endl;
408     ClusterAnalyzer* myClusterAnalyzer = new ClusterAnalyzer(producersNames_, verbosity_);
409     myClusterAnalyzer->process(iEvent, rootEvent_, lazyTools, rootBasicClusters_, "hybridSuperClusters","hybridBarrelBasicClusters", 210);
410     myClusterAnalyzer->process(iEvent, rootEvent_, lazyTools, rootBasicClusters_, "multi5x5BasicClusters", "multi5x5EndcapBasicClusters", 320);
411     myClusterAnalyzer->process(iEvent, rootEvent_, lazyTools, rootBasicClusters_, "multi5x5BasicClusters", "multi5x5BarrelBasicClusters", 310);
412     delete myClusterAnalyzer;
413    
414     if(verbosity_>1) cout << endl << "Analysing SuperClusters collection..." << endl;
415     SuperClusterAnalyzer* mySClusterAnalyzer = new SuperClusterAnalyzer(producersNames_, verbosity_);
416     mySClusterAnalyzer->process(iEvent, rootEvent_, rootSuperClusters_, "hybridSuperClusters","",210);
417     mySClusterAnalyzer->process(iEvent, rootEvent_, rootSuperClusters_, "correctedHybridSuperClusters","",211);
418     mySClusterAnalyzer->process(iEvent, rootEvent_, rootSuperClusters_, "multi5x5SuperClusters","multi5x5EndcapSuperClusters",320);
419     mySClusterAnalyzer->process(iEvent, rootEvent_, rootSuperClusters_, "multi5x5SuperClustersWithPreshower","",323);
420     mySClusterAnalyzer->process(iEvent, rootEvent_, rootSuperClusters_, "correctedMulti5x5SuperClustersWithPreshower","",322);
421     delete mySClusterAnalyzer;
422     }
423    
424    
425     if(doCluster_)
426     {
427     ClusterAssociator* myClusterAssociator = new ClusterAssociator();
428     myClusterAssociator->setVerbosity(verbosity_);
429     myClusterAssociator->process(rootSuperClusters_, rootBasicClusters_);
430     //myClusterAssociator->printSuperClusters(rootSuperClusters_, rootBasicClusters_,0); // 0 to print all types SC
431     delete myClusterAssociator;
432     }
433    
434    
435     if(doElectron_ && doCluster_)
436     {
437     ElectronAssociator* myElectronAssociator = new ElectronAssociator();
438     myElectronAssociator->setVerbosity(verbosity_);
439     myElectronAssociator->associateSuperCluster(rootElectrons_, rootSuperClusters_);
440     if(verbosity_>2) myElectronAssociator->fullPrintElectrons(rootElectrons_,rootSuperClusters_,rootBasicClusters_,0);
441     delete myElectronAssociator;
442     }
443    
444    
445     if(doPhoton_ && doCluster_)
446     {
447     PhotonAssociator* myPhotonAssociator = new PhotonAssociator();
448     myPhotonAssociator->setVerbosity(verbosity_);
449     myPhotonAssociator->associateSuperCluster(rootPhotons_, rootSuperClusters_);
450     delete myPhotonAssociator;
451     }
452    
453    
454     if(doPhoton_ && doPhotonIsolation_)
455     {
456     if(verbosity_>1) cout << endl << "Calculating Photon isolation..." << endl;
457    
458     Float_t ecalIslandIsolation;
459     Float_t ecalDoubleConeIsolation;
460     Float_t hcalRecHitIsolation;
461     Float_t isoTracks;
462     Int_t isoNTracks;
463     TRootPhoton* localPhoton;
464    
465     PhotonIsolator* myPhotonIsolator = new PhotonIsolator(&myConfig_, &producersNames_);
466     for (int iphoton=0; iphoton<rootPhotons_->GetEntriesFast(); iphoton++)
467     {
468     localPhoton = (TRootPhoton*)rootPhotons_->At(iphoton);
469    
470     ecalIslandIsolation=-1.;
471     if (doCluster_) ecalIslandIsolation = myPhotonIsolator->basicClustersIsolation(localPhoton, rootSuperClusters_, rootBasicClusters_);
472    
473     ecalDoubleConeIsolation=-1.;
474     if (doCluster_) ecalDoubleConeIsolation = myPhotonIsolator->basicClustersDoubleConeIsolation(localPhoton, rootSuperClusters_, rootBasicClusters_);
475    
476     hcalRecHitIsolation=-1.;
477     if ( myPhotonIsolator->loadHcalRecHits(iEvent, iSetup) ) hcalRecHitIsolation = myPhotonIsolator->hcalRecHitIsolation(localPhoton);
478    
479     isoTracks=-1.;
480     isoNTracks=-1;
481     if(doTrack_) isoTracks = myPhotonIsolator->trackerIsolation(localPhoton, rootTracks_, isoNTracks);
482    
483     localPhoton->setIsolationPerso(ecalIslandIsolation, ecalDoubleConeIsolation, hcalRecHitIsolation, isoTracks, isoNTracks);
484     if(verbosity_>4) cout << " After isolation - ["<< setw(3) << iphoton << "] "; localPhoton->Print(); cout << endl;
485     }
486     delete myPhotonIsolator;
487     }
488    
489    
490     if(doPhoton_ && doElectron_ && doCluster_)
491     {
492     if(verbosity_>1) cout << endl << "Runing Ambiguity Solver..." << endl;
493     AmbiguitySolver* myAmbiguitySolver = new AmbiguitySolver(myConfig_, verbosity_);
494     myAmbiguitySolver->processElectronsPhotons(rootSuperClusters_, rootElectrons_, rootPhotons_);
495     delete myAmbiguitySolver;
496     }
497    
498    
499     // Print Photons / SuperClusters / BasicClusters arborescence
500     if(doPhoton_)
501     {
502     PhotonAssociator* myPhotonAssociator = new PhotonAssociator();
503     if(verbosity_>2) myPhotonAssociator->fullPrintPhotons(rootPhotons_,rootSuperClusters_,rootBasicClusters_,0);
504     delete myPhotonAssociator;
505     }
506    
507    
508     // MET
509     if(doMET_)
510     {
511     if(verbosity_>1) cout << endl << "Analysing Missing Et..." << endl;
512     METAnalyzer* myMETAnalyzer = new METAnalyzer(producersNames_, myConfig_, verbosity_);
513     myMETAnalyzer->process(iEvent, rootMETs_);
514     delete myMETAnalyzer;
515     }
516    
517    
518     // Associate recoParticles to mcParticles
519     if(doMC_)
520     {
521     // TODO - For the moment, only for PAT format
522     if ( dataType_=="PAT")
523     {
524     MCAssociator* myMCAssociator = new MCAssociator(producersNames_, verbosity_);
525     myMCAssociator->init(iEvent, rootMCParticles_);
526     if(doPhoton_) myMCAssociator->matchGenParticlesTo(rootPhotons_);
527     if(doElectron_) myMCAssociator->matchGenParticlesTo(rootElectrons_);
528     if(doMuon_) myMCAssociator->matchGenParticlesTo(rootMuons_);
529     if(doMET_) myMCAssociator->matchGenParticlesTo(rootMETs_);
530     if(doJet_) myMCAssociator->matchGenParticlesTo(rootJets_);
531     if(doJet_ && doJetMC_) myMCAssociator->processGenJets(rootGenJets_, rootJets_);
532     if(verbosity_>4) cout << endl << "Printing recoParticles / mcParticles associations... " << endl;
533     if(verbosity_>4 && doPhoton_) myMCAssociator->printParticleAssociation(rootPhotons_);
534     if(verbosity_>4 && doElectron_) myMCAssociator->printParticleAssociation(rootElectrons_);
535     if(verbosity_>4 && doMuon_) myMCAssociator->printParticleAssociation(rootMuons_);
536     if(verbosity_>4 && doMET_) myMCAssociator->printParticleAssociation(rootMETs_);
537     if(verbosity_>4 && doJet_) myMCAssociator->printJetAssociation(rootJets_);
538     delete myMCAssociator;
539     }
540     }
541    
542    
543     if(verbosity_>1) cout << endl << "Filling rootuple..." << endl;
544     eventTree_->Fill();
545     if(verbosity_>1) cout << endl << "Start deleting objects..." << endl;
546     delete rootEvent_;
547     if(doMC_) (*rootMCParticles_).Delete();
548     if(doJetMC_) (*rootGenJets_).Delete();
549     if(doMETMC_) (*rootGenMETs_).Delete();
550     if(doSignalMuMuGamma_) delete rootMuMuGammaEvent_;
551     if(doSignalTopTop_) (*rootMCTopTop_).Delete();
552     if(doPhotonConversionMC_) (*rootMCPhotons_).Delete();
553     if(doBeamSpot_) rootBeamSpot_->clear();
554     if(doPrimaryVertex_) (*rootVertices_).Delete();
555     if(doZeePrimaryVertex_) (*rootZeeVertices_).Delete();
556     if(doTrack_) (*rootTracks_).Delete();
557     if(doJet_) (*rootJets_).Delete();
558     if(doMuon_) (*rootMuons_).Delete();
559     if(doElectron_) (*rootElectrons_).Delete();
560     if(doPhoton_) (*rootPhotons_).Delete();
561     if(doCluster_) (*rootBasicClusters_).Delete();
562     if(doCluster_) (*rootSuperClusters_).Delete();
563     if(doPhotonConversion_) (*rootConversionTracks_).Delete();
564     if(doMET_) (*rootMETs_).Delete();
565     if(verbosity_>1) cout << endl << "Objects deleted" << endl;
566     if(verbosity_>0) cout << endl;
567 mlethuil 1.1 }
568