CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Functions
zPeak.cc File Reference
#include "FWCore/Utilities/interface/InputTag.h"
#include "DataFormats/Common/interface/Handle.h"
#include "DataFormats/PatCandidates/interface/Muon.h"
#include "PhysicsTools/FWLite/interface/EventContainer.h"
#include "PhysicsTools/FWLite/interface/CommandLineParser.h"
#include "TROOT.h"

Go to the source code of this file.

Functions

int main (int argc, char *argv[])
 

Function Documentation

int main ( int  argc,
char *  argv[] 
)

Definition at line 22 of file zPeak.cc.

References fwlite::EventContainer::add(), fwlite::EventContainer::atEnd(), fwlite::EventContainer::getByLabel(), fwlite::EventContainer::hist(), edm::HandleBase::isValid(), optutl::CommandLineParser::parseArguments(), geometryXMLtoCSV::parser, optutl::VariableMapCont::stringValue(), and fwlite::EventContainer::toBegin().

23 {
25  // ////////////////////////// //
26  // // Command Line Options // //
27  // ////////////////////////// //
29 
30 
31  // Tell people what this analysis code does and setup default options.
32  optutl::CommandLineParser parser ("Plot mass of Z candidates");
33 
35  // Change any defaults or add any new command //
36  // line options you would like here. //
38  parser.stringValue ("outputFile") = "zpeak1.root";
39 
40  // Parse the command line arguments
41  parser.parseArguments (argc, argv);
42 
44  // //////////////////////////// //
45  // // Create Event Container // //
46  // //////////////////////////// //
48 
49  // This object 'event' is used both to get all information from the
50  // event as well as to store histograms, etc.
51  fwlite::EventContainer eventCont (parser);
52 
54  // ////////////////////////////////// //
55  // // Begin Run // //
56  // // (e.g., book histograms, etc) // //
57  // ////////////////////////////////// //
59 
60  // Setup a style
61  gROOT->SetStyle ("Plain");
62 
63  // Book those histograms!
64  eventCont.add( new TH1F ("Zmass", "Candidate Z mass", 50, 20, 220) );
65 
66 
68  // //////////////// //
69  // // Event Loop // //
70  // //////////////// //
72 
73  // create labels
74  edm::InputTag muonLabel ("selectedLayer1Muons");
75 
76  for (eventCont.toBegin(); ! eventCont.atEnd(); ++eventCont)
77  {
79  // Take What We Need From Event //
82  eventCont.getByLabel (muonLabel, muonHandle);
83  assert (muonHandle.isValid());
84  vector< pat::Muon > const & muonVec = *muonHandle;
85 
86  // make sure we have enough muons so that something bad doesn't
87  // happen.
88  if (muonVec.begin() == muonVec.end())
89  {
90  // If this statement is true, then we don't have at least one
91  // muon, so don't bother. The reason we need to do this is
92  // that the logic below assumes that we have at least one and
93  // we'll end up in an infinite loop without it.
94  continue;
95  }
96 
97  // O.k. Let's loop through our muons and see what we can find.
98  const vector< pat::Muon >::const_iterator kEndIter = muonVec.end();
99  const vector< pat::Muon >::const_iterator kAlmostEndIter = kEndIter - 1;
100  for (vector< pat::Muon >::const_iterator outerIter = muonVec.begin();
101  kAlmostEndIter != outerIter;
102  ++outerIter)
103  {
104  for (vector< pat::Muon >::const_iterator innerIter = outerIter + 1;
105  kEndIter != innerIter;
106  ++innerIter)
107  {
108  // make sure that we have muons of opposite charge
109  if (outerIter->charge() * innerIter->charge() >= 0) continue;
110 
111  // if we're here then we have one positively charged muon
112  // and one negatively charged muon.
113  eventCont.hist("Zmass")->Fill( (outerIter->p4() + innerIter->p4()).M() );
114  } // for innerIter
115  } // for outerIter
116  } // for eventCont
117 
118 
120  // ////////////////// //
121  // // Clean Up Job // //
122  // ////////////////// //
124 
125  // Histograms will be automatically written to the root file
126  // specificed by command line options.
127 
128  // All done! Bye bye.
129  return 0;
130 }
bool isValid() const
Definition: HandleBase.h:76