diff --git a/DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h b/DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h index 15c932e1fcd44..faedbe792ee4d 100644 --- a/DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h +++ b/DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h @@ -1,7 +1,9 @@ #ifndef DataFormats_SiStripCluster_SiStripApproximateCluster_h #define DataFormats_SiStripCluster_SiStripApproximateCluster_h - +#include "FWCore/Utilities/interface/Exception.h" #include "FWCore/Utilities/interface/typedefs.h" +#include +#include class SiStripCluster; class SiStripApproximateCluster { @@ -13,25 +15,125 @@ class SiStripApproximateCluster { cms_uint8_t avgCharge, bool filter, bool isSaturated, - bool peakFilter = false) + bool peakFilter = false, + cms_uint8_t version_ = 1) : barycenter_(barycenter), width_(width), avgCharge_(avgCharge), filter_(filter), isSaturated_(isSaturated), - peakFilter_(peakFilter) {} + peakFilter_(peakFilter), + version_(version_) {} explicit SiStripApproximateCluster(const SiStripCluster& cluster, unsigned int maxNSat, float hitPredPos, - bool peakFilter); + bool peakFilter, + cms_uint8_t version_ = 1, + float previous_barycenter = 0, + unsigned int offset_module_change = 0); + + // getBarycenter returns the barycenter as a *float* in strips (e.g. 1.0 means center of strip 1) + float getBarycenter(float previous_barycenter = 0, unsigned int offset_module_change = 0) const { + switch (version_) { + case 1: + return barycenter_ * 0.1; // in the old format barycenter_ is in tenths of strips + case 2: { + // Drop the first bit (encoding the saturation info) + // barycenterRangeMax_ is 1 in 15 bits (the 16th bit is used to encode isSaturated_ info) + double barycenter_decoded = (barycenter_ & barycenterRangeMax_); + // rescale to get the original barycenter value (float) + return barycenter_decoded / barycenterScale_ - offset_module_change + previous_barycenter; + } + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) << " of SiStripApproximateCluster not supported"; + } + } + // // kept for compatibility with v1, should be removed in the future + // cms_uint16_t barycenter(float previous_barycenter, unsigned int offset_module_change) const { + // switch (version_){ + // case 2: { + // return std::round(getBarycenter(previous_barycenter, offset_module_change)*10.); // return barycenter in tenths of strips for compatibility with v1 + // } + // default: throw cms::Exception("VersionNotSupported") << "Version " << int(version_) << " of SiStripApproximateCluster not supported for SiStripApproximateCluster::barycenter(float,unsigned int)"; + // } + // } + + // getAvgCharge returns the average charge as a *float* in ADC counts (e.g. 1.0 means 1 ADC count) + float getAvgCharge() const { + switch (version_) { + case 1: + return avgCharge_; + case 2: { + // Drop the first two bits (encoding filter_ and peakFilter_ info) + // avgChargeRangeMax_ is 1 in 6 bits (the 2 highest bits are used to encode filter_ and peakFilter_ info) + double avgCharge_decoded = (avgCharge_ & avgChargeRangeMax_); + // Rescale to get the original average charge (float) + return (avgCharge_decoded)*avgChargeScale_ + avgChargeOffset_; + } + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) << " of SiStripApproximateCluster not supported"; + } + } + + //barycenter() returns barycenter position in tenths of strip (i.e. 10 means center of strip 1) (0-1536) + //avgCharge() returns the average charge in ADC counts (0-255) + //width() returns the cluster width (0-255) + //version() returns true if the cluster is in the new format (Fall 2025) cms_uint16_t barycenter() const { return barycenter_; } cms_uint8_t width() const { return width_; } - cms_uint8_t avgCharge() const { return avgCharge_; } - bool filter() const { return filter_; } - bool isSaturated() const { return isSaturated_; } - bool peakFilter() const { return peakFilter_; } + cms_uint8_t avgCharge() const { // should be a float in the future instead of an int + switch (version_) { + case 1: + return avgCharge_; + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) + << " of SiStripApproximateCluster not supported for SiStripApproximateCluster::avgCharge()"; + } + } + + bool filter() const { + switch (version_) { + case 1: + return filter_; + // In v2, filter_ info are encoded in avgCharge_ + case 2: + return (avgCharge_ & (1 << kFilterBit)); + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) << " of SiStripApproximateCluster not supported"; + } + } + bool peakFilter() const { + switch (version_) { + case 1: + return peakFilter_; + // In v2, filter_ and peakFilter_ info are encoded in avgCharge_ + case 2: + return (avgCharge_ & (1 << kPeakFilterBit)); + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) << " of SiStripApproximateCluster not supported"; + } + } + + bool isSaturated() const { + switch (version_) { + case 1: + return isSaturated_; + // In v2, isSaturated_ info is encoded in barycenter_ + case 2: + return (barycenter_ & (1 << kSaturatedBit)); + default: + throw cms::Exception("VersionNotSupported") + << "Version " << int(version_) << " of SiStripApproximateCluster not supported"; + } + } + char version() const { return version_; } private: cms_uint16_t barycenter_ = 0; @@ -40,10 +142,43 @@ class SiStripApproximateCluster { bool filter_ = false; bool isSaturated_ = false; bool peakFilter_ = false; + // v2 --> new version + cms_uint8_t version_ = 1; static constexpr double trimMaxADC_ = 30.; static constexpr double trimMaxFracTotal_ = .15; static constexpr double trimMaxFracNeigh_ = .25; static constexpr double maxTrimmedSizeDiffNeg_ = .7; static constexpr double maxTrimmedSizeDiffPos_ = 1.; + + ////// Encoding constants for v2 /////////// + // maximum value of barycenter_ is 768 strips (128 strips/APV * 6 APVs) + // multiplied by a factor 2 as we save the distance from the previous cluster, which can be in another module + static constexpr double barycenterMax_ = 768. * 2; + // get the total number of bits in barycenter_ (16 bits for cms_uint16_t) + static constexpr int nbits_barycenter_ = sizeof(barycenter_) * CHAR_BIT; + // position of the bit used to encode isSaturated_ in barycenter_ + static constexpr int kSaturatedBit = nbits_barycenter_ - 1; + // get the largest number storable in barycenter_ with the remaining bits (2^15 -1 = 32767) + static constexpr int barycenterRangeMax_ = (1 << (nbits_barycenter_ - 1)) - 1; + + // maximum value of avgCharge_ is 255 ADC counts + static constexpr double avgChargeMax_ = 255.; + // get the number of bits in avgCharge_ (8 bits for cms_uint8_t) + static constexpr int nbits_avgCharge_ = sizeof(avgCharge_) * CHAR_BIT; + // positions of the bit used to encode filter_ and peakFilter_ in avgCharge_ + static constexpr int kPeakFilterBit = nbits_avgCharge_ - 1; + static constexpr int kFilterBit = nbits_avgCharge_ - 2; + // get the largest number storable in avgCharge_ with the remaining bits (2^6 -1 = 63) + static constexpr int avgChargeRangeMax_ = (1 << (nbits_avgCharge_ - 2)) - 1; + +public: + // constants used in compression and decompression, + // floor approximate to integer to get no approximation error for integers + static constexpr int barycenterScale_ = barycenterRangeMax_ / barycenterMax_; + static constexpr int avgChargeScale_ = avgChargeMax_ / avgChargeRangeMax_; + + static constexpr float barycenterOffset_ = +0.5; // to get no approximation error for integers + static constexpr float avgChargeOffset_ = +2.0; // to minimize bias + //////////////////////////////////////// }; #endif // DataFormats_SiStripCluster_SiStripApproximateCluster_h diff --git a/DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h b/DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h new file mode 100644 index 0000000000000..495d75cec31c0 --- /dev/null +++ b/DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h @@ -0,0 +1,111 @@ +#ifndef DataFormats_SiStripCluster_SiStripApproximateClusterCollectionV2_v1_h +#define DataFormats_SiStripCluster_SiStripApproximateClusterCollectionV2_v1_h + +#include + +#include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h" +#include +#include +/** + * This class provides a minimal interface that resembles + * edmNew::DetSetVector, but is crafted such that we are comfortable + * to provide an infinite backwards compatibility guarantee for it + * (like all RAW data). Any modifications need to be made with care. + * Please consult core software group if in doubt. +**/ +class SiStripApproximateClusterCollectionV2 { +public: + // Helper classes to make creation and iteration easier + class Filler { + public: + void push_back(SiStripApproximateCluster const& cluster) { clusters_.push_back(cluster); } + + private: + friend SiStripApproximateClusterCollectionV2; + Filler(std::vector& clusters) : clusters_(clusters) {} + + std::vector& clusters_; + }; + + class const_iterator; + class DetSet { + public: + using const_iterator = std::vector::const_iterator; + + unsigned int id() const { + return std::accumulate(coll_->detIds_.cbegin(), coll_->detIds_.cbegin() + detIndex_ + 1, 0); + } + + void move(unsigned int clusBegin) const { clusBegin_ = clusBegin; } + const_iterator begin() const { return coll_->clusters_.begin() + clusBegin_; } + const_iterator cbegin() const { return begin(); } + const_iterator end() const { return coll_->clusters_.begin() + clusEnd_; } + const_iterator cend() const { return end(); } + + private: + friend SiStripApproximateClusterCollectionV2::const_iterator; + DetSet(SiStripApproximateClusterCollectionV2 const* coll, unsigned int detIndex) + : coll_(coll), + detIndex_(detIndex), + clusEnd_(coll->clusters_.size()) + //clusBegin_(coll_->beginIndices_[detIndex]), + //clusEnd_(detIndex == coll_->beginIndices_.size() - 1 ? coll->clusters_.size() + //: coll_->beginIndices_[detIndex + 1]) {} + {} + + SiStripApproximateClusterCollectionV2 const* const coll_; + unsigned int const detIndex_; + mutable unsigned int clusBegin_ = 0; + unsigned int const clusEnd_; + }; + + class const_iterator { + public: + DetSet operator*() const { return DetSet(coll_, index_); } + + const_iterator& operator++() { + ++index_; + if (index_ == coll_->detIds_.size()) { + *this = const_iterator(); + } + return *this; + } + + const_iterator operator++(int) { + const_iterator clone = *this; + ++(*this); + return clone; + } + + bool operator==(const_iterator const& other) const { return coll_ == other.coll_ and index_ == other.index_; } + bool operator!=(const_iterator const& other) const { return not operator==(other); } + + private: + friend SiStripApproximateClusterCollectionV2; + // default-constructed object acts as the sentinel + const_iterator() = default; + const_iterator(SiStripApproximateClusterCollectionV2 const* coll) : coll_(coll) {} + + SiStripApproximateClusterCollectionV2 const* coll_ = nullptr; + unsigned int index_ = 0; + }; + + // Actual public interface + SiStripApproximateClusterCollectionV2() = default; + + void reserve(std::size_t dets, std::size_t clusters); + Filler beginDet(unsigned int detId); + + const_iterator begin() const { return clusters_.empty() ? end() : const_iterator(this); } + const_iterator cbegin() const { return begin(); } + const_iterator end() const { return const_iterator(); } + const_iterator cend() const { return end(); } + +private: + // The detIds_ and beginIndices_ have one element for each Det. An + // element of beginIndices_ points to the first cluster of the Det + // in clusters_. + std::vector detIds_; // DetId for the Det + std::vector clusters_; +}; +#endif diff --git a/DataFormats/SiStripCluster/interface/SiStripCluster.h b/DataFormats/SiStripCluster/interface/SiStripCluster.h index 83ada6f6cfc1d..72a965dce05ca 100644 --- a/DataFormats/SiStripCluster/interface/SiStripCluster.h +++ b/DataFormats/SiStripCluster/interface/SiStripCluster.h @@ -44,8 +44,10 @@ class SiStripCluster { initQB(); } - SiStripCluster(const SiStripApproximateCluster cluster, const uint16_t maxStrips); - + SiStripCluster(const SiStripApproximateCluster cluster, + const uint16_t maxStrips, + float previous_barycenter = 0, + unsigned int offset_module_change = 0); // extend the cluster template void extend(Iter begin, Iter end) { diff --git a/DataFormats/SiStripCluster/src/SiStripApproximateCluster.cc b/DataFormats/SiStripCluster/src/SiStripApproximateCluster.cc index 01b30963a083a..31f5b2ccf45c8 100644 --- a/DataFormats/SiStripCluster/src/SiStripApproximateCluster.cc +++ b/DataFormats/SiStripCluster/src/SiStripApproximateCluster.cc @@ -1,18 +1,23 @@ #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h" #include "DataFormats/SiStripCluster/interface/SiStripCluster.h" #include +#include #include SiStripApproximateCluster::SiStripApproximateCluster(const SiStripCluster& cluster, unsigned int maxNSat, float hitPredPos, - bool peakFilter) { + bool peakFilter, + unsigned char version, + float previous_barycenter, + unsigned int offset_module_change) { barycenter_ = std::round(cluster.barycenter() * 10); width_ = cluster.size(); avgCharge_ = cluster.charge() / cluster.size(); filter_ = false; isSaturated_ = false; peakFilter_ = peakFilter; + version_ = version; //mimicing the algorithm used in StripSubClusterShapeTrajectoryFilter... //Looks for 3 adjacent saturated strips (ADC>=254) @@ -54,4 +59,28 @@ SiStripApproximateCluster::SiStripApproximateCluster(const SiStripCluster& clust } else { filter_ = peakFilter_; } + + if (version_ == 2) { + // Compression of avgCharge_ to integer + avgCharge_ = floor((float(cluster.charge()) / cluster.size()) / avgChargeScale_); + // In v2, we encode the filter_ and peakFilter_ info in avgCharge_ as the two highest bits + assert(avgCharge_ <= avgChargeRangeMax_ && "Setting avgCharge > avgChargeRangeMax_"); + // filter_ and peakFilter_ are encoded in the two highest bits of avgCharge_ + avgCharge_ = (avgCharge_ | (filter_ << kFilterBit)); + assert(avgCharge_ <= ((1 << (nbits_avgCharge_ - 1)) - 1) && "Setting avgCharge with filter > max single bit range"); + avgCharge_ = (avgCharge_ | (peakFilter_ << kPeakFilterBit)); + assert(avgCharge_ <= ((1 << nbits_avgCharge_) - 1) && "Setting avgCharge with peakFilter > max full range"); + + // Compression of barycenter_ to integer [note: it contains the distance from the previous cluster] + barycenter_ = round(float(cluster.barycenter() - previous_barycenter + (offset_module_change)) * barycenterScale_); + assert(barycenter_ <= barycenterRangeMax_ && "Setting barycenter > barycenterRangeMax_"); + // isSaturated_ is encoded in the highest bit of barycenter_ + barycenter_ = (barycenter_ | (isSaturated_ << kSaturatedBit)); + assert(barycenter_ <= ((1 << nbits_barycenter_) - 1) && "Setting barycenter with isSaturated > max full range"); + + // Flags set to false to reduce event size (they should be removed when moving to v2 only) + filter_ = false; + isSaturated_ = false; + peakFilter_ = false; + } } diff --git a/DataFormats/SiStripCluster/src/SiStripApproximateClusterCollectionV2.cc b/DataFormats/SiStripCluster/src/SiStripApproximateClusterCollectionV2.cc new file mode 100644 index 0000000000000..89ba705d223d0 --- /dev/null +++ b/DataFormats/SiStripCluster/src/SiStripApproximateClusterCollectionV2.cc @@ -0,0 +1,10 @@ +#include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h" +void SiStripApproximateClusterCollectionV2::reserve(std::size_t dets, std::size_t clusters) { + detIds_.reserve(dets); + clusters_.reserve(clusters); +} + +SiStripApproximateClusterCollectionV2::Filler SiStripApproximateClusterCollectionV2::beginDet(unsigned int detId) { + detIds_.push_back((detIds_.empty()) ? detId : detId - (std::accumulate(detIds_.cbegin(), detIds_.cend(), 0))); + return Filler(clusters_); +} diff --git a/DataFormats/SiStripCluster/src/SiStripCluster.cc b/DataFormats/SiStripCluster/src/SiStripCluster.cc index 8586307a8e384..a35d7fc171ec5 100644 --- a/DataFormats/SiStripCluster/src/SiStripCluster.cc +++ b/DataFormats/SiStripCluster/src/SiStripCluster.cc @@ -1,5 +1,6 @@ #include "FWCore/Utilities/interface/Likely.h" #include "DataFormats/SiStripCluster/interface/SiStripCluster.h" +#include "FWCore/Utilities/interface/Exception.h" SiStripCluster::SiStripCluster(const SiStripDigiRange& range) : firstStrip_(range.first->strip()), error_x(-99999.9) { std::vector v; @@ -23,10 +24,51 @@ SiStripCluster::SiStripCluster(const SiStripDigiRange& range) : firstStrip_(rang initQB(); } -SiStripCluster::SiStripCluster(const SiStripApproximateCluster cluster, const uint16_t maxStrips) : error_x(-99999.9) { - barycenter_ = cluster.barycenter() / 10.0; - charge_ = cluster.width() * cluster.avgCharge(); - amplitudes_.resize(cluster.width(), cluster.avgCharge()); +SiStripCluster::SiStripCluster(const SiStripApproximateCluster cluster, + const uint16_t maxStrips, + float previous_barycenter, + unsigned int offset_module_change) + : error_x(-99999.9) { + switch (cluster.version()) { + case 1: { + barycenter_ = cluster.barycenter() / 10.0; + charge_ = cluster.width() * cluster.avgCharge(); + amplitudes_.resize(cluster.width(), cluster.avgCharge()); + break; + } + case 2: { + barycenter_ = cluster.getBarycenter(previous_barycenter, offset_module_change); + charge_ = round(cluster.getAvgCharge() * cluster.width()); + int amplitude = ceil(cluster.getAvgCharge()); + amplitudes_.resize(cluster.width(), amplitude); + int excessToBeFixed = amplitude * cluster.width() - charge_; + // Reduce the first and last amplitude to maintain the same average charge but reduce the width by 1 + if (cluster.width() > 1 && excessToBeFixed != 0) { + if (excessToBeFixed % 2 == 0) { + amplitudes_.front() = amplitude - excessToBeFixed / 2; + amplitudes_.back() = amplitude - excessToBeFixed / 2; + // assert(amplitude - excessToBeFixed / 2 >=0); + } else { + // if the excess is odd, we cannot divide it equally between the first and last strip + // we add one unit to the front and subtract one unit to the back to keep the total charge unchanged + if (int(cluster.getAvgCharge()) % 2 == 0) { // avoid bias with rounding + amplitudes_.front() = amplitude - (excessToBeFixed / 2 + 1); + amplitudes_.back() = amplitude - (excessToBeFixed / 2); + // assert(amplitude - (excessToBeFixed / 2 + 1) >=0); + } else { + amplitudes_.back() = amplitude - (excessToBeFixed / 2 + 1); + amplitudes_.front() = amplitude - (excessToBeFixed / 2); + // assert(amplitude - (excessToBeFixed / 2 + 1) >=0); + } + } + } + // assert( int(std::accumulate(amplitudes_.begin(), amplitudes_.end(), 0.0f)) == charge_); + break; + } + default: + throw cms::Exception("VersionNotSupported") + << "SiStripCluster.cc. Version " << int(cluster.version()) << " of SiStripApproximateCluster not supported"; + } filter_ = cluster.filter(); float halfwidth_ = 0.5f * float(cluster.width()); diff --git a/DataFormats/SiStripCluster/src/classes.h b/DataFormats/SiStripCluster/src/classes.h index 94eafdf2ecac7..798abef58b6f5 100644 --- a/DataFormats/SiStripCluster/src/classes.h +++ b/DataFormats/SiStripCluster/src/classes.h @@ -8,6 +8,7 @@ #include "DataFormats/SiStripCluster/interface/SiStripClustersSOA.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollection.h" +#include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h" #include "DataFormats/Common/interface/ContainerMask.h" #endif // SISTRIPCLUSTER_CLASSES_H diff --git a/DataFormats/SiStripCluster/src/classes_def.xml b/DataFormats/SiStripCluster/src/classes_def.xml index 31ff33a914296..5e65fe84a3948 100755 --- a/DataFormats/SiStripCluster/src/classes_def.xml +++ b/DataFormats/SiStripCluster/src/classes_def.xml @@ -32,7 +32,8 @@ - + + @@ -59,6 +60,12 @@ + + + + + + diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/SiStripApprox2Clusters.cc b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripApprox2Clusters.cc index 803c8949f90c6..78bc7c38ed61a 100644 --- a/RecoLocalTracker/SiStripClusterizer/plugins/SiStripApprox2Clusters.cc +++ b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripApprox2Clusters.cc @@ -1,6 +1,7 @@ #include "DataFormats/Common/interface/DetSetVectorNew.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollection.h" +#include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h" #include "DataFormats/SiStripCluster/interface/SiStripCluster.h" #include "FWCore/Framework/interface/Event.h" #include "FWCore/Framework/interface/Frameworkfwd.h" @@ -13,6 +14,9 @@ #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h" #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" +#include "CalibFormats/SiStripObjects/interface/SiStripDetInfo.h" +#include "CalibTracker/SiStripCommon/interface/SiStripDetInfoFileReader.h" +#include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" #include #include @@ -25,25 +29,60 @@ class SiStripApprox2Clusters : public edm::global::EDProducer<> { static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); private: - edm::EDGetTokenT clusterToken_; + template + void processCollection(const CollectionType& clusterCollection, + edmNew::DetSetVector& result, + const edm::EventSetup& iSetup) const; + + unsigned int collectionVersion; + edm::EDGetTokenT clusterTokenV1_; + edm::EDGetTokenT clusterTokenV2_; edm::ESGetToken tkGeomToken_; + SiStripDetInfo detInfo_; }; SiStripApprox2Clusters::SiStripApprox2Clusters(const edm::ParameterSet& conf) { - clusterToken_ = consumes(conf.getParameter("inputApproxClusters")); + const auto inputTag = conf.getParameter("inputApproxClusters"); + collectionVersion = conf.getParameter("collectionVersion"); + + clusterTokenV1_ = consumes(inputTag); + clusterTokenV2_ = consumes(inputTag); + tkGeomToken_ = esConsumes(); + detInfo_ = SiStripDetInfoFileReader::read(edm::FileInPath(SiStripDetInfoFileReader::kDefaultFile).fullPath()); produces>(); } void SiStripApprox2Clusters::produce(edm::StreamID id, edm::Event& event, const edm::EventSetup& iSetup) const { auto result = std::make_unique>(); - const auto& clusterCollection = event.get(clusterToken_); + if (collectionVersion == 1) { + const auto& clusterCollection = event.get(clusterTokenV1_); + processCollection(clusterCollection, *result, iSetup); + } else if (collectionVersion == 2) { + const auto& clusterCollection = event.get(clusterTokenV2_); + processCollection(clusterCollection, *result, iSetup); + } else { + throw cms::Exception("InvalidParameter") + << "Invalid collectionVersion: " << collectionVersion << ". Must be 1 or 2."; + } + + event.put(std::move(result)); +} + +template +void SiStripApprox2Clusters::processCollection(const CollectionType& clusterCollection, + edmNew::DetSetVector& result, + const edm::EventSetup& iSetup) const { const auto& tkGeom = &iSetup.getData(tkGeomToken_); const auto& tkDets = tkGeom->dets(); + float previous_barycenter = SiStripApproximateCluster::barycenterOffset_; + unsigned int offset_module_change = 0; + unsigned int clusBegin = 0; + for (const auto& detClusters : clusterCollection) { - edmNew::DetSetVector::FastFiller ff{*result, detClusters.id()}; + edmNew::DetSetVector::FastFiller ff{result, detClusters.id()}; unsigned int detId = detClusters.id(); uint16_t nStrips{0}; @@ -53,17 +92,36 @@ void SiStripApprox2Clusters::produce(edm::StreamID id, edm::Event& event, const const StripTopology& p = dynamic_cast(*det)->specificTopology(); nStrips = p.nstrips() - 1; + if constexpr (std::is_same::value) { + detClusters.move(clusBegin); + } for (const auto& cluster : detClusters) { - ff.push_back(SiStripCluster(cluster, nStrips)); + const auto convertedCluster = SiStripCluster(cluster, nStrips, previous_barycenter, offset_module_change); + if ((convertedCluster.barycenter()) >= nStrips + 1) { + // this should not happen for V1 + if (collectionVersion == 1) + throw cms::Exception("DataCorrupt") + << "SiStripApprox2Clusters: cluster with barycenter " << convertedCluster.barycenter() + << " out of range for module with " << nStrips + 1 << " strips."; + // in V2 is used to split clusters across modules + break; + } + previous_barycenter = convertedCluster.barycenter(); + ++clusBegin; + offset_module_change = 0; + + ff.push_back(convertedCluster); } + offset_module_change = detInfo_.getNumberOfApvsAndStripLength(detId).first * sistrip::STRIPS_PER_APV; } - - event.put(std::move(result)); } void SiStripApprox2Clusters::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { edm::ParameterSetDescription desc; desc.add("inputApproxClusters", edm::InputTag("siStripClusters")); + desc.add( + "collectionVersion", + 1); // Collection version (1 for SiStripApproximateClusterCollection, 2 for SiStripApproximateClusterCollectionV2) descriptions.add("SiStripApprox2Clusters", desc); } diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClusters2ApproxClusters.cc b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClusters2ApproxClusters.cc index e50b5886686cb..b174644b08700 100644 --- a/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClusters2ApproxClusters.cc +++ b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClusters2ApproxClusters.cc @@ -10,6 +10,7 @@ #include "DataFormats/GeometryVector/interface/LocalPoint.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateCluster.h" #include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollection.h" +#include "DataFormats/SiStripCluster/interface/SiStripApproximateClusterCollectionV2.h" #include "DataFormats/SiStripCluster/interface/SiStripCluster.h" #include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" #include "DataFormats/TrackReco/interface/Track.h" @@ -41,6 +42,12 @@ class SiStripClusters2ApproxClusters : public edm::stream::EDProducer<> { static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); private: + template + void fillCollection(CollectionType& result, + const edmNew::DetSetVector& clusterCollection, + edm::Event& event, + const edm::EventSetup& iSetup); + edm::InputTag inputClusters; edm::EDGetTokenT > clusterToken; @@ -60,6 +67,8 @@ class SiStripClusters2ApproxClusters : public edm::stream::EDProducer<> { SiStripDetInfo detInfo_; std::string csfLabel_; + unsigned int version; + unsigned int collectionVersion; edm::ESGetToken csfToken_; edm::ESGetToken stripNoiseToken_; @@ -83,15 +92,52 @@ SiStripClusters2ApproxClusters::SiStripClusters2ApproxClusters(const edm::Parame csfLabel_ = conf.getParameter("clusterShapeHitFilterLabel"); csfToken_ = esConsumes(edm::ESInputTag("", csfLabel_)); + version = conf.getParameter("version"); + collectionVersion = conf.getParameter("collectionVersion"); + + // Validate version parameter + if (version != 1 && version != 2) { + throw cms::Exception("InvalidParameter") << "Invalid version: " << version << ". Must be 1 or 2."; + } + + // Validate collectionVersion parameter + if (collectionVersion != 1 && collectionVersion != 2) { + throw cms::Exception("InvalidParameter") + << "Invalid collectionVersion: " << collectionVersion << ". Must be '1' or '2'."; + } + stripNoiseToken_ = esConsumes(); - produces(); + + if (collectionVersion == 1) { + produces(); + } else if (collectionVersion == 2) { + produces(); + } } void SiStripClusters2ApproxClusters::produce(edm::Event& event, edm::EventSetup const& iSetup) { const auto& clusterCollection = event.get(clusterToken); - auto result = std::make_unique(); - result->reserve(clusterCollection.size(), clusterCollection.dataSize()); + if (collectionVersion == 1) { + auto result = std::make_unique(); + result->reserve(clusterCollection.size(), clusterCollection.dataSize()); + + fillCollection(*result, clusterCollection, event, iSetup); + event.put(std::move(result)); + } else if (collectionVersion == 2) { + auto result = std::make_unique(); + result->reserve(clusterCollection.size(), clusterCollection.dataSize()); + + fillCollection(*result, clusterCollection, event, iSetup); + event.put(std::move(result)); + } +} + +template +void SiStripClusters2ApproxClusters::fillCollection(CollectionType& result, + const edmNew::DetSetVector& clusterCollection, + edm::Event& event, + const edm::EventSetup& iSetup) { auto const beamSpotHandle = event.getHandle(beamSpotToken_); auto const& bs = beamSpotHandle.isValid() ? *beamSpotHandle : reco::BeamSpot(); if (not beamSpotHandle.isValid()) { @@ -103,8 +149,11 @@ void SiStripClusters2ApproxClusters::produce(edm::Event& event, edm::EventSetup const auto& theFilter = &iSetup.getData(csfToken_); const auto& theNoise_ = &iSetup.getData(stripNoiseToken_); + float previous_barycenter = SiStripApproximateCluster::barycenterOffset_; + unsigned int offset_module_change = 0; + for (const auto& detClusters : clusterCollection) { - auto ff = result->beginDet(detClusters.id()); + auto ff = result.beginDet(detClusters.id()); unsigned int detId = detClusters.id(); const GeomDet* det = tkGeom->idToDet(detId); @@ -132,7 +181,10 @@ void SiStripClusters2ApproxClusters::produce(edm::Event& event, edm::EventSetup bool isTrivial = (std::abs(hitPredPos) < 2.f && hitStrips <= 2); if (!usable || isTrivial) { - ff.push_back(SiStripApproximateCluster(cluster, maxNSat, hitPredPos, true)); + SiStripApproximateCluster approxCluster( + cluster, maxNSat, hitPredPos, true, version, previous_barycenter, offset_module_change); + ff.push_back(approxCluster); + previous_barycenter = approxCluster.getBarycenter(previous_barycenter, offset_module_change); } else { bool peakFilter = false; SlidingPeakFinder pf(std::max(2, std::ceil(std::abs(hitPredPos) + subclusterWindow_))); @@ -147,12 +199,15 @@ void SiStripClusters2ApproxClusters::produce(edm::Event& event, edm::EventSetup subclusterCutSN_); peakFilter = pf.apply(cluster.amplitudes(), test); - ff.push_back(SiStripApproximateCluster(cluster, maxNSat, hitPredPos, peakFilter)); + SiStripApproximateCluster approxCluster( + cluster, maxNSat, hitPredPos, peakFilter, version, previous_barycenter, offset_module_change); + ff.push_back(approxCluster); + previous_barycenter = approxCluster.getBarycenter(previous_barycenter, offset_module_change); } + offset_module_change = 0; } + offset_module_change = nApvs * sistrip::STRIPS_PER_APV; } - - event.put(std::move(result)); } void SiStripClusters2ApproxClusters::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { @@ -161,6 +216,10 @@ void SiStripClusters2ApproxClusters::fillDescriptions(edm::ConfigurationDescript desc.add("maxSaturatedStrips", 3); desc.add("clusterShapeHitFilterLabel", "ClusterShapeHitFilter"); // add CSF label desc.add("beamSpot", edm::InputTag("offlineBeamSpot")); // add BeamSpot tag + desc.add("version", 1); // RawPrime version (1= default, 2= new v2 format) + desc.add( + "collectionVersion", + 1); // Collection version (1 for SiStripApproximateClusterCollection, 2 for SiStripApproximateClusterCollectionV2) descriptions.add("SiStripClusters2ApproxClusters", desc); }