ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/grimes/VHbbAnalysisCode/src/LightweightAnalyser.cpp

Comparing UserCode/grimes/VHbbAnalysisCode/src/LightweightAnalyser.cpp (file contents):
Revision 1.2 by grimes, Mon Feb 27 11:01:07 2012 UTC vs.
Revision 1.5 by grimes, Wed Aug 15 22:37:47 2012 UTC

#	Line 4 | Line 4
4
5		#include "TrkUpgradeAnalysis/VHbb/interface/VHbbCandidateCutSets.h"
6		#include "TrkUpgradeAnalysis/VHbb/interface/AdditionalJetStudyCutSets.h"
7	+	#include "TrkUpgradeAnalysis/VHbb/interface/tools.h" // required for trkupgradeanalysis::tools::createDirectory
8	+	#include "FWCore/Utilities/interface/EDMException.h"
9
10		#include "VHbbAnalysis/VHbbDataFormats/interface/VHbbEventAuxInfo.h"
11		#include "VHbbAnalysis/VHbbDataFormats/src/HbbCandidateFinderAlgo.cc" // Not entirely sure why this is included and not linked to
#	Line 18 | Line 20
20
21
22		const bool trkupgradeanalysis::LightweightAnalyser::verbose_=false;
23	<	const double trkupgradeanalysis::LightweightAnalyser::jetPtThresholdZ_=20;
24	<	const double trkupgradeanalysis::LightweightAnalyser::jetPtThresholdW_=20;
23	>	const double trkupgradeanalysis::LightweightAnalyser::jetPtThresholdZ_=30;
24	>	const double trkupgradeanalysis::LightweightAnalyser::jetPtThresholdW_=30;
25		const bool trkupgradeanalysis::LightweightAnalyser::useHighestPtHiggsZ_=false;
26		const bool trkupgradeanalysis::LightweightAnalyser::useHighestPtHiggsW_=true;
27
26	–	namespace // Use the unnamed namespace
27	–	{
28	–	TDirectory* createDirectory( const std::string& fullPath, TDirectory* pParent )
29	–	{
30	–	if( pParent==NULL ) throw std::runtime_error( "The parent directory is a Null pointer" );
31	–
32	–	TDirectory* pSubDirectory=pParent;
33	–	size_t currentPosition=0;
34	–	size_t nextSlash;
35	–	do
36	–	{
37	–	nextSlash=fullPath.find_first_of('/', currentPosition );
38	–	std::string directoryName=fullPath.substr(currentPosition,nextSlash-currentPosition);
39	–	currentPosition=nextSlash+1;
40	–
41	–	TDirectory* pNextSubDirectory=pSubDirectory->GetDirectory( directoryName.c_str() );
42	–	if( pNextSubDirectory==NULL ) pNextSubDirectory=pSubDirectory->mkdir( directoryName.c_str() );
43	–	if( pNextSubDirectory==NULL ) throw std::runtime_error( "Couldn't create the root directory \""+directoryName+"\"" );
44	–	pSubDirectory=pNextSubDirectory;
45	–
46	–	} while( nextSlash!=std::string::npos );
47	–
48	–	return pSubDirectory;
49	–	}
50	–	}
28
29		trkupgradeanalysis::LightweightAnalyser::LightweightAnalyser( const std::string& outputFilename )
30		: zFinder_( verbose_, jetPtThresholdZ_, useHighestPtHiggsZ_ ),
#	Line 61 | Line 38 | trkupgradeanalysis::LightweightAnalyser:
38
39
40		using namespace trkupgradeanalysis; // Don't really need this but I like to be explicit (ooh err..)
41	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZee(120)) );
42	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(120)) );
43	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionWen(120)) );
44	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionWmun(120)) );
45	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumuWithoutAdditionalJetsCut(120)) );
46	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(110)) );
47	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(115)) );
71	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(125)) );
72	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(130)) );
73	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(135)) );
74	<
75	<	//
76	<	// These are only temporary while I'm studying the number of additional jets
77	<	//
78	<	using namespace trkupgradeanalysis::additionaljetstudy;
79	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide2() ) );
80	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide3Plot1() ) );
81	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide3Plot2() ) );
82	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide3Plot3() ) );
83	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide4Plot2() ) );
84	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide5Plot1() ) );
85	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide5Plot2() ) );
86	<	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>( new WilkenSlide5Plot2AssumingTypo() ) );
87	<
88	<	//
89	<	// Don't need the following because I'm not studying the background estimation regions
90	<	//
91	<
92	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VlightRegionHWmun) );
93	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VlightRegionHWen) );
94	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VlightRegionHZmumu) );
95	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VlightRegionHZee) );
96	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new TTbarRegionHWmun) );
97	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new TTbarRegionHWen) );
98	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new TTbarRegionHZmumu) );
99	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new TTbarRegionHZee) );
100	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VbbRegionHWmun) );
101	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VbbRegionHWen) );
102	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VbbRegionHZmumu) );
103	<	//cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new VbbRegionHZee) );
41	>	//      cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZee(120)) );
42	>	//      cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionZmumu(120)) );
43	>	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new AllVariables) );
44	>	//      cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionStdGeom50PUZmumuWithoutAdditionalJetsCut(125)) );
45	>	//      cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new SignalSelectionPhase150PUZmumuWithoutAdditionalJetsCut(125)) );
46	>	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new InitialSignalSelectionZee) );
47	>	cutsToApply_.push_back( boost::shared_ptr<VHbbCandidateCutSet>(new InitialSignalSelectionZmumu) );
48		}
49
50		trkupgradeanalysis::LightweightAnalyser::~LightweightAnalyser()
#	Line 134 | Line 78 | void trkupgradeanalysis::LightweightAnal
78		void trkupgradeanalysis::LightweightAnalyser::bookHistograms( trkupgradeanalysis::LightweightAnalyser::AllPlots& plots, const std::string& directoryName )
79		{
80		TDirectory* pRootDirectory=pOutputFile_->GetDirectory("/");
81	<	TDirectory* pDirectory=createDirectory( directoryName, pRootDirectory );
81	>	TDirectory* pDirectory=trkupgradeanalysis::tools::createDirectory( directoryName, pRootDirectory );
82		//      TDirectory* pDirectory=pRootDirectory->GetDirectory( directoryName.c_str() );
83		//      if( pDirectory==NULL ) pDirectory=pRootDirectory->mkdir( directoryName.c_str() );
84
#	Line 152 | Line 96 | void trkupgradeanalysis::LightweightAnal
96		pSubDirectory=pDirectory->mkdir( "ElectronInfos" );
97		plots.allElectrons_.book( pSubDirectory );
98
99	+	pSubDirectory=pDirectory->mkdir( "HiggsCandidateJets" );
100	+	plots.higgsCandidateJets_.book( pSubDirectory );
101	+
102	+	pSubDirectory=pDirectory->mkdir( "AdditionalJets" );
103	+	plots.additionalJets_.book( pSubDirectory );
104	+
105	+	pSubDirectory=pDirectory->mkdir( "SimpleJets2" );
106	+	plots.simpleJets2_.book( pSubDirectory );
107	+
108	+	pSubDirectory=pDirectory->mkdir( "IsolationStudy" );
109	+	plots.isolationStudy_.book( pSubDirectory );
110	+
111		pSubDirectory=pDirectory->mkdir( "MonteCarloTruth" );
112		plots.monteCarloInfo_.book( pSubDirectory );
113
#	Line 168 | Line 124 | void trkupgradeanalysis::LightweightAnal
124		}
125		}
126
127	<	void trkupgradeanalysis::LightweightAnalyser::fillAllPlotsStructure( AllPlots* pPlotsToFill, const VHbbEvent& vhbbEvent, const std::vector<VHbbCandidate>& zCandidates, const std::vector<VHbbCandidate>& wCandidates )
127	>	void trkupgradeanalysis::LightweightAnalyser::fillAllPlotsStructure( AllPlots* pPlotsToFill, const VHbbEvent& vhbbEvent, const VHbbEventAuxInfo* pVHbbEventAuxInfo, const std::vector<VHbbCandidate>& zCandidates, const std::vector<VHbbCandidate>& wCandidates )
128		{
129	<	pPlotsToFill->allMuons_.fill( vhbbEvent.muInfo );
129	>	pPlotsToFill->allMuons_.fill( vhbbEvent.muInfo, pVHbbEventAuxInfo );
130		pPlotsToFill->allElectrons_.fill( vhbbEvent.eleInfo );
131	+	if( pVHbbEventAuxInfo ) pPlotsToFill->monteCarloInfo_.fill( *pVHbbEventAuxInfo );
132	+	pPlotsToFill->isolationStudy_.fill( vhbbEvent, pVHbbEventAuxInfo ); // This class checks if is null itself
133	+	pPlotsToFill->simpleJets2_.fill( vhbbEvent.simpleJets2, pVHbbEventAuxInfo );
134
135		for( std::vector<VHbbCandidate>::const_iterator iCandidate=zCandidates.begin(); iCandidate!=zCandidates.end(); ++iCandidate )
136		{
137		pPlotsToFill->candidateHistogramsForAllEvents_.fill( *iCandidate );
138	+	pPlotsToFill->higgsCandidateJets_.fill( iCandidate->H.jets, pVHbbEventAuxInfo );
139	+	pPlotsToFill->additionalJets_.fill( iCandidate->additionalJets, pVHbbEventAuxInfo );
140
141		// Loop over all of the cuts. If this candidate passes the cuts then fill the plots
142		for( std::vector<CutSetPlotSet>::iterator iCutPlotSet=pPlotsToFill->cutCollectionPlotSets_.begin();
#	Line 209 | Line 170 | void trkupgradeanalysis::LightweightAnal
170		if( pVHbbEventAuxInfo == 0 ) std::cerr << "Couldn't get the VHbbEventAuxInfo" << std::endl;
171		else
172		{
212	–	pPlotsForEverything_->allEventTypes_.monteCarloInfo_.fill( *pVHbbEventAuxInfo );
173		std::string eventType=eventTypeFromMC( *pVHbbEventAuxInfo );
174
175		std::map<std::string,AllPlots>::iterator iPlotsForThisEventType=pPlotsForEverything_->mcEventTypePlots_.find(eventType);
#	Line 226 | Line 186 | void trkupgradeanalysis::LightweightAnal
186
187		++numberOfEvents_;
188
189	+	bool applyMuonIsolation=false;
190	+	// Apply isolation to the muons. Need to modify so use a very naughty const cast
191	+	if( applyMuonIsolation )
192	+	{
193	+	std::vector<VHbbEvent::MuonInfo>& muonInfos=const_cast< std::vector<VHbbEvent::MuonInfo>& >( pVHbbEvent->muInfo );
194	+	bool allMuonsPassIsolationCut=false;
195	+	while( !allMuonsPassIsolationCut )
196	+	{
197	+	// I'm worried about looping using invalidated iterators, so whenever I erase an element I'll
198	+	// start looping from the start. Not the most efficient way but these vectors typically have
199	+	// at most 3 entries.
200	+	std::vector<VHbbEvent::MuonInfo>::iterator iMuon;
201	+	for( iMuon=muonInfos.begin(); iMuon!=muonInfos.end(); ++iMuon )
202	+	{
203	+	float isolation=MuonInfoPlotSet::thirdRho25CorrectedIsolation( *iMuon, pVHbbEventAuxInfo->puInfo.rho25 );
204	+	if( isolation>=0.15 ) break;
205	+	}
206	+	if( iMuon!=muonInfos.end() ) muonInfos.erase( iMuon );
207	+	else allMuonsPassIsolationCut=true;
208	+	}
209	+	}
210	+
211		zFinder_.run( pVHbbEvent, zCandidates );
212		wFinder_.run( pVHbbEvent, wCandidates );
213
214	<	fillAllPlotsStructure( &pPlotsForEverything_->allEventTypes_, *pVHbbEvent, zCandidates, wCandidates );
215	<	if( pCurrentPlotsForThisEventType!=NULL ) fillAllPlotsStructure( pCurrentPlotsForThisEventType, *pVHbbEvent, zCandidates, wCandidates );
214	>	fillAllPlotsStructure( &pPlotsForEverything_->allEventTypes_, *pVHbbEvent, pVHbbEventAuxInfo, zCandidates, wCandidates );
215	>	if( pCurrentPlotsForThisEventType!=NULL ) fillAllPlotsStructure( pCurrentPlotsForThisEventType, *pVHbbEvent, pVHbbEventAuxInfo, zCandidates, wCandidates );
216
217
218		if( !zCandidates.empty() || !wCandidates.empty() ) ++numberOfEventsWithAtLeastOneZOrWCandidate_;
219
220		}
221
222	<	void trkupgradeanalysis::LightweightAnalyser::processFile( TFile* pInputFile )
222	>	void trkupgradeanalysis::LightweightAnalyser::processFile( TFile* pInputFile, size_t numberOfEvents )
223		{
224		if( pPlotsForEverything_==NULL ) changeDirectory("allEvents");
225
226		fwlite::Event event( pInputFile );
227	+
228	+	size_t eventCount=0;
229		for( event.toBegin(); !event.atEnd(); ++event )
230		{
231	<	processEvent( event );
231	>	try
232	>	{
233	>	processEvent( event );
234	>	++eventCount;
235	>	}
236	>	catch( edm::Exception& error )
237	>	{
238	>	std::cerr << error.what() << std::endl;
239	>	}
240	>
241	>	if( numberOfEvents!=0 && eventCount>=numberOfEvents )
242	>	{
243	>	std::cout << "Breaking after " << eventCount << " events" << std::endl;
244	>	break;
245	>	}
246		}
247
248		std::cout << "There were " << numberOfEventsWithAtLeastOneZOrWCandidate_ << " events with a Z or W candidate out of " << numberOfEvents_ << " events." << std::endl;

Diff Legend

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