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.12 by ahart, Fri Feb 8 11:07:51 2013 UTC vs.
Revision 1.14 by lantonel, Wed Feb 13 10:05:35 2013 UTC

#	Line 20 | Line 20 | OSUAnalysis::OSUAnalysis (const edm::Par
20
21		channels_  (cfg.getParameter<vector<edm::ParameterSet> >("channels")),
22
23	<	histogramSets_ (cfg.getParameter<vector<edm::ParameterSet> >("histogramSets"))
23	>	histogramSets_ (cfg.getParameter<vector<edm::ParameterSet> >("histogramSets")),
24
25	<	{
25	>	plotAllObjectsInPassingEvents (cfg.getParameter<bool> ("plotAllObjectsInPassingEvents"))
26
27	+	{
28
29		TH1::SetDefaultSumw2 ();
30		TH2::SetDefaultSumw2 ();
#	Line 35 | Line 36 | OSUAnalysis::OSUAnalysis (const edm::Par
36
37		//always get vertex collection so we can assign the primary vertex in the event
38		allNecessaryObjects.push_back("primaryvertexs");
39	<
39	>	//always make the plot of number of primary verticex (to check pile-up reweighting)
40	>	objectsToPlot.push_back("primaryvertexs");
41
42		//parse the histogram definitions
43		for(uint currentHistogramSet = 0; currentHistogramSet != histogramSets_.size(); currentHistogramSet++){
#	Line 65 | Line 67 | OSUAnalysis::OSUAnalysis (const edm::Par
67
68
69
70	+
71		channel tempChannel;
72		//loop over all channels (event selections)
73		for(uint currentChannel = 0; currentChannel != channels_.size(); currentChannel++){
#	Line 87 | Line 90 | OSUAnalysis::OSUAnalysis (const edm::Par
90
91
92
93	+
94		//create cutFlow for this channel
95		cutFlows_.push_back (new CutFlow (fs_, channelName));
96
#	Line 99 | Line 103 | OSUAnalysis::OSUAnalysis (const edm::Par
103		std::map<std::string, TH2D*> twoDhistoMap;
104		twoDHists_.push_back(twoDhistoMap);
105
106	+
107		//book all histograms included in the configuration
108		for(uint currentHistogramIndex = 0; currentHistogramIndex != histograms.size(); currentHistogramIndex++){
109		histogram currentHistogram = histograms.at(currentHistogramIndex);
110		int numBinsElements = currentHistogram.bins.size();
111		int numInputVariables = currentHistogram.inputVariables.size();
112	+
113	+
114	+
115		if(numBinsElements != 3 && numBinsElements !=6) {
116		std::cout << "Error: Didn't find correct number of bin specifications for histogram named '" << currentHistogram.name << "'\n";
117		exit(0);
#	Line 112 | Line 120 | OSUAnalysis::OSUAnalysis (const edm::Par
120		std::cout << "Error: Didn't find correct number of input variables for histogram named '" << currentHistogram.name << "'\n";
121		exit(0);
122		}
123	<	else if(numBinsElements == 3)
123	>	else if(numBinsElements == 3){
124		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));
125	<	else if(numBinsElements == 6)
125	>	}
126	>	else if(numBinsElements == 6){
127		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));
128
129	+	}
130	+
131
132		}
133		//book a histogram for the number of each object type to be plotted
#	Line 124 | Line 135 | OSUAnalysis::OSUAnalysis (const edm::Par
135		string currentObject = objectsToPlot.at(currentObjectIndex);
136		int maxNum = 10;
137		if(currentObject == "mcparticles") maxNum = 50;
138	<	currentObject[0] = toupper(currentObject[0]);
138	>	currentObject.at(0) = toupper(currentObject.at(0));
139		string histoName = "num" + currentObject;
140		oneDHists_.at(currentChannel)[histoName] = directories.at(currentChannel).make<TH1D> (TString(histoName),channelLabel+" channel: Number of Selected "+currentObject+"; # "+currentObject, maxNum, 0, maxNum);
141		}
#	Line 201 | Line 212 | OSUAnalysis::OSUAnalysis (const edm::Par
212		}//end loop over channels
213
214
204	–
205	–
215		//make unique vector of objects we need to cut on (so we make sure to get them from the event)
216		sort( allNecessaryObjects.begin(), allNecessaryObjects.end() );
217		allNecessaryObjects.erase( unique( allNecessaryObjects.begin(), allNecessaryObjects.end() ), allNecessaryObjects.end() );
#	Line 282 | Line 291 | OSUAnalysis::analyze (const edm::Event &
291		event.getByLabel (superclusters_, superclusters);
292
293
294	+
295		//loop over all channels
296
297		for(uint currentChannelIndex = 0; currentChannelIndex != channels.size(); currentChannelIndex++){
#	Line 307 | Line 317 | OSUAnalysis::analyze (const edm::Event &
317		individualFlags[currentObject].push_back (vector<bool> ());
318		cumulativeFlags[currentObject].push_back (vector<bool> ());
319
320	+
321		if(currentObject == "jets") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,jets.product(),"jets");
322		else if(currentObject == "muons") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,muons.product(),"muons");
323		else if(currentObject == "electrons") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,electrons.product(),"electrons");
#	Line 322 | Line 333 | OSUAnalysis::analyze (const edm::Event &
333		else if(currentObject == "superclusters") setObjectFlags(currentCut,currentCutIndex,individualFlags,cumulativeFlags,superclusters.product(),"superclusters");
334
335
336	+
337		}
338
339
#	Line 345 | Line 357 | OSUAnalysis::analyze (const edm::Event &
357		cut currentCut = currentChannel.cuts.at(currentCutIndex);
358		int numberPassing = 0;
359
360	<	for (uint object = 0; object != cumulativeFlags[currentCut.inputCollection].at(currentCutIndex).size() ; object++)
361	<	if(cumulativeFlags[currentCut.inputCollection].at(currentCutIndex).at(object)) numberPassing++;
360	>	for (uint object = 0; object != cumulativeFlags.at(currentCut.inputCollection).at(currentCutIndex).size() ; object++)
361	>	if(cumulativeFlags.at(currentCut.inputCollection).at(currentCutIndex).at(object)) numberPassing++;
362
363
364		bool cutDecision = evaluateComparison(numberPassing,currentCut.eventComparativeOperator,currentCut.numberRequired);
#	Line 362 | Line 374 | OSUAnalysis::analyze (const edm::Event &
374		if(!eventPassedAllCuts)continue;
375
376
377	+
378	+
379		//set position of primary vertex in event, in order to calculate quantities relative to it
380		primaryVertex_ = 0;
381	<	vector<bool> vertexFlags = cumulativeFlags["primaryvertexs"].back();
381	>	vector<bool> vertexFlags = cumulativeFlags.at("primaryvertexs").back();
382		for (uint vertexIndex = 0; vertexIndex != vertexFlags.size(); vertexIndex++){
383		if(!vertexFlags.at(vertexIndex)) continue;
384		primaryVertex_ = new BNprimaryvertex (primaryvertexs->at (vertexIndex));
#	Line 380 | Line 394 | OSUAnalysis::analyze (const edm::Event &
394		if(currentHistogram.inputVariables.size() == 1){
395		TH1D* histo;
396		histo = oneDHists_.at(currentChannelIndex).at(currentHistogram.name);
397	<	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product());
398	<	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product());
399	<	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product());
400	<	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product());
401	<	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product());
402	<	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product());
403	<	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product());
404	<	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product());
405	<	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product());
406	<	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product());
407	<	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product());
408	<	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product());
409	<	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product());
397	>	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product(),cumulativeFlags.at("jets").back());
398	>	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product(),cumulativeFlags.at("muons").back());
399	>	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product(),cumulativeFlags.at("electrons").back());
400	>	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product(),cumulativeFlags.at("events").back());
401	>	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product(),cumulativeFlags.at("taus").back());
402	>	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product(),cumulativeFlags.at("mets").back());
403	>	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product(),cumulativeFlags.at("tracks").back());
404	>	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product(),cumulativeFlags.at("genjets").back());
405	>	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product(),cumulativeFlags.at("mcparticles").back());
406	>	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product(),cumulativeFlags.at("primaryvertexs").back());
407	>	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product(),cumulativeFlags.at("bxlumis").back());
408	>	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product(),cumulativeFlags.at("photons").back());
409	>	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product(),cumulativeFlags.at("superclusters").back());
410		}
411		else if(currentHistogram.inputVariables.size() == 2){
412		TH2D* histo;
413		histo = twoDHists_.at(currentChannelIndex).at(currentHistogram.name);
414	<	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product());
415	<	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product());
416	<	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product());
417	<	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product());
418	<	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product());
419	<	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product());
420	<	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product());
421	<	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product());
422	<	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product());
423	<	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product());
424	<	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product());
425	<	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product());
426	<	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product());
414	>	if(currentHistogram.inputCollection == "jets") fillHistogram(histo,currentHistogram,jets.product(),cumulativeFlags.at("jets").back());
415	>	else if(currentHistogram.inputCollection == "muons") fillHistogram(histo,currentHistogram,muons.product(),cumulativeFlags.at("muons").back());
416	>	else if(currentHistogram.inputCollection == "electrons") fillHistogram(histo,currentHistogram,electrons.product(),cumulativeFlags.at("electrons").back());
417	>	else if(currentHistogram.inputCollection == "events") fillHistogram(histo,currentHistogram,events.product(),cumulativeFlags.at("events").back());
418	>	else if(currentHistogram.inputCollection == "taus") fillHistogram(histo,currentHistogram,taus.product(),cumulativeFlags.at("taus").back());
419	>	else if(currentHistogram.inputCollection == "mets") fillHistogram(histo,currentHistogram,mets.product(),cumulativeFlags.at("mets").back());
420	>	else if(currentHistogram.inputCollection == "tracks") fillHistogram(histo,currentHistogram,tracks.product(),cumulativeFlags.at("tracks").back());
421	>	else if(currentHistogram.inputCollection == "genjets") fillHistogram(histo,currentHistogram,genjets.product(),cumulativeFlags.at("genjets").back());
422	>	else if(currentHistogram.inputCollection == "mcparticles") fillHistogram(histo,currentHistogram,mcparticles.product(),cumulativeFlags.at("mcparticles").back());
423	>	else if(currentHistogram.inputCollection == "primaryvertexs") fillHistogram(histo,currentHistogram,primaryvertexs.product(),cumulativeFlags.at("primaryvertexs").back());
424	>	else if(currentHistogram.inputCollection == "bxlumis") fillHistogram(histo,currentHistogram,bxlumis.product(),cumulativeFlags.at("bxlumis").back());
425	>	else if(currentHistogram.inputCollection == "photons") fillHistogram(histo,currentHistogram,photons.product(),cumulativeFlags.at("photons").back());
426	>	else if(currentHistogram.inputCollection == "superclusters") fillHistogram(histo,currentHistogram,superclusters.product(),cumulativeFlags.at("superclusters").back());
427		}
428		}
429
#	Line 418 | Line 432 | OSUAnalysis::analyze (const edm::Event &
432		for (uint currentObjectIndex = 0; currentObjectIndex != objectsToPlot.size(); currentObjectIndex++){
433		string currentObject = objectsToPlot.at(currentObjectIndex);
434		string tempCurrentObject = currentObject;
435	<	tempCurrentObject[0] = toupper(tempCurrentObject[0]);
435	>	tempCurrentObject.at(0) = toupper(tempCurrentObject.at(0));
436		string histoName = "num" + tempCurrentObject;
437
438	<	if(currentObject == "jets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(jets->size());
439	<	else if(currentObject == "muons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(muons->size());
440	<	else if(currentObject == "electrons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(electrons->size());
441	<	else if(currentObject == "events") oneDHists_.at(currentChannelIndex)[histoName]->Fill(events->size());
442	<	else if(currentObject == "taus") oneDHists_.at(currentChannelIndex)[histoName]->Fill(taus->size());
443	<	else if(currentObject == "mets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(mets->size());
444	<	else if(currentObject == "tracks") oneDHists_.at(currentChannelIndex)[histoName]->Fill(tracks->size());
445	<	else if(currentObject == "genjets") oneDHists_.at(currentChannelIndex)[histoName]->Fill(genjets->size());
446	<	else if(currentObject == "mcparticles") oneDHists_.at(currentChannelIndex)[histoName]->Fill(mcparticles->size());
447	<	else if(currentObject == "primaryvertexs") oneDHists_.at(currentChannelIndex)[histoName]->Fill(primaryvertexs->size());
448	<	else if(currentObject == "bxlumis") oneDHists_.at(currentChannelIndex)[histoName]->Fill(bxlumis->size());
449	<	else if(currentObject == "photons") oneDHists_.at(currentChannelIndex)[histoName]->Fill(photons->size());
450	<	else if(currentObject == "superclusters") oneDHists_.at(currentChannelIndex)[histoName]->Fill(superclusters->size());
438	>	if(currentObject == "jets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(jets->size());
439	>	else if(currentObject == "muons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(muons->size());
440	>	else if(currentObject == "electrons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(electrons->size());
441	>	else if(currentObject == "events") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(events->size());
442	>	else if(currentObject == "taus") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(taus->size());
443	>	else if(currentObject == "mets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mets->size());
444	>	else if(currentObject == "tracks") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(tracks->size());
445	>	else if(currentObject == "genjets") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(genjets->size());
446	>	else if(currentObject == "mcparticles") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(mcparticles->size());
447	>	else if(currentObject == "primaryvertexs") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(primaryvertexs->size());
448	>	else if(currentObject == "bxlumis") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(bxlumis->size());
449	>	else if(currentObject == "photons") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(photons->size());
450	>	else if(currentObject == "superclusters") oneDHists_.at(currentChannelIndex).at(histoName)->Fill(superclusters->size());
451
452		}
453		} //end loop over channel
#	Line 441 | Line 455 | OSUAnalysis::analyze (const edm::Event &
455		masterCutFlow_->fillCutFlow();
456
457
458	+
459		}
460
461
#	Line 776 | Line 791 | OSUAnalysis::valueLookup (const BNmuon*
791		//user-defined variables
792		else if(variable == "correctedD0VertexErr") value =  hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
793		else if(variable == "correctedD0VertexSig") value =  object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
794	<	else if(variable == "detIso") value = (object->trackIso + object->caloIso) / object->pt;
794	>	else if(variable == "detIso") value = (object->trackIso) / object->pt;
795		else if(variable == "relPFdBetaIso") value = (object->pfIsoR04SumChargedHadronPt + max(0.0, object->pfIsoR04SumNeutralHadronEt + object->pfIsoR04SumPhotonEt - 0.5*object->pfIsoR04SumPUPt)) / object->pt;
796		else if(variable == "relPFrhoIso") value = ( object->chargedHadronIso + max(0.0, object->neutralHadronIso + object->photonIso - object->AEffDr03*object->rhoPrime) ) / object->pt;
797		else if(variable == "tightID") {
#	Line 974 | Line 989 | OSUAnalysis::valueLookup (const BNelectr
989		//user-defined variables
990		else if(variable == "correctedD0VertexErr") value =  hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
991		else if(variable == "correctedD0VertexSig") value =  object->correctedD0Vertex / hypot (object->tkD0err, hypot (primaryVertex_->xError, primaryVertex_->yError));
992	<	else if(variable == "detIso") value = (object->trackIso + object->caloIso) / object->pt;
992	>	else if(variable == "detIso") value = (object->trackIso) / object->pt;
993		else if(variable == "relPFrhoIso") value = ( object->chargedHadronIsoDR03 + max(0.0, object->neutralHadronIsoDR03 + object->photonIsoDR03 - object->AEffDr03*object->rhoPrime) ) / object->pt;
994
995
#	Line 1541 | Line 1556 | void OSUAnalysis::setObjectFlags(cut &cu
1556
1557		decision = evaluateComparison(value,currentCut.comparativeOperator,currentCut.cutValue);
1558		}
1559	<	individualFlags[inputType].at(currentCutIndex).push_back(decision);
1559	>	individualFlags.at(inputType).at(currentCutIndex).push_back(decision);
1560
1561
1562		//set flags for objects that pass each cut AND all the previous cuts
1563		bool previousCumulativeFlag = true;
1564		for(uint previousCutIndex = 0; previousCutIndex != currentCutIndex; previousCutIndex++){
1565	<	if(previousCumulativeFlag && individualFlags[inputType].at(previousCutIndex).at(object)) previousCumulativeFlag = true;
1565	>	if(previousCumulativeFlag && individualFlags.at(inputType).at(previousCutIndex).at(object)) previousCumulativeFlag = true;
1566		else{ previousCumulativeFlag = false; break;}
1567		}
1568	<	cumulativeFlags[inputType].at(currentCutIndex).push_back(previousCumulativeFlag && decision);
1568	>	cumulativeFlags.at(inputType).at(currentCutIndex).push_back(previousCumulativeFlag && decision);
1569
1570		}
1571
#	Line 1558 | Line 1573 | void OSUAnalysis::setObjectFlags(cut &cu
1573
1574
1575		template <class InputCollection>
1576	<	void OSUAnalysis::fillHistogram(TH1* histo, histogram parameters, InputCollection inputCollection){
1562	<
1576	>	void OSUAnalysis::fillHistogram(TH1* histo, histogram parameters, InputCollection inputCollection,vector<bool> flags){
1577
1578		for (uint object = 0; object != inputCollection->size(); object++){
1579	+
1580	+	if(!plotAllObjectsInPassingEvents && !flags.at(object)) continue;
1581	+
1582		if(parameters.inputVariables.size() == 1){
1583		double value = valueLookup(&inputCollection->at(object), parameters.inputVariables.at(0), parameters.function);
1584		histo->Fill(value);
#	Line 1571 | Line 1588 | void OSUAnalysis::fillHistogram(TH1* his
1588		double valueY = valueLookup(&inputCollection->at(object), parameters.inputVariables.at(1), parameters.function);
1589		histo->Fill(valueX,valueY);
1590		}
1591	+
1592		}
1593
1594		}

Diff Legend

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