Skip to content

Commit

Permalink
Merge pull request #177 from DUNE-DAQ/aeo/threshold_by_plane
Browse files Browse the repository at this point in the history 
TPG Threshold By Plane
  • Loading branch information
@aeoranday
aeoranday committed May 14, 2024
2 parents aaa56d9 + 119fc48 commit cb2144d
Show file tree
Hide file tree
Showing 10 changed files with 66 additions and 36 deletions.
16 changes: 12 additions & 4 deletions include/fdreadoutlibs/wibeth/WIBEthFrameProcessor.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,7 @@ class WIBEthFrameHandler {

void reset();

void initialize(uint16_t threshold_value, uint16_t memory_factor, uint16_t scale_factor, int16_t frug_streaming_acclimt);
void initialize(uint16_t memory_factor, uint16_t scale_factor, int16_t frug_streaming_acclimt);

uint16_t* get_hits_dest();

Expand Down Expand Up @@ -146,7 +146,7 @@ class WIBEthFrameProcessor : public readoutlibs::TaskRawDataProcessorModel<types
private:
bool m_tpg_enabled;

bool m_enable_simple_threshold_on_collection = false;
bool m_enable_simple_threshold_on_plane2 = false;
// Selected TPG algorithm properties from configuration
std::string m_tpg_algorithm;
uint16_t m_tpg_rs_memory_factor;
Expand All @@ -157,7 +157,10 @@ class WIBEthFrameProcessor : public readoutlibs::TaskRawDataProcessorModel<types
uint32_t m_tp_max_width;
std::vector<int> m_channel_mask_vec;
std::set<uint> m_channel_mask_set;
uint16_t m_tpg_threshold;

uint16_t m_tpg_threshold_plane0;
uint16_t m_tpg_threshold_plane1;
uint16_t m_tpg_threshold_plane2;

// Algorithm used to form a trigger primitive
dunedaq::trgdataformats::TriggerPrimitive::Algorithm m_tp_algo = trgdataformats::TriggerPrimitive::Algorithm::kUnknown;
Expand All @@ -181,10 +184,15 @@ class WIBEthFrameProcessor : public readoutlibs::TaskRawDataProcessorModel<types


// Create an array to store the values of the memory factor
// AAA: silver bullet to be able to use SimpleThreshold on collection and RS on induction planes
// AAA: silver bullet to be able to use SimpleThreshold on plane 2 and RS on planes 0 and 1
// By default set all the values to the selected memory factor
std::array<uint16_t, swtpg_wibeth::NUM_REGISTERS_PER_FRAME * swtpg_wibeth::SAMPLES_PER_REGISTER> m_register_memory_factor = {0};

// Create an array to store the values of the TPG threshold
// This is to be used for setting a different value by plane
// By default set it to 150
std::array<uint16_t, swtpg_wibeth::NUM_REGISTERS_PER_FRAME * swtpg_wibeth::SAMPLES_PER_REGISTER> m_tpg_threshold = {150};


std::function<void(swtpg_wibeth::ProcessingInfo<swtpg_wibeth::NUM_REGISTERS_PER_FRAME>& info)> m_assigned_tpg_algorithm_function;

Expand Down
7 changes: 4 additions & 3 deletions include/fdreadoutlibs/wibeth/tpg/ProcessAVX2.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,7 @@ process_window_avx2(ProcessingInfo<NREGISTERS>& info)
// Pointer to keep track of where we'll write the next output hit
__m256i* output_loc = (__m256i*)(info.output); // NOLINT(readability/casting)

const __m256i iota = _mm256_set_epi16(15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0);
const __m256i iota = _mm256_lddqu_si256((__m256i*) IOTA);

int nhits = 0;

Expand All @@ -42,6 +42,8 @@ process_window_avx2(ProcessingInfo<NREGISTERS>& info)
// from the previous go-around.

ChanState<NREGISTERS>& state = info.chanState;
__m256i threshold = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.threshold) + ireg); // NOLINT

