diff --git a/DataFormats/SiStripCluster/interface/SiStripCluster.h b/DataFormats/SiStripCluster/interface/SiStripCluster.h index 245fcdaa25969..7ac6cb96a354e 100644 --- a/DataFormats/SiStripCluster/interface/SiStripCluster.h +++ b/DataFormats/SiStripCluster/interface/SiStripCluster.h @@ -30,6 +30,9 @@ class SiStripCluster { initQB(); } + SiStripCluster(uint16_t firstStrip, std::vector&& data, float barycenter, float charge) + : amplitudes_(std::move(data)), firstStrip_(firstStrip), barycenter_(barycenter), charge_(charge) {} + template SiStripCluster(const uint16_t& firstStrip, Iter begin, Iter end) : amplitudes_(begin, end), firstStrip_(firstStrip) { initQB(); diff --git a/DataFormats/SiStripClusterSoA/BuildFile.xml b/DataFormats/SiStripClusterSoA/BuildFile.xml new file mode 100644 index 0000000000000..0fa8c933eeecd --- /dev/null +++ b/DataFormats/SiStripClusterSoA/BuildFile.xml @@ -0,0 +1,8 @@ + + + + + + + + diff --git a/DataFormats/SiStripClusterSoA/README.md b/DataFormats/SiStripClusterSoA/README.md new file mode 100644 index 0000000000000..f1fff393d9d05 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/README.md @@ -0,0 +1,26 @@ +# SiStripClusterSoA +The `SiStripClusterHost`/`SiStripClusterDevice` is a portable collection based on the `SiStripClusterSoALayout` (`DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h`). + +It is used to store the collection of SiStrip cluster candidates from the heterogeneous SiStripClusterizer module (`RecoLocalTracker/SiStripClusterizer`). + +## Data members +The fields in the structure have the following meaning: + +| SoA type | C-type | Name | Description | +| --- | --- | --- | --- | +| column | uint32_t | clusterIndex | Index for the first strip amplitude in this cluster candidate, to be fetched by the Digi collection | +| column | uint16_t | clusterSize | Cluster cand. strip count | +| column | uint32_t | clusterDetId | Cluster cand. detector ID | +| column | uint16_t | firstStrip | First strip ID of the cluster cand. | +| column | bool | candidateAccepted | Is the cluster candidate a good one? [^ThreeThresholdAlgorithm.cc#L103-L107] | +| column | float | barycenter | Cluster cand. barycenter [^SiStripCluster.h#L164] | +| column | float | charge | Cluster cand. charge [^SiStripCluster.h#L152-L159] | +| column | uint32_t | candidateAcceptedPrefix | Prefix sum of the candidateAccepted colum, used to index the good clusters | +| scalar | uint32_t | nClusterCandidates | Number of cluster candidates in the collection [^noteCand] | +| scalar | uint32_t | maxClusterSize | Max number of contiguous strips for clustering (setup) | + + +[^ThreeThresholdAlgorithm.cc#L103-L107]: [RecoLocalTracker/SiStripClusterizer/src/ThreeThresholdAlgorithm.cc](https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/RecoLocalTracker/SiStripClusterizer/src/ThreeThresholdAlgorithm.cc#L103-L107) +[^SiStripCluster.h#L164]: [DataFormats/SiStripCluster/interface/SiStripCluster.h#L164](https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/DataFormats/SiStripCluster/interface/SiStripCluster.h#L164) +[^SiStripCluster.h#L152-L159]: [DataFormats/SiStripCluster/interface/SiStripCluster.h#L152-L159](https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/DataFormats/SiStripCluster/interface/SiStripCluster.h#L152-L159) +[^noteCand]: this is typically lower than the collection pre-allocated size (i.e., from `collection->metadata().size()`). It indicates the number of non-contiguous strip seeds over which the clustering is perfomed. It can be used while loop over the collection to break earlier than the collection size. diff --git a/DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h b/DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h new file mode 100644 index 0000000000000..dac2e12c57940 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h @@ -0,0 +1,12 @@ +#ifndef DataFormats_SiStripClusterSoA_interface_SiStripClusterHost_h +#define DataFormats_SiStripClusterSoA_interface_SiStripClusterHost_h + +#include "DataFormats/Portable/interface/PortableHostCollection.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h" + +namespace sistrip { + // SoA with SiStripCluster fields in host memory + using SiStripClusterHost = PortableHostCollection; +} // namespace sistrip + +#endif // DataFormats_SiStripClusterSoA_interface_SiStripClusterHost_h diff --git a/DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h b/DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h new file mode 100644 index 0000000000000..7016df42353a0 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h @@ -0,0 +1,26 @@ + +#ifndef DataFormats_SiStripClusterSoA_interface_SiStripClusterSoA_h +#define DataFormats_SiStripClusterSoA_interface_SiStripClusterSoA_h + +#include "DataFormats/SoATemplate/interface/SoALayout.h" + +namespace sistrip { + + GENERATE_SOA_LAYOUT(SiStripClusterSoALayout, + SOA_COLUMN(uint32_t, clusterIndex), + SOA_COLUMN(uint16_t, clusterSize), + SOA_COLUMN(uint32_t, clusterDetId), + SOA_COLUMN(uint16_t, firstStrip), + SOA_COLUMN(bool, candidateAccepted), + SOA_COLUMN(float, barycenter), + SOA_COLUMN(float, charge), + SOA_COLUMN(uint32_t, candidateAcceptedPrefix), + SOA_SCALAR(uint32_t, nClusterCandidates), + SOA_SCALAR(uint32_t, maxClusterSize)) + + using SiStripClusterSoA = SiStripClusterSoALayout<>; + using SiStripClusterView = SiStripClusterSoA::View; + using SiStripClusterConstView = SiStripClusterSoA::ConstView; +} // namespace sistrip + +#endif // DataFormats_SiStripClusterSoA_interface_SiStripClustersSoA_h diff --git a/DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h b/DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h new file mode 100644 index 0000000000000..0f8c4d300ca53 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h @@ -0,0 +1,19 @@ +#ifndef DataFormats_SiStripClusterSoA_interface_alpaka_SiStripClustersDevice_h +#define DataFormats_SiStripClusterSoA_interface_alpaka_SiStripClustersDevice_h + +#include "DataFormats/Portable/interface/alpaka/PortableCollection.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h" +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + // make the names from the top-level sistrip namespace visible for unqualified lookup + // inside the ALPAKA_ACCELERATOR_NAMESPACE::sistrip namespace + using namespace ::sistrip; + using SiStripClusterDevice = PortableCollection; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +// check that the sistrip device collection for the host device is the same as the sistrip host collection +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::SiStripClusterDevice, sistrip::SiStripClusterHost); + +#endif // DataFormats_SiStripClusterSoA_interface_alpaka_SiStripClustersDevice_h diff --git a/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda.h b/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda.h new file mode 100644 index 0000000000000..625b171abf832 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda.h @@ -0,0 +1,9 @@ +#ifndef DataFormats_SiStripClusterSoA_src_alpaka_classes_cuda_h +#define DataFormats_SiStripClusterSoA_src_alpaka_classes_cuda_h + +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/Common/interface/DeviceProduct.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h" +#include "DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h" + +#endif diff --git a/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda_def.xml b/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda_def.xml new file mode 100644 index 0000000000000..db173765b65f0 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/alpaka/classes_cuda_def.xml @@ -0,0 +1,5 @@ + + + + + diff --git a/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm.h b/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm.h new file mode 100644 index 0000000000000..899b184219c5a --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm.h @@ -0,0 +1,9 @@ +#ifndef DataFormats_SiStripClusterSoA_src_alpaka_classes_rocm_h +#define DataFormats_SiStripClusterSoA_src_alpaka_classes_rocm_h + +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/Common/interface/DeviceProduct.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterSoA.h" +#include "DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h" + +#endif diff --git a/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm_def.xml b/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm_def.xml new file mode 100644 index 0000000000000..2a0dfaa604c48 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/alpaka/classes_rocm_def.xml @@ -0,0 +1,5 @@ + + + + + diff --git a/DataFormats/SiStripClusterSoA/src/classes.cc b/DataFormats/SiStripClusterSoA/src/classes.cc new file mode 100644 index 0000000000000..7e0fdf398a2b4 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/classes.cc @@ -0,0 +1,4 @@ +#include "DataFormats/Portable/interface/PortableHostCollectionReadRules.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h" + +SET_PORTABLEHOSTCOLLECTION_READ_RULES(sistrip::SiStripClusterHost); diff --git a/DataFormats/SiStripClusterSoA/src/classes.h b/DataFormats/SiStripClusterSoA/src/classes.h new file mode 100644 index 0000000000000..f66185247bcdd --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/classes.h @@ -0,0 +1,7 @@ +#ifndef DataFormats_SiStripClusterSoA_src_classes_h +#define DataFormats_SiStripClusterSoA_src_classes_h + +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h" + +#endif // DataFormats_SiStripClusterSoA_src_classes_h diff --git a/DataFormats/SiStripClusterSoA/src/classes_def.xml b/DataFormats/SiStripClusterSoA/src/classes_def.xml new file mode 100644 index 0000000000000..a5fee8d19ae11 --- /dev/null +++ b/DataFormats/SiStripClusterSoA/src/classes_def.xml @@ -0,0 +1,7 @@ + + + + + + + diff --git a/DataFormats/SiStripDigiSoA/BuildFile.xml b/DataFormats/SiStripDigiSoA/BuildFile.xml new file mode 100644 index 0000000000000..0fa8c933eeecd --- /dev/null +++ b/DataFormats/SiStripDigiSoA/BuildFile.xml @@ -0,0 +1,8 @@ + + + + + + + + diff --git a/DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h b/DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h new file mode 100644 index 0000000000000..198a821e278ca --- /dev/null +++ b/DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h @@ -0,0 +1,12 @@ +#ifndef DataFormats_SiStripDigiSoA_interface_SiStripDigiHost_h +#define DataFormats_SiStripDigiSoA_interface_SiStripDigiHost_h + +#include "DataFormats/Portable/interface/PortableHostCollection.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h" + +namespace sistrip { + // SoA with SiStripClusters fields in host memory + using SiStripDigiHost = PortableHostCollection; +} // namespace sistrip + +#endif diff --git a/DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h b/DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h new file mode 100644 index 0000000000000..fa15cecbf7025 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h @@ -0,0 +1,20 @@ + +#ifndef DataFormats_SiStripDigiSoA_interface_SiStripDigiSoA_h +#define DataFormats_SiStripDigiSoA_interface_SiStripDigiSoA_h + +#include "DataFormats/SoATemplate/interface/SoALayout.h" + +namespace sistrip { + GENERATE_SOA_LAYOUT(SiStripDigiSoALayout, + SOA_COLUMN(uint8_t, adc), + SOA_COLUMN(uint16_t, channel), + SOA_COLUMN(uint16_t, stripId), + SOA_SCALAR(uint32_t, nbGoodCandidates), + SOA_SCALAR(uint32_t, nbCandidates)) + + using SiStripDigiSoA = SiStripDigiSoALayout<>; + using SiStripDigiView = SiStripDigiSoA::View; + using SiStripDigiConstView = SiStripDigiSoA::ConstView; +} // namespace sistrip + +#endif diff --git a/DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h b/DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h new file mode 100644 index 0000000000000..a91dcb3929d17 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h @@ -0,0 +1,18 @@ +#ifndef DataFormats_SiStripDigiSoA_interface_alpaka_SiStripDigiDevice_h +#define DataFormats_SiStripDigiSoA_interface_alpaka_SiStripDigiDevice_h + +#include "DataFormats/Portable/interface/alpaka/PortableCollection.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h" +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + // make the names from the top-level sistrip namespace visible for unqualified lookup + // inside the ALPAKA_ACCELERATOR_NAMESPACE::sistrip namespace + using namespace ::sistrip; + using SiStripDigiDevice = PortableCollection; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::SiStripDigiDevice, sistrip::SiStripDigiHost); + +#endif diff --git a/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda.h b/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda.h new file mode 100644 index 0000000000000..cdf51b439c92d --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda.h @@ -0,0 +1,8 @@ +#ifndef DataFormats_SiStripDigiSoA_src_alpaka_classes_cuda_h +#define DataFormats_SiStripDigiSoA_src_alpaka_classes_cuda_h + +#include "DataFormats/Common/interface/DeviceProduct.h" +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h" + +#endif diff --git a/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda_def.xml b/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda_def.xml new file mode 100644 index 0000000000000..b11906ad6a4b7 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/alpaka/classes_cuda_def.xml @@ -0,0 +1,5 @@ + + + + + diff --git a/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm.h b/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm.h new file mode 100644 index 0000000000000..078c441d6d3e4 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm.h @@ -0,0 +1,8 @@ +#ifndef DataFormats_SiStripDigiSoA_src_alpaka_classes_rocm_h +#define DataFormats_SiStripDigiSoA_src_alpaka_classes_rocm_h + +#include "DataFormats/Common/interface/DeviceProduct.h" +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h" + +#endif diff --git a/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm_def.xml b/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm_def.xml new file mode 100644 index 0000000000000..84133e9a43b4c --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/alpaka/classes_rocm_def.xml @@ -0,0 +1,5 @@ + + + + + diff --git a/DataFormats/SiStripDigiSoA/src/classes.cc b/DataFormats/SiStripDigiSoA/src/classes.cc new file mode 100644 index 0000000000000..b0e725d5f1fd9 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/classes.cc @@ -0,0 +1,4 @@ +#include "DataFormats/Portable/interface/PortableHostCollectionReadRules.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h" + +SET_PORTABLEHOSTCOLLECTION_READ_RULES(sistrip::SiStripDigiHost); diff --git a/DataFormats/SiStripDigiSoA/src/classes.h b/DataFormats/SiStripDigiSoA/src/classes.h new file mode 100644 index 0000000000000..08daeb82015cc --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/classes.h @@ -0,0 +1,8 @@ +#ifndef DataFormats_SiStripDigiSoA_src_classes_h +#define DataFormats_SiStripDigiSoA_src_classes_h + +#include "DataFormats/Common/interface/Wrapper.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiSoA.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h" + +#endif diff --git a/DataFormats/SiStripDigiSoA/src/classes_def.xml b/DataFormats/SiStripDigiSoA/src/classes_def.xml new file mode 100644 index 0000000000000..e8534e51ee9c7 --- /dev/null +++ b/DataFormats/SiStripDigiSoA/src/classes_def.xml @@ -0,0 +1,7 @@ + + + + + + + diff --git a/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h b/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h index 0c73e848fb9b7..88e5525f641f5 100644 --- a/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h +++ b/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h @@ -194,18 +194,18 @@ namespace sistrip { namespace detail { template - uint16_t getADC_W(const uint8_t* data, uint_fast16_t offset, uint8_t bits_shift) { + constexpr uint16_t getADC_W(const uint8_t* data, uint_fast16_t offset, uint8_t bits_shift) { // get ADC from one or two bytes (at most 10 bits), and shift if needed return (data[offset ^ 7] + (num_words == 2 ? ((data[(offset + 1) ^ 7] & 0x03) << 8) : 0)) << bits_shift; } template - uint16_t getADC_B2(const uint8_t* data, uint_fast16_t wOffset, uint_fast8_t bOffset) { + constexpr uint16_t getADC_B2(const uint8_t* data, uint_fast16_t wOffset, uint_fast8_t bOffset) { // get ADC from two bytes, from wOffset until bOffset bits from the next byte (maximum decided by mask) return (((data[wOffset ^ 7]) << bOffset) + (data[(wOffset + 1) ^ 7] >> (BITS_PER_BYTE - bOffset))) & mask; } template - uint16_t getADC_B1(const uint8_t* data, uint_fast16_t wOffset, uint_fast8_t bOffset) { + constexpr uint16_t getADC_B1(const uint8_t* data, uint_fast16_t wOffset, uint_fast8_t bOffset) { // get ADC from one byte, until bOffset into the byte at wOffset (maximum decided by mask) return (data[wOffset ^ 7] >> (BITS_PER_BYTE - bOffset)) & mask; } diff --git a/EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h b/EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h index a30349e6506cb..6446d8765b87a 100644 --- a/EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h +++ b/EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h @@ -615,12 +615,15 @@ namespace sistrip { //holds information about position of a channel in the buffer for use by unpacker class FEDChannel { public: - FEDChannel(const uint8_t* const data, const uint32_t offset, const uint16_t length); + enum class ZSROMode { nonLite = 7, lite = 2 }; + constexpr FEDChannel(const uint8_t* const data, const uint32_t offset, const uint16_t length); + constexpr FEDChannel(const uint8_t* const data, const uint32_t offset, const ZSROMode zsRoMode); //gets length from first 2 bytes (assuming normal FED channel) - FEDChannel(const uint8_t* const data, const uint32_t offset); - uint16_t length() const; - const uint8_t* data() const; - uint32_t offset() const; + constexpr FEDChannel(const uint8_t* const data, const uint32_t offset); + constexpr uint16_t stripsInCh(uint8_t num_bits) const; + constexpr uint16_t length() const; + constexpr const uint8_t* data() const; + constexpr uint32_t offset() const; /** * Retrieve the APV CM median for a non-lite zero-suppressed channel * @@ -628,15 +631,16 @@ namespace sistrip { * No additional checks are done here, so the caller should check * the readout mode and/or packet code. */ - uint16_t cmMedian(const uint8_t apvIndex) const; + constexpr uint16_t cmMedian(const uint8_t apvIndex) const; //third byte of channel data for normal FED channels - uint8_t packetCode() const; + constexpr uint8_t packetCode() const; private: friend class FEDBuffer; const uint8_t* data_; uint32_t offset_; uint16_t length_; + uint8_t headerLen_ = 0; }; //base class for sistrip FED buffers which have a DAQ header/trailer and tracker special header @@ -1557,18 +1561,38 @@ namespace sistrip { //FEDChannel - inline FEDChannel::FEDChannel(const uint8_t* const data, const uint32_t offset) : data_(data), offset_(offset) { + constexpr FEDChannel::FEDChannel(const uint8_t* const data, const uint32_t offset, const ZSROMode zsRoMode) + : data_(data), offset_(offset), headerLen_(zsRoMode == ZSROMode::nonLite ? 7 : 2) { + length_ = (data_[(offset_) ^ 7] + (data_[(offset_ + 1) ^ 7] << 8)); + } + + constexpr FEDChannel::FEDChannel(const uint8_t* const data, const uint32_t offset) : data_(data), offset_(offset) { length_ = (data_[(offset_) ^ 7] + (data_[(offset_ + 1) ^ 7] << 8)); } - inline FEDChannel::FEDChannel(const uint8_t* const data, const uint32_t offset, const uint16_t length) + constexpr FEDChannel::FEDChannel(const uint8_t* const data, const uint32_t offset, const uint16_t length) : data_(data), offset_(offset), length_(length) {} - inline uint16_t FEDChannel::length() const { return length_; } + constexpr uint16_t FEDChannel::stripsInCh(uint8_t num_bits) const { + const bool emptyCh = (headerLen_ + 2) >= (length_); + const uint16_t start = offset_ + headerLen_; + const uint16_t end = offset_ + length_; + uint16_t stripN = 0; + if (!emptyCh) { + for (uint16_t nStrip_wOfs = start + 1; nStrip_wOfs < end;) { + const uint8_t clustStripN = data_[(nStrip_wOfs) ^ 7]; + nStrip_wOfs += ((uint32_t)clustStripN) * num_bits / 8 + 2; + stripN += clustStripN; + } + } + return stripN; + } + + constexpr uint16_t FEDChannel::length() const { return length_; } - inline uint8_t FEDChannel::packetCode() const { return data_[(offset_ + 2) ^ 7]; } + constexpr uint8_t FEDChannel::packetCode() const { return data_[(offset_ + 2) ^ 7]; } - inline uint16_t FEDChannel::cmMedian(const uint8_t apvIndex) const { + constexpr uint16_t FEDChannel::cmMedian(const uint8_t apvIndex) const { uint16_t result = 0; //CM median is 10 bits with lowest order byte first. First APV CM median starts in 4th byte of channel data result |= data_[(offset_ + 3 + 2 * apvIndex) ^ 7]; @@ -1576,9 +1600,9 @@ namespace sistrip { return result; } - inline const uint8_t* FEDChannel::data() const { return data_; } + constexpr const uint8_t* FEDChannel::data() const { return data_; } - inline uint32_t FEDChannel::offset() const { return offset_; } + constexpr uint32_t FEDChannel::offset() const { return offset_; } } // namespace sistrip #endif //ndef EventFilter_SiStripRawToDigi_FEDBufferComponents_H diff --git a/RecoLocalTracker/SiStripClusterizer/BuildFile.xml b/RecoLocalTracker/SiStripClusterizer/BuildFile.xml index e7fdcd6cbf34a..ba5eb6f7adab5 100644 --- a/RecoLocalTracker/SiStripClusterizer/BuildFile.xml +++ b/RecoLocalTracker/SiStripClusterizer/BuildFile.xml @@ -9,6 +9,12 @@ + + + + + + diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h new file mode 100644 index 0000000000000..e6d2ffd1745b3 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h @@ -0,0 +1,13 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsHostObject_h +#define RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsHostObject_h + +#include "DataFormats/Portable/interface/PortableHostObject.h" +#include "DataFormats/Portable/interface/PortableObject.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h" + +namespace sistrip { + using SiStripClusterizerConditionsDetToFedsHostObject = PortableHostObject; + using SiStripClusterizerConditionsGainNoiseCalsHostObject = PortableHostObject; +} // namespace sistrip + +#endif diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h new file mode 100644 index 0000000000000..cef7a537a7efa --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h @@ -0,0 +1,19 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsRecord_h +#define RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsRecord_h + +#include "CalibTracker/Records/interface/SiStripDependentRecords.h" +#include "FWCore/Framework/interface/DependentRecordImplementation.h" +#include "FWCore/Utilities/interface/mplVector.h" + +namespace sistrip { + class SiStripClusterizerConditionsDetToFedsRecord + : public edm::eventsetup::DependentRecordImplementation> {}; + + class SiStripClusterizerConditionsGainNoiseCalsRecord + : public edm::eventsetup::DependentRecordImplementation< + SiStripClusterizerConditionsGainNoiseCalsRecord, + edm::mpl::Vector> {}; +} // namespace sistrip + +#endif // RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsRecord_h diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h new file mode 100644 index 0000000000000..34176db643c72 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h @@ -0,0 +1,22 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsStruct_h +#define RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsStruct_h + +#include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" + +namespace sistrip { + // Quality, detID, invThick and iPair are indexed by the `channelIndex` map + // Gain and noise are indexed by the stripIndex and apvIndex, respectively + struct DetToFeds { + std::array qualityOk; + }; + + struct GainNoiseCals { + std::array detID; + std::array invthick; + std::array gain; + std::array iPair; + std::array noise; + }; +} // namespace sistrip + +#endif // RecoLocalTracker_SiStripClusterizer_interface_SiStripClusterizerConditionsStruct_h diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingHost.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingHost.h new file mode 100644 index 0000000000000..9e0c493ce17c1 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingHost.h @@ -0,0 +1,11 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingHost_h +#define RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingHost_h + +#include "DataFormats/Portable/interface/PortableHostCollection.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingSoA.h" + +namespace sistrip { + using SiStripMappingHost = PortableHostCollection; +} + +#endif // RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingHost_h diff --git a/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingSoA.h b/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingSoA.h new file mode 100644 index 0000000000000..d887ad3838b5a --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingSoA.h @@ -0,0 +1,29 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingSoA_h +#define RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingSoA_h + +#include "DataFormats/SoATemplate/interface/SoACommon.h" +#include "DataFormats/SoATemplate/interface/SoALayout.h" + +namespace sistrip { + GENERATE_SOA_LAYOUT(SiStripMappingSoALayout, + // Detector ID, FED ID and FED channel for indexing + SOA_COLUMN(uint32_t, detID), + SOA_COLUMN(uint16_t, fedID), + SOA_COLUMN(uint8_t, fedCh), + // FEDChannel.offset() + SOA_COLUMN(uint16_t, fedChOfs), + // Offset of the FEDChannel.data() in the flatten raw buffer + SOA_COLUMN(uint32_t, fedChDataOfsBuf), + // Number of strips in the fedCh + SOA_COLUMN(uint32_t, fedChStripsN), + // readout mode of the buffer the FED channels are taken + SOA_COLUMN(uint8_t, readoutMode), + // packet code of the buffer the FED channels are taken + SOA_COLUMN(uint8_t, packetCode)) + + using SiStripMappingSoA = SiStripMappingSoALayout<>; + using SiStripMappingView = SiStripMappingSoA::View; + using SiStripMappingConstView = SiStripMappingSoA::ConstView; +} // namespace sistrip + +#endif // RecoLocalTracker_SiStripClusterizer_interface_SiStripMappingSoA_h diff --git a/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h b/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h new file mode 100644 index 0000000000000..d307d2a243e83 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h @@ -0,0 +1,27 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_alpaka_SiStripClusterizerConditionsDeviceObject_h +#define RecoLocalTracker_SiStripClusterizer_interface_alpaka_SiStripClusterizerConditionsDeviceObject_h + +#include "DataFormats/Portable/interface/alpaka/PortableObject.h" + +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" + +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE { + namespace sistrip { + // make the names from the top-level sistrip namespace visible for unqualified lookup + using namespace ::sistrip; + using SiStripClusterizerConditionsDetToFedsDeviceObject = PortableObject; + using SiStripClusterizerConditionsGainNoiseCalsDeviceObject = PortableObject; + } // namespace sistrip + +} // namespace ALPAKA_ACCELERATOR_NAMESPACE + +// check that the portable device collection for the host device is the same as the portable host collection +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::SiStripClusterizerConditionsDetToFedsDeviceObject, + sistrip::SiStripClusterizerConditionsDetToFedsHostObject); +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::SiStripClusterizerConditionsGainNoiseCalsDeviceObject, + sistrip::SiStripClusterizerConditionsGainNoiseCalsHostObject); + +#endif diff --git a/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripMappingDevice.h b/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripMappingDevice.h new file mode 100644 index 0000000000000..b5c6e30871254 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripMappingDevice.h @@ -0,0 +1,17 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_interface_alpaka_SiStripMappingDevice_h +#define RecoLocalTracker_SiStripClusterizer_interface_alpaka_SiStripMappingDevice_h + +#include "DataFormats/Portable/interface/alpaka/PortableCollection.h" +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingHost.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripMappingSoA.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + using namespace ::sistrip; + using SiStripMappingDevice = PortableCollection; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +// check that the sistrip device collection for the host device is the same as the sistrip host collection +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::SiStripMappingDevice, sistrip::SiStripMappingHost); + +#endif // RecoLocalTracker_SiStripClusterizer_interface_alpaka_SiStripMappingDevice_h diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/BuildFile.xml b/RecoLocalTracker/SiStripClusterizer/plugins/BuildFile.xml index 8933e253d7be6..993015c581591 100644 --- a/RecoLocalTracker/SiStripClusterizer/plugins/BuildFile.xml +++ b/RecoLocalTracker/SiStripClusterizer/plugins/BuildFile.xml @@ -11,3 +11,24 @@ + + + + + + + + + + + + + + + + + + + + + diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClustersToLegacy.cc b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClustersToLegacy.cc new file mode 100644 index 0000000000000..610b7466cea26 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/plugins/SiStripClustersToLegacy.cc @@ -0,0 +1,151 @@ +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" + +#include "FWCore/Framework/interface/global/EDProducer.h" +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/Utilities/interface/StreamID.h" +#include "FWCore/Utilities/interface/InputTag.h" + +#include "DataFormats/Common/interface/DetSetVectorNew.h" +#include "DataFormats/SiStripCluster/interface/SiStripCluster.h" +#include "DataFormats/SiStripClusterSoA/interface/SiStripClusterHost.h" +#include "DataFormats/SiStripDigiSoA/interface/SiStripDigiHost.h" + +// using namespace ::sistrip; +namespace sistrip { + class SiStripClustersToLegacy : public edm::global::EDProducer<> { + public: + SiStripClustersToLegacy(edm::ParameterSet const& iConfig) + : siStripClustersToken_(consumes(iConfig.getParameter("source"))), + siStripDigiToken_(consumes(iConfig.getParameter("source"))), + siStripClustersSetVecPutToken_(produces()) {} + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("source"); + descriptions.addWithDefaultLabel(desc); + } + + void produce(edm::StreamID, edm::Event& iEvent, edm::EventSetup const&) const override { + using out_t = edmNew::DetSetVector; + auto output = std::make_unique(); + + // Get clusters and amplitudes + auto const& clusters_onHost = iEvent.get(siStripClustersToken_); + auto const& amplitudes_onHost = iEvent.get(siStripDigiToken_); + + const uint32_t nClusterCandidates = clusters_onHost->metadata().size(); + const auto& clusterSizeArr = clusters_onHost->clusterSize(); + const auto& detIdArr = clusters_onHost->clusterDetId(); + const auto& firstStripArr = clusters_onHost->firstStrip(); + const auto& candidateAcceptedArr = clusters_onHost->candidateAccepted(); + const auto& clusterIndexArr = clusters_onHost->clusterIndex(); + const auto& barycenterArr = clusters_onHost->barycenter(); + const auto& chargeArr = clusters_onHost->charge(); + + const auto& clusterAmplArr = amplitudes_onHost->adc(); + + // Return immediately if there are no clusters + if (nClusterCandidates == 0) { + iEvent.put(siStripClustersSetVecPutToken_, std::move(output)); + return; + } + + // Educated guess for the total number of detector IDs, + // based on Run: 386593 Event: 536278171 with 13883 detectors. + // const uint32_t nModulesWithClustersGuess = 15000; + // The number of clusters from x->nClusterCandidates() is an upper limit, + // the flag trueCluster then mask the real clusters. + // From Run: 386593 Event: 536278171 there are nClusterCandidates=112735 with + // 99863 real clusters (so 112735-99863 = 12872 clusters are masked out ) + const uint32_t clusterCandidatesNb = clusters_onHost->nClusterCandidates(); + const uint32_t goodClustersNb = clusters_onHost->candidateAcceptedPrefix(nClusterCandidates - 1); + // output->reserve(nModulesWithClustersGuess, goodClustersNb); + + uint32_t clusterN = 0; + for (uint32_t i = 0; i < clusterCandidatesNb && (clusterN < goodClustersNb);) { + const auto detid = detIdArr[i]; + out_t::FastFiller record(*output, detid); + + while (i < clusterCandidatesNb && detIdArr[i] == detid) { + if (candidateAcceptedArr[i]) { + const auto size = clusterSizeArr[i]; + const auto firstStrip = firstStripArr[i]; + + const float barycenter = barycenterArr[i]; + const float charge = chargeArr[i]; + + const auto index = clusterIndexArr[i]; + auto clusterAdc = clusterAmplArr.subspan(index, size); + std::vector adcs(clusterAdc.begin(), clusterAdc.end()); + + record.push_back(SiStripCluster(firstStrip, std::move(adcs), barycenter, charge)); + clusterN++; + } + i++; + } + } + + if (edm::isDebugEnabled()) { + dumpClusters(output.get()); + } + iEvent.put(siStripClustersSetVecPutToken_, std::move(output)); + } + + void dumpClusters(edmNew::DetSetVector* detSetClusters) const; + + private: + const edm::EDGetTokenT siStripClustersToken_; + const edm::EDGetTokenT siStripDigiToken_; + const edm::EDPutTokenT> siStripClustersSetVecPutToken_; + }; + + void SiStripClustersToLegacy::dumpClusters(edmNew::DetSetVector* detSetClusters) const { + int clustersAlloc = detSetClusters->dataSize(); + + std::ostringstream dumpMsg; + dumpMsg << "#clDump,Produced:" << clustersAlloc << "\n"; + dumpMsg << "i,cIdx,cSz,cDetId,chg,1st,tCl,bary,|clusterADCs|\n"; + + int i = 0; + int cIdx = 0; + const int trueCluster = 1; + for (auto it = detSetClusters->begin(); it != detSetClusters->end(); ++cIdx, it++) { + if (true || i < 100 || i > (clustersAlloc - 100)) { + auto detSet = *it; + auto cDetId = detSet.detId(); + // int clNum = detSet.size(); + for (auto j = detSet.begin(); j != detSet.end(); ++j, ++i) { + auto cluster = *j; + // + auto cSz = cluster.size(); + auto chg = cluster.charge(); + auto firstStrip = cluster.firstStrip(); + + auto bary = cluster.barycenter(); + + dumpMsg << i << "," << cIdx << "," << cSz << "," << cDetId << "," << chg << "," << firstStrip << "," + << trueCluster << "," << bary << ",|"; + for (int k = 0; k < (int)cluster.amplitudes().size(); ++k) { + auto adc = cluster.amplitudes()[k]; + dumpMsg << (int)(adc); + if (k != ((int)cluster.amplitudes().size() - 1)) { + dumpMsg << "/"; + } + } + dumpMsg << "|\n"; + } + } + } + dumpMsg << "#zClDump\n"; + std::cout << dumpMsg.str(); + } + +} // namespace sistrip + +#include "FWCore/Framework/interface/MakerMacros.h" +DEFINE_FWK_MODULE(sistrip::SiStripClustersToLegacy); diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripClusterizerConditionsESProducerAlpaka.cc b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripClusterizerConditionsESProducerAlpaka.cc new file mode 100644 index 0000000000000..c6756c483a56a --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripClusterizerConditionsESProducerAlpaka.cc @@ -0,0 +1,213 @@ +#include "CalibFormats/SiStripObjects/interface/SiStripGain.h" +#include "CalibFormats/SiStripObjects/interface/SiStripQuality.h" +#include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h" +#include "CalibTracker/Records/interface/SiStripDependentRecords.h" + +#include "CondFormats/SiStripObjects/interface/SiStripNoises.h" + +#include "DataFormats/SiStripCommon/interface/ConstantsForHardwareSystems.h" +#include "DataFormats/SiStripCluster/interface/SiStripClusterTools.h" + +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/Utilities/interface/ESGetToken.h" + +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/ESProducer.h" +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/ModuleFactory.h" +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" +#include "HeterogeneousCore/AlpakaInterface/interface/host.h" + +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + // Unqualified lookup of the top level namespace + using namespace ::sistrip; + + class DetToFed { + public: + DetToFed(uint32_t detid, uint16_t fedid, uint16_t fedch, uint16_t ipair) + : detID_(detid), fedID_(fedid), fedCh_(fedch), iPair_(ipair) {} + + inline uint32_t detID() const { return detID_; } + inline uint16_t fedID() const { return fedID_; } + inline uint16_t fedCh() const { return fedCh_; } + inline uint16_t iPair() const { return iPair_; } + + private: + uint32_t detID_; + uint16_t fedID_; + uint16_t fedCh_; + uint16_t iPair_; + }; + + class SiStripClusterizerConditionsESProducerAlpaka : public ESProducer { + public: + SiStripClusterizerConditionsESProducerAlpaka(edm::ParameterSet const& iConfig) : ESProducer(iConfig) { + // Two tokens for quality because of the exc message: + // You may need multiple tokens if you want to get the same data in multiple transitions. + auto listenA = setWhatProduced(this, &SiStripClusterizerConditionsESProducerAlpaka::produceDetToFeds); + qualityTokenA_ = listenA.consumes(iConfig.getParameter("QualityLabel")); + + auto listenB = setWhatProduced(this, &SiStripClusterizerConditionsESProducerAlpaka::produceData); + qualityTokenB_ = listenB.consumes(iConfig.getParameter("QualityLabel")); + gainsToken_ = listenB.consumes(); + noisesToken_ = listenB.consumes(); + } + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + edm::ParameterSetDescription desc; + desc.add("QualityLabel"); + desc.add("Label"); + descriptions.addWithDefaultLabel(desc); + } + + std::unique_ptr produceDetToFeds( + SiStripClusterizerConditionsDetToFedsRecord const& iRecord) { + const auto& quality = iRecord.get(qualityTokenA_); + + auto product = std::make_unique(cms::alpakatools::host()); + product->zeroInitialise(); + + auto& detToFeds_qualityFlags = (*product).data()->qualityOk; + + // connected: map> + // map of KEY=detid DATA=vector of apvs, maximum 6 APVs per detector module : + const auto& connected = quality.cabling()->connected(); + // detCabling: map + // map of KEY=detid DATA=vector + const auto& detCabling = quality.cabling()->getDetCabling(); + + for (const auto& conn : connected) { + const auto det = conn.first; + if (!quality.IsModuleBad(det)) { + const auto detConn_it = detCabling.find(det); + + if (detCabling.end() != detConn_it) { + for (const auto& conn : (*detConn_it).second) { + if (conn && conn->fedId() && conn->isConnected()) { + const auto conn_fedID = conn->fedId(); + const auto conn_fedCh = conn->fedCh(); + + uint32_t idx = channelIndex(conn_fedID, conn_fedCh); + detToFeds_qualityFlags[idx] = true; + // ss << idx << ","; + } + } + } + } + } + + return product; + } + + std::unique_ptr produceData( + SiStripClusterizerConditionsGainNoiseCalsRecord const& iRecord) { + auto gains = iRecord.getTransientHandle(gainsToken_); + auto noises = iRecord.getTransientHandle(noisesToken_); + auto quality = iRecord.getTransientHandle(qualityTokenB_); + + // Prepare the conditions on the host + auto product = std::make_unique(cms::alpakatools::host()); + // Fill the collections + fillSiStripClusterizerConditions(quality.product(), gains.product(), noises.product(), *product->data()); + // + return product; + } + + // Auxiliary functions to translate indexes on the arrays + static constexpr uint16_t badBit = 1 << 15; + inline uint16_t fedIndex(uint16_t fed) { return (fed - FED_ID_MIN); }; + inline uint32_t stripIndex(uint16_t fed, uint16_t channel, uint16_t strip) { + return (fedIndex(fed) * FEDCH_PER_FED * STRIPS_PER_FEDCH + channel * STRIPS_PER_FEDCH + + (strip % STRIPS_PER_FEDCH)); + }; + inline uint32_t apvIndex(uint16_t fed, uint16_t channel, uint16_t strip) { + return (fedIndex(fed) * APVS_PER_FEDCH * FEDCH_PER_FED + APVS_PER_CHAN * channel + + (strip % STRIPS_PER_FEDCH) / STRIPS_PER_APV); + }; + inline uint32_t channelIndex(uint16_t fedId, uint16_t fedCh) { return (fedIndex(fedId) * FEDCH_PER_FED + fedCh); }; + + private: + edm::ESGetToken gainsToken_; + edm::ESGetToken noisesToken_; + // Two tokens for quality because of the exc message: + // You may need multiple tokens if you want to get the same data in multiple transitions. + edm::ESGetToken qualityTokenA_; + edm::ESGetToken qualityTokenB_; + + // Make conditions as in the RecoLocalTracker/SiStripClusterizer/plugins/ClustersFromRawProducer.cc module + void fillSiStripClusterizerConditions(const SiStripQuality* quality, + const SiStripGain* gains, + const SiStripNoises* noises, + GainNoiseCals& calibs); + }; + + void SiStripClusterizerConditionsESProducerAlpaka::fillSiStripClusterizerConditions(const SiStripQuality* quality, + const SiStripGain* gains, + const SiStripNoises* noises, + GainNoiseCals& calibs) { + // Alias the members + auto& invthick = calibs.invthick; + auto& detID = calibs.detID; + auto& iPair = calibs.iPair; + auto& noise = calibs.noise; + auto& gain = calibs.gain; + + // connected: map> + // map of KEY=detid DATA=vector of apvs, maximum 6 APVs per detector module : + const auto& connected = quality->cabling()->connected(); + // detCabling: map + // map of KEY=detid DATA=vector + const auto& detCabling = quality->cabling()->getDetCabling(); + + for (const auto& conn : connected) { + const auto det = conn.first; + + if (!quality->IsModuleBad(det)) { + const auto detConn_it = detCabling.find(det); + + const auto gainRange = gains->getRange(det); + if (detCabling.end() != detConn_it) { + for (const auto& chan : (*detConn_it).second) { + if (chan && chan->fedId() && chan->isConnected()) { + const uint32_t chan_detID = chan->detId(); + const uint16_t chan_fedID = chan->fedId(); + const uint16_t chan_fedCh = chan->fedCh(); + const uint16_t chan_apvPairNumber = chan->apvPairNumber(); + + // Fill the data structures + // Note: the channelIndex is used to access the data structures + detID[channelIndex(chan_fedID, chan_fedCh)] = chan_detID; + iPair[channelIndex(chan_fedID, chan_fedCh)] = chan_apvPairNumber; + invthick[channelIndex(chan_fedID, chan_fedCh)] = siStripClusterTools::sensorThicknessInverse(chan_detID); + + auto offset = STRIPS_PER_FEDCH * chan_apvPairNumber; + for (uint16_t strip = 0; strip < STRIPS_PER_FEDCH; ++strip) { + const auto detstrip = strip + offset; + const uint16_t strip_noise = SiStripNoises::getRawNoise(detstrip, noises->getRange(det)); + const auto bad = quality->IsStripBad(quality->getRange(det), detstrip); + + // setStrip_(chan_fedID, chan_fedCh, detstrip, noise, gain, bad); + if (bad) [[unlikely]] { + noise[stripIndex(chan_fedID, chan_fedCh, detstrip)] = badBit; + } else { + noise[stripIndex(chan_fedID, chan_fedCh, detstrip)] = strip_noise; + } + } + + // gain is actually stored per-APV, not per-strip (so stored for the strpis 0-127 and 128-255) + gain[apvIndex(chan_fedID, chan_fedCh, 0)] = SiStripGain::getApvGain(0, gainRange); + gain[apvIndex(chan_fedID, chan_fedCh, 255)] = SiStripGain::getApvGain(1, gainRange); + } + } + } + } + } + } +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/ModuleFactory.h" +DEFINE_FWK_EVENTSETUP_ALPAKA_MODULE(sistrip::SiStripClusterizerConditionsESProducerAlpaka); diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToCluster.cc b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToCluster.cc new file mode 100644 index 0000000000000..8d113583122ed --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToCluster.cc @@ -0,0 +1,298 @@ +#include "FWCore/Framework/interface/ESWatcher.h" +#include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/ParameterSet/interface/ParameterSetDescription.h" +#include "FWCore/Utilities/interface/InputTag.h" + +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/stream/SynchronizingEDProducer.h" +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/EDPutToken.h" +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/ESGetToken.h" +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" + +#include "CalibFormats/SiStripObjects/interface/SiStripClusterizerConditions.h" + +#include "DataFormats/FEDRawData/interface/FEDRawDataCollection.h" +#include "DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h" +#include "DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h" + +#include "EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h" +#include "EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h" + +#include "RecoLocalTracker/SiStripClusterizer/interface/StripClusterizerAlgorithmFactory.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h" +#include "RecoLocalTracker/Records/interface/SiStripClusterizerConditionsRcd.h" + +#include "SiStripRawToClusterAlgo.h" + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + using namespace ::sistrip; + + class SiStripRawToCluster : public stream::SynchronizingEDProducer<> { + public: + SiStripRawToCluster(edm::ParameterSet const& iConfig); + ~SiStripRawToCluster() override = default; + + static void fillDescriptions(edm::ConfigurationDescriptions& descriptions); + + private: + void acquire(device::Event const& iEvent, device::EventSetup const& iSetup) override; + void produce(device::Event& iEvent, device::EventSetup const& iSetup) override; + + std::unique_ptr fillFedIdFedChBuffer(Queue& queue, + const SiStripClusterizerConditions& stripCond, + const FEDRawDataCollection& rawColl); + + // RAW unpacking and clustering algorithm + SiStripRawToClusterAlgo algo_; + + edm::EDGetTokenT fedRawGetToken_; + edm::ESGetToken stripCondGetToken_; + device::ESGetToken + stripCablCondGetToken_; + device::ESGetToken + stripDataCondGetToken_; + device::EDPutToken stripClustPutToken_; + device::EDPutToken stripDigiPutToken_; + + // Setup + bool doAPVEmulatorCheck_; + + // Helper functions to fill valid, condition-passing raw/buffers + std::unique_ptr fillBuffer(int fedId, const FEDRawData& rawData); + }; + + SiStripRawToCluster::SiStripRawToCluster(const edm::ParameterSet& iConfig) + : SynchronizingEDProducer(iConfig), + algo_(iConfig.getParameter("Clusterizer")), + fedRawGetToken_(consumes(iConfig.getParameter("ProductLabel"))), + stripCondGetToken_(esConsumes(edm::ESInputTag{"", iConfig.getParameter("ConditionsLabel")})), + stripCablCondGetToken_( + esConsumes(edm::ESInputTag{"", iConfig.getParameter("CablingConditionsLabel")})), + stripDataCondGetToken_(esConsumes()), + stripClustPutToken_(produces()), + stripDigiPutToken_(produces()), + doAPVEmulatorCheck_(iConfig.getParameter("DoAPVEmulatorCheck")) {} + + void SiStripRawToCluster::fillDescriptions(edm::ConfigurationDescriptions& descriptions) { + // Add some custom parameter description to the automatically-created ones by the addWithDefaultLabel methos + edm::ParameterSetDescription desc; + desc.add("ProductLabel", edm::InputTag("rawDataCollector")); + desc.add("ConditionsLabel", ""); + desc.add("CablingConditionsLabel", ""); + + // Setup parameters (from SiStripClusterizerFromRaw) + desc.add("DoAPVEmulatorCheck", false); + + // Inherit all the parameters from the clusterizers (all var.s from the Clusterizer PSet) + // 1:Algorithm, ConditionsLabel, ChannelThreshold, SeedThreshold, ClusterThreshold, + // MaxSequentialHoles, MaxSequentialBad, MaxAdjacentBad, MaxClusterSize, RemoveApvShots, + // setDetId, 12:clusterChargeCut + edm::ParameterSetDescription clusterizer; + StripClusterizerAlgorithmFactory::fillDescriptions(clusterizer); + // The parameter MaxSeedStrips determines the size for pre-allocation of the clusters collection SoA + clusterizer.add("MaxSeedStrips", 200000u); + desc.add("Clusterizer", clusterizer); + + descriptions.addWithDefaultLabel(desc); + } + + void SiStripRawToCluster::acquire(device::Event const& iEvent, device::EventSetup const& iSetup) { + // If this is the first time the module is called / the record signalled a change in conditions, + // then load/update the conditions + cabling map + + // Host conditions (for buffer preconstruct + fedCh mapping preparation) + const auto& stripCond = iSetup.getData(stripCondGetToken_); + + // Device conditions (cabling map) + const auto& stripCablCond = iSetup.getData(stripCablCondGetToken_); + LogDebug("SiStripRawToCluster") << "#sizeB,stripCablCond_," << alpaka::getExtentProduct(stripCablCond.buffer()); + + // Get FED raw data collection + const auto& rawCollection = iEvent.get(fedRawGetToken_); + + // Prepare the LUT for FEDCh mapping in the raw[] on the device + std::unique_ptr fedChMover = fillFedIdFedChBuffer(iEvent.queue(), stripCond, rawCollection); + + // Move PortableFEDMover class to algo + algo_.prepareUnpackCluster(iEvent.queue(), stripCablCond.const_data(), std::move(fedChMover)); + } + + void SiStripRawToCluster::produce(device::Event& iEvent, device::EventSetup const& iSetup) { + // Device conditions (strip noise) + const auto& stripDataCond = iSetup.getData(stripDataCondGetToken_); + LogDebug("SiStripRawToCluster") << "#sizeB,stripDataCond_," << alpaka::getExtentProduct(stripDataCond.buffer()); + + // Unpack the raw FED data into strip digi + auto nStrips = algo_.unpackStrips(iEvent.queue(), stripDataCond.const_data()); + if (nStrips == 0) { + // No strips to unpack, empty cluster collection + iEvent.emplace(stripClustPutToken_, 0, iEvent.queue()); + iEvent.emplace(stripDigiPutToken_, 0, iEvent.queue()); + } else { + // Run the clusterization algorithm (ThreeThresholdAlgorithm) + auto cluster_d = algo_.makeClusters(iEvent.queue(), stripDataCond.const_data()); + + // Get the clusters amplitudes + auto clusterAmpls_d = algo_.releaseDigiAmplitudes(); + + iEvent.put(stripClustPutToken_, std::move(cluster_d)); + iEvent.put(stripDigiPutToken_, std::move(clusterAmpls_d)); + } + } + + std::unique_ptr SiStripRawToCluster::fillFedIdFedChBuffer( + Queue& queue, const SiStripClusterizerConditions& stripCond, const FEDRawDataCollection& rawColl) { + // Containers for the condition-passing raw data + const auto fedID_maxIdx = sistrip::FED_ID_MAX + 1; + std::vector raw(fedID_maxIdx, nullptr); + std::vector> buffers(fedID_maxIdx); + + std::vector fedChOfs_wrt_rawFedId_; + uint32_t rawBuffFlattenSize = 0; + + // loop over good det in cabling + for (auto idet : stripCond.allDetIds()) { + // it populates raw, buffers with only connected fed + auto const& det = stripCond.findDetId(idet); + if (!det.valid()) { + continue; // idet loop + } + + // Loop over apv-pairs of det + // int connIdx = 0; + for (auto const conn : stripCond.currentConnection(det)) { + // connIdx++; + if (!conn) { + continue; // conn loop + } + + const uint16_t fedId = conn->fedId(); + + // If fed id is null or connection is invalid continue + if (!fedId || !conn->isConnected()) { + continue; // conn loop + } + + if (!raw[fedId]) { + // pointer to the raw data in the collection + raw[fedId] = &rawColl.FEDData(fedId); + rawBuffFlattenSize += raw[fedId]->size(); + } + + // If Fed hasnt already been initialised, extract data and initialise + sistrip::FEDBuffer* buffer = buffers[fedId].get(); + if (!buffer) { + buffers[fedId] = fillBuffer(fedId, *raw[fedId]); + if (!buffers[fedId]) { + continue; + } else { + buffer = buffers[fedId].get(); + } + } + assert(buffer); + + // Set legacy unpacking + // buffers[fedidx]->setLegacyMode(legacyUnpacker_); + + // Check the readout mode of the buffer + const FEDReadoutMode buffROMode = buffer->readoutMode(); + // Make sure EACH buffer has a readout mode supported by the module + if (!(buffROMode >= READOUT_MODE_ZERO_SUPPRESSED_LITE10 && + buffROMode <= READOUT_MODE_ZERO_SUPPRESSED_LITE8_BOTBOT_CMOVERRIDE && + buffROMode != READOUT_MODE_PROC_RAW)) { + if (edm::isDebugEnabled()) { + std::ostringstream ss; + ss << "Unsupported buffer ROmode=" << buffROMode << " fedID= " << fedId; + edm::LogWarning("fillFedIdFedChBuffer") << ss.str(); + } + continue; // conn loop + } + + // check channel + const uint8_t fedCh = conn->fedCh(); + if (!buffer->channelGood(fedCh, doAPVEmulatorCheck_)) { + if (edm::isDebugEnabled()) { + std::ostringstream ss; + ss << "Problem unpacking channel " << (int)fedCh << " on FED " << fedId; + edm::LogWarning("fillFedIdFedChBuffer") << ss.str(); + } + continue; // conn loop + } + + // Good fedId and fedCh in the raw data + const FEDChannel& channel = buffer->channel(fedCh); + const uint32_t fedChOfs = channel.offset(); + + auto diff = channel.data() - raw[fedId]->data(); + if (diff < 0 || diff > std::numeric_limits::max()) { + if (edm::isDebugEnabled()) { + std::ostringstream ss; + ss << "Large diff " << diff << " for fedId " << fedId << " fedCh " << fedCh; + edm::LogWarning("fillFedIdFedChBuffer") << ss.str(); + } + continue; + } + const uint32_t fedChOfs_wrt_rawFedId = static_cast(diff); + + fedChOfs_wrt_rawFedId_.push_back( + FEDChMetadata(idet, fedId, fedCh, fedChOfs, fedChOfs_wrt_rawFedId, buffROMode)); + } + } + + // Do I have any channel to unpack ( FEDRaw data empty || mapping empty) ? + if (fedChOfs_wrt_rawFedId_.empty()) { + return nullptr; + } + + // Create container for the buffer and mapping + auto fedMover = std::make_unique(queue, rawBuffFlattenSize, fedChOfs_wrt_rawFedId_.size()); + + // Fill buffer + fedMover->fillBuffer(raw); + + // Fill the mapping + fedMover->fillMapping(fedChOfs_wrt_rawFedId_); + + return fedMover; + } + + std::unique_ptr SiStripRawToCluster::fillBuffer(int fedId, const FEDRawData& rawData) { + std::unique_ptr buffer; + + // Check on FEDRawData pointer + const auto st_buffer = sistrip::preconstructCheckFEDBuffer(rawData); + if UNLIKELY (sistrip::FEDBufferStatusCode::SUCCESS != st_buffer) { + if (edm::isDebugEnabled()) { + edm::LogWarning(sistrip::mlRawToCluster_) + << "[ClustersFromRawProducer::" << __func__ << "]" << st_buffer << " for FED ID " << fedId; + } + return buffer; + } + buffer = std::make_unique(rawData); + const auto st_chan = buffer->findChannels(); + if UNLIKELY (sistrip::FEDBufferStatusCode::SUCCESS != st_chan) { + if (edm::isDebugEnabled()) { + edm::LogWarning(sistrip::mlRawToCluster_) + << "Exception caught when creating FEDBuffer object for FED " << fedId << ": " << st_chan; + } + buffer.reset(); + return buffer; + } + if UNLIKELY (!buffer->doChecks(false)) { + if (edm::isDebugEnabled()) { + edm::LogWarning(sistrip::mlRawToCluster_) + << "Exception caught when creating FEDBuffer object for FED " << fedId << ": FED Buffer check fails"; + } + buffer.reset(); + return buffer; + } + + return buffer; + } + +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/MakerMacros.h" +DEFINE_FWK_ALPAKA_MODULE(sistrip::SiStripRawToCluster); diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.dev.cc b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.dev.cc new file mode 100644 index 0000000000000..11e77bcf1edec --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.dev.cc @@ -0,0 +1,1020 @@ +#include + +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" +#include "HeterogeneousCore/AlpakaInterface/interface/warpsize.h" +#include "HeterogeneousCore/AlpakaInterface/interface/workdivision.h" +#include "HeterogeneousCore/AlpakaInterface/interface/prefixScan.h" +#include "HeterogeneousCore/AlpakaInterface/interface/moveToDeviceAsync.h" + +#include "RecoLocalTracker/SiStripClusterizer/interface/ClusterChargeCut.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsStruct.h" + +#include "EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h" +#include "EventFilter/SiStripRawToDigi/interface/SiStripFEDBufferComponents.h" + +#include "SiStripRawToClusterAlgo.h" + +// Generic raw unpackers + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip::fedchannelunpacker { + enum class StatusCode { SUCCESS = 0, BAD_CHANNEL_LENGTH, UNORDERED_DATA, BAD_PACKET_CODE, ZERO_PACKET_CODE }; + namespace detail { + // (adapted from https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h#L263 + template + ALPAKA_FN_HOST_ACC StatusCode unpackZSW(uint32_t chan, + const uint8_t* channel_data, + uint16_t channel_length, + uint32_t channel_offset, + SiStripDigiView out, + uint32_t* aoffIdx, + uint8_t headerLength, + uint16_t stripStart, + uint8_t bits_shift = 0) { + constexpr auto num_words = num_bits / 8; + static_assert(((num_bits % 8) == 0) && (num_words > 0) && (num_words < 3)); + if (channel_length & 0xF000) { + return StatusCode::BAD_CHANNEL_LENGTH; + } + const uint8_t* const data = channel_data; + uint_fast16_t offset = channel_offset + headerLength; // header is 2 (lite) or 7 + uint_fast8_t firstStrip{0}, nInCluster{0}, inCluster{0}; + const uint_fast16_t end = channel_offset + channel_length; + while (offset != end) { + if (inCluster == nInCluster) { + if (offset + 2 >= end) { + // offset should already be at end then (empty cluster) + // printf("offset should already set"); + break; + } + const uint_fast8_t newFirstStrip = data[(offset++) ^ 7]; + if (newFirstStrip < (firstStrip + inCluster)) { + // LogDebug("FEDBuffer") << "First strip of new cluster is not greater than last strip of previous cluster. " + // << "Last strip of previous cluster is " << uint16_t(firstStrip + inCluster) << ". " + // << "First strip of new cluster is " << uint16_t(newFirstStrip) << "."; + return StatusCode::UNORDERED_DATA; + } + firstStrip = newFirstStrip; + nInCluster = data[(offset++) ^ 7]; + inCluster = 0; + } + out.channel(*aoffIdx) = chan; + out.stripId(*aoffIdx) = stripStart + firstStrip + inCluster; + out.adc(*aoffIdx) = ::sistrip::fedchannelunpacker::detail::getADC_W(data, offset, bits_shift); + (*aoffIdx)++; + offset += num_words; + ++inCluster; + } + return StatusCode::SUCCESS; + } + + // (adapted from https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h#L301) + template + ALPAKA_FN_HOST_ACC StatusCode unpackZSB(uint32_t chan, + const uint8_t* channel_data, + uint16_t channel_length, + uint32_t channel_offset, + SiStripDigiView out, + uint32_t* idx, + uint8_t headerLength, + uint16_t stripStart) { + constexpr uint16_t mask = (1 << num_bits) - 1; + if (channel_length & 0xF000) { + return StatusCode::BAD_CHANNEL_LENGTH; + } + uint_fast16_t wOffset = channel_offset + headerLength; // header is 2 (lite) or 7 + uint_fast8_t bOffset{0}, firstStrip{0}, nInCluster{0}, inCluster{0}; + const uint_fast16_t chEnd = channel_offset + channel_length; + while (((wOffset + 1) < chEnd) || + ((inCluster != nInCluster) && ((chEnd - wOffset) * BITS_PER_BYTE - bOffset >= num_bits))) { + if (inCluster == nInCluster) { + if (wOffset + 2 >= chEnd) { + // offset should already be at end then (empty cluster) + break; + } + if (bOffset) { + ++wOffset; + bOffset = 0; + } + const uint_fast8_t newFirstStrip = channel_data[(wOffset++) ^ 7]; + if (newFirstStrip < (firstStrip + inCluster)) { + // LogDebug("FEDBuffer") << "First strip of new cluster is not greater than last strip of previous cluster. " + // << "Last strip of previous cluster is " << uint16_t(firstStrip + inCluster) << ". " + // << "First strip of new cluster is " << uint16_t(newFirstStrip) << "."; + // printf("[%i] | UNORDERED_DATA - %i\t%i\t%i\n", *idx, newFirstStrip, firstStrip, inCluster); + return StatusCode::UNORDERED_DATA; + } + firstStrip = newFirstStrip; + nInCluster = channel_data[(wOffset++) ^ 7]; + inCluster = 0; + bOffset = 0; + } + bOffset += num_bits; + if ((num_bits > BITS_PER_BYTE) || (bOffset > BITS_PER_BYTE)) { + bOffset -= BITS_PER_BYTE; + out.adc(*idx) = ::sistrip::fedchannelunpacker::detail::getADC_B2(channel_data, wOffset, bOffset); + (*idx)++; + ++wOffset; + } else { + out.adc(*idx) = ::sistrip::fedchannelunpacker::detail::getADC_B1(channel_data, wOffset, bOffset); + (*idx)++; + } + out.channel(*idx) = chan; + out.stripId(*idx) = stripStart + firstStrip + inCluster; + ++inCluster; + if (bOffset == BITS_PER_BYTE) { + bOffset = 0; + ++wOffset; + } + } + return StatusCode::SUCCESS; + } + } // namespace detail + + namespace checks { + // (adapted from https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h#L391) + ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE bool isNonLiteZS(uint8_t mode) { return (mode == 10 || mode == 11); } + } // namespace checks + + namespace unpackers { + // (adapted from https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/EventFilter/SiStripRawToDigi/interface/SiStripFEDBuffer.h#L451) + ALPAKA_FN_HOST_ACC StatusCode unpackZeroSuppressed2(uint32_t chanIdx, + const uint8_t* channel_data, + uint16_t channel_length, + uint32_t channel_offset, + SiStripDigiView stripDigis, + uint32_t* aoff, + uint16_t stripStart, + bool isNonLite, + uint8_t mode, + uint8_t packetCode = 0) { + if ((isNonLite && packetCode == PACKET_CODE_ZERO_SUPPRESSED10) || + (mode == READOUT_MODE_ZERO_SUPPRESSED_LITE10 || mode == READOUT_MODE_ZERO_SUPPRESSED_LITE10_CMOVERRIDE)) { + // printf("[%i] unpackZSB<10>\n", *aoff); + return detail::unpackZSB<10>( + chanIdx, channel_data, channel_length, channel_offset, stripDigis, aoff, (isNonLite ? 7 : 2), stripStart); + } else if (mode == READOUT_MODE_PREMIX_RAW) { + // printf("[%i] unpackZSB<16>\n", *aoff); + return detail::unpackZSW<16>( + chanIdx, channel_data, channel_length, channel_offset, stripDigis, aoff, 7, stripStart); + } else { // 8bit + uint8_t bits_shift = 0; + if (isNonLite) { + if (packetCode == PACKET_CODE_ZERO_SUPPRESSED8_TOPBOT) + bits_shift = 1; + else if (packetCode == PACKET_CODE_ZERO_SUPPRESSED8_BOTBOT) + bits_shift = 2; + } else { // lite + if (mode == READOUT_MODE_ZERO_SUPPRESSED_LITE8_TOPBOT || + mode == READOUT_MODE_ZERO_SUPPRESSED_LITE8_TOPBOT_CMOVERRIDE) + bits_shift = 1; + else if (mode == READOUT_MODE_ZERO_SUPPRESSED_LITE8_BOTBOT || + mode == READOUT_MODE_ZERO_SUPPRESSED_LITE8_BOTBOT_CMOVERRIDE) + bits_shift = 2; + } + // printf("[%i] unpackZSB<8> [isNonLite %i] [packetCode %i] [mode %i] [channel_length %i] [channel_offset %i]\n", *aoff, isNonLite, packetCode, mode, channel_length, channel_offset); + auto st = detail::unpackZSW<8>(chanIdx, + channel_data, + channel_length, + channel_offset, + stripDigis, + aoff, + (isNonLite ? 7 : 2), + stripStart, + bits_shift); + if (isNonLite && packetCode == 0 && StatusCode::SUCCESS == st) { + // workaround for a pre-2015 bug in the packer: assume default ZS packing + // printf("[%i] ZERO_PACKET_CODE\n", *aoff); + return StatusCode::ZERO_PACKET_CODE; + } + return st; + } + } + } // namespace unpackers +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip::fedchannelunpacker + +// kernels and related objects +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + using namespace cms::alpakatools; + + constexpr uint16_t invalidStrip = std::numeric_limits::max(); + constexpr uint16_t invalidFed = std::numeric_limits::max(); + constexpr uint16_t badBit = (1 << 15); + constexpr uint16_t stripIndexMask = 0x7FFF; + + ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr uint16_t fedIndex(uint16_t fedId) { return fedId - FED_ID_MIN; } + ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr uint32_t stripIndex(uint16_t fedID, uint8_t fedCH, uint16_t strip) { + return fedIndex(fedID) * FEDCH_PER_FED * STRIPS_PER_FEDCH + fedCH * STRIPS_PER_FEDCH + (strip % STRIPS_PER_FEDCH); + } + ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr uint32_t apvIndex(uint16_t fed, uint8_t channel, uint16_t strip) { + return fedIndex(fed) * APVS_PER_FEDCH * FEDCH_PER_FED + APVS_PER_CHAN * channel + + (strip % STRIPS_PER_FEDCH) / STRIPS_PER_APV; + } + ALPAKA_FN_HOST_ACC ALPAKA_FN_INLINE constexpr uint32_t channelIndex(uint16_t fedId, uint8_t fedCh) { + return fedIndex(fedId) * FEDCH_PER_FED + fedCh; + } + + class SiStripKer_applyQualConds { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + uint32_t* blockStripN, + uint32_t* invalidFedCh, + uint8_t* fedChannelsData, + SiStripMappingView mapping, + const DetToFeds* qualityConditions) const { + // + for (auto chan : uniform_elements(acc, mapping.metadata().size())) { + const auto fedId = mapping.fedID(chan); + const auto fedCh = mapping.fedCh(chan); + + uint32_t index = channelIndex(fedId, fedCh); + + if (qualityConditions->qualityOk[index]) { + const uint32_t fedChOfs = mapping.fedChOfs(chan); + const uint32_t fedChDataOfsBuf = mapping.fedChDataOfsBuf(chan); + + bool isNonLite = fedchannelunpacker::checks::isNonLiteZS(mapping.readoutMode(chan)); + + FEDChannel fedChan(fedChannelsData + fedChDataOfsBuf, + fedChOfs, + isNonLite ? FEDChannel::ZSROMode::nonLite : FEDChannel::ZSROMode::lite); + + // Calculate the number of strips in the channel + uint8_t numBits = 8; + const uint8_t packetCode = fedChan.packetCode(); + const auto roMode = mapping.readoutMode(chan); + if ((isNonLite && packetCode == PACKET_CODE_ZERO_SUPPRESSED10) || + (roMode == READOUT_MODE_ZERO_SUPPRESSED_LITE10 || + roMode == READOUT_MODE_ZERO_SUPPRESSED_LITE10_CMOVERRIDE)) { + // unpackZSB<10> + numBits = 10; + } else if (roMode == READOUT_MODE_PREMIX_RAW) { + numBits = 16; + } + uint16_t chStripNb = fedChan.stripsInCh(numBits); + mapping.fedChStripsN(chan) = chStripNb; + + // uint8_t chLength = fedChan.length(); + // printf("#chLen,%i,%i,%i,%i,%i\n", fedId, fedCh, fedChOfs, chStripNb, blockStripN[0]); + // alpaka::atomicAdd(acc, blockStripN, static_cast(chStripNb), alpaka::hierarchy::Blocks{}); + } else { + // Atomic add the total number of strips which cannot be unpacked (unlikely case) + mapping.fedID(chan) = invalidFed; + mapping.fedChStripsN(chan) = 0; + // printf("#invFed,%i,%i,%i,%i\n", index, fedId, fedCh, invalidFedCh[0]); + alpaka::atomicAdd(acc, invalidFedCh, 1u, alpaka::hierarchy::Blocks{}); + } + } + } + }; + + class SiStripKer_init { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + const float r_channelThreshold, + const float r_seedThreshold, + const float r_clusterThresholdSquared, + const uint8_t r_maxSequentialHoles, + const uint8_t r_maxSequentialBad, + const uint8_t r_maxAdjacentBad, + const float r_minGoodCharge, + const uint32_t r_clusterSizeLimit, + StripClustersAuxView clusterDataObj) const { + if (once_per_grid(acc)) { + // Initialize the members of the clusterizer + clusterDataObj.channelThreshold() = r_channelThreshold; + clusterDataObj.seedThreshold() = r_seedThreshold; + clusterDataObj.clusterThresholdSquared() = r_clusterThresholdSquared; + clusterDataObj.maxSequentialHoles() = r_maxSequentialHoles; + clusterDataObj.maxSequentialBad() = r_maxSequentialBad; + clusterDataObj.maxAdjacentBad() = r_maxAdjacentBad; + clusterDataObj.minGoodCharge() = r_minGoodCharge; + clusterDataObj.clusterSizeLimit() = r_clusterSizeLimit; + } + } + }; + + class SiStripKer_unpackZS2 { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + uint8_t* fedChannelsData, + SiStripDigiView stripDigis, + SiStripMappingConstView mapping, + const GainNoiseCals* calibs) const { + // Loop over the FEDChannel collection to be digitized + for (auto chan : uniform_elements(acc, mapping.metadata().size())) { + // for (uint32_t chan=0; chan<(uint32_t)mapping.metadata().size(); chan++) { + const auto fedId = mapping.fedID(chan); + const auto fedCh = mapping.fedCh(chan); + // const auto detId = mapping.detID(chan); + + // reject strips which are in the conditions but not in the fed data collection + if (fedId == invalidFed) { + continue; + } + + const auto ipair = calibs->iPair[channelIndex(fedId, fedCh)]; + int ipoff = STRIPS_PER_FEDCH * ipair; + + const uint32_t fedChDataOfsBuf = mapping.fedChDataOfsBuf(chan); + uint16_t fedChOfs = mapping.fedChOfs(chan); + const uint8_t mode = mapping.readoutMode(chan); + const bool isNonLite = fedchannelunpacker::checks::isNonLiteZS(mode); + + const FEDChannel fedChan(fedChannelsData + fedChDataOfsBuf, fedChOfs, FEDChannel::ZSROMode(isNonLite)); + uint8_t packetCode = fedChan.packetCode(); + // packetCode for non-lite generic READOUT_MODE_ZERO_SUPPRESSED + // (https://github.com/cms-sw/cmssw/blob/CMSSW_16_0_X/RecoLocalTracker/SiStripClusterizer/plugins/ClustersFromRawProducer.cc#L398) + if (mode != READOUT_MODE_ZERO_SUPPRESSED) { + packetCode = 0; + } + packetCode = isNonLite ? packetCode : 0; + + uint32_t absoluteOffset = 0; + if (chan > 0) [[likely]] { + absoluteOffset = mapping.fedChStripsN(chan - 1); + } + + auto retCode = fedchannelunpacker::unpackers::unpackZeroSuppressed2(chan, + fedChan.data(), + fedChan.length(), + fedChan.offset(), + stripDigis, + &absoluteOffset, + ipoff, + isNonLite, + mode, + packetCode); + if (retCode != fedchannelunpacker::StatusCode::SUCCESS) { + printf("[%i] [fedID %hu] [fedCH %hhu] - Returned %i \n", chan, fedId, fedCh, (int)retCode); + } + } // data != nullptr && len > 0 + } + }; + + class SiStripKer_setSeedStrips { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + SiStripDigiConstView stripDigi, + StripClustersAuxView clusterDataObj, + SiStripMappingConstView mapping, + const GainNoiseCals* calibs) const { + auto nStrips = stripDigi.metadata().size(); + const float seedThreshold = clusterDataObj.seedThreshold(); + for (auto i : uniform_elements(acc, nStrips)) { + clusterDataObj.seedStripsMask(i) = 0; + clusterDataObj.seedStripsNCMask(i) = 0; + clusterDataObj.prefixSeedStripsNCMask(i) = 0; + + const auto ch = stripDigi.channel(i); + const auto fedId = mapping.fedID(ch); + const auto fedCh = mapping.fedCh(ch); + const auto stripID = stripDigi.stripId(i); + + const uint32_t idx = stripIndex(fedId, fedCh, stripID); + const uint16_t noise_tmp = calibs->noise[idx]; + const bool isBad = (noise_tmp & badBit) > 0; + + if (!isBad) { + const float noise_i = 0.1f * (noise_tmp & ~badBit); + const uint8_t adc_i = stripDigi.adc(i); + clusterDataObj.seedStripsMask(i) = (adc_i >= static_cast(noise_i * seedThreshold)) ? 1 : 0; + clusterDataObj.seedStripsNCMask(i) = clusterDataObj.seedStripsMask(i); + } + } + } + }; + + class SiStripKer_setNCSeedStrips { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + SiStripDigiConstView stripDigi, + StripClustersAuxView clusterDataObj, + SiStripMappingConstView mapping) const { + auto channels = stripDigi.channel(); + // Loop over the strips + for (auto stripIdx : uniform_elements(acc, 1, stripDigi.metadata().size())) { + const auto detid = mapping.detID(channels[stripIdx]); + const auto detid1 = mapping.detID(channels[stripIdx - 1]); + + if (clusterDataObj.seedStripsMask(stripIdx) && clusterDataObj.seedStripsMask(stripIdx - 1) && + (stripDigi.stripId(stripIdx) - stripDigi.stripId(stripIdx - 1)) == 1 && (detid == detid1)) { + clusterDataObj.seedStripsNCMask(stripIdx) = 0; + } + } + } + }; + + class SiStripKer_setNCStripIndex { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, StripClustersAuxView clusterDataObj) const { + // Loop over the strips + for (auto stripIdx : uniform_elements(acc, clusterDataObj.metadata().size())) { + if (clusterDataObj.seedStripsNCMask(stripIdx) == 1) { + const int index = (clusterDataObj.prefixSeedStripsNCMask(stripIdx) - 1); + clusterDataObj.seedStripsNCIndex(index) = stripIdx; + } + } + } + }; + + class SiStripKer_makeCandidates { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + const uint32_t kMaxSeedStrips, + SiStripDigiConstView stripDataObj, + StripClustersAuxConstView clusterDataObj, + SiStripClusterView clusters, + SiStripMappingConstView mapping, + const GainNoiseCals* calibs) const { + // + const int32_t nStrips = stripDataObj.metadata().size(); + const uint32_t nSeedStripsNC = (kMaxSeedStrips < clusterDataObj.prefixSeedStripsNCMask(nStrips - 1)) + ? kMaxSeedStrips + : clusterDataObj.prefixSeedStripsNCMask(nStrips - 1); + + const auto& chanArr = stripDataObj.channel(); + const auto& stripIdArr = stripDataObj.stripId(); + const auto& adcArr = stripDataObj.adc(); + + const float channelThreshold = clusterDataObj.channelThreshold(); + const float clusterThresholdSquared = clusterDataObj.clusterThresholdSquared(); + const uint8_t maxSequentialHoles = clusterDataObj.maxSequentialHoles(); + const uint16_t clusterSizeLimit = clusterDataObj.clusterSizeLimit(); + + // set this only once in the whole kernel grid + if (once_per_grid(acc)) { + clusters.nClusterCandidates() = nSeedStripsNC; + clusters.maxClusterSize() = clusterSizeLimit; + } + + // Loop over only the non-contiguous strips (flagged in setStripIndex) + for (auto i : uniform_elements(acc, nSeedStripsNC)) { + const uint32_t ch = clusterDataObj.seedStripsNCIndex(i); + + const auto chan = chanArr[ch]; + const auto fedId = mapping.fedID(chan); + const auto fedCh = mapping.fedCh(chan); + const auto detId = mapping.detID(chan); + const auto stripId = stripIdArr[ch]; + // + const uint32_t idx = stripIndex(fedId, fedCh, stripId); + const uint16_t noise_tmp = calibs->noise[idx]; + const float noise_i = 0.1f * (noise_tmp & ~badBit); + + // Calculate the accumulated noise2 and ADC + uint16_t clSize = 1; + float noiseSquared_i = noise_i * noise_i; + float adcSum_i = adcArr[ch]; + int32_t testIndex = ch - 1; + + auto addtocluster = [&](int& indexLR) { + const auto test_chan = chanArr[testIndex]; + const auto test_fedId = mapping.fedID(test_chan); + const auto test_fedCh = mapping.fedCh(test_chan); + const auto test_stripId = stripIdArr[testIndex]; + + const uint32_t test_idx = stripIndex(test_fedId, test_fedCh, test_stripId); + const uint16_t test_noise_tmp = calibs->noise[test_idx]; + const bool test_isBad = (test_noise_tmp & badBit) > 0; + + const float test_noise_i = 0.1f * (test_noise_tmp & ~badBit); + const uint8_t testADC = adcArr[testIndex]; + + if (!test_isBad && (testADC >= static_cast(test_noise_i * channelThreshold))) { + ++clSize; + indexLR = testIndex; + noiseSquared_i += test_noise_i * test_noise_i; + adcSum_i += testADC; + } + }; + + // find left boundary + int32_t indexLeft = ch; + + // if (stripIdArr[testIndex] == invalidStrip && testIndex >= 0) { + // testIndex -= 2; + // } + + if (testIndex >= 0) { + const auto testchan = chanArr[testIndex]; + const auto testDetId = mapping.detID(testchan); + + int16_t rangeLeft = stripIdArr[indexLeft] - stripIdArr[testIndex] - 1; + bool isSameDetLeft = (detId == testDetId); + + while (isSameDetLeft && (rangeLeft >= 0) && (rangeLeft <= maxSequentialHoles) && + (clSize < (clusterSizeLimit + 1))) { + addtocluster(indexLeft); + --testIndex; + // if (testIndex >= 0 && stripIdArr[testIndex] == invalidStrip) { + // testIndex -= 2; + // } + if (testIndex >= 0) { + rangeLeft = stripIdArr[indexLeft] - stripIdArr[testIndex] - 1; + const auto newchan = chanArr[testIndex]; + const auto newdet = mapping.detID(newchan); + isSameDetLeft = (detId == newdet); + } else { + isSameDetLeft = false; + } + } // while loop + } // testIndex >= 0 + + // find right boundary + int32_t indexRight = ch; + testIndex = ch + 1; + + // if (stripIdArr[testIndex] == invalidStrip && testIndex < nStrips) { + // testIndex += 2; + // } + + if (testIndex < nStrips) { + const auto testchan = chanArr[testIndex]; + const auto testDetId = mapping.detID(testchan); + + int16_t rangeRight = stripIdArr[testIndex] - stripIdArr[indexRight] - 1; + bool isSameDetRight = (detId == testDetId); + + while (isSameDetRight && (rangeRight >= 0) && (rangeRight <= maxSequentialHoles) && + (clSize < (clusterSizeLimit + 1))) { + addtocluster(indexRight); + ++testIndex; + // if (testIndex < nStrips && stripIdArr[testIndex] == invalidStrip) { + // testIndex += 2; + // } + if (testIndex < nStrips) { + rangeRight = stripIdArr[testIndex] - stripIdArr[indexRight] - 1; + const auto newchan = chanArr[testIndex]; + const auto newdet = mapping.detID(newchan); + isSameDetRight = (detId == newdet); + } else { + isSameDetRight = false; + } + } // while loop + } // testIndex < nStrips + + clusters.clusterIndex(i) = indexLeft; + clusters.clusterSize(i) = indexRight - indexLeft + 1; + clusters.clusterDetId(i) = detId; + clusters.firstStrip(i) = stripIdArr[indexLeft]; + // Flag candidates which do not pass (cluster noise threshold) && (cluster size) conditions. Max number of holes already accounted in candidate finder + // Floating point approximation leads sensitivity on the above condition for O(12/500) events, with NCluster deviation w.r.t. legacy module of O(1). + clusters.candidateAccepted(i) = (noiseSquared_i * clusterThresholdSquared <= adcSum_i * adcSum_i) && + (clusters.clusterSize(i) <= clusterSizeLimit); + clusters.candidateAcceptedPrefix(i) = static_cast(clusters.candidateAccepted(i)); + } // i < nSeedStripsNC + } + }; + + class SiStripKer_endCandidates { + public: + ALPAKA_FN_ACC void operator()(Acc1D const& acc, + SiStripDigiView stripDataObj, + StripClustersAuxConstView clusterDataObj, + SiStripClusterView clusters, + SiStripMappingConstView mapping, + const GainNoiseCals* calibs) const { + // + constexpr uint8_t adc_low_saturation = 254; + constexpr uint8_t adc_high_saturation = 255; + constexpr int charge_low_saturation = 253; + constexpr int charge_high_saturation = 1022; + // + const auto& clusterIndexArr = clusters.clusterIndex(); + + for (auto i : uniform_elements(acc, clusters.nClusterCandidates())) { + if (clusters.candidateAccepted(i)) { + const uint32_t indexLeft = clusterIndexArr[i]; + const uint16_t size = clusters.clusterSize(i); + + if (i > 0 && clusterIndexArr[i - 1] == indexLeft) { + clusters.candidateAccepted(i) = false; // ignore duplicates + clusters.candidateAcceptedPrefix(i) = 0; + } else { + float chargeSum = 0.0f; + uint16_t sumx = 0; + uint16_t suma = 0; + + uint16_t j = 0; + for (uint16_t k = 0; k < size; k++) { + uint32_t index = indexLeft + k; + const auto chan = stripDataObj.channel(index); + const auto fedId = mapping.fedID(chan); + const auto fedCh = mapping.fedCh(chan); + const auto stripId = stripDataObj.stripId(index); + + // ThreeThresholdAlgorithm::applyGains + const float gain_j = calibs->gain[apvIndex(fedId, fedCh, stripId)]; + uint8_t amplitudes_j = stripDataObj.adc(index); + const uint16_t charge = static_cast(static_cast(amplitudes_j) / gain_j + 0.5f); + + if (amplitudes_j < adc_low_saturation) { + amplitudes_j = ((charge > charge_high_saturation) + ? adc_high_saturation + : (charge > charge_low_saturation ? adc_low_saturation : charge)); + } + + // Overrides the ADC value in the StripDigi with the corrected amplitude + stripDataObj.adc(index) = amplitudes_j; + // clusters.clusterADCs(i)[j] = amplitudes_j; + + // SiStripCluster::initQB() + chargeSum += static_cast(amplitudes_j); + sumx += j * amplitudes_j; + suma += amplitudes_j; + j++; + } + + clusters.charge(i) = chargeSum; + + const auto chanL = stripDataObj.channel(indexLeft); + const auto fedIdL = mapping.fedID(chanL); + const auto fedChL = mapping.fedCh(chanL); + + // Flags for siStripClusterTools::chargePerCM in current candidate + clusters.candidateAccepted(i) = + (clusterDataObj.minGoodCharge() <= 0 || + ((chargeSum * calibs->invthick[channelIndex(fedIdL, fedChL)]) > clusterDataObj.minGoodCharge())); + + // SiStripCluster::initQB() -> barycenter_ + const float bary_i = static_cast(sumx) / static_cast(suma); + clusters.barycenter(i) = + static_cast(stripDataObj.stripId(indexLeft) & stripIndexMask) + bary_i + 0.5f; + + clusters.clusterSize(i) = j; + } // not a duplicate cluster + } + } // i < nSeedStripsNC + } + }; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +// kernels launchers +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + using namespace ::sistrip; + using namespace cms::alpakatools; + + SiStripRawToClusterAlgo::SiStripRawToClusterAlgo(const edm::ParameterSet& clustPar) + : channelThreshold_(clustPar.getParameter("ChannelThreshold")), + seedThreshold_(clustPar.getParameter("SeedThreshold")), + clusterThresholdSquared_(std::pow(clustPar.getParameter("ClusterThreshold"), 2.0f)), + maxSequentialHoles_(clustPar.getParameter("MaxSequentialHoles")), + maxSequentialBad_(clustPar.getParameter("MaxSequentialBad")), + maxAdjacentBad_(clustPar.getParameter("MaxAdjacentBad")), + maxClusterSize_(clustPar.getParameter("MaxClusterSize")), + minGoodCharge_(clusterChargeCut(clustPar)), + kMaxSeedStrips_(clustPar.getParameter("MaxSeedStrips")) { + // Checks non-sensical parameters + if (maxClusterSize_ > 768) { + throw cms::Exception("SiStripRawToClstAlg", "MaxClusterSize must be <= 768"); + } + } + + void SiStripRawToClusterAlgo::prepareUnpackCluster(Queue& queue, + const DetToFeds* conditions_DetToFeds, + std::unique_ptr rFEDChMover) { + // Move ownership of the host-data container this class + fedChMover_ = std::move(rFEDChMover); + + // If not, allocate nStrips_h_ for the first time + if (!nStrips_h_) { + nStrips_h_ = cms::alpakatools::make_host_buffer(queue); + } + + // There are no channels (=>strips) to unpack, set to strip-to-unpack to zero and return + if (!fedChMover_ || fedChMover_->channelNb() == 0) { + alpaka::memset(queue, *nStrips_h_, 0); + return; + } + + // Move the data to the device + fedBuffer_d_.emplace(cms::alpakatools::make_device_buffer(queue, fedChMover_->bufferSize())); + alpaka::memcpy(queue, *fedBuffer_d_, fedChMover_->buffer(), fedChMover_->bufferSize()); + + // Move fedchannels to the device + stripMapping_d_ = cms::alpakatools::moveToDeviceAsync(queue, fedChMover_->mapping()); + + // Apply quality conditions to mapping and calculate the number of strips to unpack + uint32_t divider = 256u; + uint32_t groups = divide_up_by(fedChMover_->channelNb(), divider); + auto workDiv = make_workdiv(groups, divider); + + auto nstrips_d = make_device_buffer(queue); + alpaka::memset(queue, nstrips_d, 0); + + auto invalidFedChN_d = make_device_buffer(queue); + alpaka::memset(queue, invalidFedChN_d, 0); + + alpaka::exec(queue, + workDiv, + SiStripKer_applyQualConds{}, + nstrips_d.data(), + invalidFedChN_d.data(), + fedBuffer_d_->data(), + stripMapping_d_->view(), + conditions_DetToFeds); + + const uint32_t threads = 1024u; + const uint32_t nBlocks = divide_up_by(stripMapping_d_->view().metadata().size(), threads); + const auto workDivMultiBlock = make_workdiv(nBlocks, threads); + auto blockCounter_d = make_device_buffer(queue); + alpaka::memset(queue, blockCounter_d, 0); + alpaka::exec(queue, + workDivMultiBlock, + multiBlockPrefixScan(), + stripMapping_d_->const_view().fedChStripsN().data(), + stripMapping_d_->view().fedChStripsN().data(), + (uint32_t)stripMapping_d_->view().metadata().size(), + (int32_t)nBlocks, + blockCounter_d.data(), + alpaka::getPreferredWarpSize(alpaka::getDev(queue))); + + // To do: merge the two kernels for applyQualConds and sum into a single one. + + // Get the number of strips to unpack + nStrips_h_ = cms::alpakatools::make_host_buffer(queue); + const uint32_t size = stripMapping_d_->view().metadata().size(); + auto viewSrc = make_device_view(queue, stripMapping_d_->view().fedChStripsN(size - 1)); + alpaka::memcpy(queue, *nStrips_h_, viewSrc); + } + + uint32_t SiStripRawToClusterAlgo::unpackStrips(Queue& queue, const GainNoiseCals* calibs) { + // Allocate the SiStripDigi collection on device + const uint32_t nStrips = *nStrips_h_->data(); + // Skip all if there are no strips to be unpacked + if (nStrips == 0) { + return nStrips; + } + digis_d_ = std::make_unique(nStrips, queue); + + // Run the unpacking kernel + uint32_t divider = 256u; + uint32_t nChannels = (*stripMapping_d_)->metadata().size(); + uint32_t groups = divide_up_by(nChannels, divider); + auto workDiv = make_workdiv(groups, divider); + alpaka::exec(queue, + workDiv, + SiStripKer_unpackZS2{}, + fedBuffer_d_->data(), + digis_d_->view(), + stripMapping_d_->const_view(), + calibs); + + // dumpUnpackedStrips(queue, digis_d_.get()); // (for debugging) + + // Allocate and initialize the StripClustersAux collection + sClustersAux_d_ = std::make_unique(nStrips, queue); + // LogDebug("sClustersAux") << "Size of StripClustersAuxDevice (bytes): " << alpaka::getExtentProduct(sClustersAux_d_->buffer()) * sizeof(std::byte); + alpaka::exec(queue, + make_workdiv(1u, 1u), + SiStripKer_init{}, + channelThreshold_, + seedThreshold_, + clusterThresholdSquared_, + maxSequentialHoles_, + maxSequentialBad_, + maxAdjacentBad_, + minGoodCharge_, + maxClusterSize_, + sClustersAux_d_->view()); + + // Cluster seeding + alpaka::exec(queue, + workDiv, + SiStripKer_setSeedStrips{}, + digis_d_->const_view(), + sClustersAux_d_->view(), + stripMapping_d_->const_view(), + calibs); + + // Un-seed any contiguous strips in the same detector + alpaka::exec(queue, + workDiv, + SiStripKer_setNCSeedStrips{}, + digis_d_->const_view(), + sClustersAux_d_->view(), + stripMapping_d_->const_view()); + + // Calculate the discrete integral (prefix sum) of seedStripsNCMask. + // When the integral increase AND I am at a non-contigous strip, the beginning of new cluster is marked. + const uint32_t nThreads = 1024u; + const int32_t nBlocks = divide_up_by(nStrips, nThreads); + auto blockCounter_d = make_device_buffer(queue); + alpaka::memset(queue, blockCounter_d, 0); + alpaka::exec(queue, + make_workdiv(nBlocks, nThreads), + multiBlockPrefixScan(), + sClustersAux_d_->const_view().seedStripsNCMask().data(), + sClustersAux_d_->view().prefixSeedStripsNCMask().data(), + nStrips, + nBlocks, + blockCounter_d.data(), + alpaka::getPreferredWarpSize(alpaka::getDev(queue))); + + // Get the total number of non-contiguous seeds (ready on the host in produce) + nSeeds_h_ = cms::alpakatools::make_host_buffer(queue); + auto viewSrc = make_device_view(queue, sClustersAux_d_->view().prefixSeedStripsNCMask(nStrips - 1)); + alpaka::memcpy(queue, *nSeeds_h_, viewSrc); + + // Find index of the non-contiguous strip seeds + alpaka::exec(queue, workDiv, SiStripKer_setNCStripIndex{}, sClustersAux_d_->view()); + + // dumpSeeds(queue, digis_d_.get(), sClustersAux_d_.get(), &stripMapping_d_.value()); // (for debugging) + return nStrips; + } + + std::unique_ptr SiStripRawToClusterAlgo::makeClusters(Queue& queue, + const GainNoiseCals* calibs) { + // The maximum number of clusters is set to kMaxSeedStrips + auto clusters_d = std::make_unique(kMaxSeedStrips_, queue); + clusters_d->zeroInitialise(queue); + + // The number of seed over which to loop for clusters is the min between the number of strips and the kMaxSeeds + const uint32_t nStrips = *nStrips_h_->data(); + const uint32_t nSeeds = std::min(kMaxSeedStrips_, nStrips); + + uint32_t divider = 256u; + uint32_t groups = divide_up_by(nSeeds, divider); + auto workDiv = make_workdiv(groups, divider); + + // Three-threshold clusterization algo + alpaka::exec(queue, + workDiv, + SiStripKer_makeCandidates{}, + kMaxSeedStrips_, + digis_d_->const_view(), + sClustersAux_d_->const_view(), + clusters_d->view(), + stripMapping_d_->const_view(), + calibs); + + // dumpClusters(queue, clusters_d.get(), digis_d_.get()); + + // Apply the conditions + alpaka::exec(queue, + workDiv, + SiStripKer_endCandidates{}, + digis_d_->view(), + sClustersAux_d_->const_view(), + clusters_d->view(), + stripMapping_d_->const_view(), + calibs); + + // Fill the prefix indexes for the candidateAccepted + const uint32_t nThreads = 1024u; + const int32_t nBlocks = divide_up_by(kMaxSeedStrips_, nThreads); + auto blockCounter_d = make_device_buffer(queue); + alpaka::memset(queue, blockCounter_d, 0); + alpaka::exec(queue, + make_workdiv(nBlocks, nThreads), + multiBlockPrefixScan(), + clusters_d->const_view().candidateAcceptedPrefix().data(), + clusters_d->view().candidateAcceptedPrefix().data(), + kMaxSeedStrips_, + nBlocks, + blockCounter_d.data(), + alpaka::getPreferredWarpSize(alpaka::getDev(queue))); + + // Store the total number of good cluster candidates into the scalar of the StripDigi SoA + auto viewSrc_realClusters = + make_device_view(queue, clusters_d->view().candidateAcceptedPrefix(kMaxSeedStrips_ - 1)); + auto viewDst_realClusters = make_device_view(queue, digis_d_->view().nbGoodCandidates()); + alpaka::memcpy(queue, viewDst_realClusters, viewSrc_realClusters); + + // Store also the number of cluster candidates, in order to reduce as much as possible the loop for the slimming in the legacy converter + auto viewSrc_candidatesN = make_device_view(queue, clusters_d->view().nClusterCandidates()); + auto viewDst_candidatesN = make_device_view(queue, digis_d_->view().nbCandidates()); + alpaka::memcpy(queue, viewDst_candidatesN, viewSrc_candidatesN); + + // dumpClusters(queue, clusters_d.get(), digis_d_.get()); + return clusters_d; + } + + std::unique_ptr SiStripRawToClusterAlgo::releaseDigiAmplitudes() { + // Allow the release of device memory (when queue completes) + sClustersAux_d_.reset(); + return std::move(digis_d_); + } +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +// Debugging functions +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + void SiStripRawToClusterAlgo::dumpUnpackedStrips(Queue& queue, SiStripDigiDevice* digis_d) { + const int digisSize = digis_d->const_view().metadata().size(); + auto digis_h = SiStripDigiHost(digisSize, queue); + alpaka::memcpy(queue, digis_h.buffer(), digis_d->const_buffer()); + alpaka::wait(queue); + std::ostringstream dumpMsg(""); + dumpMsg << "[SiStripRawToClusterAlgo::unpackStrips] Dumping unpacked strips\n"; + dumpMsg << "Allocated " << digisSize << " strips\n"; + dumpMsg << "i,chan,stripId,adc\n"; + for (int i = 0; i < digisSize; ++i) { + if (i < std::min(digisSize, 100) || i > (digisSize - 100)) { + dumpMsg << i << "," << (int)(digis_h->channel(i)) << "," << (int)(digis_h->stripId(i)) << "," + << (int)digis_h->adc(i) << "\n"; + } + } + std::cout << dumpMsg.str(); + } + + void SiStripRawToClusterAlgo::dumpSeeds(Queue& queue, + SiStripDigiDevice* digis_d, + StripClustersAuxDevice* sClustersAux_d, + SiStripMappingDevice* mapping_d) { + // Store the size of the digi to avoid repetitions + const int digisSize = digis_d->const_view().metadata().size(); + auto digis_h = SiStripDigiHost(digisSize, queue); + alpaka::memcpy(queue, digis_h.buffer(), digis_d->const_buffer()); + // Seed table and digis have the same size + auto sClustersAux_h = StripClustersAuxHost(digisSize, queue); + alpaka::memcpy(queue, sClustersAux_h.buffer(), sClustersAux_d->const_buffer()); + // Mapping table + auto mapping_h = SiStripMappingHost(mapping_d->const_view().metadata().size(), queue); + alpaka::memcpy(queue, mapping_h.buffer(), mapping_d->const_buffer()); + alpaka::wait(queue); + + std::ostringstream dumpMsg; + dumpMsg << "[SiStripRawToClusterAlgo::setSeedsAndMakeIndexes] Dumping seeds table\n"; + dumpMsg << "i,adc,chan,stripId,detId,seed,seedNC,pfxNCMask,seedNCIndex\n"; + for (int i = 0; i < digisSize; ++i) { + if (i < 600 || i > (digisSize - 600) || i % 10000 == 0) { + if (digis_h->stripId(i) != invalidStrip) { + const auto chan = digis_h->channel(i); + dumpMsg << i << "," << (int)(digis_h->adc(i)) << "," << digis_h->channel(i) << "," << digis_h->stripId(i) + << "," << mapping_h->detID(chan) << "," << sClustersAux_h->seedStripsMask(i) << "," + << sClustersAux_h->seedStripsNCMask(i) << "," << sClustersAux_h->prefixSeedStripsNCMask(i) << "," + << sClustersAux_h->seedStripsNCIndex(i) << "\n"; + } + } + } + std::cout << dumpMsg.str(); + } + + void SiStripRawToClusterAlgo::dumpClusters(Queue& queue, + SiStripClusterDevice* clusters_d, + SiStripDigiDevice* digis_d, + bool fullDump) { + // Store the size of the digi to avoid repetitions + const int clustersPrealloc = clusters_d->view().metadata().size(); + auto clusters_h = SiStripClusterHost(clustersPrealloc, queue); + alpaka::memcpy(queue, clusters_h.buffer(), clusters_d->const_buffer()); + + const uint32_t nStrips = digis_d->view().metadata().size(); + auto digis_h = SiStripDigiHost(nStrips, queue); + alpaka::memcpy(queue, digis_h.buffer(), digis_d->const_buffer()); + + alpaka::wait(queue); + + const int clustersN = clusters_h->nClusterCandidates(); + + std::ostringstream dumpMsg(""); + dumpMsg << "#clDump,Pre-allocated:" << clustersPrealloc << ",Candidates:" << clustersN << "\n"; + dumpMsg << "i,cIdx,cSz,cDetId,chg,1st,tCl,tClIdx,bary,|clusterADCs|\n"; + + for (int i = 0; i < clustersN; ++i) { + if (fullDump || i < 100 || i > (clustersN - 100)) { + dumpMsg << i << "," << clusters_h->clusterIndex(i) << "," << clusters_h->clusterSize(i) << "," + << clusters_h->clusterDetId(i) << "," << clusters_h->charge(i) << "," << clusters_h->firstStrip(i) + << "," << clusters_h->candidateAccepted(i) << "," << clusters_h->candidateAcceptedPrefix(i) << "," + << clusters_h->barycenter(i) << ",|"; + if (clusters_h->candidateAccepted(i)) { + for (int j = 0; j < clusters_h->clusterSize(i); ++j) { + uint32_t index = clusters_h->clusterIndex(i) + j; + dumpMsg << (int)(digis_h->adc(index)); + if (j != (clusters_h->clusterSize(i) - 1)) { + dumpMsg << "/"; + } + } + } else { + dumpMsg << "-"; + } + dumpMsg << "|\n"; + } + } + dumpMsg << "#zClDump\n"; + dumpMsg << "#goodCandidates," << digis_h.view().nbGoodCandidates() << ":" + << clusters_h->candidateAcceptedPrefix(clustersPrealloc - 1) << "\n"; + + dumpMsg << "#last\n"; + dumpMsg << "i,cIdx,cSz,cDetId,chg,1st,tCl,tClIdx,bary,|clusterADCs|\n"; + for (int i = clustersPrealloc - 100; i > 0 && i < clustersPrealloc; ++i) { + dumpMsg << i << "," << clusters_h->clusterIndex(i) << "," << clusters_h->clusterSize(i) << "," + << clusters_h->clusterDetId(i) << "," << clusters_h->charge(i) << "," << clusters_h->firstStrip(i) << "," + << clusters_h->candidateAccepted(i) << "," << clusters_h->candidateAcceptedPrefix(i) << "," + << clusters_h->barycenter(i) << ",|"; + if (clusters_h->candidateAccepted(i)) { + for (int j = 0; j < clusters_h->clusterSize(i); ++j) { + uint32_t index = clusters_h->clusterIndex(i) + j; + dumpMsg << (int)(digis_h->adc(index)); + if (j != (clusters_h->clusterSize(i) - 1)) { + dumpMsg << "/"; + } + } + } else { + dumpMsg << "-"; + } + dumpMsg << "|\n"; + } + std::cout << dumpMsg.str(); + } +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip diff --git a/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.h b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.h new file mode 100644 index 0000000000000..4d735c0cc70e6 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/plugins/alpaka/SiStripRawToClusterAlgo.h @@ -0,0 +1,154 @@ +#ifndef RecoLocalTracker_SiStripClusterizer_plugins_alpaka_SiStripRawToClusterAlgo_h +#define RecoLocalTracker_SiStripClusterizer_plugins_alpaka_SiStripRawToClusterAlgo_h + +#include "DataFormats/Portable/interface/PortableHostCollection.h" +#include "DataFormats/Portable/interface/alpaka/PortableCollection.h" + +#include "DataFormats/FEDRawData/interface/FEDRawData.h" + +#include "DataFormats/SiStripClusterSoA/interface/alpaka/SiStripClusterDevice.h" +#include "DataFormats/SiStripDigiSoA/interface/alpaka/SiStripDigiDevice.h" + +#include "HeterogeneousCore/AlpakaInterface/interface/config.h" +#include "HeterogeneousCore/AlpakaInterface/interface/memory.h" + +#include "RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripMappingDevice.h" + +namespace edm { + class ParameterSet; +} + +///////////////////////////////////////////////////////////////////////////////////////////////////////////////////// +// The StripClustersAuxSoALayout is an auxiliary SoA, with the same size of the SiStripDigiSoA, which helps +// in the clusterization process. +namespace sistrip { + GENERATE_SOA_LAYOUT(StripClustersAuxSoALayout, + SOA_COLUMN(uint32_t, seedStripsMask), + SOA_COLUMN(uint32_t, seedStripsNCMask), + SOA_COLUMN(uint32_t, seedStripsNCIndex), + SOA_COLUMN(uint32_t, prefixSeedStripsNCMask), + // + SOA_SCALAR(float, channelThreshold), + SOA_SCALAR(float, seedThreshold), + SOA_SCALAR(float, clusterThresholdSquared), + SOA_SCALAR(uint8_t, maxSequentialHoles), + SOA_SCALAR(uint8_t, maxSequentialBad), + SOA_SCALAR(uint8_t, maxAdjacentBad), + SOA_SCALAR(float, minGoodCharge), + SOA_SCALAR(uint16_t, clusterSizeLimit)) + + using StripClustersAuxSoA = StripClustersAuxSoALayout<>; + using StripClustersAuxView = StripClustersAuxSoA::View; + using StripClustersAuxConstView = StripClustersAuxSoA::ConstView; + using StripClustersAuxHost = PortableHostCollection; +} // namespace sistrip + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + using namespace ::sistrip; + using StripClustersAuxDevice = PortableCollection; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +ASSERT_DEVICE_MATCHES_HOST_COLLECTION(sistrip::StripClustersAuxDevice, sistrip::StripClustersAuxHost); +///////////////////////////////////////////////////////////////////////////////////////////////////////////////////// + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + struct FEDChMetadata { + uint32_t detId; + uint16_t fedId; + uint8_t fedCh; + // + uint16_t fedChOfs; + uint32_t fedChOfs_wrt_rawFedId; + uint8_t fedBufferRO; + }; + + class PortableFEDMover { + public: + PortableFEDMover(Queue& queue, uint32_t rawBufferSize, uint32_t fedChannelsNb) + : buffer_(cms::alpakatools::make_host_buffer(queue, rawBufferSize)), + mapping_(fedChannelsNb, queue), + bufferSize_(0), + channelNb_(fedChannelsNb), + offset4FedId_(sistrip::FED_ID_MAX + 1) {} + + void fillBuffer(const std::vector& raw) { + for (uint16_t fedId = sistrip::FED_ID_MIN; fedId <= sistrip::FED_ID_MAX; ++fedId) { + if (raw[fedId]) { + std::memcpy(buffer_.data() + bufferSize_, raw[fedId]->data(), raw[fedId]->size()); + offset4FedId_[fedId] = bufferSize_; + bufferSize_ += raw[fedId]->size(); + } + } + } + + void fillMapping(const std::vector& channelsMeta) { + for (uint32_t i = 0; i < channelNb_; i++) { + const auto& channel = channelsMeta[i]; + mapping_->detID(i) = channel.detId; + mapping_->fedID(i) = channel.fedId; + mapping_->fedCh(i) = channel.fedCh; + mapping_->fedChOfs(i) = channel.fedChOfs; + mapping_->fedChDataOfsBuf(i) = channel.fedChOfs_wrt_rawFedId + offset4FedId_[channel.fedId]; + mapping_->readoutMode(i) = channel.fedBufferRO; + } + } + + auto buffer() { return buffer_; } + uint32_t bufferSize() const { return bufferSize_; } + + auto mapping() { return std::move(mapping_); } + uint32_t channelNb() const { return channelNb_; } + + private: + cms::alpakatools::host_buffer buffer_; + SiStripMappingHost mapping_; + + uint32_t bufferSize_; + uint32_t channelNb_; + std::vector offset4FedId_; + }; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip { + class SiStripRawToClusterAlgo { + public: + SiStripRawToClusterAlgo(const edm::ParameterSet& clustPar); + void prepareUnpackCluster(Queue& queue, + const DetToFeds* conditions_DetToFeds, + std::unique_ptr FEDChMover); + uint32_t unpackStrips(Queue& queue, const GainNoiseCals* calibs); + + std::unique_ptr makeClusters(Queue& queue, const GainNoiseCals* calibs); + std::unique_ptr releaseDigiAmplitudes(); + + private: + // Clusterizer parameters + const float channelThreshold_, seedThreshold_, clusterThresholdSquared_; + const uint8_t maxSequentialHoles_, maxSequentialBad_, maxAdjacentBad_; + const uint32_t maxClusterSize_; + const float minGoodCharge_; + const uint32_t kMaxSeedStrips_; + + std::unique_ptr fedChMover_; + + std::optional> fedBuffer_d_; + std::optional stripMapping_d_; + std::optional> nStrips_h_; + std::optional> nSeeds_h_; + + int nStripsBytes_ = 0; + std::unique_ptr digis_d_; + std::unique_ptr sClustersAux_d_; + + // Debug functions + void dumpUnpackedStrips(Queue& queue, SiStripDigiDevice* digis_d); + void dumpSeeds(Queue& queue, + SiStripDigiDevice* digis_d, + StripClustersAuxDevice* sClustersAux_d, + SiStripMappingDevice* mapping_d); + void dumpClusters(Queue& queue, SiStripClusterDevice* clusters_d, SiStripDigiDevice* digis_d, bool fullDump = false); + }; +} // namespace ALPAKA_ACCELERATOR_NAMESPACE::sistrip + +#endif // RecoLocalTracker_SiStripClusterizer_plugins_alpaka_SiStripRawToClusterAlgo_h diff --git a/RecoLocalTracker/SiStripClusterizer/python/customizeStripClustersFromRaw.py b/RecoLocalTracker/SiStripClusterizer/python/customizeStripClustersFromRaw.py index 0d5e53c56539c..e2afa2d575146 100644 --- a/RecoLocalTracker/SiStripClusterizer/python/customizeStripClustersFromRaw.py +++ b/RecoLocalTracker/SiStripClusterizer/python/customizeStripClustersFromRaw.py @@ -22,3 +22,68 @@ def customizeHLTStripClustersFromRaw(process): cms.Sequence(process.hltSiStripRawToClustersFacility, process.siStripClustersTaskCUDA)) return process + +def customizeHLTStripClustersFromRaw_alpaka(process: cms.Process, MaxClusterSize:int = 768, doNotReplaceInPath = []): + if hasattr(process, 'hltSiStripRawToClustersFacility'): + from RecoLocalTracker.SiStripZeroSuppression.DefaultAlgorithms_cff import DefaultAlgorithms + from RecoLocalTracker.SiStripClusterizer.DefaultClusterizer_cff import DefaultClusterizer + + # Take the parameters from the original producer + # (parameters_() makes a copy, this could cause issues if the parameters are changed after this customizer is applied) + initialPars = process.hltSiStripRawToClustersFacility.parameters_() + + # Make sure the module configuration is compatible with the heterogeneous producer + if 'LegacyUnpacker' in initialPars and initialPars['LegacyUnpacker'] == True: + raise RuntimeError("The LegacyUnpacker is not supported by the alpaka producer. Please set LegacyUnpacker to False.") + if 'HybridZeroSuppressed' in initialPars and initialPars['HybridZeroSuppressed'] == True: + # The alpaka version supports READOUT_MODE_ZERO_SUPPRESSED* modes + raise RuntimeError("The HybridZeroSuppressed is not supported by the alpaka producer. Please set HybridZeroSuppressed to False.") + + # Create the alpaka-producer and set its parameters from the original one + hltSiStripRawToClustersFacilityAlpaka = cms.EDProducer("sistrip::SiStripRawToCluster@alpaka", **initialPars) + + # Add the extra pars (if not present) + if not hasattr(hltSiStripRawToClustersFacilityAlpaka, "CablingConditionsLabel"): hltSiStripRawToClustersFacilityAlpaka.CablingConditionsLabel = cms.string("") + ## Make sure the Clusterizer PSet has the MaxClusterSize argument + if not hasattr(hltSiStripRawToClustersFacilityAlpaka.Clusterizer, "MaxClusterSize"): + # print(f"[hltSiStripRawToClustersFacility] No Clusterizer.MaxClusterSize defined. Defaulting to {MaxClusterSize}") + hltSiStripRawToClustersFacilityAlpaka.Clusterizer.MaxClusterSize = cms.uint32(MaxClusterSize) + ## Add the MaxSeedStrips to the clusterizer PSet if not present + if not hasattr(hltSiStripRawToClustersFacilityAlpaka.Clusterizer, "MaxSeedStrips"): + hltSiStripRawToClustersFacilityAlpaka.Clusterizer.MaxSeedStrips = cms.uint32(200000) + + # Remove illegal parameters from the configuration from comparison of + # the legacy `edmPluginHelp -b -p SiStripClusterizerFromRaw` and + # the heterogeneous `edmPluginHelp -b -p sistrip::SiStripRawToCluster@alpaka` parameters. + for par in ['onDemand', 'LegacyUnpacker', 'HybridZeroSuppressed', 'Algorithms']: + if hasattr(hltSiStripRawToClustersFacilityAlpaka, par): delattr(hltSiStripRawToClustersFacilityAlpaka, par) + + # Create the converter bringing the alpaka-made cluster into legacy objects + hltSiStripClustersToLegacy = cms.EDProducer("sistrip::SiStripClustersToLegacy", + source = cms.InputTag("hltSiStripRawToClustersFacilityAlpaka") + ) + + ####### Produce ES for the alpaka clusterizer ####### + # Produce the SiStripClusterizerConditionsHost object + hltSiStripClusterizerConditionsESProducerAlpaka = cms.ESProducer("sistrip::SiStripClusterizerConditionsESProducerAlpaka@alpaka", + QualityLabel = cms.ESInputTag(""), + Label = cms.ESInputTag(""), + ) + + # Add to process, if not already present + if not hasattr(process, "hltSiStripClusterizerConditionsESProducerAlpaka"): + process.hltSiStripClusterizerConditionsESProducerAlpaka = hltSiStripClusterizerConditionsESProducerAlpaka + process.hltSiStripRawToClustersFacilityAlpaka = hltSiStripRawToClustersFacilityAlpaka + process.hltSiStripRawToClustersFacility = hltSiStripClustersToLegacy + + sequencesToFix = [] + for path_name, path in process.paths_().items(): + if path.contains(process.hltSiStripRawToClustersFacility): + # If Path is in doNotReplaceInPath, then don't replace there + if path_name in doNotReplaceInPath: continue + + sequencesToFix.append(path_name) + pathToFix_ptr = getattr(process, path_name) + pathToFix_ptr.replace(process.hltSiStripRawToClustersFacility, process.hltSiStripRawToClustersFacilityAlpaka + process.hltSiStripRawToClustersFacility) + + return process diff --git a/RecoLocalTracker/SiStripClusterizer/src/ES_SiStripClusterizerConditionsHost.cc b/RecoLocalTracker/SiStripClusterizer/src/ES_SiStripClusterizerConditionsHost.cc new file mode 100644 index 0000000000000..015ba84c53d80 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/src/ES_SiStripClusterizerConditionsHost.cc @@ -0,0 +1,5 @@ +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsHostObject.h" +#include "FWCore/Utilities/interface/typelookup.h" + +TYPELOOKUP_DATA_REG(sistrip::SiStripClusterizerConditionsDetToFedsHostObject); +TYPELOOKUP_DATA_REG(sistrip::SiStripClusterizerConditionsGainNoiseCalsHostObject); diff --git a/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterizerConditionsRecord.cc b/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterizerConditionsRecord.cc new file mode 100644 index 0000000000000..79d6d44ccd7c7 --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/src/SiStripClusterizerConditionsRecord.cc @@ -0,0 +1,5 @@ +#include "RecoLocalTracker/SiStripClusterizer/interface/SiStripClusterizerConditionsRecord.h" +#include "FWCore/Framework/interface/eventsetuprecord_registration_macro.h" + +EVENTSETUP_RECORD_REG(sistrip::SiStripClusterizerConditionsDetToFedsRecord); +EVENTSETUP_RECORD_REG(sistrip::SiStripClusterizerConditionsGainNoiseCalsRecord); diff --git a/RecoLocalTracker/SiStripClusterizer/src/alpaka/ES_SiStripClusterizerConditionsDevice.cc b/RecoLocalTracker/SiStripClusterizer/src/alpaka/ES_SiStripClusterizerConditionsDevice.cc new file mode 100644 index 0000000000000..b7b69cbc5ef1f --- /dev/null +++ b/RecoLocalTracker/SiStripClusterizer/src/alpaka/ES_SiStripClusterizerConditionsDevice.cc @@ -0,0 +1,5 @@ +#include "HeterogeneousCore/AlpakaCore/interface/alpaka/typelookup.h" +#include "RecoLocalTracker/SiStripClusterizer/interface/alpaka/SiStripClusterizerConditionsDeviceObject.h" + +TYPELOOKUP_ALPAKA_DATA_REG(sistrip::SiStripClusterizerConditionsDetToFedsDeviceObject); +TYPELOOKUP_ALPAKA_DATA_REG(sistrip::SiStripClusterizerConditionsGainNoiseCalsDeviceObject);