ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/OSUT3Analysis/AnaTools/plugins/OSUAnalysis.cc

Comparing UserCode/OSUT3Analysis/AnaTools/plugins/OSUAnalysis.cc (file contents):
Revision 1.8 by ahart, Wed Jan 30 22:29:59 2013 UTC vs.
Revision 1.12 by ahart, Fri Feb 8 11:07:51 2013 UTC

#	Line 18 | Line 18 | OSUAnalysis::OSUAnalysis (const edm::Par
18		triggers_ (cfg.getParameter<edm::InputTag> ("triggers")),
19
20
21	<	channels_  (cfg.getParameter<vector<edm::ParameterSet> >("channels"))
21	>	channels_  (cfg.getParameter<vector<edm::ParameterSet> >("channels")),
22	>
23	>	histogramSets_ (cfg.getParameter<vector<edm::ParameterSet> >("histogramSets"))
24
25		{
26	+
27	+
28		TH1::SetDefaultSumw2 ();
29	+	TH2::SetDefaultSumw2 ();
30
31		// Construct Cutflow Objects. These store the results of cut decisions and
32		// handle filling cut flow histograms.
33		masterCutFlow_ = new CutFlow (fs_);
34		std::vector<TFileDirectory> directories;
35
36	+	//always get vertex collection so we can assign the primary vertex in the event
37	+	allNecessaryObjects.push_back("primaryvertexs");
38	+
39	+
40	+	//parse the histogram definitions
41	+	for(uint currentHistogramSet = 0; currentHistogramSet != histogramSets_.size(); currentHistogramSet++){
42	+
43	+	string tempInputCollection = histogramSets_.at(currentHistogramSet).getParameter<string>("inputCollection");
44	+	objectsToPlot.push_back(tempInputCollection);
45	+	allNecessaryObjects.push_back(tempInputCollection);
46	+	vector<edm::ParameterSet> histogramList_  (histogramSets_.at(currentHistogramSet).getParameter<vector<edm::ParameterSet> >("histograms"));
47	+
48	+	for(uint currentHistogram = 0; currentHistogram != histogramList_.size(); currentHistogram++){
49	+
50	+	histogram tempHistogram;
51	+	tempHistogram.inputCollection = tempInputCollection;
52	+	tempHistogram.name = histogramList_.at(currentHistogram).getParameter<string>("name");
53	+	tempHistogram.title = histogramList_.at(currentHistogram).getParameter<string>("title");
54	+	tempHistogram.bins = histogramList_.at(currentHistogram).getParameter<vector<double> >("bins");
55	+	tempHistogram.inputVariables = histogramList_.at(currentHistogram).getParameter<vector<string> >("inputVariables");
56	+	if(histogramList_.at(currentHistogram).exists("function"))
57	+	tempHistogram.function = histogramList_.at(currentHistogram).getParameter<string>("function");
58	+	else
59	+	tempHistogram.function = "";
60	+
61	+	histograms.push_back(tempHistogram);
62	+
63	+	}
64	+	}
65	+
66	+
67
68		channel tempChannel;
69		//loop over all channels (event selections)
#	Line 50 | Line 86 | OSUAnalysis::OSUAnalysis (const edm::Par
86		}
87
88
89	+
90		//create cutFlow for this channel
91		cutFlows_.push_back (new CutFlow (fs_, channelName));
92
93		//book a directory in the output file with the name of the channel
94		TFileDirectory subDir = fs_->mkdir( channelName );
95		directories.push_back(subDir);
59	–	std::map<std::string, TH1D*> histoMap;
60	–	oneDHists_.push_back(histoMap);
96
97	<	//gen histograms
98	<	oneDHists_.at(currentChannel)["genD0"] = directories.at(currentChannel).make<TH1D> ("genD0",channelLabel+" channel: Gen d_{0}; d_{0} [cm] ", 1000, -1, 1);
99	<	oneDHists_.at(currentChannel)["genAbsD0"] = directories.at(currentChannel).make<TH1D> ("genAbsD0",channelLabel+" channel: Gen d_{0}; |d_{0}| [cm] ", 1000, 0, 1);
100	<	oneDHists_.at(currentChannel)["genDZ"] = directories.at(currentChannel).make<TH1D> ("genDZ",channelLabel+" channel: Gen d_{z}; d_{z} [cm] ", 1000, -10, 10);
101	<	oneDHists_.at(currentChannel)["genAbsDZ"] = directories.at(currentChannel).make<TH1D> ("genAbsDZ",channelLabel+" channel: Gen d_{z}; |d_{z}| [cm] ", 1000, 0, 10);
102	<
103	<	//muon histograms
104	<	allNecessaryObjects.push_back("muons");
105	<	oneDHists_.at(currentChannel)["muonPt"]  = directories.at(currentChannel).make<TH1D> ("muonPt",channelLabel+" channel: Muon Transverse Momentum; p_{T} [GeV]", 100, 0, 500);
106	<	oneDHists_.at(currentChannel)["muonEta"] = directories.at(currentChannel).make<TH1D> ("muonEta",channelLabel+" channel: Muon Eta; #eta", 100, -5, 5);
107	<	oneDHists_.at(currentChannel)["muonD0"] = directories.at(currentChannel).make<TH1D> ("muonD0",channelLabel+" channel: Muon d_{0}; d_{0} [cm] ", 1000, -1, 1);
108	<	oneDHists_.at(currentChannel)["muonDz"] = directories.at(currentChannel).make<TH1D> ("muonDz",channelLabel+" channel: Muon d_{z}; d_{z} [cm] ", 1000, -20, 20);
109	<	oneDHists_.at(currentChannel)["muonAbsD0"] = directories.at(currentChannel).make<TH1D> ("muonAbsD0",channelLabel+" channel: Muon d_{0}; |d_{0}| [cm] ", 1000, 0, 1);
110	<	oneDHists_.at(currentChannel)["muonDZ"] = directories.at(currentChannel).make<TH1D> ("muonDZ",channelLabel+" channel: Muon d_{z}; d_{z} [cm] ", 1000, -10, 10);
111	<	oneDHists_.at(currentChannel)["muonAbsDZ"] = directories.at(currentChannel).make<TH1D> ("muonAbsDZ",channelLabel+" channel: Muon d_{z}; |d_{z}| [cm] ", 1000, 0, 10);
112	<	oneDHists_.at(currentChannel)["muonD0Sig"] = directories.at(currentChannel).make<TH1D> ("muonD0Sig",channelLabel+" channel: Muon d_{0} Significance; d_{0} / #sigma_{d_{0}} ", 1000, -10.0, 10.0);
113	<	oneDHists_.at(currentChannel)["muonAbsD0Sig"] = directories.at(currentChannel).make<TH1D> ("muonAbsD0Sig",channelLabel+" channel: Muon d_{0} Significance; |d_{0}| / #sigma_{d_{0}} ", 1000, 0, 10.0);
114	<	oneDHists_.at(currentChannel)["muonIso"] = directories.at(currentChannel).make<TH1D> ("muonIso",channelLabel+" channel: Muon Combined Relative Isolation; rel. iso. ", 1000, 0, 1);
115	<	oneDHists_.at(currentChannel)["numMuons"] = directories.at(currentChannel).make<TH1D> ("numMuons",channelLabel+" channel: Number of Selected Muons; # muons", 10, 0, 10);
116	<
117	<
118	<	//electron histograms
84	<	allNecessaryObjects.push_back("electrons");
85	<	oneDHists_.at(currentChannel)["electronPt"]  = directories.at(currentChannel).make<TH1D> ("electronPt",channelLabel+" channel: Electron Transverse Momentum; p_{T} [GeV]", 100, 0, 500);
86	<	oneDHists_.at(currentChannel)["electronEta"] = directories.at(currentChannel).make<TH1D> ("electronEta",channelLabel+" channel: Electron Eta; #eta", 50, -5, 5);
87	<	oneDHists_.at(currentChannel)["electronD0"] = directories.at(currentChannel).make<TH1D> ("electronD0",channelLabel+" channel: Electron d_{0}; d_{0} [cm] ", 1000, -1, 1);
88	<	oneDHists_.at(currentChannel)["electronDz"] = directories.at(currentChannel).make<TH1D> ("electronDz",channelLabel+" channel: Electron d_{z}; d_{z} [cm] ", 1000, -20, 20);
89	<	oneDHists_.at(currentChannel)["electronAbsD0"] = directories.at(currentChannel).make<TH1D> ("electronAbsD0",channelLabel+" channel: Electron d_{0}; |d_{0}| [cm] ", 1000, 0, 1);
90	<	oneDHists_.at(currentChannel)["electronDZ"] = directories.at(currentChannel).make<TH1D> ("electronDZ",channelLabel+" channel: Electron d_{z}; d_{z} [cm] ", 1000, -10, 10);
91	<	oneDHists_.at(currentChannel)["electronAbsDZ"] = directories.at(currentChannel).make<TH1D> ("electronAbsDZ",channelLabel+" channel: Electron d_{z}; |d_{z}| [cm] ", 1000, 0, 10);
92	<	oneDHists_.at(currentChannel)["electronD0Sig"] = directories.at(currentChannel).make<TH1D> ("electronD0Sig",channelLabel+" channel: Electron d_{0} Significance; d_{0} / #sigma_{d_{0}} ", 1000, -10.0, 10.0);
93	<	oneDHists_.at(currentChannel)["electronAbsD0Sig"] = directories.at(currentChannel).make<TH1D> ("electronAbsD0Sig",channelLabel+" channel: Electron d_{0} Significance; |d_{0}| / #sigma_{d_{0}} ", 1000, 0, 10.0);
94	<	oneDHists_.at(currentChannel)["electronIso"] = directories.at(currentChannel).make<TH1D> ("electronIso",channelLabel+" channel: Electron Combined Relative Isolation; rel. iso. ", 1000, 0, 1);
95	<	oneDHists_.at(currentChannel)["electronFbrem"] = directories.at(currentChannel).make<TH1D> ("electronFbrem",channelLabel+" channel: Electron Brem. Energy Fraction; fbrem", 1000, 0, 2);
96	<	oneDHists_.at(currentChannel)["electronMVATrig"] = directories.at(currentChannel).make<TH1D> ("electronMVATrig",channelLabel+" channel: Electron ID MVA Output", 100, -1, 1);
97	<	oneDHists_.at(currentChannel)["electronMVANonTrig"] = directories.at(currentChannel).make<TH1D> ("electronMVANonTrig",channelLabel+" channel: Electron ID MVA Output", 100, -1, 1);
98	<	oneDHists_.at(currentChannel)["numElectrons"] = directories.at(currentChannel).make<TH1D> ("numElectrons",channelLabel+" channel: Number of Selected Electrons; # electrons", 10, 0, 10);
97	>	std::map<std::string, TH1D*> oneDhistoMap;
98	>	oneDHists_.push_back(oneDhistoMap);
99	>	std::map<std::string, TH2D*> twoDhistoMap;
100	>	twoDHists_.push_back(twoDhistoMap);
101	>
102	>	//book all histograms included in the configuration
103	>	for(uint currentHistogramIndex = 0; currentHistogramIndex != histograms.size(); currentHistogramIndex++){
104	>	histogram currentHistogram = histograms.at(currentHistogramIndex);
105	>	int numBinsElements = currentHistogram.bins.size();
106	>	int numInputVariables = currentHistogram.inputVariables.size();
107	>	if(numBinsElements != 3 && numBinsElements !=6) {
108	>	std::cout << "Error: Didn't find correct number of bin specifications for histogram named '" << currentHistogram.name << "'\n";
109	>	exit(0);
110	>	}
111	>	else if((numBinsElements == 3 && numInputVariables !=1) || (numBinsElements == 6 && numInputVariables!=2)){
112	>	std::cout << "Error: Didn't find correct number of input variables for histogram named '" << currentHistogram.name << "'\n";
113	>	exit(0);
114	>	}
115	>	else if(numBinsElements == 3)
116	>	oneDHists_.at(currentChannel)[currentHistogram.name] = directories.at(currentChannel).make<TH1D> (TString(currentHistogram.name),channelLabel+" channel: "+currentHistogram.title, currentHistogram.bins.at(0), currentHistogram.bins.at(1), currentHistogram.bins.at(2));
117	>	else if(numBinsElements == 6)
118	>	twoDHists_.at(currentChannel)[currentHistogram.name] = directories.at(currentChannel).make<TH2D> (TString(currentHistogram.name),channelLabel+" channel: "+currentHistogram.title, currentHistogram.bins.at(0), currentHistogram.bins.at(1), currentHistogram.bins.at(2),currentHistogram.bins.at(3),currentHistogram.bins.at(4),currentHistogram.bins.at(5));
119
120
121	+	}
122	+	//book a histogram for the number of each object type to be plotted
123	+	for (uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){
124	+	string currentObject = objectsToPlot.at(currentObjectIndex);
125	+	int maxNum = 10;
126	+	if(currentObject == "mcparticles") maxNum = 50;
127	+	currentObject[0] = toupper(currentObject[0]);
128	+	string histoName = "num" + currentObject;
129	+	oneDHists_.at(currentChannel)[histoName] = directories.at(currentChannel).make<TH1D> (TString(histoName),channelLabel+" channel: Number of Selected "+currentObject+"; # "+currentObject, maxNum, 0, maxNum);
130	+	}
131	+
132
133		//get list of cuts for this channel
134		vector<edm::ParameterSet> cuts_  (channels_.at(currentChannel).getParameter<vector<edm::ParameterSet> >("cuts"));
135
105	–	vector<string> tempInputCollections;
106	–
107	–
136		//loop over and parse all cuts
137		for(uint currentCut = 0; currentCut != cuts_.size(); currentCut++){
138
#	Line 112 | Line 140 | OSUAnalysis::OSUAnalysis (const edm::Par
140		//store input collection for cut
141		string inputCollection = cuts_.at(currentCut).getParameter<string> ("inputCollection");
142		tempCut.inputCollection = inputCollection;
115	–
116	–	tempInputCollections.push_back(inputCollection);
143		allNecessaryObjects.push_back(inputCollection);
144
145		//split cut string into parts and store them
146		string cutString = cuts_.at(currentCut).getParameter<string> ("cutString");
147		std::vector<string> cutStringVector = splitString(cutString);
148	+	if(cutStringVector.size()!=3){
149	+	std::cout << "Error: Didn't find three elements in the following cut string: '" <<cutString << "'\n";
150	+	exit(0);
151	+	}
152	+
153		tempCut.variable = cutStringVector.at(0);// variable to cut on
154		tempCut.comparativeOperator = cutStringVector.at(1);// comparison to make
155		tempCut.cutValue = atof(cutStringVector.at(2).c_str());// threshold value to pass cut
#	Line 126 | Line 157 | OSUAnalysis::OSUAnalysis (const edm::Par
157		//get number of objects required to pass cut for event to pass
158		string numberRequiredString = cuts_.at(currentCut).getParameter<string> ("numberRequired");
159		std::vector<string> numberRequiredVector = splitString(numberRequiredString);
160	+	if(numberRequiredVector.size()!=2){
161	+	std::cout << "Error: Didn't find two elements in the following number requirement string: '" << numberRequiredString << "'\n";
162	+	exit(0);
163	+	}
164
165		// determine number required if comparison contains "="
166		int numberRequiredInt = atoi(numberRequiredVector.at(1).c_str());
#	Line 156 | Line 191 | OSUAnalysis::OSUAnalysis (const edm::Par
191		tempChannel.cuts.push_back(tempCut);
192
193
194	+
195		}//end loop over cuts
196
161	–	//make unique vector of all objects that this channel cuts on (so we know to set flags for those)
162	–	sort( tempInputCollections.begin(), tempInputCollections.end() );
163	–	tempInputCollections.erase( unique( tempInputCollections.begin(), tempInputCollections.end() ), tempInputCollections.end() );
164	–	tempChannel.inputCollections = tempInputCollections;
197
198		channels.push_back(tempChannel);
199		tempChannel.cuts.clear();
200
201		}//end loop over channels
202
203	+
204	+
205	+
206		//make unique vector of objects we need to cut on (so we make sure to get them from the event)
207		sort( allNecessaryObjects.begin(), allNecessaryObjects.end() );
208		allNecessaryObjects.erase( unique( allNecessaryObjects.begin(), allNecessaryObjects.end() ), allNecessaryObjects.end() );
#	Line 324 | Line 359 | OSUAnalysis::analyze (const edm::Event &
359		cutFlows_.at(currentChannelIndex)->fillCutFlow();
360
361
327	–
362		if(!eventPassedAllCuts)continue;
363
364
365	<	TVector3 vertexPosition;
365	>	//set position of primary vertex in event, in order to calculate quantities relative to it
366	>	primaryVertex_ = 0;
367		vector<bool> vertexFlags = cumulativeFlags["primaryvertexs"].back();
333	–
368		for (uint vertexIndex = 0; vertexIndex != vertexFlags.size(); vertexIndex++){
369		if(!vertexFlags.at(vertexIndex)) continue;
370	<
337	<	vertexPosition.SetXYZ (primaryvertexs->at(vertexIndex).x,primaryvertexs->at(vertexIndex).y,primaryvertexs->at(vertexIndex).z);
370	>	primaryVertex_ = new BNprimaryvertex (primaryvertexs->at (vertexIndex));
371		break;
372		}
373
341	–	vector<bool> mcparticleFlags = cumulativeFlags["mcparticles"].back();
374
343	–	for (uint mcparticleIndex = 0; mcparticleIndex != mcparticleFlags.size(); mcparticleIndex++){
344	–	if(!mcparticleFlags.at(mcparticleIndex)) continue;
345	–
346	–	double vx = mcparticles->at(mcparticleIndex).vx - vertexPosition.X (),
347	–	vy = mcparticles->at(mcparticleIndex).vy - vertexPosition.Y (),
348	–	vz = mcparticles->at(mcparticleIndex).vz - vertexPosition.Z (),
349	–	px = mcparticles->at(mcparticleIndex).px,
350	–	py = mcparticles->at(mcparticleIndex).py,
351	–	pt = mcparticles->at(mcparticleIndex).pt,
352	–	d0;
353	–
354	–	d0 = (-vx * py + vy * px) / pt;
355	–	oneDHists_.at(currentChannelIndex)["genD0"]->Fill (d0);
356	–	oneDHists_.at(currentChannelIndex)["genAbsD0"]->Fill (fabs (d0));
357	–	oneDHists_.at(currentChannelIndex)["genDZ"]->Fill (vz);
358	–	oneDHists_.at(currentChannelIndex)["genAbsDZ"]->Fill (fabs (vz));
359	–	}
360	–
361	–	vector<bool> electronFlags = cumulativeFlags["electrons"].back();
362	–	int electronCounter = 0;
363	–
364	–	for (uint electronIndex = 0; electronIndex != electronFlags.size(); electronIndex++){
365	–	if(!electronFlags.at(electronIndex)) continue;
366	–
367	–
368	–	electronCounter++;
369	–	oneDHists_.at(currentChannelIndex)["electronPt"]->Fill (electrons->at(electronIndex).pt);
370	–	oneDHists_.at(currentChannelIndex)["electronEta"]->Fill (electrons->at(electronIndex).eta);
371	–	oneDHists_.at(currentChannelIndex)["electronD0"]->Fill (electrons->at(electronIndex).correctedD0Vertex);
372	–	oneDHists_.at(currentChannelIndex)["electronDz"]->Fill (electrons->at(electronIndex).correctedDZ);
373	–	oneDHists_.at(currentChannelIndex)["electronAbsD0"]->Fill (fabs(electrons->at(electronIndex).correctedD0Vertex));
374	–	oneDHists_.at(currentChannelIndex)["electronDZ"]->Fill (electrons->at(electronIndex).correctedDZ);
375	–	oneDHists_.at(currentChannelIndex)["electronAbsDZ"]->Fill (fabs(electrons->at(electronIndex).correctedDZ));
376	–	oneDHists_.at(currentChannelIndex)["electronD0Sig"]->Fill (electrons->at(electronIndex).correctedD0Vertex / electrons->at(electronIndex).tkD0err);
377	–	oneDHists_.at(currentChannelIndex)["electronAbsD0Sig"]->Fill (fabs(electrons->at(electronIndex).correctedD0Vertex) / electrons->at(electronIndex).tkD0err);
378	–	oneDHists_.at(currentChannelIndex)["electronIso"]->Fill (electrons->at(electronIndex).puChargedHadronIso / electrons->at(electronIndex).pt);
379	–	oneDHists_.at(currentChannelIndex)["electronFbrem"]->Fill (electrons->at(electronIndex).fbrem);
380	–	oneDHists_.at(currentChannelIndex)["electronMVATrig"]->Fill (electrons->at(electronIndex).mvaTrigV0);
381	–	oneDHists_.at(currentChannelIndex)["electronMVANonTrig"]->Fill (electrons->at(electronIndex).mvaNonTrigV0);
375
376	+	//filling histograms
377	+	for (uint histogramIndex = 0; histogramIndex != histograms.size(); histogramIndex++){
378	+	histogram currentHistogram = histograms.at(histogramIndex);
379	+
380	+	if(currentHistogram.inputVariables.size() == 1){
381	+	TH1D* histo;
382	+	histo = oneDHists_.at(currentChannelIndex).at(currentHistogram.name);
383	+	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product());
384	+	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product());
385	+	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product());
386	+	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product());
387	+	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product());
388	+	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product());
389	+	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product());
390	+	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product());
391	+	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product());
392	+	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product());
393	+	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product());
394	+	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product());
395	+	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product());
396	+	}
397	+	else if(currentHistogram.inputVariables.size() == 2){
398	+	TH2D* histo;
399	+	histo = twoDHists_.at(currentChannelIndex).at(currentHistogram.name);
400	+	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product());
401	+	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product());
402	+	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product());
403	+	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product());
404	+	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product());
405	+	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product());
406	+	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product());
407	+	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product());
408	+	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product());
409	+	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product());
410	+	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product());
411	+	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product());
412	+	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product());
413	+	}
414		}
384	–	oneDHists_.at(currentChannelIndex)["numElectrons"]->Fill (electronCounter);
385	–
415
387	–	vector<bool> muonFlags = cumulativeFlags["muons"].back();
388	–	int muonCounter = 0;
416
417	+	//fills histograms with the sizes of collections
418	+	for (uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){
419	+	string currentObject = objectsToPlot.at(currentObjectIndex);
420	+	string tempCurrentObject = currentObject;
421	+	tempCurrentObject[0] = toupper(tempCurrentObject[0]);
422	+	string histoName = "num" + tempCurrentObject;
423	+
424	+	if(currentObject == "jets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(jets->size());
425	+	else if(currentObject == "muons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(muons->size());
426	+	else if(currentObject == "electrons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(electrons->size());
427	+	else if(currentObject == "events") oneDHists_.at(currentChannelIndex)[histoName]->Fill(events->size());
428	+	else if(currentObject == "taus") oneDHists_.at(currentChannelIndex)[histoName]->Fill(taus->size());
429	+	else if(currentObject == "mets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(mets->size());
430	+	else if(currentObject == "tracks") oneDHists_.at(currentChannelIndex)[histoName]->Fill(tracks->size());
431	+	else if(currentObject == "genjets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(genjets->size());
432	+	else if(currentObject == "mcparticles") oneDHists_.at(currentChannelIndex)[histoName]->Fill(mcparticles->size());
433	+	else if(currentObject == "primaryvertexs") oneDHists_.at(currentChannelIndex)[histoName]->Fill(primaryvertexs->size());
434	+	else if(currentObject == "bxlumis") oneDHists_.at(currentChannelIndex)[histoName]->Fill(bxlumis->size());
435	+	else if(currentObject == "photons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(photons->size());
436	+	else if(currentObject == "superclusters") oneDHists_.at(currentChannelIndex)[histoName]->Fill(superclusters->size());
437
391	–	for (uint muonIndex = 0; muonIndex != muonFlags.size(); muonIndex++){
392	–	if(!muonFlags.at(muonIndex)) continue;
393	–	muonCounter++;
394	–
395	–	oneDHists_.at(currentChannelIndex)["muonPt"]->Fill (muons->at(muonIndex).pt);
396	–	oneDHists_.at(currentChannelIndex)["muonEta"]->Fill (muons->at(muonIndex).eta);
397	–	oneDHists_.at(currentChannelIndex)["muonD0"]->Fill (muons->at(muonIndex).correctedD0Vertex);
398	–	oneDHists_.at(currentChannelIndex)["muonDz"]->Fill (muons->at(muonIndex).correctedDZ);
399	–	oneDHists_.at(currentChannelIndex)["muonAbsD0"]->Fill (fabs(muons->at(muonIndex).correctedD0Vertex));
400	–	oneDHists_.at(currentChannelIndex)["muonDZ"]->Fill (muons->at(muonIndex).correctedDZ);
401	–	oneDHists_.at(currentChannelIndex)["muonAbsDZ"]->Fill (fabs(muons->at(muonIndex).correctedDZ));
402	–	oneDHists_.at(currentChannelIndex)["muonD0Sig"]->Fill (muons->at(muonIndex).correctedD0Vertex / muons->at(muonIndex).tkD0err);
403	–	oneDHists_.at(currentChannelIndex)["muonAbsD0Sig"]->Fill (fabs(muons->at(muonIndex).correctedD0Vertex) / muons->at(muonIndex).tkD0err);
404	–	oneDHists_.at(currentChannelIndex)["muonIso"]->Fill (muons->at(muonIndex).puChargedHadronIso / muons->at(muonIndex).pt);
438		}
406	–	oneDHists_.at(currentChannelIndex)["numMuons"]->Fill (muonCounter);
407	–
408	–
409	–
410	–
411	–
439		} //end loop over channel
440
441		masterCutFlow_->fillCutFlow();
#	Line 746 | Line 773 | OSUAnalysis::valueLookup (const BNmuon*
773		else if(variable == "numberOfMatchedStations") value = object->numberOfMatchedStations;
774		else if(variable == "time_ndof") value = object->time_ndof;
775
776	+	//user-defined variables
777	+	else if(variable == "correctedD0VertexErr") value =  hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
778	+	else if(variable == "correctedD0VertexSig") value =  object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
779	+	else if(variable == "detIso") value = (object->trackIso + object->caloIso) / object->pt;
780	+	else if(variable == "relPFdBetaIso") value = (object->pfIsoR04SumChargedHadronPt + max(0.0, object->pfIsoR04SumNeutralHadronEt + object->pfIsoR04SumPhotonEt - 0.5*object->pfIsoR04SumPUPt)) / object->pt;
781	+	else if(variable == "relPFrhoIso") value = ( object->chargedHadronIso + max(0.0, object->neutralHadronIso + object->photonIso - object->AEffDr03*object->rhoPrime) ) / object->pt;
782	+	else if(variable == "tightID") {
783	+	value = object->isGlobalMuon > 0            \
784	+	&& object->isPFMuon > 0                   \
785	+	&& object->normalizedChi2 < 10            \
786	+	&& object->numberOfValidMuonHits > 0      \
787	+	&& object->numberOfMatchedStations > 1    \
788	+	&& fabs(object->correctedD0Vertex) < 0.2  \
789	+	&& fabs(object->correctedDZ) < 0.5        \
790	+	&& object->numberOfValidPixelHits > 0     \
791	+	&& object->numberOfLayersWithMeasurement > 5;
792	+
793	+
794	+	}
795	+	else if(variable == "tightIDdisplaced"){
796	+	value = object->isGlobalMuon > 0            \
797	+	&& object->isPFMuon > 0                   \
798	+	&& object->normalizedChi2 < 10            \
799	+	&& object->numberOfValidMuonHits > 0      \
800	+	&& object->numberOfMatchedStations > 1    \
801	+	&& object->numberOfValidPixelHits > 0     \
802	+	&& object->numberOfLayersWithMeasurement > 5;
803	+	}
804	+
805		else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;}
806
807		value = applyFunction(function, value);
#	Line 915 | Line 971 | OSUAnalysis::valueLookup (const BNelectr
971		else if(variable == "eidHyperTight4MC") value = object->eidHyperTight4MC;
972		else if(variable == "passConvVeto") value = object->passConvVeto;
973
974	+	//user-defined variables
975	+	else if(variable == "correctedD0VertexErr") value =  hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
976	+	else if(variable == "correctedD0VertexSig") value =  object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
977	+	else if(variable == "detIso") value = (object->trackIso + object->caloIso) / object->pt;
978	+	else if(variable == "relPFrhoIso") value = ( object->chargedHadronIsoDR03 + max(0.0, object->neutralHadronIsoDR03 + object->photonIsoDR03 - object->AEffDr03*object->rhoPrime) ) / object->pt;
979	+
980
981		else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;}
982
#	Line 1275 | Line 1337 | OSUAnalysis::valueLookup (const BNmcpart
1337		else if(variable == "grandMother11Status") value = object->grandMother11Status;
1338		else if(variable == "grandMother11Charge") value = object->grandMother11Charge;
1339
1340	+	//user-defined variables
1341	+	else if (variable == "d0"){
1342	+	double vx = object->vx - primaryVertex_->x,
1343	+	vy = object->vy - primaryVertex_->y,
1344	+	px = object->px,
1345	+	py = object->py,
1346	+	pt = object->pt;
1347	+	value = (-vx * py + vy * px) / pt;
1348	+	}
1349	+	else if (variable == "dz"){
1350	+	double vx = object->vx - primaryVertex_->x,
1351	+	vy = object->vy - primaryVertex_->y,
1352	+	vz = object->vz - primaryVertex_->z,
1353	+	px = object->px,
1354	+	py = object->py,
1355	+	pz = object->pz,
1356	+	pt = object->pt;
1357	+	value = vz - (vx * px + vy * py)/pt * (pz/pt);
1358	+	}
1359	+
1360	+
1361	+
1362		else{std::cout << "WARNING: invalid variable '" << variable << "'\n"; value = -999;}
1363
1364		value = applyFunction(function, value);
#	Line 1435 | Line 1519 | OSUAnalysis::valueLookup (const BNsuperc
1519		double
1520		OSUAnalysis::applyFunction(string function, double value){
1521
1522	<	if(function == "abs") value = abs(value);
1522	>	if(function == "abs") value = fabs(value);
1523
1524
1525		return value;
#	Line 1473 | Line 1557 | void OSUAnalysis::setObjectFlags(cut &cu
1557		}
1558
1559
1560	+	template <class InputCollection>
1561	+	void OSUAnalysis::fillHistogram(TH1* histo, histogram parameters, InputCollection inputCollection){
1562	+
1563	+
1564	+	for (uint object = 0; object != inputCollection->size(); object++){
1565	+	if(parameters.inputVariables.size() == 1){
1566	+	double value = valueLookup(&inputCollection->at(object), parameters.inputVariables.at(0), parameters.function);
1567	+	histo->Fill(value);
1568	+	}
1569	+	else if(parameters.inputVariables.size() == 2){
1570	+	double valueX = valueLookup(&inputCollection->at(object), parameters.inputVariables.at(0), parameters.function);
1571	+	double valueY = valueLookup(&inputCollection->at(object), parameters.inputVariables.at(1), parameters.function);
1572	+	histo->Fill(valueX,valueY);
1573	+	}
1574	+	}
1575	+
1576	+	}
1577	+
1578
1579
1580

Diff Legend

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