__m256i median = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.pedestals) + ireg); // NOLINT
// The accumulator that we increase/decrease when the current
// sample is greater/less than the median
Expand Down Expand Up @@ -93,8 +95,7 @@ process_window_avx2(ProcessingInfo<NREGISTERS>& info)
// --------------------------------------------------------------
// Mask for channels that are over the threshold in this step

// FIXED THRESHOLD
__m256i threshold = _mm256_set1_epi16(info.threshold);
// Define a register for elements above the threhsold
__m256i is_over = _mm256_cmpgt_epi16(s, threshold);

// Mask for channels that left "over threshold" state this step
Expand Down
7 changes: 4 additions & 3 deletions include/fdreadoutlibs/wibeth/tpg/ProcessAbsRSAVX2.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ process_window_rs_avx2(ProcessingInfo<NREGISTERS>& info)
// Pointer to keep track of where we'll write the next output hit
__m256i* output_loc = (__m256i*)(info.output); // NOLINT(readability/casting)

const __m256i iota = _mm256_set_epi16(15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0);
const __m256i iota = _mm256_lddqu_si256((__m256i*) IOTA);

int nhits = 0;

Expand All @@ -55,6 +55,8 @@ process_window_rs_avx2(ProcessingInfo<NREGISTERS>& info)
// from the previous go-around.

ChanState<NREGISTERS>& state = info.chanState;
__m256i threshold = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.threshold) + ireg); // NOLINT

__m256i median = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.pedestals) + ireg); // NOLINT
//__m256i quantile25 = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.quantile25) + ireg); // NOLINT
//__m256i quantile75 = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.quantile75) + ireg); // NOLINT
Expand Down Expand Up @@ -183,8 +185,7 @@ process_window_rs_avx2(ProcessingInfo<NREGISTERS>& info)
// IQR-based THRESHOLD
//__m256i is_over = _mm256_cmpgt_epi16(RS, sigma * info.threshold);

// FIXED THRESHOLD
__m256i threshold = _mm256_set1_epi16(info.threshold);
// Define a register for elements above the threhsold
__m256i is_over = _mm256_cmpgt_epi16(RS, threshold);


Expand Down
3 changes: 2 additions & 1 deletion include/fdreadoutlibs/wibeth/tpg/ProcessNaive.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,7 @@ process_window_naive(ProcessingInfo<NREGISTERS>& info)
int16_t& accum = state.accum[ichan];

// Variables for hit finding
int16_t& threshold = state.threshold[ichan]; // Threshold for this channel.
uint16_t& prev_was_over = state.prev_was_over[ichan]; // was the previous sample over threshold?
uint16_t& hit_charge = state.hit_charge[ichan];
uint16_t& hit_tover = state.hit_tover[ichan]; // time over threshold
Expand Down Expand Up @@ -90,7 +91,7 @@ process_window_naive(ProcessingInfo<NREGISTERS>& info)
// --------------------------------------------------------------
// Hit finding
// --------------------------------------------------------------
bool is_over = sample > info.threshold;
bool is_over = sample > threshold;
//printf("% 5d % 5d % 5d % 5d\n", (uint16_t)ichan, (uint16_t)itime, sample, info.threshold); // NOLINT
if (is_over) {
// Simulate saturated add
Expand Down
3 changes: 2 additions & 1 deletion include/fdreadoutlibs/wibeth/tpg/ProcessNaiveRS.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,7 @@ process_window_naive_RS(ProcessingInfo<NREGISTERS>& info)
//int16_t& quantile75 = state.quantile75[ichan];

// Variables for hit finding
int16_t& threshold = state.threshold[ichan]; // Threshold for this channel.
uint16_t& prev_was_over = state.prev_was_over[ichan]; // was the previous sample over threshold?
uint16_t& hit_charge = state.hit_charge[ichan];
uint16_t& hit_tover = state.hit_tover[ichan]; // time over threshold
Expand Down Expand Up @@ -116,7 +117,7 @@ process_window_naive_RS(ProcessingInfo<NREGISTERS>& info)
// --------------------------------------------------------------
// Hit finding
// --------------------------------------------------------------
bool is_over = RS > info.threshold;
bool is_over = RS > threshold;
if (is_over) {
// Simulate saturated add
int32_t tmp_charge = hit_charge;
Expand Down
9 changes: 5 additions & 4 deletions include/fdreadoutlibs/wibeth/tpg/ProcessStandardRSAVX2.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -31,15 +31,15 @@ process_window_standard_rs_avx2(ProcessingInfo<NREGISTERS>& info)

// Scaling factor to stop the ADCs from overflowing
// (may not needs this, depends on magnitude of FIR output)
const __m256i scale_factor = _mm256_set1_epi16(info.rs_scale_factor);
//const __m256i scale_factor = _mm256_set1_epi16(info.rs_scale_factor);

// The maximum value that sigma can have before the threshold overflows a 16-bit signed integer
//const __m256i sigmaMax = _mm256_set1_epi16((1 << 15) / (info.multiplier * info.threshold));

// Pointer to keep track of where we'll write the next output hit
__m256i* output_loc = (__m256i*)(info.output); // NOLINT(readability/casting)

const __m256i iota = _mm256_set_epi16(15, 14, 13, 12, 11, 10, 9, 8, 7, 6, 5, 4, 3, 2, 1, 0);
const __m256i iota = _mm256_lddqu_si256((__m256i*) IOTA);

int nhits = 0;

Expand All @@ -55,6 +55,8 @@ process_window_standard_rs_avx2(ProcessingInfo<NREGISTERS>& info)
// from the previous go-around.

ChanState<NREGISTERS>& state = info.chanState;
__m256i threshold = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.threshold) + ireg); // NOLINT

