Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-25 23:59:50

0001 /** \class HLTPMMassFilter
0002  *
0003  * Original Author: Jeremy Werner
0004  * Institution: Princeton University, USA
0005  * Contact: Jeremy.Werner@cern.ch
0006  * Date: February 21, 2007
0007  */
0008 #include "HLTPMMassFilter.h"
0009 
0010 //
0011 // constructors and destructor
0012 //
0013 HLTPMMassFilter::HLTPMMassFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig), magFieldToken_(esConsumes()) {
0014   candTag_ = iConfig.getParameter<edm::InputTag>("candTag");
0015   beamSpot_ = iConfig.getParameter<edm::InputTag>("beamSpot");
0016   l1EGTag_ = iConfig.getParameter<edm::InputTag>("l1EGCand");
0017 
0018   lowerMassCut_ = iConfig.getParameter<double>("lowerMassCut");
0019   upperMassCut_ = iConfig.getParameter<double>("upperMassCut");
0020   nZcandcut_ = iConfig.getParameter<int>("nZcandcut");
0021   reqOppCharge_ = iConfig.getUntrackedParameter<bool>("reqOppCharge", false);
0022   isElectron1_ = iConfig.getUntrackedParameter<bool>("isElectron1", true);
0023   isElectron2_ = iConfig.getUntrackedParameter<bool>("isElectron2", true);
0024 
0025   candToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(candTag_);
0026   beamSpotToken_ = consumes<reco::BeamSpot>(beamSpot_);
0027 }
0028 
0029 HLTPMMassFilter::~HLTPMMassFilter() = default;
0030 
0031 void HLTPMMassFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0032   edm::ParameterSetDescription desc;
0033   makeHLTFilterDescription(desc);
0034   desc.add<edm::InputTag>("candTag", edm::InputTag("hltL1NonIsoDoublePhotonEt5UpsHcalIsolFilter"));
0035   desc.add<edm::InputTag>("beamSpot", edm::InputTag("hltOfflineBeamSpot"));
0036   desc.add<double>("lowerMassCut", 8.0);
0037   desc.add<double>("upperMassCut", 11.0);
0038   desc.add<int>("nZcandcut", 1);
0039   desc.addUntracked<bool>("reqOppCharge", true);
0040   desc.addUntracked<bool>("isElectron1", false);
0041   desc.addUntracked<bool>("isElectron2", false);
0042   desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0043   descriptions.add("hltPMMassFilter", desc);
0044 }
0045 
0046 // ------------ method called to produce the data  ------------
0047 bool HLTPMMassFilter::hltFilter(edm::Event& iEvent,
0048                                 const edm::EventSetup& iSetup,
0049                                 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0050   using namespace std;
0051   using namespace edm;
0052   using namespace reco;
0053   using namespace trigger;
0054 
0055   // The filter object
0056   if (saveTags()) {
0057     filterproduct.addCollectionTag(l1EGTag_);
0058   }
0059 
0060   auto const& theMagField = iSetup.getData(magFieldToken_);
0061 
0062   edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0063   iEvent.getByToken(candToken_, PrevFilterOutput);
0064 
0065   // beam spot
0066   edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
0067   iEvent.getByToken(beamSpotToken_, recoBeamSpotHandle);
0068   // gets its position
0069   const GlobalPoint vertexPos(
0070       recoBeamSpotHandle->position().x(), recoBeamSpotHandle->position().y(), recoBeamSpotHandle->position().z());
0071 
0072   int n = 0;
0073 
0074   // REMOVED USAGE OF STATIC ARRAYS
0075   // double px[66];
0076   // double py[66];
0077   // double pz[66];
0078   // double energy[66];
0079   std::vector<TLorentzVector> pEleCh1;
0080   std::vector<TLorentzVector> pEleCh2;
0081   std::vector<double> charge;
0082   std::vector<double> etaOrig;
0083 
0084   if (isElectron1_ && isElectron2_) {
0085     Ref<ElectronCollection> refele;
0086 
0087     vector<Ref<ElectronCollection> > electrons;
0088     PrevFilterOutput->getObjects(TriggerElectron, electrons);
0089 
0090     for (auto& electron : electrons) {
0091       refele = electron;
0092 
0093       TLorentzVector pThisEle(refele->px(), refele->py(), refele->pz(), refele->energy());
0094       pEleCh1.push_back(pThisEle);
0095       charge.push_back(refele->charge());
0096     }
0097 
0098     for (unsigned int jj = 0; jj < electrons.size(); jj++) {
0099       TLorentzVector p1Ele = pEleCh1.at(jj);
0100       for (unsigned int ii = jj + 1; ii < electrons.size(); ii++) {
0101         TLorentzVector p2Ele = pEleCh1.at(ii);
0102 
0103         if (fabs(p1Ele.E() - p2Ele.E()) < 0.00001)
0104           continue;
0105         if (reqOppCharge_ && charge[jj] * charge[ii] > 0)
0106           continue;
0107 
0108         TLorentzVector pTot = p1Ele + p2Ele;
0109         double mass = pTot.M();
0110 
0111         if (mass >= lowerMassCut_ && mass <= upperMassCut_) {
0112           n++;
0113           refele = electrons[ii];
0114           filterproduct.addObject(TriggerElectron, refele);
0115           refele = electrons[jj];
0116           filterproduct.addObject(TriggerElectron, refele);
0117         }
0118       }
0119     }
0120 
0121   } else {
0122     Ref<RecoEcalCandidateCollection> refsc;
0123 
0124     vector<Ref<RecoEcalCandidateCollection> > scs;
0125     PrevFilterOutput->getObjects(TriggerCluster, scs);
0126     if (scs.empty())
0127       PrevFilterOutput->getObjects(TriggerPhoton, scs);  //we dont know if its type trigger cluster or trigger photon
0128 
0129     for (auto& i : scs) {
0130       refsc = i;
0131       const reco::SuperClusterRef sc = refsc->superCluster();
0132       TLorentzVector pscPos = approxMomAtVtx(theMagField, vertexPos, sc, 1);
0133       pEleCh1.push_back(pscPos);
0134 
0135       TLorentzVector pscEle = approxMomAtVtx(theMagField, vertexPos, sc, -1);
0136       pEleCh2.push_back(pscEle);
0137       etaOrig.push_back(sc->eta());
0138     }
0139 
0140     for (unsigned int jj = 0; jj < scs.size(); jj++) {
0141       TLorentzVector p1Ele = pEleCh1.at(jj);
0142       for (unsigned int ii = 0; ii < scs.size(); ii++) {
0143         TLorentzVector p2Ele = pEleCh2.at(ii);
0144 
0145         if (fabs(p1Ele.E() - p2Ele.E()) < 0.00001)
0146           continue;
0147 
0148         TLorentzVector pTot = p1Ele + p2Ele;
0149         double mass = pTot.M();
0150 
0151         if (mass >= lowerMassCut_ && mass <= upperMassCut_) {
0152           n++;
0153           refsc = scs[ii];
0154           filterproduct.addObject(TriggerCluster, refsc);
0155           refsc = scs[jj];
0156           filterproduct.addObject(TriggerCluster, refsc);
0157         }
0158       }
0159     }
0160   }
0161 
0162   // filter decision
0163   bool accept(n >= nZcandcut_);
0164   // if (accept) std::cout << "n size = " << n << std::endl;
0165 
0166   return accept;
0167 }
0168 
0169 TLorentzVector HLTPMMassFilter::approxMomAtVtx(const MagneticField& magField,
0170                                                const GlobalPoint& xvert,
0171                                                const reco::SuperClusterRef sc,
0172                                                int charge) const {
0173   GlobalPoint xsc(sc->position().x(), sc->position().y(), sc->position().z());
0174   float energy = sc->energy();
0175   auto theFTS = trackingTools::ftsFromVertexToPoint(magField, xsc, xvert, energy, charge);
0176   float theApproxMomMod = theFTS.momentum().x() * theFTS.momentum().x() +
0177                           theFTS.momentum().y() * theFTS.momentum().y() + theFTS.momentum().z() * theFTS.momentum().z();
0178   TLorentzVector theApproxMom(
0179       theFTS.momentum().x(), theFTS.momentum().y(), theFTS.momentum().z(), sqrt(theApproxMomMod + 2.61121E-7));
0180   return theApproxMom;
0181 }
0182 
0183 // declare this class as a framework plugin
0184 #include "FWCore/Framework/interface/MakerMacros.h"
0185 DEFINE_FWK_MODULE(HLTPMMassFilter);