emulate_from_tpstream.cxx

Go to the documentation of this file.

#include "trgtools/TCEmulationUnit.hpp"
#include "trgtools/TAEmulationWorker.hpp"
 
#include "CLI/App.hpp"
#include "CLI/Config.hpp"
#include "CLI/Formatter.hpp"
 
#include <fmt/core.h>
#include <fmt/format.h>
#include <fmt/chrono.h>
#include <filesystem>
#include <optional>
 
#include "hdf5libs/HDF5RawDataFile.hpp"
#include "hdf5libs/HDF5SourceIDHandler.hpp"
 
#include "triggeralgs/TriggerCandidateFactory.hpp"
 
//-----------------------------------------------------------------------------
 
using namespace dunedaq;
using namespace trgtools;
 
void save_fragments(const std::string& _outputfilename,
                   const hdf5libs::HDF5SourceIDHandler::source_id_geo_id_map_t& _sourceid_geoid_map,
                   std::map<uint64_t, std::vector<std::unique_ptr<daqdataformats::Fragment>>> _frags,
                   bool _quiet)
{
  std::string output_filename = _outputfilename;
  if (output_filename.size() < 5 ||
      output_filename.compare(output_filename.size() - 5, 5, ".hdf5") != 0) {
    output_filename += ".hdf5";
  }
 
  // Create layout parameter object required for HDF5 creation
  hdf5libs::HDF5FileLayoutParameters layout_params;
 
  // Create HDF5 parameter path required for the layout
  hdf5libs::HDF5PathParameters params_trigger;
  params_trigger.detector_group_type = "Detector_Readout";
  params_trigger.detector_group_name = "TPC";
  params_trigger.element_name_prefix = "Link";
  params_trigger.digits_for_element_number = 5;
 
  // Fill the HDF5 layout
  std::vector<hdf5libs::HDF5PathParameters> params;
  params.push_back(params_trigger);
  layout_params.record_name_prefix = "TimeSlice";
  layout_params.digits_for_record_number = 6;
  layout_params.digits_for_sequence_number = 0;
  layout_params.record_header_dataset_name = "TimeSliceHeader";
  layout_params.raw_data_group_name = "RawData";
  layout_params.view_group_name = "Views";
  layout_params.path_params_list = {params};
 
  // Create pointer to a new output HDF5 file
  std::unique_ptr<hdf5libs::HDF5RawDataFile> output_file = std::make_unique<hdf5libs::HDF5RawDataFile>(
      output_filename,
      _frags.begin()->second[0]->get_run_number(),
      0,
      "emulate_from_tpstream",
      layout_params,
      _sourceid_geoid_map);
 
  // Iterate over the time slices & save all the fragments
  for (auto& [slice_id, vec_frags]: _frags) {
    // Create a new timeslice header
    daqdataformats::TimeSliceHeader tsh;
    tsh.timeslice_number = slice_id;
    tsh.run_number = vec_frags[0]->get_run_number();
    tsh.element_id = dunedaq::daqdataformats::SourceID(dunedaq::daqdataformats::SourceID::Subsystem::kTRBuilder, 0);
 
    // Create a new timeslice
    dunedaq::daqdataformats::TimeSlice ts(tsh);
    if (!_quiet) {
      std::cout << "Time slice number: " << slice_id << std::endl;
    }
    // Add the fragments to the timeslices
    for (std::unique_ptr<daqdataformats::Fragment>& frag_ptr: vec_frags) {
      if (!_quiet) {
        std::cout << "  Writing elementid: " << frag_ptr->get_element_id()  << " trigger number: " << frag_ptr->get_trigger_number() << " trigger_timestamp: " << frag_ptr->get_trigger_timestamp() << " window_begin: " << frag_ptr->get_window_begin() << " sequence_no: " << frag_ptr->get_sequence_number() << std::endl;
      }
      ts.add_fragment(std::move(frag_ptr));
    }
 
    // Write the timeslice to output file
    output_file->write(ts);
  }
}
 