__m256i median = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.pedestals) + ireg); // NOLINT
//__m256i quantile25 = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.quantile25) + ireg); // NOLINT
//__m256i quantile75 = _mm256_lddqu_si256(reinterpret_cast<__m256i*>(state.quantile75) + ireg); // NOLINT
Expand Down Expand Up @@ -178,8 +180,7 @@ process_window_standard_rs_avx2(ProcessingInfo<NREGISTERS>& info)
// IQR-based THRESHOLD
//__m256i is_over = _mm256_cmpgt_epi16(RS, sigma * info.threshold);

// FIXED THRESHOLD
__m256i threshold = _mm256_set1_epi16(info.threshold);
// Define a register for elements above the threhsold
__m256i is_over = _mm256_cmpgt_epi16(RS, threshold);


Expand Down
17 changes: 14 additions & 3 deletions include/fdreadoutlibs/wibeth/tpg/ProcessingInfo.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@ struct ChanState
ChanState()
{
for (size_t i = 0; i < NREGISTERS * SAMPLES_PER_REGISTER; ++i) {
threshold[i] = 0;
pedestals[i] = 0;
accum[i] = 0;
RS[i] = 0;
Expand All @@ -39,6 +40,7 @@ struct ChanState
}
}

alignas(32) int16_t __restrict__ threshold[NREGISTERS * SAMPLES_PER_REGISTER];
alignas(32) int16_t __restrict__ pedestals[NREGISTERS * SAMPLES_PER_REGISTER];
alignas(32) int16_t __restrict__ accum[NREGISTERS * SAMPLES_PER_REGISTER];

Expand Down Expand Up @@ -74,7 +76,6 @@ struct ProcessingInfo
uint8_t last_register_, // NOLINT
uint16_t* __restrict__ output_, // NOLINT
const uint8_t exponent_, // NOLINT
uint16_t threshold_, // NOLINT
uint16_t rs_memory_factor_, // NOLINT
uint16_t rs_scale_factor_, // NOLINT
int16_t frugal_streaming_accumulator_limit_, // NOLINT
Expand All @@ -86,7 +87,6 @@ struct ProcessingInfo
, last_register(last_register_)
, output(output_)
, exponent(exponent_)
, threshold(threshold_)
, rs_memory_factor(rs_memory_factor_)
, rs_scale_factor(rs_scale_factor_)
, frugal_streaming_accumulator_limit(frugal_streaming_accumulator_limit_)
Expand All @@ -95,6 +95,18 @@ struct ProcessingInfo
, nhits(nhits_)
{}

