Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-10-25 09:38:54

0001 #ifndef ECALTRIGPRIMCOMPACTCOLL_H
0002 #define ECALTRIGPRIMCOMPACTCOLL_H
0003 
0004 #include <vector>
0005 #include <cinttypes>
0006 #include "FWCore/Utilities/interface/Exception.h"
0007 
0008 #include "DataFormats/EcalDigi/interface/EcalTriggerPrimitiveDigi.h"
0009 #include "DataFormats/Common/interface/SortedCollection.h"
0010 typedef edm::SortedCollection<EcalTriggerPrimitiveDigi> EcalTrigPrimDigiCollection;
0011 
0012 /** \class EcalTrigPrimCompactColl
0013 
0014 This collection is used to store ECAL trigger primitive with a low footpring in ROOT EDM file.
0015 
0016 Note that only the time sample of interest time is stored. In nornal operation data contains only this time sample.
0017 
0018 The interface is similar to EcalTriggerPrimitiveSample one. The collection can also be converted to an EcalTrigPrimDigiCollection
0019 with the method toEcalTrigPrimDigiCollection().
0020 
0021 The collection is generated with EcalCompactTrigPrimProducer module from package RecoLocalCalo/EcalRecProducers.
0022 
0023 author: Ph. Gras CEA/IRFU Saclay
0024 
0025 */
0026 
0027 class EcalTrigPrimCompactColl {
0028 private:
0029   static const int nPhiBins = 72;
0030   static const int nEtaBins = 56;
0031   static const int nBins = nPhiBins * nEtaBins;
0032 
0033 private:
0034   static size_t index(int ieta, int iphi) {
0035     size_t r = unsigned(((ieta < 0) ? (ieta + 28) : (ieta + 27)) * nPhiBins + (iphi - 1));
0036     if (r >= (unsigned)nBins)
0037       throw cms::Exception("Invalid argument")
0038           << "Trigger tower index (" << ieta << "," << iphi << ") are out of range";
0039     return r;
0040   }
0041 
0042 public:
0043   EcalTrigPrimCompactColl() : formatVers_(0), data_(nBins){};
0044 
0045   ///Set data
0046   void setValue(int ieta, int iphi, uint16_t sample) { data_[index(ieta, iphi)] = sample; }
0047 
0048   //@{
0049   /// get the raw word
0050   uint16_t raw(int ieta, int iphi) const { return data_[index(ieta, iphi)]; }
0051   uint16_t raw(const EcalTrigTowerDetId& ttId) const { return raw(ttId.ieta(), ttId.iphi()); }
0052   //@}
0053 
0054   //@{
0055   /// get the encoded/compressed Et (8 bits)
0056   int compressedEt(int ieta, int iphi) const { return raw(ieta, iphi) & 0xFF; }
0057   int compressedEt(const EcalTrigTowerDetId& ttId) const { return compressedEt(ttId.ieta(), ttId.iphi()); }
0058   //@}
0059 
0060   //@{
0061   /// get the fine-grain bit (1 bit)
0062   bool fineGrain(int ieta, int iphi) const { return (raw(ieta, iphi) & 0x100) != 0; }
0063   bool fineGrain(const EcalTrigTowerDetId& ttId) const { return fineGrain(ttId.ieta(), ttId.iphi()); }
0064   //@}
0065 
0066   //@{
0067   /// get the Trigger tower Flag (3 bits)
0068   int ttFlag(int ieta, int iphi) const { return (raw(ieta, iphi) >> 9) & 0x7; }
0069   int ttFlag(const EcalTrigTowerDetId& ttId) const { return ttFlag(ttId.ieta(), ttId.iphi()); }
0070   //@}
0071 
0072   //@{
0073   /// Gets the "strip fine grain veto bit" (sFGVB)  used as L1A spike detection
0074   /// @return 0 spike like pattern
0075   ///         1 EM shower like pattern
0076   int sFGVB(int ieta, int iphi) const { return (raw(ieta, iphi) >> 12) & 0x1; }
0077   int sFGVB(const EcalTrigTowerDetId& ttId) const { return sFGVB(ttId.ieta(), ttId.iphi()); }
0078   //@}
0079 
0080   //@{
0081   /// Gets the "strip fine grain veto bit" (sFGVB)  used as L1A spike detection
0082   /// Deprecated. Use instead sFGVB() method. Indeed the name of the method being missleading,
0083   /// since it returns 0 for spike-compatible deposit.
0084   /// @return 0 spike-like pattern
0085   ///         1 EM-shower-like pattern
0086   int l1aSpike(int ieta, int iphi) const { return sFGVB(ieta, iphi); }
0087   int l1aSpike(const EcalTrigTowerDetId& ttId) const { return sFGVB(ttId); }
0088   //@}
0089 
0090   void toEcalTrigPrimDigiCollection(EcalTrigPrimDigiCollection& dest) const;
0091 
0092 private:
0093   int16_t formatVers_;
0094   std::vector<uint16_t> data_;
0095 };
0096 
0097 #endif  //ECALTRIGPRIMCOMPACTCOLL_H not defined