Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:20:31

0001 // -*- C++ -*-
0002 //
0003 // Class:      HLTHcalLaserFilter
0004 //
0005 /**\class HLTHcalLaserFilter
0006 
0007  Description: HLT filter module for rejecting events with HCAL laser firing
0008 
0009  Implementation:
0010      <Notes on implementation>
0011 */
0012 //
0013 // Original Author:  Alex Mott
0014 //         Created:  Wed Aug 15 10:37:03 EST 2012
0015 //
0016 //
0017 
0018 #include "HLTrigger/JetMET/interface/HLTHcalLaserFilter.h"
0019 
0020 #include "DataFormats/Common/interface/Handle.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 #include "FWCore/Framework/interface/ESHandle.h"
0023 #include "FWCore/Framework/interface/Event.h"
0024 #include "FWCore/Framework/interface/EventSetup.h"
0025 
0026 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0027 
0028 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0029 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0030 #include "FWCore/Utilities/interface/InputTag.h"
0031 
0032 #include <iostream>
0033 
0034 HLTHcalLaserFilter::HLTHcalLaserFilter(const edm::ParameterSet& iConfig)
0035     : hcalDigiCollection_(iConfig.getParameter<edm::InputTag>("hcalDigiCollection")),
0036       timeSlices_(iConfig.getParameter<std::vector<int> >("timeSlices")),
0037       thresholdsfC_(iConfig.getParameter<std::vector<double> >("thresholdsfC")),
0038       CalibCountFilterValues_(iConfig.getParameter<std::vector<int> >("CalibCountFilterValues")),
0039       CalibChargeFilterValues_(iConfig.getParameter<std::vector<double> >("CalibChargeFilterValues")),
0040       maxTotalCalibCharge_(iConfig.getParameter<double>("maxTotalCalibCharge")),
0041       maxAllowedHFcalib_(iConfig.getParameter<int>("maxAllowedHFcalib"))
0042 
0043 {
0044   //maxAllowedHFcalib_=10;
0045 
0046   m_theCalibToken = consumes<HcalCalibDigiCollection>(hcalDigiCollection_);
0047 }
0048 
0049 HLTHcalLaserFilter::~HLTHcalLaserFilter() = default;
0050 
0051 void HLTHcalLaserFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0052   edm::ParameterSetDescription desc;
0053   desc.add<edm::InputTag>("hcalDigiCollection", edm::InputTag("hltHcalDigis"));
0054   desc.add<double>("maxTotalCalibCharge", -1);
0055 
0056   std::vector<int> dummy_vint;
0057   std::vector<double> dummy_vdouble;
0058 
0059   desc.add<std::vector<int> >("timeSlices", dummy_vint);
0060   desc.add<std::vector<double> >("thresholdsfC", dummy_vdouble);
0061   desc.add<std::vector<int> >("CalibCountFilterValues", dummy_vint);
0062   desc.add<std::vector<double> >("CalibChargeFilterValues", dummy_vdouble);
0063   desc.add<int>("maxAllowedHFcalib", -1);
0064   descriptions.add("hltHcalLaserFilter", desc);
0065 }
0066 
0067 //
0068 // member functions
0069 //
0070 
0071 bool HLTHcalLaserFilter::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0072   edm::Handle<HcalCalibDigiCollection> hCalib;
0073   iEvent.getByToken(m_theCalibToken, hCalib);
0074 
0075   int numHFcalib = 0;
0076 
0077   // Set up potential filter variables
0078   double totalCalibCharge = 0;
0079 
0080   // Track multiplicity and total charge for each fC threshold
0081   std::vector<int> CalibCount;
0082   std::vector<double> CalibCharge;
0083   for (unsigned int i = 0; i < thresholdsfC_.size(); ++i) {
0084     CalibCount.push_back(0);
0085     CalibCharge.push_back(0);
0086   }
0087 
0088   const float adc2fC[128] = {
0089       -0.5,   0.5,    1.5,    2.5,    3.5,    4.5,    5.5,    6.5,    7.5,    8.5,    9.5,    10.5,   11.5,
0090       12.5,   13.5,   15.,    17.,    19.,    21.,    23.,    25.,    27.,    29.5,   32.5,   35.5,   38.5,
0091       42.,    46.,    50.,    54.5,   59.5,   64.5,   59.5,   64.5,   69.5,   74.5,   79.5,   84.5,   89.5,
0092       94.5,   99.5,   104.5,  109.5,  114.5,  119.5,  124.5,  129.5,  137.,   147.,   157.,   167.,   177.,
0093       187.,   197.,   209.5,  224.5,  239.5,  254.5,  272.,   292.,   312.,   334.5,  359.5,  384.5,  359.5,
0094       384.5,  409.5,  434.5,  459.5,  484.5,  509.5,  534.5,  559.5,  584.5,  609.5,  634.5,  659.5,  684.5,
0095       709.5,  747.,   797.,   847.,   897.,   947.,   997.,   1047.,  1109.5, 1184.5, 1259.5, 1334.5, 1422.,
0096       1522.,  1622.,  1734.5, 1859.5, 1984.5, 1859.5, 1984.5, 2109.5, 2234.5, 2359.5, 2484.5, 2609.5, 2734.5,
0097       2859.5, 2984.5, 3109.5, 3234.5, 3359.5, 3484.5, 3609.5, 3797.,  4047.,  4297.,  4547.,  4797.,  5047.,
0098       5297.,  5609.5, 5984.5, 6359.5, 6734.5, 7172.,  7672.,  8172.,  8734.5, 9359.5, 9984.5};
0099 
0100   if (hCalib.isValid() == true) {
0101     // loop over calibration channels
0102 
0103     //for timing reasons, we abort within the loop if a field ever goes out of bounds
0104     for (auto const& digi : *hCalib) {
0105       if (digi.id().hcalSubdet() == 0)
0106         continue;
0107 
0108       HcalCalibDetId myid = (HcalCalibDetId)digi.id();
0109       if (myid.hcalSubdet() == HcalBarrel || myid.hcalSubdet() == HcalEndcap) {
0110         if (myid.calibFlavor() == HcalCalibDetId::HOCrosstalk)
0111           continue;  // ignore HOCrosstalk channels
0112 
0113         // Add this digi to total calibration charge
0114         for (int i = 0; i < (int)digi.size(); i++)
0115           totalCalibCharge = totalCalibCharge + adc2fC[digi.sample(i).adc() & 0xff];
0116 
0117         if (maxTotalCalibCharge_ >= 0 && totalCalibCharge > maxTotalCalibCharge_)
0118           return false;
0119 
0120         // Compute total charge found in the provided subset of timeslices
0121         double sumCharge = 0;
0122         unsigned int NTS = timeSlices_.size();
0123         int digisize = (int)digi.size();  // gives value of largest time slice
0124 
0125         for (unsigned int ts = 0; ts < NTS; ++ts)  // loop over provided timeslices
0126         {
0127           if (timeSlices_[ts] < 0 || timeSlices_[ts] > digisize)
0128             continue;
0129           sumCharge += adc2fC[digi.sample(timeSlices_[ts]).adc() & 0xff];
0130         }
0131 
0132         // Check multiplicity and charge against filter settings for each charge threshold
0133         for (unsigned int thresh = 0; thresh < thresholdsfC_.size(); ++thresh) {
0134           if (sumCharge > thresholdsfC_[thresh]) {
0135             ++CalibCount[thresh];
0136             CalibCharge[thresh] += sumCharge;
0137             // FilterValues must be >=0 in order for filter to be applied
0138             if (CalibCount[thresh] >= CalibCountFilterValues_[thresh] && CalibCountFilterValues_[thresh] >= 0) {
0139               //std::cout <<"Number of channels > "<<thresholdsfC_[thresh]<<" = "<<CalibCount[thresh]<<"; vetoing!"<<std::endl;
0140               return false;
0141             }
0142             if (CalibCharge[thresh] >= CalibChargeFilterValues_[thresh] && CalibChargeFilterValues_[thresh] >= 0) {
0143               //std::cout <<"FILTERED BY HBHE"<<std::endl;
0144               return false;
0145             }
0146           }  //if (sumCharge > thresholdsfC_[thresh])
0147         }    //for (unsigned int thresh=0;thresh<thresholdsfC_.size();++thresh)
0148       }      // if HB or HE Calib
0149       else if (myid.hcalSubdet() == HcalForward && maxAllowedHFcalib_ >= 0) {
0150         ++numHFcalib;
0151         //std::cout <<"numHFcalib = "<<numHFcalib<<"  Max allowed = "<<maxAllowedHFcalib_<<std::endl;
0152         if (numHFcalib > maxAllowedHFcalib_) {
0153           //std::cout <<"FILTERED BY HF; "<<maxAllowedHFcalib_<<std::endl;
0154           return false;
0155         }
0156       }
0157     }  // loop on calibration digis:  for (HcalCalibDigiCollection::...)
0158 
0159     /*
0160       for (unsigned int thresh=0;thresh<thresholdsfC_.size();++thresh)
0161     {
0162       std::cout <<"Thresh = "<<thresholdsfC_[thresh]<<"  Num channels found = "<<CalibCount[thresh]<<std::endl;
0163     }
0164       */
0165   }  // if (hCalib.isValid()==true)
0166   //std::cout <<"UNFILTERED"<<std::endl;
0167   return true;
0168 }