// Set TPG threshold by plane
void setThresholdState(std::array<uint16_t, swtpg_wibeth::NUM_REGISTERS_PER_FRAME * swtpg_wibeth::SAMPLES_PER_REGISTER>& register_threshold
)
{
// Set the threshold values
for (size_t j = 0; j < NREGISTERS * SAMPLES_PER_REGISTER; ++j) {
const size_t i = IOTA[j % SAMPLES_PER_REGISTER];
const size_t reg_idx = j / SAMPLES_PER_REGISTER;
chanState.threshold[j] = register_threshold[i + reg_idx * SAMPLES_PER_REGISTER];
}

}

// Set the initial state from the window starting at first_msg_p
template<size_t N>
Expand Down Expand Up @@ -153,7 +165,6 @@ struct ProcessingInfo
uint8_t last_register; // NOLINT
uint16_t* __restrict__ output; // NOLINT
uint8_t exponent; // NOLINT
uint16_t threshold; // NOLINT
uint16_t rs_memory_factor; // NOLINT
uint16_t rs_scale_factor; // NOLINT
int16_t frugal_streaming_accumulator_limit; // NOLINT
Expand Down
3 changes: 3 additions & 0 deletions include/fdreadoutlibs/wibeth/tpg/TPGConstants_wibeth.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,9 @@ const constexpr std::size_t BYTES_PER_REGISTER = 32;
// How many samples are in a register
const constexpr std::size_t SAMPLES_PER_REGISTER = 16;

// Order of the channels after frame expansion.
const constexpr std::uint16_t IOTA[SAMPLES_PER_REGISTER] = {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15};

// One netio message's worth of channel ADCs after
// expansion: 64 frames per message times 4 registers per frame times
// 32 bytes (256 bits) per register
Expand Down
2 changes: 1 addition & 1 deletion src/wib2/WIB2FrameProcessor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -209,7 +209,7 @@ WIB2FrameProcessor::conf(const nlohmann::json& cfg)
// AAA: The set provides faster look up than a std::vector
m_channel_mask_set.insert(m_channel_mask_vec.begin(), m_channel_mask_vec.end());

m_tpg_threshold_selected = config.tpg_threshold;
m_tpg_threshold_selected = config.tpg_threshold_default;

m_crate_no = config.crate_id;
m_slot_no = config.slot_id;
Expand Down
35 changes: 19 additions & 16 deletions src/wibeth/WIBEthFrameProcessor.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,7 @@ WIBEthFrameHandler::reset()
}

