Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:25:11

0001 // Authors: Marco Rovere - marco.rovere@cern.ch, Felice Pantaleo - felice.pantaleo@cern.ch
0002 // Date: 09/2020
0003 
0004 #ifndef RecoHGCal_TICL_ClusterFilterByAlgoAndSizeAndLayerRange_H__
0005 #define RecoHGCal_TICL_ClusterFilterByAlgoAndSizeAndLayerRange_H__
0006 
0007 #include "DataFormats/CaloRecHit/interface/CaloCluster.h"
0008 #include "ClusterFilterBase.h"
0009 
0010 #include <memory>
0011 #include <utility>
0012 
0013 // Filter clusters that belong to a specific algorithm
0014 namespace ticl {
0015   class ClusterFilterByAlgoAndSizeAndLayerRange final : public ClusterFilterBase {
0016   public:
0017     ClusterFilterByAlgoAndSizeAndLayerRange(const edm::ParameterSet& ps)
0018         : ClusterFilterBase(ps),
0019           algo_number_(ps.getParameter<std::vector<int>>("algo_number")),
0020           min_cluster_size_(ps.getParameter<int>("min_cluster_size")),
0021           max_cluster_size_(ps.getParameter<int>("max_cluster_size")),
0022           min_layerId_(ps.getParameter<int>("min_layerId")),
0023           max_layerId_(ps.getParameter<int>("max_layerId")) {}
0024     ~ClusterFilterByAlgoAndSizeAndLayerRange() override{};
0025 
0026     void filter(const std::vector<reco::CaloCluster>& layerClusters,
0027                 const TICLClusterFilterMask& availableLayerClusters,
0028                 std::vector<float>& layerClustersMask,
0029                 hgcal::RecHitTools& rhtools) const override {
0030       auto filteredLayerClusters = std::make_unique<TICLClusterFilterMask>();
0031       for (auto const& cl : availableLayerClusters) {
0032         auto const& layerCluster = layerClusters[cl.first];
0033         auto const& haf = layerCluster.hitsAndFractions();
0034         auto layerId = rhtools.getLayerWithOffset(haf[0].first);
0035         if (find(algo_number_.begin(), algo_number_.end(), layerCluster.algo()) != algo_number_.end() and
0036             layerId <= max_layerId_ and layerId >= min_layerId_ and haf.size() <= max_cluster_size_ and
0037             (haf.size() >= min_cluster_size_ or !(rhtools.isSilicon(haf[0].first)))) {
0038           filteredLayerClusters->emplace_back(cl);
0039         } else {
0040           layerClustersMask[cl.first] = 0.;
0041         }
0042       }
0043     }
0044 
0045   private:
0046     std::vector<int> algo_number_;
0047     unsigned int min_cluster_size_;
0048     unsigned int max_cluster_size_;
0049     unsigned int min_layerId_;
0050     unsigned int max_layerId_;
0051   };
0052 }  // namespace ticl
0053 
0054 #endif