CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_5_3_13_patch3/src/ElectroWeakAnalysis/ZMuMu/plugins/ZToLLEdmNtupleDumper.cc

Go to the documentation of this file.
00001 /* \class ZToLLEdmNtupleDumper
00002  *
00003  * \author Luca Lista, INFN
00004  *
00005  */
00006 #include "FWCore/Framework/interface/EDProducer.h"
00007 #include "FWCore/Utilities/interface/InputTag.h"
00008 #include "DataFormats/Common/interface/Handle.h"
00009 #include "FWCore/Framework/interface/Event.h"
00010 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00011 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00012 #include "DataFormats/Candidate/interface/CandMatchMap.h"
00013 #include "DataFormats/Candidate/interface/Candidate.h"
00014 //#include "DataFormats/Common/interface/AssociationVector.h"
00015 #include "DataFormats/PatCandidates/interface/Muon.h"
00016 #include "DataFormats/PatCandidates/interface/GenericParticle.h"
00017 #include "FWCore/Utilities/interface/EDMException.h"
00018 #include "DataFormats/PatCandidates/interface/TriggerObjectStandAlone.h"
00019 
00020 #include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
00021 #include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
00022 #include "DataFormats/PatCandidates/interface/Isolation.h"
00023 #include "DataFormats/Common/interface/ValueMap.h"
00024 #include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
00025 #include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
00026 
00027 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00028 
00029 #include <vector>
00030 
00031 using namespace edm;
00032 using namespace std;
00033 using namespace reco;
00034 using namespace isodeposit;
00035 //using namespace pat;
00036 
00037 class ZToLLEdmNtupleDumper : public edm::EDProducer {
00038 public:
00039   typedef math::XYZVector Vector;
00040   ZToLLEdmNtupleDumper( const edm::ParameterSet & );
00041    
00042 private:
00043   void produce( edm::Event &, const edm::EventSetup & );
00044     std::vector<std::string> zName_;
00045   std::vector<edm::InputTag> z_, zGenParticlesMatch_ ;
00046   edm::InputTag beamSpot_,  primaryVertices_;
00047 
00048   std::vector<double> ptThreshold_, etEcalThreshold_, etHcalThreshold_ ,dRVetoTrk_, dRTrk_, dREcal_ , dRHcal_,  alpha_, beta_; 
00049   std::vector<double> relativeIsolation_;
00050   std::vector<string> hltPath_;
00051   int counter;
00052   
00053   
00054 };
00055 
00056 template<typename T>
00057   double isolation(const T * t, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta, bool relativeIsolation) {
00058    // on 34X: 
00059 const pat::IsoDeposit * trkIso = t->isoDeposit(pat::TrackIso);
00060 //  const pat::IsoDeposit * trkIso = t->trackerIsoDeposit();
00061    // on 34X 
00062 const pat::IsoDeposit * ecalIso = t->isoDeposit(pat::EcalIso);
00063 //  const pat::IsoDeposit * ecalIso = t->ecalIsoDeposit();
00064 //    on 34X 
00065 const pat::IsoDeposit * hcalIso = t->isoDeposit(pat::HcalIso);   
00066 //    const pat::IsoDeposit * hcalIso = t->hcalIsoDeposit();
00067 
00068     Direction dir = Direction(t->eta(), t->phi());
00069    
00070     pat::IsoDeposit::AbsVetos vetosTrk;
00071     vetosTrk.push_back(new ConeVeto( dir, dRVetoTrk ));
00072     vetosTrk.push_back(new ThresholdVeto( ptThreshold ));
00073     
00074     pat::IsoDeposit::AbsVetos vetosEcal;
00075     vetosEcal.push_back(new ConeVeto( dir, 0.));
00076     vetosEcal.push_back(new ThresholdVeto( etEcalThreshold ));
00077     
00078     pat::IsoDeposit::AbsVetos vetosHcal;
00079     vetosHcal.push_back(new ConeVeto( dir, 0. ));
00080     vetosHcal.push_back(new ThresholdVeto( etHcalThreshold ));
00081 
00082     double isovalueTrk = (trkIso->sumWithin(dRTrk,vetosTrk));
00083     double isovalueEcal = (ecalIso->sumWithin(dREcal,vetosEcal));
00084     double isovalueHcal = (hcalIso->sumWithin(dRHcal,vetosHcal));
00085     
00086 
00087     double iso = alpha*( ((1+beta)/2*isovalueEcal) + ((1-beta)/2*isovalueHcal) ) + ((1-alpha)*isovalueTrk) ;
00088     if(relativeIsolation) iso /= t->pt();
00089     return iso;
00090   }
00091 
00092 
00093 double candIsolation( const reco::Candidate* c, double ptThreshold, double etEcalThreshold, double etHcalThreshold , double dRVetoTrk, double dRTrk, double dREcal , double dRHcal,  double alpha, double beta, bool relativeIsolation) {
00094     const pat::Muon * mu = dynamic_cast<const pat::Muon *>(c);
00095     if(mu != 0) return isolation(mu, ptThreshold, etEcalThreshold, etHcalThreshold ,dRVetoTrk, dRTrk, dREcal , dRHcal,  alpha, beta, relativeIsolation);
00096     const pat::GenericParticle * trk = dynamic_cast<const pat::GenericParticle*>(c);
00097     if(trk != 0) return isolation(trk,  ptThreshold, etEcalThreshold, etHcalThreshold ,dRVetoTrk, dRTrk, dREcal , dRHcal,  alpha, beta, relativeIsolation);
00098     throw edm::Exception(edm::errors::InvalidReference) 
00099       << "Candidate daughter #0 is neither pat::Muons nor pat::GenericParticle\n";      
00100     return -1;
00101   }
00102 
00103 
00104 
00105 ZToLLEdmNtupleDumper::ZToLLEdmNtupleDumper( const ParameterSet & cfg ) {
00106   string alias;
00107   vector<ParameterSet> psets = cfg.getParameter<vector<ParameterSet> > ( "zBlocks" );
00108   for( std::vector<edm::ParameterSet>::const_iterator i = psets.begin(); i != psets.end(); ++ i ) {
00109     string zName = i->getParameter<string>( "zName" );
00110     InputTag z =  i->getParameter<InputTag>( "z" );
00111     InputTag zGenParticlesMatch = i->getParameter<InputTag>( "zGenParticlesMatch" );
00112     beamSpot_ =  i->getParameter<InputTag>( "beamSpot" );
00113     primaryVertices_= i->getParameter<InputTag>( "primaryVertices" ) ;
00114       double ptThreshold = i->getParameter<double>("ptThreshold");
00115     double etEcalThreshold = i->getParameter<double>("etEcalThreshold");
00116     double etHcalThreshold= i->getParameter<double>("etHcalThreshold");
00117     double dRVetoTrk=i->getParameter<double>("deltaRVetoTrk");
00118     double dRTrk=i->getParameter<double>("deltaRTrk");
00119     double dREcal=i->getParameter<double>("deltaREcal");
00120     double dRHcal=i->getParameter<double>("deltaRHcal");
00121     double alpha=i->getParameter<double>("alpha");
00122     double beta=i->getParameter<double>("beta");
00123     bool relativeIsolation = i->getParameter<bool>("relativeIsolation");
00124     string hltPath = i ->getParameter<std::string >("hltPath");
00125     zName_.push_back( zName );
00126     z_.push_back( z );
00127     zGenParticlesMatch_.push_back( zGenParticlesMatch );
00128       ptThreshold_.push_back( ptThreshold );
00129     etEcalThreshold_.push_back(etEcalThreshold);  
00130     etHcalThreshold_.push_back(etHcalThreshold);
00131     dRVetoTrk_.push_back(dRVetoTrk);
00132     dRTrk_.push_back(dRTrk);
00133     dREcal_.push_back(dREcal);
00134     dRHcal_.push_back(dRHcal);
00135     alpha_.push_back(alpha);
00136     beta_.push_back(beta);
00137     relativeIsolation_.push_back(relativeIsolation);
00138     hltPath_.push_back(hltPath);
00139     produces<vector<unsigned int> >( alias = zName + "EventNumber" ).setBranchAlias( alias );
00140     produces<vector<unsigned int> >( alias = zName + "RunNumber" ).setBranchAlias( alias );
00141     produces<vector<unsigned int> >( alias = zName + "LumiBlock" ).setBranchAlias( alias );
00142     produces<vector<float> >( alias = zName + "Mass" ).setBranchAlias( alias );
00143     produces<vector<float> >( alias = zName + "MassSa" ).setBranchAlias( alias );
00144     produces<vector<float> >( alias = zName + "Pt" ).setBranchAlias( alias );
00145     produces<vector<float> >( alias = zName + "Eta" ).setBranchAlias( alias );
00146     produces<vector<float> >( alias = zName + "Phi" ).setBranchAlias( alias );
00147     produces<vector<float> >( alias = zName + "Y" ).setBranchAlias( alias );
00148     produces<vector<float> >( alias = zName + "Dau1Pt" ).setBranchAlias( alias );
00149     produces<vector<float> >( alias = zName + "Dau2Pt" ).setBranchAlias( alias );
00150     produces<vector<float> >( alias = zName + "Dau1SaPt" ).setBranchAlias( alias );
00151     produces<vector<float> >( alias = zName + "Dau2SaPt" ).setBranchAlias( alias );
00152     produces<vector<unsigned int> >( alias = zName + "Dau1HLTBit" ).setBranchAlias( alias );
00153     produces<vector<unsigned int> >( alias = zName + "Dau2HLTBit" ).setBranchAlias( alias );
00154     produces<vector<int> >( alias = zName + "Dau1Q" ).setBranchAlias( alias );
00155     produces<vector<int> >( alias = zName + "Dau2Q" ).setBranchAlias( alias );
00156     produces<vector<float> >( alias = zName + "Dau1Eta" ).setBranchAlias( alias );
00157     produces<vector<float> >( alias = zName + "Dau2Eta" ).setBranchAlias( alias );
00158     produces<vector<float> >( alias = zName + "Dau1SaEta" ).setBranchAlias( alias );
00159     produces<vector<float> >( alias = zName + "Dau2SaEta" ).setBranchAlias( alias );
00160     produces<vector<float> >( alias = zName + "Dau1Phi" ).setBranchAlias( alias );
00161     produces<vector<float> >( alias = zName + "Dau2Phi" ).setBranchAlias( alias );
00162     produces<vector<float> >( alias = zName + "Dau1SaPhi" ).setBranchAlias( alias );
00163     produces<vector<float> >( alias = zName + "Dau2SaPhi" ).setBranchAlias( alias );
00164     produces<vector<float> >( alias = zName + "Dau1Iso" ).setBranchAlias( alias );
00165     produces<vector<float> >( alias = zName + "Dau2Iso" ).setBranchAlias( alias );
00166     produces<vector<float> >( alias = zName + "Dau1TrkIso" ).setBranchAlias( alias );
00167     produces<vector<float> >( alias = zName + "Dau2TrkIso" ).setBranchAlias( alias );
00168     produces<vector<float> >( alias = zName + "Dau1EcalIso" ).setBranchAlias( alias );
00169     produces<vector<float> >( alias = zName + "Dau2EcalIso" ).setBranchAlias( alias );
00170     produces<vector<float> >( alias = zName + "Dau1HcalIso" ).setBranchAlias( alias );
00171     produces<vector<float> >( alias = zName + "Dau2HcalIso" ).setBranchAlias( alias );
00172     produces<vector<float> >( alias = zName + "Dau1MuEnergyEm" ).setBranchAlias( alias );
00173     produces<vector<float> >( alias = zName + "Dau1MuEnergyHad" ).setBranchAlias( alias );
00174     produces<vector<float> >( alias = zName + "Dau2MuEnergyEm" ).setBranchAlias( alias );
00175     produces<vector<float> >( alias = zName + "Dau2MuEnergyHad" ).setBranchAlias( alias );
00176 
00177     produces<vector<float> >( alias = zName + "VtxNormChi2" ).setBranchAlias( alias );
00178     produces<vector<unsigned int> >( alias = zName + "Dau1NofHit" ).setBranchAlias( alias );
00179     produces<vector<unsigned int> >( alias = zName + "Dau2NofHit" ).setBranchAlias( alias );
00180     produces<vector<unsigned int> >( alias = zName + "Dau1NofHitTk" ).setBranchAlias( alias );
00181     produces<vector<unsigned int> >( alias = zName + "Dau1NofHitSta" ).setBranchAlias( alias );
00182     produces<vector<unsigned int> >( alias = zName + "Dau2NofHitTk" ).setBranchAlias( alias );
00183     produces<vector<unsigned int> >( alias = zName + "Dau2NofHitSta" ).setBranchAlias( alias );
00184     produces<vector<unsigned int> >( alias = zName + "Dau1NofMuChambers" ).setBranchAlias( alias );
00185     produces<vector<unsigned int> >( alias = zName + "Dau2NofMuChambers" ).setBranchAlias( alias );
00186     produces<vector<unsigned int> >( alias = zName + "Dau1NofMuMatches" ).setBranchAlias( alias );
00187     produces<vector<unsigned int> >( alias = zName + "Dau2NofMuMatches" ).setBranchAlias( alias );
00188     produces<vector<float> >( alias = zName + "Dau1Chi2" ).setBranchAlias( alias );
00189     produces<vector<float> >( alias = zName + "Dau2Chi2" ).setBranchAlias( alias );
00190     produces<vector<float> >( alias = zName + "Dau1TrkChi2" ).setBranchAlias( alias );
00191     produces<vector<float> >( alias = zName + "Dau2TrkChi2" ).setBranchAlias( alias );
00192     produces<vector<float> >( alias = zName + "Dau1dxyFromBS" ).setBranchAlias( alias );
00193     produces<vector<float> >( alias = zName + "Dau2dxyFromBS" ).setBranchAlias( alias );
00194     produces<vector<float> >( alias = zName + "Dau1dzFromBS" ).setBranchAlias( alias );
00195     produces<vector<float> >( alias = zName + "Dau2dzFromBS" ).setBranchAlias( alias );
00196     produces<vector<float> >( alias = zName + "Dau1dxyFromPV" ).setBranchAlias( alias );
00197     produces<vector<float> >( alias = zName + "Dau2dxyFromPV" ).setBranchAlias( alias );
00198     produces<vector<float> >( alias = zName + "Dau1dzFromPV" ).setBranchAlias( alias );
00199     produces<vector<float> >( alias = zName + "Dau2dzFromPV" ).setBranchAlias( alias );
00200     produces<vector<float> >( alias = zName + "TrueMass" ).setBranchAlias( alias );
00201     produces<vector<float> >( alias = zName + "TruePt" ).setBranchAlias( alias );
00202     produces<vector<float> >( alias = zName + "TrueEta" ).setBranchAlias( alias );
00203     produces<vector<float> >( alias = zName + "TruePhi" ).setBranchAlias( alias );
00204     produces<vector<float> >( alias = zName + "TrueY" ).setBranchAlias( alias );
00205   }
00206 }
00207 
00208 
00209 
00210 void ZToLLEdmNtupleDumper::produce( Event & evt, const EventSetup & ) {
00211     Handle<reco::BeamSpot> beamSpotHandle;
00212   if (!evt.getByLabel(beamSpot_, beamSpotHandle)) {
00213     std::cout << ">>> No beam spot found !!!"<<std::endl;
00214   }
00215   Handle<reco::VertexCollection> primaryVertices;  // Collection of primary Vertices
00216   if (!evt.getByLabel(primaryVertices_, primaryVertices)){
00217     std::cout << ">>> No primary verteces  found !!!"<<std::endl;
00218   }
00219     
00220   unsigned int size = z_.size();
00221   for( unsigned int c = 0; c < size; ++ c ) {
00222     Handle<CandidateView> zColl;
00223     evt.getByLabel( z_[c], zColl );
00224     bool isMCMatchTrue = false;  
00225     //if (zGenParticlesMatch_[c] != "")  isMCMatchTrue = true;     
00226     Handle<GenParticleMatch> zGenParticlesMatch;
00227     if (evt.getByLabel( zGenParticlesMatch_[c], zGenParticlesMatch )) {
00228       isMCMatchTrue=true;
00229     }
00230     unsigned int zSize = zColl->size();
00231     auto_ptr<vector<unsigned int> > event( new vector<unsigned int> );
00232     auto_ptr<vector<unsigned int> > run( new vector<unsigned int> );
00233     auto_ptr<vector<unsigned int> > lumi( new vector<unsigned int > );
00234     auto_ptr<vector<float> > zMass( new vector<float> );
00235     auto_ptr<vector<float> > zMassSa( new vector<float> );
00236     auto_ptr<vector<float> > zPt( new vector<float> );
00237     auto_ptr<vector<float> > zEta( new vector<float> );
00238     auto_ptr<vector<float> > zPhi( new vector<float> );
00239     auto_ptr<vector<float> > zY( new vector<float> );
00240     auto_ptr<vector<float> > zDau1Pt( new vector<float> );
00241     auto_ptr<vector<float> > zDau2Pt( new vector<float> );
00242     auto_ptr<vector<float> > zDau1SaPt( new vector<float> );
00243     auto_ptr<vector<float> > zDau2SaPt( new vector<float> );
00244     auto_ptr<vector<unsigned int> > zDau1HLTBit( new vector<unsigned int> );
00245     auto_ptr<vector<unsigned int> > zDau2HLTBit( new vector<unsigned int> );
00246     auto_ptr<vector<int> > zDau1Q( new vector<int> );
00247     auto_ptr<vector<int> > zDau2Q( new vector<int> );
00248     auto_ptr<vector<float> > zDau1Eta( new vector<float> );
00249     auto_ptr<vector<float> > zDau2Eta( new vector<float> );
00250     auto_ptr<vector<float> > zDau1SaEta( new vector<float> );
00251     auto_ptr<vector<float> > zDau2SaEta( new vector<float> );
00252     auto_ptr<vector<float> > zDau1Phi( new vector<float> );
00253     auto_ptr<vector<float> > zDau2Phi( new vector<float> );
00254     auto_ptr<vector<float> > zDau1SaPhi( new vector<float> );
00255     auto_ptr<vector<float> > zDau2SaPhi( new vector<float> );
00256     auto_ptr<vector<float> > zDau1Iso( new vector<float> );
00257     auto_ptr<vector<float> > zDau2Iso( new vector<float> );
00258     auto_ptr<vector<float> > zDau1TrkIso( new vector<float> );
00259     auto_ptr<vector<float> > zDau2TrkIso( new vector<float> );
00260     auto_ptr<vector<float> > zDau1EcalIso( new vector<float> );
00261     auto_ptr<vector<float> > zDau2EcalIso( new vector<float> );
00262     auto_ptr<vector<float> > zDau1HcalIso( new vector<float> );
00263     auto_ptr<vector<float> > zDau2HcalIso( new vector<float> );
00264     auto_ptr<vector<float> > zDau1MuEnergyEm( new vector<float> );
00265     auto_ptr<vector<float> > zDau2MuEnergyEm( new vector<float> );
00266     auto_ptr<vector<float> > zDau1MuEnergyHad( new vector<float> );
00267     auto_ptr<vector<float> > zDau2MuEnergyHad( new vector<float> );
00268     auto_ptr<vector<float> > vtxNormChi2( new vector<float> );
00269     auto_ptr<vector<unsigned int> > zDau1NofHit( new vector<unsigned int> );
00270     auto_ptr<vector<unsigned int> > zDau2NofHit( new vector<unsigned int> );
00271     auto_ptr<vector<unsigned int> > zDau1NofHitTk( new vector<unsigned int> );
00272     auto_ptr<vector<unsigned int> > zDau2NofHitTk( new vector<unsigned int> ); 
00273     auto_ptr<vector<unsigned int> > zDau1NofHitSta( new vector<unsigned int> );
00274     auto_ptr<vector<unsigned int> > zDau2NofHitSta( new vector<unsigned int> );
00275     auto_ptr<vector<unsigned int> > zDau1NofMuChambers( new vector<unsigned int> );
00276     auto_ptr<vector<unsigned int> > zDau2NofMuChambers( new vector<unsigned int> );
00277     auto_ptr<vector<unsigned int> > zDau1NofMuMatches( new vector<unsigned int> );
00278     auto_ptr<vector<unsigned int> > zDau2NofMuMatches( new vector<unsigned int> );
00279     auto_ptr<vector<float> > zDau1Chi2( new vector<float> );
00280     auto_ptr<vector<float> > zDau2Chi2( new vector<float> );
00281     auto_ptr<vector<float> > zDau1TrkChi2( new vector<float> );
00282     auto_ptr<vector<float> > zDau2TrkChi2( new vector<float> );
00283     auto_ptr<vector<float> > zDau1dxyFromBS( new vector<float> );
00284     auto_ptr<vector<float> > zDau2dxyFromBS( new vector<float> );
00285     auto_ptr<vector<float> > zDau1dzFromBS( new vector<float> );
00286     auto_ptr<vector<float> > zDau2dzFromBS( new vector<float> );
00287     auto_ptr<vector<float> > zDau1dxyFromPV( new vector<float> );
00288     auto_ptr<vector<float> > zDau2dxyFromPV( new vector<float> );
00289     auto_ptr<vector<float> > zDau1dzFromPV( new vector<float> );
00290     auto_ptr<vector<float> > zDau2dzFromPV( new vector<float> );
00291     auto_ptr<vector<float> > trueZMass( new vector<float> );
00292     auto_ptr<vector<float> > trueZPt( new vector<float> );
00293     auto_ptr<vector<float> > trueZEta( new vector<float> );
00294     auto_ptr<vector<float> > trueZPhi( new vector<float> );
00295     auto_ptr<vector<float> > trueZY( new vector<float> );
00296     event -> push_back(evt.id().event());
00297     run -> push_back(evt.id().run());
00298     lumi -> push_back(evt.luminosityBlock());
00299     for( unsigned int i = 0; i < zSize; ++ i ) {
00300       const Candidate & z = (*zColl)[ i ];
00301       CandidateBaseRef zRef = zColl->refAt(i);
00302       zMass->push_back( z.mass() );
00303       zPt->push_back( z.pt() );
00304       zEta->push_back( z.eta() );
00305       zPhi->push_back( z.phi() );
00306       zY->push_back( z.rapidity() );
00307       vtxNormChi2->push_back(z.vertexNormalizedChi2() );
00308       const Candidate * dau1 = z.daughter(0); 
00309       const Candidate * dau2 = z.daughter(1); 
00310       zDau1Pt->push_back( dau1->pt() );
00311       zDau2Pt->push_back( dau2->pt() );
00312       zDau1Q->push_back( dau1->charge() );
00313       zDau2Q->push_back( dau2->charge() );
00314       zDau1Eta->push_back( dau1->eta() );
00315       zDau2Eta->push_back( dau2->eta() );
00316       zDau1Phi->push_back( dau1->phi() );
00317       zDau2Phi->push_back( dau2->phi() );      
00318       if(!(dau1->hasMasterClone()&&dau2->hasMasterClone()))
00319         throw edm::Exception(edm::errors::InvalidReference) 
00320           << "Candidate daughters have no master clone\n"; 
00321       const CandidateBaseRef & mr1 = dau1->masterClone(), & mr2 = dau2->masterClone();
00322        
00323        const Candidate * m1 = &*mr1, * m2 = &*mr2;
00324 
00325       // isolation as defined by us into the analyzer
00326       double iso1 = candIsolation(m1,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c],  alpha_[c], beta_[c], relativeIsolation_[c]);
00327       double iso2 = candIsolation(m2,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c],  alpha_[c], beta_[c], relativeIsolation_[c] );
00328       // tracker isolation : alpha =0
00329       double trkIso1 = candIsolation(m1,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c], 0.0, beta_[c], relativeIsolation_[c]);
00330       double trkIso2 = candIsolation(m2,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c],  0.0, beta_[c], relativeIsolation_[c] );
00331       // ecal isolation : alpha =1, beta =1
00332       double ecalIso1 = candIsolation(m1,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c], 1.0, 1.0, relativeIsolation_[c]);
00333       double ecalIso2 = candIsolation(m2,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c],  1.0, 1.0, relativeIsolation_[c] );
00334       // hcal isolation : alpha =1, beta =-1
00335       double hcalIso1 = candIsolation(m1,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c], 1.0, -1.0, relativeIsolation_[c]);
00336       double hcalIso2 = candIsolation(m2,ptThreshold_[c], etEcalThreshold_[c], etHcalThreshold_[c] ,dRVetoTrk_[c], dRTrk_[c], dREcal_[c] , dRHcal_[c],  1.0, -1.0, relativeIsolation_[c] );
00337 
00338       zDau1Iso->push_back( iso1 );
00339       zDau2Iso->push_back( iso2 );
00340       zDau1TrkIso->push_back( trkIso1 );
00341       zDau2TrkIso->push_back( trkIso2 );
00342       zDau1EcalIso->push_back( ecalIso1 );
00343       zDau2EcalIso->push_back( ecalIso2 );
00344       zDau1HcalIso->push_back( hcalIso1 );
00345       zDau2HcalIso->push_back( hcalIso2 );
00346       if (isMCMatchTrue){
00347         GenParticleRef trueZRef  = (*zGenParticlesMatch)[zRef]; 
00348         //CandidateRef trueZRef = trueZIter->val;
00349         if( trueZRef.isNonnull() ) {
00350           const Candidate & z = * trueZRef;
00351           trueZMass->push_back( z.mass() );
00352           trueZPt->push_back( z.pt() );
00353           trueZEta->push_back( z.eta() );
00354           trueZPhi->push_back( z.phi() );
00355           trueZY->push_back( z.rapidity() );
00356         } else {
00357           trueZMass->push_back( -100 ); 
00358           trueZPt->push_back( -100 );
00359           trueZEta->push_back( -100 );
00360           trueZPhi->push_back( -100 );
00361           trueZY->push_back( -100 );  
00362         }
00363       }      
00364       // quality variables 
00365       const pat::Muon * mu1 = dynamic_cast<const pat::Muon*>(m1);
00366         // protection for standalone and trackerMuon
00367       if (mu1->isGlobalMuon() == true){
00368         zDau1NofHit->push_back(mu1->numberOfValidHits());
00369         zDau1NofHitTk->push_back(mu1->innerTrack()->numberOfValidHits());
00370         zDau1NofHitSta->push_back(mu1->outerTrack()->numberOfValidHits());
00371         zDau1Chi2->push_back(mu1->normChi2());
00372         TrackRef mu1TrkRef = mu1->innerTrack();
00373         zDau1TrkChi2->push_back( mu1TrkRef->normalizedChi2());
00374         zDau1dxyFromBS->push_back(mu1TrkRef->dxy(beamSpotHandle->position()));
00375         zDau1dzFromBS->push_back(mu1TrkRef->dz(beamSpotHandle->position()));
00376         zDau1dxyFromPV->push_back(mu1TrkRef->dxy(primaryVertices->begin()->position() ));
00377         zDau1dzFromPV->push_back(mu1TrkRef->dz(primaryVertices->begin()->position() ));     
00378         zDau1MuEnergyEm->push_back( mu1->calEnergy().em);
00379         zDau1MuEnergyHad->push_back( mu1->calEnergy().had);
00380 
00381       } else  if (mu1->isStandAloneMuon() == true) {
00382         // the muon is a standalone
00383         TrackRef mu1StaRef = mu1->outerTrack(); 
00384         zDau1NofHit->push_back(mu1StaRef->numberOfValidHits());
00385         zDau1NofHitTk->push_back(0);
00386         zDau1NofHitSta->push_back(mu1StaRef->numberOfValidHits());
00387         zDau1Chi2->push_back(mu1StaRef->normalizedChi2()); 
00388         zDau1TrkChi2->push_back(0);
00389         zDau1dxyFromBS->push_back(mu1StaRef->dxy(beamSpotHandle->position()));
00390         zDau1dzFromBS->push_back(mu1StaRef->dz(beamSpotHandle->position()));
00391         zDau1dxyFromPV->push_back(mu1StaRef->dxy(primaryVertices->begin()->position() ));
00392         zDau1dzFromPV->push_back(mu1StaRef->dz(primaryVertices->begin()->position() ));     
00393         zDau1MuEnergyEm->push_back( -1);
00394         zDau1MuEnergyHad->push_back( -1);
00395       } else  if (mu1->isTrackerMuon() == true) {
00396         // the muon is a trackerMuon
00397         TrackRef mu1TrkRef = mu1->innerTrack(); 
00398         zDau1NofHit->push_back(mu1TrkRef->numberOfValidHits());
00399         zDau1NofHitTk->push_back(mu1TrkRef->numberOfValidHits());
00400         zDau1NofHitSta->push_back(0);
00401         zDau1Chi2->push_back(mu1TrkRef->normalizedChi2()); 
00402         zDau1TrkChi2->push_back(mu1TrkRef->normalizedChi2());
00403         zDau1dxyFromBS->push_back(mu1TrkRef->dxy(beamSpotHandle->position()));
00404         zDau1dzFromBS->push_back(mu1TrkRef->dz(beamSpotHandle->position()));
00405         zDau1dxyFromPV->push_back(mu1TrkRef->dxy(primaryVertices->begin()->position() ));
00406         zDau1dzFromPV->push_back(mu1TrkRef->dz(primaryVertices->begin()->position() ));     
00407         zDau1MuEnergyEm->push_back( mu1->calEnergy().em);
00408         zDau1MuEnergyHad->push_back( mu1->calEnergy().had);
00409       }
00410       zDau1NofMuChambers->push_back(mu1->numberOfChambers());
00411       zDau1NofMuMatches->push_back(mu1->numberOfMatches());
00412 
00413       // would we like to add another variables??? 
00414       // HLT trigger  bit
00415       const pat::TriggerObjectStandAloneCollection mu1HLTMatches =  mu1->triggerObjectMatchesByPath( hltPath_[c] );
00416 
00417       int dimTrig1 = mu1HLTMatches.size();
00418       if(dimTrig1 !=0 ){
00419         zDau1HLTBit->push_back(1);
00420       } else {
00421         zDau1HLTBit->push_back(0); 
00422       }
00423       const pat::Muon * mu2 = dynamic_cast<const pat::Muon*>(m2);
00424       if (mu2!=0 ) {
00425         if (mu2->isGlobalMuon() == true) {
00426           zDau2NofHit->push_back(mu2->numberOfValidHits());
00427           zDau2NofHitTk->push_back(mu2->innerTrack()->numberOfValidHits());
00428           zDau2NofHitSta->push_back(mu2->outerTrack()->numberOfValidHits());
00429           zDau2Chi2->push_back(mu2->normChi2());
00430           TrackRef mu2TrkRef = mu2->innerTrack();
00431           zDau1TrkChi2->push_back( mu2TrkRef->normalizedChi2());
00432           zDau2dxyFromBS->push_back(mu2TrkRef->dxy(beamSpotHandle->position()));
00433           zDau2dzFromBS->push_back(mu2TrkRef->dz(beamSpotHandle->position()));
00434           zDau2dxyFromPV->push_back(mu2TrkRef->dxy(primaryVertices->begin()->position() ));
00435           zDau2dzFromPV->push_back(mu2TrkRef->dz(primaryVertices->begin()->position() ));  
00436           zDau2MuEnergyEm->push_back( mu2->calEnergy().em);
00437           zDau2MuEnergyHad->push_back( mu2->calEnergy().had); 
00438         } else if (mu2->isStandAloneMuon() == true){
00439           // its' a standalone
00440           zDau2HLTBit->push_back(0);
00441           TrackRef mu2StaRef = mu2->outerTrack(); 
00442           zDau2NofHit->push_back(mu2StaRef->numberOfValidHits());
00443           zDau2NofHitTk->push_back(0);
00444           zDau2NofHitSta->push_back(mu2StaRef->numberOfValidHits());
00445           zDau2Chi2->push_back(mu2StaRef->normalizedChi2()); 
00446           zDau2TrkChi2->push_back(0);
00447           zDau2dxyFromBS->push_back(mu2StaRef->dxy(beamSpotHandle->position()));
00448           zDau2dzFromBS->push_back(mu2StaRef->dz(beamSpotHandle->position()));
00449           zDau2dxyFromPV->push_back(mu2StaRef->dxy(primaryVertices->begin()->position() ));
00450           zDau2dzFromPV->push_back(mu2StaRef->dz(primaryVertices->begin()->position() ));     
00451           zDau1MuEnergyEm->push_back( -1);
00452           zDau1MuEnergyHad->push_back(  -1);
00453         } else  if (mu2->isTrackerMuon() == true) {
00454           // the muon is a trackerMuon
00455           TrackRef mu2TrkRef = mu2->innerTrack();
00456           zDau2NofHit->push_back(mu2TrkRef->numberOfValidHits());
00457           zDau2NofHitSta->push_back(0);
00458           zDau2NofHitTk->push_back(mu2TrkRef->numberOfValidHits());
00459           zDau2Chi2->push_back(mu2TrkRef->normalizedChi2()); 
00460           zDau2TrkChi2->push_back(mu2TrkRef->normalizedChi2());
00461           zDau2dxyFromBS->push_back(mu2TrkRef->dxy(beamSpotHandle->position()));
00462           zDau2dzFromBS->push_back(mu2TrkRef->dz(beamSpotHandle->position()));
00463           zDau2dxyFromPV->push_back(mu2TrkRef->dxy(primaryVertices->begin()->position() ));
00464           zDau2dzFromPV->push_back(mu2TrkRef->dz(primaryVertices->begin()->position() ));     
00465           zDau2MuEnergyEm->push_back( mu2->calEnergy().em);
00466           zDau2MuEnergyHad->push_back( mu2->calEnergy().had);
00467         }
00468 
00469         // HLT trigger  bit
00470           const pat::TriggerObjectStandAloneCollection mu2HLTMatches = mu2->triggerObjectMatchesByPath( hltPath_[c] );
00471           int dimTrig2 = mu2HLTMatches.size();
00472           if(dimTrig2 !=0 ){
00473             zDau2HLTBit->push_back(1);
00474           } 
00475           else {
00476             zDau2HLTBit->push_back(0);
00477           }
00479           if ( mu1->isGlobalMuon() && mu2->isGlobalMuon() ) {
00480           TrackRef stAloneTrack1;
00481           TrackRef stAloneTrack2;
00482           Vector momentum;
00483           Candidate::PolarLorentzVector p4_1;
00484           double mu_mass;
00485           stAloneTrack1 = dau1->get<TrackRef,reco::StandAloneMuonTag>();
00486           stAloneTrack2 = dau2->get<TrackRef,reco::StandAloneMuonTag>();
00487           zDau1SaEta->push_back(stAloneTrack1->eta());
00488           zDau2SaEta->push_back(stAloneTrack2->eta());
00489           zDau1SaPhi->push_back(stAloneTrack1->phi());
00490           zDau2SaPhi->push_back(stAloneTrack2->phi());
00491           if(counter % 2 == 0) {
00492             momentum = stAloneTrack1->momentum();
00493             p4_1 = dau2->polarP4();
00494             mu_mass = dau1->mass();
00496             zDau1SaPt->push_back(stAloneTrack1 ->pt());
00497             zDau2SaPt->push_back(- stAloneTrack2->pt());
00498           }else{
00499             momentum = stAloneTrack2->momentum();
00500             p4_1= dau1->polarP4();
00501             mu_mass = dau2->mass();
00503             zDau1SaPt->push_back( - stAloneTrack1->pt());
00504             zDau2SaPt->push_back( stAloneTrack2->pt());
00505           }
00506 
00507           Candidate::PolarLorentzVector p4_2(momentum.rho(), momentum.eta(),momentum.phi(), mu_mass);
00508           double mass = (p4_1+p4_2).mass();
00509           zMassSa->push_back(mass);     
00510           ++counter;
00511           }
00512 
00513 
00514           zDau2NofMuChambers->push_back(mu2->numberOfChambers());
00515           zDau2NofMuMatches->push_back(mu2->numberOfMatches());
00516       } else{
00517         // for ZMuTk case...
00518         // it's a track......
00519       const pat::GenericParticle * trk2 = dynamic_cast<const pat::GenericParticle*>(m2);
00520       TrackRef mu2TrkRef = trk2->track(); 
00521       zDau2NofHit->push_back(mu2TrkRef->numberOfValidHits());
00522       zDau2NofHitTk->push_back( mu2TrkRef->numberOfValidHits());
00523       zDau2NofHitSta->push_back( 0);
00524       zDau2NofMuChambers->push_back(0);
00525       zDau2NofMuMatches->push_back(0);
00526       zDau2Chi2->push_back(mu2TrkRef->normalizedChi2());
00527       zDau2dxyFromBS->push_back(mu2TrkRef->dxy(beamSpotHandle->position()));
00528       zDau2dzFromBS->push_back(mu2TrkRef->dz(beamSpotHandle->position()));
00529       zDau2dxyFromPV->push_back(mu2TrkRef->dxy(primaryVertices->begin()->position() ));
00530       zDau2dzFromPV->push_back(mu2TrkRef->dz(primaryVertices->begin()->position() ));   
00531         zDau1MuEnergyEm->push_back( -1);
00532         zDau1MuEnergyHad->push_back( -1);
00533     }
00534   }
00535   const string & zName = zName_[c];
00536   evt.put( event,zName + "EventNumber" );
00537   evt.put( run, zName + "RunNumber" );
00538   evt.put( lumi,zName + "LumiBlock" );
00539   evt.put( zMass, zName +  "Mass" );
00540   evt.put( zMassSa, zName +  "MassSa" );
00541   evt.put( zPt, zName + "Pt" );
00542   evt.put( zEta, zName + "Eta" );
00543   evt.put( zPhi, zName + "Phi" );
00544   evt.put( zY, zName + "Y" );
00545   evt.put( zDau1Pt, zName + "Dau1Pt" );
00546   evt.put( zDau2Pt, zName + "Dau2Pt" );
00547   evt.put( zDau1SaPt, zName + "Dau1SaPt" );
00548   evt.put( zDau2SaPt, zName + "Dau2SaPt" );
00549   evt.put( zDau1HLTBit, zName + "Dau1HLTBit" );
00550   evt.put( zDau2HLTBit, zName + "Dau2HLTBit" );
00551   evt.put( zDau1Q, zName + "Dau1Q" );
00552   evt.put( zDau2Q, zName + "Dau2Q" );
00553   evt.put( zDau1Eta, zName + "Dau1Eta" );
00554   evt.put( zDau2Eta, zName + "Dau2Eta" );
00555   evt.put( zDau1SaEta, zName + "Dau1SaEta" );
00556   evt.put( zDau2SaEta, zName + "Dau2SaEta" );
00557   evt.put( zDau1Phi, zName + "Dau1Phi" );
00558   evt.put( zDau2Phi, zName + "Dau2Phi" );
00559   evt.put( zDau1SaPhi, zName + "Dau1SaPhi" );
00560   evt.put( zDau2SaPhi, zName + "Dau2SaPhi" );
00561   evt.put( zDau1Iso, zName + "Dau1Iso" );
00562   evt.put( zDau2Iso, zName + "Dau2Iso" );
00563   evt.put( zDau1TrkIso, zName + "Dau1TrkIso" );
00564   evt.put( zDau2TrkIso, zName + "Dau2TrkIso" );
00565   evt.put( zDau1EcalIso, zName + "Dau1EcalIso" );
00566   evt.put( zDau2EcalIso, zName + "Dau2EcalIso" );
00567   evt.put( zDau1HcalIso, zName + "Dau1HcalIso" );
00568   evt.put( zDau2HcalIso, zName + "Dau2HcalIso" );
00569   evt.put( zDau1MuEnergyEm, zName + "Dau1MuEnergyEm" );
00570   evt.put( zDau2MuEnergyEm, zName + "Dau2MuEnergyEm" );
00571   evt.put( zDau1MuEnergyHad, zName + "Dau1MuEnergyHad" );
00572   evt.put( zDau2MuEnergyHad, zName + "Dau2MuEnergyHad" );
00573   evt.put( vtxNormChi2, zName + "VtxNormChi2" );
00574   evt.put( zDau1NofHit, zName + "Dau1NofHit" );
00575   evt.put( zDau2NofHit, zName + "Dau2NofHit" );
00576   evt.put( zDau1NofHitTk, zName + "Dau1NofHitTk" );
00577   evt.put( zDau2NofHitTk, zName + "Dau2NofHitTk" );
00578   evt.put( zDau1NofHitSta, zName + "Dau1NofHitSta" );
00579   evt.put( zDau2NofHitSta, zName + "Dau2NofHitSta" );
00580   evt.put( zDau1NofMuChambers, zName + "Dau1NofMuChambers" );
00581   evt.put( zDau1NofMuMatches, zName + "Dau1NofMuMatches" );
00582   evt.put( zDau2NofMuChambers, zName + "Dau2NofMuChambers" );
00583   evt.put( zDau2NofMuMatches, zName + "Dau2NofMuMatches" );
00584   evt.put( zDau1Chi2, zName + "Dau1Chi2" );
00585   evt.put( zDau2Chi2, zName + "Dau2Chi2" );
00586   evt.put( zDau1TrkChi2, zName + "Dau1TrkChi2" );
00587   evt.put( zDau2TrkChi2, zName + "Dau2TrkChi2" );
00588   evt.put( zDau1dxyFromBS, zName + "Dau1dxyFromBS" );
00589   evt.put( zDau2dxyFromBS, zName + "Dau2dxyFromBS" );
00590   evt.put( zDau1dxyFromPV, zName + "Dau1dxyFromPV" );
00591   evt.put( zDau2dxyFromPV, zName + "Dau2dxyFromPV" );
00592   evt.put( zDau1dzFromBS, zName + "Dau1dzFromBS" );
00593   evt.put( zDau2dzFromBS, zName + "Dau2dzFromBS" );
00594   evt.put( zDau1dzFromPV, zName + "Dau1dzFromPV" );
00595   evt.put( zDau2dzFromPV, zName + "Dau2dzFromPV" );
00596   evt.put( trueZMass, zName +  "TrueMass" );
00597   evt.put( trueZPt, zName + "TruePt" );
00598   evt.put( trueZEta, zName + "TrueEta" );
00599   evt.put( trueZPhi, zName + "TruePhi" );
00600   evt.put( trueZY, zName + "TrueY" );
00601 }
00602 }
00603 
00604 #include "FWCore/Framework/interface/MakerMacros.h"
00605 
00606 DEFINE_FWK_MODULE( ZToLLEdmNtupleDumper );
00607