std::map<std::string, std::vector<std::shared_ptr<hdf5libs::HDF5RawDataFile>>>
sort_files_per_writer(const std::vector<std::string>& _files)
{
  std::map<std::string, std::vector<std::shared_ptr<hdf5libs::HDF5RawDataFile>>> files_sorted;
 
  // Put files into per-writer app groups
  for (const std::string& file : _files) {
    std::shared_ptr<hdf5libs::HDF5RawDataFile> fl = std::make_shared<hdf5libs::HDF5RawDataFile>(file);
    if (!fl->is_timeslice_type()) {
      throw std::runtime_error(fmt::format("ERROR: input file '{}' not of type 'TimeSlice'", file));
    }
 
    std::string application_name = fl->get_attribute<std::string>("application_name");
    files_sorted[application_name].push_back(fl);
  }
 
  // Sort files for each writer application individually
  for (auto& [app_name, vec_files]: files_sorted) {
    // Don't sort if we have 0 or 1 files in the application...
    if (vec_files.size() <= 1) {
      continue;
    }
 
    // Sort w.r.t. file index attribute
    std::sort(vec_files.begin(), vec_files.end(),
        [](const std::shared_ptr<hdf5libs::HDF5RawDataFile>& a, const std::shared_ptr<hdf5libs::HDF5RawDataFile>& b) {
        return a->get_attribute<size_t>("file_index") <
               b->get_attribute<size_t>("file_index");
        });
  }
 
  return files_sorted;
};
 
std::pair<uint64_t, uint64_t>
get_available_slice_id_range(const std::map<std::string, std::vector<std::shared_ptr<hdf5libs::HDF5RawDataFile>>>& _files,
                         bool _quiet)
{
  if (_files.empty()) {
    throw std::runtime_error("No files provided");
  }
 
  uint64_t global_start = std::numeric_limits<uint64_t>::min();
  uint64_t global_end = std::numeric_limits<uint64_t>::max();
 
  // Get the min & max record id per application, and the global
  for (auto& [appid, vec_files]: _files) {
    uint64_t app_start = std::numeric_limits<uint64_t>::max();
    uint64_t app_end = std::numeric_limits<uint64_t>::min();
 
    // Find min / max record id for this application
    for (auto& file : vec_files) {
      auto record_ids = file->get_all_record_ids();
      if (record_ids.empty()) {
        throw std::runtime_error(fmt::format("File from application {} contains no records.", appid));
      }
      app_start = std::min(app_start, record_ids.begin()->first);
      app_end = std::max(app_end, record_ids.rbegin()->first);
    }
 
    // Update the global min / max record id
    global_start = std::max(global_start, app_start);
    global_end = std::min(global_end, app_end);
 
    if (!_quiet) {
      std::cout << "Application: " << appid << " " << " TimeSliceID start: " << app_start << " end: " << app_end << std::endl;
    }
  }
 
  if (!_quiet) {
    std::cout << "Global start: " << global_start << " Global end: " << global_end << std::endl;
  }
  if (global_start > global_end) {
    throw std::runtime_error("One of the provided files' id range did not overlap with the rest. Please select files with overlapping TimeSlice IDs");
  }
 
  // Extra validation / error handling
  for (auto& [appid, vec_files]: _files) {
    for (auto& file : vec_files) {
      auto record_ids = file->get_all_record_ids();
 
      uint64_t file_start = record_ids.begin()->first;
      uint64_t file_end = record_ids.rbegin()->first;
      if ((file_start > global_end || file_end < global_start)) {
        uint64_t file_index = file->get_attribute<size_t>("file_index");
        throw std::runtime_error(fmt::format(
          "File from TPStreamWrite application '{}' (index '{}') has record id range [{}, {}], which does not overlap with global range [{}, {}].",
           appid, file_index, file_start, file_end, global_start, global_end
          ));
      }
    }
  }
 
  return {global_start, global_end};
}
 
struct Options
{
  std::vector<std::string> input_files;
  std::string output_filename;
  std::string config_name;
  bool quiet = false;
  bool latencies = false;
  bool run_parallel = false;
  std::optional<uint64_t> slice_start = std::nullopt;
  std::optional<uint64_t> num_slices = std::nullopt;
};
 
