File indexing completed on 2024-04-06 12:18:41
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "HLTPixelAsymmetryFilter.h"
0012
0013
0014
0015
0016
0017 HLTPixelAsymmetryFilter::HLTPixelAsymmetryFilter(const edm::ParameterSet& config)
0018 : HLTFilter(config),
0019 inputTag_(config.getParameter<edm::InputTag>("inputTag")),
0020 min_asym_(config.getParameter<double>("MinAsym")),
0021 max_asym_(config.getParameter<double>("MaxAsym")),
0022 clus_thresh_(config.getParameter<double>("MinCharge")),
0023 bmincharge_(config.getParameter<double>("MinBarrel")) {
0024 inputToken_ = consumes<edmNew::DetSetVector<SiPixelCluster> >(inputTag_);
0025 LogDebug("") << "Using the " << inputTag_ << " input collection";
0026 LogDebug("") << "Requesting events with a charge repartition asymmetry between " << min_asym_ << " and " << max_asym_;
0027 LogDebug("") << "Mean cluster charge in the barrel should be higher than" << bmincharge_ << " electrons ";
0028 LogDebug("") << "Only clusters with a charge larger than " << clus_thresh_ << " electrons will be used ";
0029 }
0030
0031 HLTPixelAsymmetryFilter::~HLTPixelAsymmetryFilter() = default;
0032
0033 void HLTPixelAsymmetryFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0034 edm::ParameterSetDescription desc;
0035 makeHLTFilterDescription(desc);
0036 desc.add<edm::InputTag>("inputTag", edm::InputTag("hltSiPixelClusters"));
0037 desc.add<double>("MinAsym", 0.);
0038 desc.add<double>("MaxAsym", 1.);
0039 desc.add<double>("MinCharge", 4000.);
0040 desc.add<double>("MinBarrel", 10000.);
0041 descriptions.add("hltPixelAsymmetryFilter", desc);
0042 }
0043
0044
0045
0046
0047
0048
0049 bool HLTPixelAsymmetryFilter::hltFilter(edm::Event& event,
0050 const edm::EventSetup& iSetup,
0051 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0052
0053
0054
0055
0056
0057 if (saveTags())
0058 filterproduct.addCollectionTag(inputTag_);
0059
0060
0061 edm::Handle<edmNew::DetSetVector<SiPixelCluster> > clusterColl;
0062 event.getByToken(inputToken_, clusterColl);
0063
0064 unsigned int clusterSize = clusterColl->dataSize();
0065 LogDebug("") << "Number of clusters accepted: " << clusterSize;
0066
0067 bool accept = (clusterSize >= 2);
0068
0069 if (!accept)
0070 return accept;
0071
0072 double asym_pix_1 = -1;
0073 double asym_pix_2 = -1;
0074
0075 int n_clus[3] = {0, 0, 0};
0076 double e_clus[3] = {0., 0., 0.};
0077
0078 for (edmNew::DetSetVector<SiPixelCluster>::const_iterator DSViter = clusterColl->begin();
0079 DSViter != clusterColl->end();
0080 DSViter++) {
0081 edmNew::DetSet<SiPixelCluster>::const_iterator begin = DSViter->begin();
0082 edmNew::DetSet<SiPixelCluster>::const_iterator end = DSViter->end();
0083 uint32_t detid = DSViter->id();
0084
0085 bool barrel = DetId(detid).subdetId() == static_cast<int>(PixelSubdetector::PixelBarrel);
0086 bool endcap = DetId(detid).subdetId() == static_cast<int>(PixelSubdetector::PixelEndcap);
0087
0088 int detpos = -1;
0089
0090
0091
0092
0093 if (endcap) {
0094 PixelEndcapName::HalfCylinder position = PixelEndcapName(DetId(detid)).halfCylinder();
0095
0096 if (position == PixelEndcapName::mI || position == PixelEndcapName::mO)
0097 detpos = 0;
0098
0099 if (position == PixelEndcapName::pI || position == PixelEndcapName::pO)
0100 detpos = 2;
0101 }
0102
0103 if (barrel)
0104 detpos = 1;
0105
0106 if (detpos < 0)
0107 continue;
0108
0109 for (edmNew::DetSet<SiPixelCluster>::const_iterator iter = begin; iter != end; ++iter) {
0110 if (iter->charge() > clus_thresh_) {
0111 ++n_clus[detpos];
0112 e_clus[detpos] += iter->charge();
0113 }
0114 }
0115 }
0116
0117 for (int i = 0; i < 3; ++i) {
0118 if (n_clus[i])
0119 e_clus[i] /= n_clus[i];
0120 }
0121
0122 if (e_clus[1] < bmincharge_)
0123 return false;
0124
0125 if (e_clus[0] + e_clus[2] != 0)
0126 {
0127 asym_pix_1 = e_clus[0] / (e_clus[0] + e_clus[2]);
0128 asym_pix_2 = e_clus[2] / (e_clus[0] + e_clus[2]);
0129 } else
0130 {
0131 return false;
0132 }
0133
0134 bool pkam_1 = (asym_pix_1 <= max_asym_ && asym_pix_1 >= min_asym_);
0135 bool pkam_2 = (asym_pix_2 <= max_asym_ && asym_pix_2 >= min_asym_);
0136
0137 if (pkam_1 || pkam_2)
0138 return accept;
0139
0140 return false;
0141 }
0142
0143
0144 #include "FWCore/Framework/interface/MakerMacros.h"
0145 DEFINE_FWK_MODULE(HLTPixelAsymmetryFilter);