void
WIBEthFrameHandler::initialize(uint16_t threshold_value, uint16_t memory_factor, uint16_t scale_factor, int16_t frug_streaming_acclimt)
WIBEthFrameHandler::initialize(uint16_t memory_factor, uint16_t scale_factor, int16_t frug_streaming_acclimt)
{

if(m_hits_dest == nullptr) {m_hits_dest = new uint16_t[100000];}
Expand All @@ -83,7 +83,6 @@ WIBEthFrameHandler::initialize(uint16_t threshold_value, uint16_t memory_factor,
swtpg_wibeth::NUM_REGISTERS_PER_FRAME,
m_hits_dest,
m_tpg_exponent,
threshold_value,
memory_factor,
scale_factor,
frug_streaming_acclimt,
Expand Down Expand Up @@ -116,8 +115,7 @@ WIBEthFrameProcessor::start(const nlohmann::json& args)
m_tps_suppressed_too_long = 0;
m_tps_send_failed = 0;

m_wibeth_frame_handler->initialize(m_tpg_threshold,
m_tpg_rs_memory_factor,
m_wibeth_frame_handler->initialize(m_tpg_rs_memory_factor,
m_tpg_rs_scale_factor,
m_tpg_frugal_streaming_accumulator_limit
);
Expand Down Expand Up @@ -185,13 +183,13 @@ WIBEthFrameProcessor::conf(const nlohmann::json& cfg)
} else if (m_tpg_algorithm == "AbsRS" ) {
m_tp_algo = trgdataformats::TriggerPrimitive::Algorithm::kAbsRunningSum;
m_assigned_tpg_algorithm_function = &swtpg_wibeth::process_window_rs_avx2<swtpg_wibeth::NUM_REGISTERS_PER_FRAME>;
// Enable simple threshold on collection is used only if you are using a Running Sum algorithm
m_enable_simple_threshold_on_collection = config.enable_simple_threshold_on_collection;
// Enable simple threshold on plane 2 is used only if you are using a Running Sum algorithm
m_enable_simple_threshold_on_plane2 = config.enable_simple_threshold_on_plane2;
} else if (m_tpg_algorithm == "StandardRS" ) {
m_tp_algo = trgdataformats::TriggerPrimitive::Algorithm::kRunningSum;
m_assigned_tpg_algorithm_function = &swtpg_wibeth::process_window_standard_rs_avx2<swtpg_wibeth::NUM_REGISTERS_PER_FRAME>;
// Enable simple threshold on collection is used only if you are using a Running Sum algorithm
m_enable_simple_threshold_on_collection = config.enable_simple_threshold_on_collection;
// Enable simple threshold on plane 2 is used only if you are using a Running Sum algorithm
m_enable_simple_threshold_on_plane2 = config.enable_simple_threshold_on_plane2;
} else {
throw TPGAlgorithmInexistent(ERS_HERE, m_tpg_algorithm);
}
Expand All @@ -217,7 +215,10 @@ WIBEthFrameProcessor::conf(const nlohmann::json& cfg)
// AAA: The set provides faster look up than a std::vector
m_channel_mask_set.insert(m_channel_mask_vec.begin(), m_channel_mask_vec.end());

m_tpg_threshold = config.tpg_threshold;
// Use config.<plane>_threshold if it is non-zero, else default to config.tpg_threshold.
m_tpg_threshold_plane2 = config.tpg_threshold_plane2 ? config.tpg_threshold_plane2 : config.tpg_threshold_default;
m_tpg_threshold_plane1 = config.tpg_threshold_plane1 ? config.tpg_threshold_plane1 : config.tpg_threshold_default;
m_tpg_threshold_plane0 = config.tpg_threshold_plane0 ? config.tpg_threshold_plane0 : config.tpg_threshold_default;

m_crate_no = config.crate_id;
m_slot_no = config.slot_id;
Expand Down Expand Up @@ -438,15 +439,16 @@ WIBEthFrameProcessor::find_hits(constframeptr fp, WIBEthFrameHandler* frame_hand
auto chan_value = frame_handler->register_channel_map.channel[i];
m_register_channels[i] = chan_value;

if (m_enable_simple_threshold_on_collection) {
// If the given channel is a collection then set R (memory factor) to zero
if (m_channel_map->get_plane_from_offline_channel(chan_value) == 0 ) {
m_register_memory_factor[i] = 0;
} else {
m_register_memory_factor[i] = m_tpg_rs_memory_factor;
}
if (m_channel_map->get_plane_from_offline_channel(chan_value) == 2 ) {
// If SimpleThreshold on plane 2, then set the memory factor to 0, else use the common memory factor.
m_register_memory_factor[i] = m_enable_simple_threshold_on_plane2 ? 0 : m_tpg_rs_memory_factor;
m_tpg_threshold[i] = m_tpg_threshold_plane2;
} else if (m_channel_map->get_plane_from_offline_channel(chan_value) == 1) {
m_register_memory_factor[i] = m_tpg_rs_memory_factor;
m_tpg_threshold[i] = m_tpg_threshold_plane1;
} else {
m_register_memory_factor[i] = m_tpg_rs_memory_factor;
m_tpg_threshold[i] = m_tpg_threshold_plane0;
}

//TLOG () << "Index number " << i << " offline channel " << frame_handler->register_channel_map.channel[i];
Expand All @@ -455,6 +457,7 @@ WIBEthFrameProcessor::find_hits(constframeptr fp, WIBEthFrameHandler* frame_hand
m_tp_channel_rate_map[frame_handler->register_channel_map.channel[i]] = 0;
}

frame_handler->m_tpg_processing_info->setThresholdState(m_tpg_threshold);
frame_handler->m_tpg_processing_info->setState(registers_array, m_register_memory_factor);


Expand Down

0 comments on commit cb2144d

Please sign in to comment.