ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/kiesel/TreeWriter/treeWriter.cc

Comparing UserCode/kiesel/TreeWriter/treeWriter.cc (file contents):
Revision 1.8 by kiesel, Wed Apr 10 20:46:33 2013 UTC vs.
Revision 1.9 by kiesel, Sat Apr 13 09:53:02 2013 UTC

#	Line 8 | Line 8
8
9		using namespace std;
10
11	<	TreeWriter::TreeWriter(TString inputName, TString outputName, int loggingVerbosity_ ) {
11	>	TreeWriter::TreeWriter( TString inputName, TString outputName, int loggingVerbosity_ ) {
12		// read the input file
13		inputTree = new TChain("susyTree");
14		if (loggingVerbosity_ > 0)
#	Line 20 | Line 20 | TreeWriter::TreeWriter(TString inputName
20		event = new susy::Event;
21		inputTree->SetBranchAddress("susyEvent", &event);
22
23	<	if (outputName == "") {
24	<	// Here the output filename is generated automaticly
25	<	try{
26	<	string fName = (string) inputName;
27	<	// Search for 'susyEvents_', and get the unique characters afterwards.
28	<	// eg /path/susyEvents_123_Myl.root -> susyTree_123_Mly.root
29	<	outputName = "susyTree_" + fName.substr(fName.find("susyEvents_")+11,-1);
30	<	} catch (int e ) {
31	<	outputName = "susyTree.root";
32	<	}
33	<	}
34	<
23	>	// open the output file
24		if (loggingVerbosity_>0)
25		std::cout << "Open file " << outputName << " for writing." << std::endl;
37	–	// open the output file
26		outFile = new TFile( outputName, "recreate" );
27		tree = new TTree("susyTree","Tree for single photon analysis");
28
#	Line 43 | Line 31 | TreeWriter::TreeWriter(TString inputName
31		reportEvery = 1000;
32		loggingVerbosity = loggingVerbosity_;
33		skim = true;
34	+	pileupHisto = 0;
35	+	}
36
37	+	void TreeWriter::PileUpWeightFile( string pileupFileName ) {
38	+	TFile *puFile = new TFile( pileupFileName.c_str() );
39	+	pileupHisto = (TH1F*) puFile->Get("pileup");
40		}
41
42	+
43	+
44		TreeWriter::~TreeWriter() {
45	<	if (!inputTree) return;
46	<	delete inputTree->GetCurrentFile();
45	>	if (pileupHisto != 0 )
46	>	delete pileupHisto;
47	>	inputTree->GetCurrentFile()->Close();
48		}
49
50		// useful functions
#	Line 111 | Line 107 | float photonIso_corrected(susy::Photon g
107		return iso;
108		}
109
114	–
115	–
116	–
110		float TreeWriter::getPtFromMatchedJet( susy::Photon myPhoton, susy::Event myEvent ) {
111		/**
112		* \brief Takes jet p_T as photon p_T
#	Line 132 | Line 125 | float TreeWriter::getPtFromMatchedJet( s
125		} else {
126		susy::PFJetCollection& jetColl = pfJets_it->second;
127		for(std::vector<susy::PFJet>::iterator it = jetColl.begin();
128	<	it != jetColl.end(); it++) {
128	>	it != jetColl.end(); ++it) {
129		std::map<TString,Float_t>::iterator s_it = it->jecScaleFactors.find("L2L3");
130		if (s_it == it->jecScaleFactors.end()) {
131		std::cout << "JEC is not available for this jet!!!" << std::endl;
#	Line 156 | Line 149 | float TreeWriter::getPtFromMatchedJet( s
149		float pt = 0;
150		float minPtDifferenz = 1E20; // should be very high
151		for( std::vector<susy::PFJet>::iterator it = nearJets.begin(), jetEnd = nearJets.end();
152	<	it != jetEnd; it++ ) {
152	>	it != jetEnd; ++it ) {
153		float ptDiff = fabs(myPhoton.momentum.Et() - it->momentum.Et());
154		if (  ptDiff < minPtDifferenz ) {
155		minPtDifferenz = ptDiff;
#	Line 190 | Line 183 | void TreeWriter::Loop() {
183
184		tree->Branch("photon", &photon);
185		tree->Branch("jet", &jet);
186	+	tree->Branch("electron", &electron);
187	+	tree->Branch("muon", &muon);
188		tree->Branch("met", &met, "met/F");
189		tree->Branch("ht", &ht, "ht/F");
190		tree->Branch("nVertex", &nVertex, "nVertex/I");
191	<	tree->Branch("nElectron", &nElectron, "nElectron/I");
191	>	tree->Branch("pu_weight", &pu_weight, "pu_weight/F");
192
193
194	<	for (long jentry=0; jentry < processNEvents; jentry++) {
194	>	for (long jentry=0; jentry < processNEvents; ++jentry) {
195		if ( loggingVerbosity>0 && jentry%reportEvery==0 )
196		std::cout << jentry << " / " << processNEvents << std :: endl;
197		inputTree->GetEntry(jentry);
198
199		photon.clear();
200		jet.clear();
201	+	electron.clear();
202	+	muon.clear();
203		ht = 0;
204
205	+	// weights
206	+	if (pileupHisto == 0) {
207	+	pu_weight = 1.;
208	+	} else {
209	+	float trueNumInteractions = -1;
210	+	for( susy::PUSummaryInfoCollection::const_iterator iBX = event->pu.begin();
211	+	iBX != event->pu.end() && trueNumInteractions < 0; ++iBX) {
212	+	if (iBX->BX == 0)
213	+	trueNumInteractions = iBX->trueNumInteractions;
214	+	}
215	+	pu_weight = pileupHisto->GetBinContent( pileupHisto->FindBin( trueNumInteractions ) );
216	+	}
217	+
218	+
219		// photons
220		if( loggingVerbosity > 1 )
221		std::cout << "Process photons" << std::endl;
222		std::map<TString, std::vector<susy::Photon> >::iterator phoMap = event->photons.find("photons");
223		for(std::vector<susy::Photon>::iterator it = phoMap->second.begin();
224	<	it != phoMap->second.end() && phoMap != event->photons.end(); it++ ) {
224	>	it != phoMap->second.end() && phoMap != event->photons.end(); ++it ) {
225		if( ! it->isEB() && skim )
226		continue;
227		thisphoton->pt = getPtFromMatchedJet( *it, *event );
#	Line 228 | Line 239 | void TreeWriter::Loop() {
239		thisphoton->hadTowOverEm = it->hadTowOverEm;
240		thisphoton->pixelseed = it->nPixelSeeds;
241		photon.push_back( *thisphoton );
231	–	ht += thisphoton->pt;
242		if( loggingVerbosity > 2 )
243		std::cout << " p_T, gamma = " << thisphoton->pt << std::endl;
244		}
#	Line 248 | Line 258 | void TreeWriter::Loop() {
258		susy::PFJetCollection& jetColl = pfJets_it->second;
259
260		for(std::vector<susy::PFJet>::iterator it = jetColl.begin();
261	<	it != jetColl.end(); it++) {
261	>	it != jetColl.end(); ++it) {
262		std::map<TString,Float_t>::iterator s_it = it->jecScaleFactors.find("L2L3");
263		if (s_it == it->jecScaleFactors.end()) {
264		std::cout << "JEC is not available for this jet!!!" << std::endl;
#	Line 262 | Line 272 | void TreeWriter::Loop() {
272		thisjet->pt = corrP4.Et();
273		thisjet->eta = corrP4.Eta();
274		thisjet->phi = corrP4.Phi();
275	+	thisjet->bCSV = it->bTagDiscriminators[susy::kCSV];
276	+	if( loggingVerbosity > 2 )
277	+	std::cout << " p_T, jet = " << thisjet->pt << std::endl;
278	+
279		jet.push_back( *thisjet );
280		ht += thisjet->pt;
281		}// for jet
#	Line 269 | Line 283 | void TreeWriter::Loop() {
283		if( jet.size() < 2 && skim )
284		continue;
285		std::sort( jet.begin(), jet.end(), tree::EtGreater);
286	+	if( loggingVerbosity > 1 )
287	+	std::cout << "Found " << jet.size() << " jets" << std::endl;
288	+
289
290		// met
291		std::map<TString, susy::MET>::iterator met_it = event->metMap.find("pfMet");
292		susy::MET* metobj = &(met_it->second);
293		met = metobj->met();
294	+	met_phi = metobj->mEt.Phi();
295	+	if( loggingVerbosity > 2 )
296	+	std::cout << " met = " << met << std::endl;
297	+
298	+	std::map<TString, susy::MET>::iterator type1met_it = event->metMap.find("pfType1CorrectedMet");
299	+	susy::MET* type1metobj = &(type1met_it->second);
300	+	type1met = type1metobj->met();
301	+	type1met_phi = type1metobj->mEt.Phi();
302	+	if( loggingVerbosity > 2 )
303	+	std::cout << " type1met = " << type1met << std::endl;
304
305		// electrons
306	<	//std::vector<susy::Electron> eVector = event->electrons["gsfElectrons"];
280	<	//nElectron = eVector.size();
306	>	tree::Particle* thiselectron = new tree::Particle();
307
308	<	// vertices
308	>	map<TString, vector<susy::Electron> >::iterator eleMap = event->electrons.find("gsfElectrons");
309	>	if(eleMap == event->electrons.end() && loggingVerbosity > 0) {
310	>	cout << "gsfElectrons not found!" << endl;
311	>	} else {
312	>	for(vector<susy::Electron>::iterator it = eleMap->second.begin(); it < eleMap->second.end(); ++it) {
313	>	thiselectron->pt = it->momentum.Et();
314	>	if( thiselectron->pt < 20 )
315	>	continue;
316	>	if( loggingVerbosity > 2 )
317	>	std::cout << " p_T, electron = " << it->momentum.Et() << std::endl;
318	>	thiselectron->eta = it->momentum.Eta();
319	>	thiselectron->phi = it->momentum.Phi();
320	>	electron.push_back( *thiselectron );
321	>	}
322	>	}
323	>	if( loggingVerbosity > 1 )
324	>	std::cout << "Found " << electron.size() << " electrons" << std::endl;
325	>
326	>	// this seems not to work yet, where is the bug?
327	>	/*
328	>	std::vector<susy::Electron> eVector = event->electrons["gsfElectronsx"];
329	>	for( std::vector<susy::Electron>::iterator it = eVector.begin(); it != eVector.end(); ++it) {
330	>	thiselectron->pt = it->momentum.Et();
331	>	if( thiselectron->pt < 20 )
332	>	continue;
333	>	if( loggingVerbosity > 2 )
334	>	std::cout << " p_T, electron = " << it->momentum.Et() << std::endl;
335	>	thiselectron->eta = it->momentum.Eta();
336	>	thiselectron->phi = it->momentum.Phi();
337	>	electron.push_back( *thiselectron );
338	>	}
339	>	*/
340
341	+	// muons
342	+	std::vector<susy::Muon> mVector = event->muons["muons"];
343	+	tree::Particle* thismuon = new tree::Particle();
344	+	for( std::vector<susy::Muon>::iterator it = mVector.begin(); it != mVector.end(); ++it) {
345	+	thismuon->pt = it->momentum.Et();
346	+	if( thismuon->pt < 20 )
347	+	continue;
348	+	thismuon->eta = it->momentum.Eta();
349	+	thismuon->phi = it->momentum.Phi();
350	+	muon.push_back( *thismuon );
351	+	}
352	+	if( loggingVerbosity > 1 )
353	+	std::cout << "Found " << muon.size() << " muons" << std::endl;
354	+
355	+
356	+	// vertices
357		nVertex = event->vertices.size();
358
359		if( ht < 450 && skim)
#	Line 290 | Line 363 | void TreeWriter::Loop() {
363		} // for jentry
364
365
293	–
294	–	tree->Write();
366		outFile->cd();
367	+	tree->Write();
368		outFile->Write();
369		outFile->Close();
370		}

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines