Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:20:26 

0001 // JetFinderMethods.cc
0002 // Author: Alex Barbieri
0003 //
0004 // This file should contain the different algorithms used to find jets.
0005 // Currently the standard is the sliding window method, used by both
0006 // HI and PP.
0007 
0008 #include "DataFormats/L1TCalorimeter/interface/CaloRegion.h"
0009 #include "DataFormats/L1Trigger/interface/Jet.h"
0010 #include "L1Trigger/L1TCalorimeter/interface/JetFinderMethods.h"
0011 
0012 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegionDetId.h"
0013 
0014 #include <vector>
0015 
0016 namespace l1t {
0017 
0018   int deltaGctPhi(const CaloRegion& region, const CaloRegion& neighbor) {
0019     int phi1 = region.hwPhi();
0020     int phi2 = neighbor.hwPhi();
0021     int diff = phi1 - phi2;
0022     if (std::abs(phi1 - phi2) == L1CaloRegionDetId::N_PHI - 1) {  //18 regions in phi
0023       diff = -diff / std::abs(diff);
0024     }
0025     return diff;
0026   }
0027 
0028   // turn each central region into a jet
0029   void passThroughJets(const std::vector<l1t::CaloRegion>* regions, std::vector<l1t::Jet>* uncalibjets) {
0030     for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
0031       int jetQual = 0;
0032       if (region->hwEta() < 4 || region->hwEta() > 17)
0033         jetQual = 2;
0034       int jetET = region->hwPt();
0035       int jetEta = region->hwEta();
0036       int jetPhi = region->hwPhi();
0037 
0038       ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
0039       l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
0040       uncalibjets->push_back(theJet);
0041     }
0042   }
0043 
0044   void slidingWindowJetFinder(const int jetSeedThreshold,
0045                               const std::vector<l1t::CaloRegion>* regions,
0046                               std::vector<l1t::Jet>* uncalibjets) {
0047     for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
0048       int regionET = region->hwPt();  //regionPhysicalEt(*region);
0049       if (regionET <= jetSeedThreshold)
0050         continue;
0051       int neighborN_et = 0;
0052       int neighborS_et = 0;
0053       int neighborE_et = 0;
0054       int neighborW_et = 0;
0055       int neighborNE_et = 0;
0056       int neighborSW_et = 0;
0057       int neighborNW_et = 0;
0058       int neighborSE_et = 0;
0059       unsigned int nNeighbors = 0;
0060       for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
0061            neighbor++) {
0062         int neighborET = neighbor->hwPt();  //regionPhysicalEt(*neighbor);
0063         if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
0064           neighborN_et = neighborET;
0065           nNeighbors++;
0066           continue;
0067         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
0068           neighborS_et = neighborET;
0069           nNeighbors++;
0070           continue;
0071         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
0072           neighborE_et = neighborET;
0073           nNeighbors++;
0074           continue;
0075         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
0076           neighborW_et = neighborET;
0077           nNeighbors++;
0078           continue;
0079         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0080           neighborNE_et = neighborET;
0081           nNeighbors++;
0082           continue;
0083         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0084           neighborSW_et = neighborET;
0085           nNeighbors++;
0086           continue;
0087         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0088           neighborNW_et = neighborET;
0089           nNeighbors++;
0090           continue;
0091         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0092           neighborSE_et = neighborET;
0093           nNeighbors++;
0094           continue;
0095         }
0096       }
0097       if (regionET > neighborN_et && regionET > neighborNW_et && regionET > neighborW_et && regionET > neighborSW_et &&
0098           regionET >= neighborNE_et && regionET >= neighborE_et && regionET >= neighborSE_et &&
0099           regionET >= neighborS_et) {
0100         unsigned int jetET = regionET + neighborN_et + neighborS_et + neighborE_et + neighborW_et + neighborNE_et +
0101                              neighborSW_et + neighborSE_et + neighborNW_et;
0102 
0103         int jetPhi = region->hwPhi();
0104         int jetEta = region->hwEta();
0105 
0106         bool neighborCheck = (nNeighbors == 8);
0107         // On the eta edge we only expect 5 neighbors
0108         if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
0109           neighborCheck = true;
0110 
0111         if (!neighborCheck) {
0112           std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
0113           assert(false);
0114         }
0115 
0116         //first iteration, eta cut defines forward
0117         //const bool forward = (jetEta <= 4 || jetEta >= 17);
0118         const bool forward = (jetEta < 4 || jetEta > 17);
0119         int jetQual = 0;
0120         if (forward)
0121           jetQual |= 0x2;
0122 
0123         // check for input overflow regions
0124         if (forward && regionET == 255) {
0125           jetET = 1023;  // 10 bit max
0126         } else if (!forward && regionET == 1023) {
0127           jetET = 1023;  // 10 bit max
0128         } else if (region->hwEta() == 17) {
0129           if (neighborNE_et == 255 || neighborE_et == 255 || neighborSE_et == 255)
0130             jetET = 1023;  // 10 bit max
0131         } else if (region->hwEta() == 4) {
0132           if (neighborNW_et == 255 || neighborW_et == 255 || neighborSW_et == 255)
0133             jetET = 1023;  // 10 bit max
0134         }
0135 
0136         ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
0137         l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
0138         //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
0139 
0140         uncalibjets->push_back(theJet);
0141       }
0142     }
0143   }
0144 
0145   void TwelveByTwelveFinder(const int jetSeedThreshold,
0146                             const std::vector<l1t::CaloRegion>* regions,
0147                             std::vector<l1t::Jet>* uncalibjets) {
0148     for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
0149       int regionET = region->hwPt();  //regionPhysicalEt(*region);
0150       if (regionET < jetSeedThreshold)
0151         continue;
0152       int neighborN_et = 0;
0153       int neighborS_et = 0;
0154       int neighborE_et = 0;
0155       int neighborW_et = 0;
0156       int neighborNE_et = 0;
0157       int neighborSW_et = 0;
0158       int neighborNW_et = 0;
0159       int neighborSE_et = 0;
0160       unsigned int nNeighbors = 0;
0161       for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
0162            neighbor++) {
0163         int neighborET = neighbor->hwPt();  //regionPhysicalEt(*neighbor);
0164         if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
0165           neighborN_et = neighborET;
0166           nNeighbors++;
0167           continue;
0168         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
0169           neighborS_et = neighborET;
0170           nNeighbors++;
0171           continue;
0172         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
0173           neighborE_et = neighborET;
0174           nNeighbors++;
0175           continue;
0176         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
0177           neighborW_et = neighborET;
0178           nNeighbors++;
0179           continue;
0180         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0181           neighborNE_et = neighborET;
0182           nNeighbors++;
0183           continue;
0184         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0185           neighborSW_et = neighborET;
0186           nNeighbors++;
0187           continue;
0188         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0189           neighborNW_et = neighborET;
0190           nNeighbors++;
0191           continue;
0192         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0193           neighborSE_et = neighborET;
0194           nNeighbors++;
0195           continue;
0196         }
0197       }
0198       unsigned int jetET = regionET + neighborN_et + neighborS_et + neighborE_et + neighborW_et + neighborNE_et +
0199                            neighborSW_et + neighborSE_et + neighborNW_et;
0200 
0201       int jetPhi = region->hwPhi();
0202       int jetEta = region->hwEta();
0203 
0204       bool neighborCheck = (nNeighbors == 8);
0205       // On the eta edge we only expect 5 neighbors
0206       if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
0207         neighborCheck = true;
0208 
0209       if (!neighborCheck) {
0210         std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
0211         assert(false);
0212       }
0213 
0214       //first iteration, eta cut defines forward
0215       //const bool forward = (jetEta <= 4 || jetEta >= 17);
0216       const bool forward = (jetEta < 4 || jetEta > 17);
0217       int jetQual = 0;
0218       if (forward)
0219         jetQual |= 0x2;
0220 
0221       ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
0222       l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
0223       //l1t::Jet theJet(0, jetET, jetEta, jetPhi);
0224 
0225       uncalibjets->push_back(theJet);
0226     }
0227   }
0228 
0229   void TwoByTwoFinder(const int jetSeedThreshold,
0230                       const int etaMask,
0231                       const std::vector<l1t::CaloRegion>* regions,
0232                       std::vector<l1t::Jet>* uncalibjets) {
0233     for (std::vector<CaloRegion>::const_iterator region = regions->begin(); region != regions->end(); region++) {
0234       int regionET = region->hwPt();
0235       if (regionET <= jetSeedThreshold)
0236         continue;
0237       int subEta = region->hwEta();
0238       if ((etaMask & (1 << subEta)) >> subEta)
0239         regionET = 0;
0240       int neighborN_et = 0;
0241       int neighborS_et = 0;
0242       int neighborE_et = 0;
0243       int neighborW_et = 0;
0244       int neighborNE_et = 0;
0245       int neighborSW_et = 0;
0246       int neighborNW_et = 0;
0247       int neighborSE_et = 0;
0248       unsigned int nNeighbors = 0;
0249       for (std::vector<CaloRegion>::const_iterator neighbor = regions->begin(); neighbor != regions->end();
0250            neighbor++) {
0251         int neighborET = neighbor->hwPt();
0252         int subEta2 = neighbor->hwEta();
0253         if ((etaMask & (1 << subEta2)) >> subEta2)
0254           neighborET = 0;
0255 
0256         if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta()) == neighbor->hwEta()) {
0257           neighborN_et = neighborET;
0258           nNeighbors++;
0259           continue;
0260         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta()) == neighbor->hwEta()) {
0261           neighborS_et = neighborET;
0262           nNeighbors++;
0263           continue;
0264         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() + 1) == neighbor->hwEta()) {
0265           neighborE_et = neighborET;
0266           nNeighbors++;
0267           continue;
0268         } else if (deltaGctPhi(*region, *neighbor) == 0 && (region->hwEta() - 1) == neighbor->hwEta()) {
0269           neighborW_et = neighborET;
0270           nNeighbors++;
0271           continue;
0272         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0273           neighborNE_et = neighborET;
0274           nNeighbors++;
0275           continue;
0276         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0277           neighborSW_et = neighborET;
0278           nNeighbors++;
0279           continue;
0280         } else if (deltaGctPhi(*region, *neighbor) == 1 && (region->hwEta() - 1) == neighbor->hwEta()) {
0281           neighborNW_et = neighborET;
0282           nNeighbors++;
0283           continue;
0284         } else if (deltaGctPhi(*region, *neighbor) == -1 && (region->hwEta() + 1) == neighbor->hwEta()) {
0285           neighborSE_et = neighborET;
0286           nNeighbors++;
0287           continue;
0288         }
0289       }
0290       if (regionET > neighborN_et && regionET > neighborNW_et && regionET > neighborW_et && regionET > neighborSW_et &&
0291           regionET >= neighborNE_et && regionET >= neighborE_et && regionET >= neighborSE_et &&
0292           regionET >= neighborS_et) {
0293         // use the highest-pT 2x2 jet inside this 3x3
0294         unsigned int jetET_NW;
0295         unsigned int jetET_NE;
0296         unsigned int jetET_SW;
0297         unsigned int jetET_SE;
0298 
0299         jetET_NW = regionET + neighborW_et + neighborNW_et + neighborN_et;
0300         jetET_NE = regionET + neighborE_et + neighborNE_et + neighborN_et;
0301         jetET_SW = regionET + neighborS_et + neighborSW_et + neighborW_et;
0302         jetET_SE = regionET + neighborS_et + neighborSE_et + neighborE_et;
0303 
0304         unsigned int jetET = std::max(jetET_NW, jetET_NE);
0305         jetET = std::max(jetET, jetET_SW);
0306         jetET = std::max(jetET, jetET_SE);
0307 
0308         int jetPhi = region->hwPhi();
0309         int jetEta = region->hwEta();
0310 
0311         bool neighborCheck = (nNeighbors == 8);
0312         // On the eta edge we only expect 5 neighbor
0313         if (!neighborCheck && (jetEta == 0 || jetEta == 21) && nNeighbors == 5)
0314           neighborCheck = true;
0315 
0316         if (!neighborCheck) {
0317           std::cout << "phi: " << jetPhi << " eta: " << jetEta << " n: " << nNeighbors << std::endl;
0318           assert(false);
0319         }
0320 
0321         //first iteration, eta cut defines forward
0322         const bool forward = (jetEta < 4 || jetEta > 17);
0323         int jetQual = 0;
0324         if (forward)
0325           jetQual |= 0x2;
0326 
0327         ROOT::Math::LorentzVector<ROOT::Math::PxPyPzE4D<double> > jetLorentz(0, 0, 0, 0);
0328         l1t::Jet theJet(*&jetLorentz, jetET, jetEta, jetPhi, jetQual);
0329         uncalibjets->push_back(theJet);
0330       }
0331     }
0332   }
0333 }  // namespace l1t
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!