File indexing completed on 2023-03-17 11:09:46
0001 #include "HLTrigger/HLTcore/interface/HLTFilter.h"
0002
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/EventSetup.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0007 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0008 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0009 #include "DataFormats/Common/interface/Handle.h"
0010 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0011 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0012 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0013 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0014 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0015
0016
0017
0018
0019
0020 class HLTPixelActivityFilter : public HLTFilter {
0021 public:
0022 explicit HLTPixelActivityFilter(const edm::ParameterSet&);
0023 ~HLTPixelActivityFilter() override;
0024 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0025
0026 private:
0027 bool hltFilter(edm::Event&,
0028 const edm::EventSetup&,
0029 trigger::TriggerFilterObjectWithRefs& filterproduct) const override;
0030
0031 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> const trackerTopologyRcdToken_;
0032 edm::InputTag inputTag_;
0033 unsigned int min_clusters_;
0034 unsigned int max_clusters_;
0035 unsigned int min_clustersBPix_;
0036 unsigned int max_clustersBPix_;
0037 unsigned int min_clustersFPix_;
0038 unsigned int max_clustersFPix_;
0039 unsigned int min_layersBPix_;
0040 unsigned int max_layersBPix_;
0041 unsigned int min_layersFPix_;
0042 unsigned int max_layersFPix_;
0043 edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster> > inputToken_;
0044 };
0045
0046
0047
0048
0049
0050 HLTPixelActivityFilter::HLTPixelActivityFilter(const edm::ParameterSet& config)
0051 : HLTFilter(config),
0052 trackerTopologyRcdToken_(esConsumes()),
0053 inputTag_(config.getParameter<edm::InputTag>("inputTag")),
0054 min_clusters_(config.getParameter<unsigned int>("minClusters")),
0055 max_clusters_(config.getParameter<unsigned int>("maxClusters")),
0056 min_clustersBPix_(config.getParameter<unsigned int>("minClustersBPix")),
0057 max_clustersBPix_(config.getParameter<unsigned int>("maxClustersBPix")),
0058 min_clustersFPix_(config.getParameter<unsigned int>("minClustersFPix")),
0059 max_clustersFPix_(config.getParameter<unsigned int>("maxClustersFPix")),
0060 min_layersBPix_(config.getParameter<unsigned int>("minLayersBPix")),
0061 max_layersBPix_(config.getParameter<unsigned int>("maxLayersBPix")),
0062 min_layersFPix_(config.getParameter<unsigned int>("minLayersFPix")),
0063 max_layersFPix_(config.getParameter<unsigned int>("maxLayersFPix")) {
0064 inputToken_ = consumes<edmNew::DetSetVector<SiPixelCluster> >(inputTag_);
0065 LogDebug("") << "Using the " << inputTag_ << " input collection";
0066 LogDebug("") << "Requesting at least " << min_clusters_ << " clusters";
0067 if (max_clusters_ > 0)
0068 LogDebug("") << "...but no more than " << max_clusters_ << " clusters";
0069 if (min_layersBPix_ > 0)
0070 LogDebug("") << "Also requesting at least " << min_layersBPix_ << " layers with clusters in BPix";
0071 if (max_layersBPix_ > 0)
0072 LogDebug("") << "...but no more than " << max_layersBPix_ << " layers with clusters in BPix";
0073 if (min_layersFPix_ > 0)
0074 LogDebug("") << "Also requesting at least " << min_layersFPix_ << " layers with clusters in FPix";
0075 if (max_layersFPix_ > 0)
0076 LogDebug("") << "...but no more than " << max_layersFPix_ << " layers with clusters in FPix";
0077 }
0078
0079 HLTPixelActivityFilter::~HLTPixelActivityFilter() = default;
0080
0081 void HLTPixelActivityFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0082 edm::ParameterSetDescription desc;
0083 makeHLTFilterDescription(desc);
0084 desc.add<edm::InputTag>("inputTag", edm::InputTag("hltSiPixelClusters"));
0085 desc.add<unsigned int>("minClusters", 3);
0086 desc.add<unsigned int>("maxClusters", 0);
0087 desc.add<unsigned int>("minClustersBPix", 0);
0088 desc.add<unsigned int>("maxClustersBPix", 0);
0089 desc.add<unsigned int>("minClustersFPix", 0);
0090 desc.add<unsigned int>("maxClustersFPix", 0);
0091 desc.add<unsigned int>("minLayersBPix", 0);
0092 desc.add<unsigned int>("maxLayersBPix", 0);
0093 desc.add<unsigned int>("minLayersFPix", 0);
0094 desc.add<unsigned int>("maxLayersFPix", 0);
0095 descriptions.add("hltPixelActivityFilter", desc);
0096 }
0097
0098
0099
0100
0101
0102 bool HLTPixelActivityFilter::hltFilter(edm::Event& event,
0103 const edm::EventSetup& iSetup,
0104 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0105
0106
0107
0108
0109
0110 if (saveTags())
0111 filterproduct.addCollectionTag(inputTag_);
0112
0113
0114 edm::Handle<edmNew::DetSetVector<SiPixelCluster> > clusterColl;
0115 event.getByToken(inputToken_, clusterColl);
0116
0117 unsigned int clusterSize = clusterColl->dataSize();
0118 LogDebug("") << "Number of clusters accepted: " << clusterSize;
0119 bool accept = (clusterSize >= min_clusters_);
0120 if (max_clusters_ > 0)
0121 accept &= (clusterSize <= max_clusters_);
0122 if (min_layersBPix_ > 0 || max_layersBPix_ > 0 || min_layersFPix_ > 0 || max_layersFPix_ > 0 ||
0123 min_clustersBPix_ > 0 || max_clustersBPix_ > 0 || min_clustersFPix_ > 0 || max_clustersFPix_ > 0) {
0124 auto const& tTopoHandle = iSetup.getHandle(trackerTopologyRcdToken_);
0125
0126 const TrackerTopology& tTopo = *tTopoHandle;
0127 unsigned int layerCountBPix = 0;
0128 unsigned int layerCountFPix = 0;
0129 unsigned int clusterCountBPix = 0;
0130 unsigned int clusterCountFPix = 0;
0131 const edmNew::DetSetVector<SiPixelCluster>& clusters = *clusterColl;
0132
0133 edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = clusters.begin();
0134
0135 std::vector<int> foundLayersB;
0136 std::vector<int> foundLayersEp;
0137 std::vector<int> foundLayersEm;
0138 for (; DSViter != clusters.end(); DSViter++) {
0139 unsigned int detid = DSViter->detId();
0140 DetId detIdObject(detid);
0141 const auto nCluster = DSViter->size();
0142 const auto subdet = detIdObject.subdetId();
0143 if (subdet == PixelSubdetector::PixelBarrel) {
0144 if (!(std::find(foundLayersB.begin(), foundLayersB.end(), tTopo.layer(detIdObject)) != foundLayersB.end()) &&
0145 nCluster > 0)
0146 foundLayersB.push_back(tTopo.layer(detIdObject));
0147 clusterCountBPix += nCluster;
0148 } else if (subdet == PixelSubdetector::PixelEndcap) {
0149 clusterCountFPix += nCluster;
0150 if (tTopo.side(detIdObject) == 2) {
0151 if (!(std::find(foundLayersEp.begin(), foundLayersEp.end(), tTopo.layer(detIdObject)) !=
0152 foundLayersEp.end()) &&
0153 nCluster > 0)
0154 foundLayersEp.push_back(tTopo.layer(detIdObject));
0155 } else if (tTopo.side(detIdObject) == 1) {
0156 if (!(std::find(foundLayersEm.begin(), foundLayersEm.end(), tTopo.layer(detIdObject)) !=
0157 foundLayersEm.end()) &&
0158 nCluster > 0)
0159 foundLayersEm.push_back(tTopo.layer(detIdObject));
0160 }
0161 }
0162 }
0163 layerCountBPix = foundLayersB.size();
0164 layerCountFPix = foundLayersEp.size() + foundLayersEm.size();
0165 if (max_layersBPix_ > 0)
0166 accept &= (layerCountBPix <= max_layersBPix_);
0167 if (min_layersBPix_ > 0)
0168 accept &= (layerCountBPix >= min_layersBPix_);
0169 if (max_layersFPix_ > 0)
0170 accept &= (layerCountFPix <= max_layersFPix_);
0171 if (min_layersFPix_ > 0)
0172 accept &= (layerCountFPix >= min_layersFPix_);
0173 if (min_clustersBPix_ > 0)
0174 accept &= (clusterCountBPix >= min_clustersBPix_);
0175 if (max_clustersBPix_ > 0)
0176 accept &= (clusterCountBPix <= max_clustersBPix_);
0177 if (min_clustersFPix_ > 0)
0178 accept &= (clusterCountFPix >= min_clustersFPix_);
0179 if (max_clustersFPix_ > 0)
0180 accept &= (clusterCountFPix <= max_clustersFPix_);
0181 }
0182
0183 return accept;
0184 }
0185
0186
0187 #include "FWCore/Framework/interface/MakerMacros.h"
0188 DEFINE_FWK_MODULE(HLTPixelActivityFilter);