void parse_app(CLI::App& _app, Options& _opts)
{
  _app.add_option("-i,--input-files", _opts.input_files, "List of input files (required)")
    ->required()
    ->check(CLI::ExistingFile); // Validate that each file exists
 
  _app.add_option("-o,--output-file", _opts.output_filename, "Output file (required)")
    ->required(); // make the argument required
 
  _app.add_option("-j,--json-config", _opts.config_name, "Trigger Activity and Candidate config JSON to use (required)")
    ->required()
    ->check(CLI::ExistingFile);
  
  _app.add_flag("--parallel", _opts.run_parallel, "Run the TAMakers in parallel");
 
  _app.add_flag("--quiet", _opts.quiet, "Quiet outputs.");
 
  _app.add_flag("--latencies", _opts.latencies, "Saves latencies per TP into csv");
 
  _app.add_option("-s,--slice-start", _opts.slice_start, "Inclusive lower bound for TimeSlice ID to process");
  _app.add_option("-n,--num-slices", _opts.num_slices, "Number of TimeSlices to process");
}
 
int main(int argc, char const *argv[])
{
  // Do all the CLI processing first
  CLI::App app{"Offline trigger TriggerActivity & TriggerCandidate emulatior"};
  Options opts{};
 
  parse_app(app, opts);
 
  try {
    app.parse(argc, argv);
  }
  catch (const CLI::ParseError &e) {
    return app.exit(e);
  }
 
  // Get the configuration file
  std::ifstream config_stream(opts.config_name);
  nlohmann::json config = nlohmann::json::parse(config_stream);
 
  if (!opts.quiet) {
    std::cout << "Files to process:\n";
    for (const std::string& file : opts.input_files) {
      std::cout << "- " << file << "\n";
    }
  }
 
  // Sort the files into a map writer_id::vector<HDF5>
  std::map<std::string, std::vector<std::shared_ptr<hdf5libs::HDF5RawDataFile>>> sorted_files =
    sort_files_per_writer(opts.input_files);
 
  // Get the available record_id range
  std::pair<uint64_t, uint64_t> recordid_range = get_available_slice_id_range(sorted_files, opts.quiet);
  std::pair<uint64_t, uint64_t> processing_range = recordid_range;
 
  if (opts.slice_start.has_value()) {
    processing_range.first = std::max(processing_range.first, opts.slice_start.value());
  }
 
  if (opts.num_slices.has_value()) {
    if (opts.num_slices.value() == 0) {
      throw std::runtime_error("Invalid --num-slices value: must be greater than 0");
    }
 
    uint64_t requested_end = processing_range.first + (opts.num_slices.value() - 1);
    if (requested_end < processing_range.first) {
      requested_end = std::numeric_limits<uint64_t>::max();
    }
 
    processing_range.second = std::min(processing_range.second, requested_end);
  }
 
  if (processing_range.first > processing_range.second) {
    throw std::runtime_error(fmt::format(
      "Requested timeslice subset (slice-start={}, num-slices={}) does not overlap available range [{}, {}]",
      opts.slice_start.value_or(0),
      opts.num_slices.value_or(0),
      recordid_range.first,
      recordid_range.second));
  }
 
  if (!opts.quiet) {
    std::cout << "Processing TimeSliceID range: [" << processing_range.first << ", " << processing_range.second << "]" << std::endl;
  }
 
  // Create the file handlers
  std::vector<std::unique_ptr<TAEmulationWorker>> ta_emu_workers;
  for (auto [name, files] : sorted_files) {
    ta_emu_workers.push_back(std::make_unique<TAEmulationWorker>(files, config, processing_range, opts.run_parallel, opts.quiet));
  }
 
  // Start each file handler
  for (const auto& handler : ta_emu_workers) {
    handler->start_processing();
  }
 
  // Output map of TA vectors & function that appends TAs to that vector
  std::map<uint64_t, std::vector<triggeralgs::TriggerActivity>> tas;
  auto append_tas = [&tas](std::map<uint64_t, std::vector<triggeralgs::TriggerActivity>>&& _tas) {
    for (auto& [sliceid, src_vec] : _tas) {
      auto& dest_vec = tas[sliceid];
      dest_vec.reserve(dest_vec.size() + src_vec.size());
 
      dest_vec.insert(dest_vec.end(), std::make_move_iterator(src_vec.begin()), std::make_move_iterator(src_vec.end()));
      src_vec.clear();
    }
  };
 
  // Output map of Fragment vectors & function that appends TAs to that vector
  std::map<uint64_t, std::vector<std::unique_ptr<daqdataformats::Fragment>>> frags;
  auto append_frags = [&frags](std::map<uint64_t, std::vector<std::unique_ptr<daqdataformats::Fragment>>>&& _frags) {
    for (auto& [sliceid, src_vec] : _frags) {
      auto& dest_vec = frags[sliceid];
 
      dest_vec.reserve(dest_vec.size() + src_vec.size());
 
      dest_vec.insert(dest_vec.end(), std::make_move_iterator(src_vec.begin()), std::make_move_iterator(src_vec.end()));
      src_vec.clear();
    }
  };
 
  // Iterate over the handlers, wait for them to complete their job & append
  // their TAs to our vector when ready.
  hdf5libs::HDF5SourceIDHandler::source_id_geo_id_map_t sourceid_geoid_map;
  for (const auto& handler : ta_emu_workers) {
    // Wait for all TAs to be made
    handler->wait_to_complete_work();
 
    // Append output TA/fragments
    append_tas(std::move(handler->get_tas()));
    append_frags(std::move(handler->get_frags()));
 
    // Get and merge the source ID map
    hdf5libs::HDF5SourceIDHandler::source_id_geo_id_map_t map = handler->get_sourceid_geoid_map();
    sourceid_geoid_map.insert(map.begin(), map.end());
  }
 
  // Sort the TAs in each slice before pushing them into the TCMaker
  size_t n_tas = 0;
  for (auto& [sliceid, vec_tas]: tas) {
    std::sort(vec_tas.begin(), vec_tas.end(),
        [](const triggeralgs::TriggerActivity& a, const triggeralgs::TriggerActivity& b) {
        return std::tie(a.time_start, a.channel_start, a.time_end) <
               std::tie(b.time_start, b.channel_start, b.time_end);
        });
    n_tas += vec_tas.size();
  }
 
  if (!opts.quiet) {
    std::cout << "Total number of TAs made: " << n_tas << std::endl;
    std::cout << "Creating a TCMaker..." << std::endl;
  }
  // Create the TC emulator
  std::string algo_name = config["trigger_candidate_plugin"][0];
  nlohmann::json algo_config = config["trigger_candidate_config"][0];
 
  std::unique_ptr<triggeralgs::TriggerCandidateMaker> tc_maker =
    triggeralgs::TriggerCandidateFactory::get_instance()->build_maker(algo_name);
  tc_maker->configure(algo_config);
 
  trgtools::TCEmulationUnit tc_emulator;
  tc_emulator.set_maker(tc_maker);
 
  // Emulate the TriggerCandidates
  std::vector<triggeralgs::TriggerCandidate> tcs;
  for (auto& [sliceid, vec_tas]: tas) {
    std::unique_ptr<daqdataformats::Fragment> tc_frag = tc_emulator.emulate_vector(vec_tas);
    if (!tc_frag) {
      continue;
    }
 
    // Manipulate the fragment header
    daqdataformats::FragmentHeader frag_hdr = tc_frag->get_header();
    frag_hdr.element_id = daqdataformats::SourceID{daqdataformats::SourceID::Subsystem::kTrigger, tc_frag->get_element_id().id+10000};
 
    tc_frag->set_header_fields(frag_hdr);
    tc_frag->set_type(daqdataformats::FragmentType::kTriggerCandidate);
    tc_frag->set_trigger_number(sliceid);
    tc_frag->set_window_begin(frags[sliceid][0]->get_window_begin());
    tc_frag->set_window_end(frags[sliceid][0]->get_window_end());
 
    // Push the fragment to our list
    frags[sliceid].push_back(std::move(tc_frag));
 
    std::vector<triggeralgs::TriggerCandidate> tmp_tcs = tc_emulator.get_last_output_buffer();
    tcs.reserve(tcs.size() + tmp_tcs.size());
    tcs.insert(tcs.end(), std::make_move_iterator(tmp_tcs.begin()), std::make_move_iterator(tmp_tcs.end()));
  }
  if (!opts.quiet) {
    std::cout << "Total number of TCs made: " << tcs.size() << std::endl;
  }
 
  if (!frags.empty()) {
    save_fragments(opts.output_filename, sourceid_geoid_map, std::move(frags), opts.quiet);
  } else {
    std::cout << "No TA/TC fragments generated. Output file will not be generated" << std::endl;
  }
 
  return 0;
}