26 static inline int quant(
double min,
double max,
unsigned nSteps,
double val);
27 static inline double unquant(
double min,
double max,
unsigned nSteps,
int step);
29 static inline std::string
to_string(
const std::vector<T> &
v);
47 ATH_MSG_INFO(
"Convolution: " << to_string(
const_cast<std::vector<int>&
>(
m_conv.value())));
59 ATH_MSG_FATAL(
"initialize() Image size must be greater than 0");
65 ATH_MSG_FATAL(
"initialize() Hit extentsion list must have size % nLayers");
67 ATH_MSG_FATAL(
"initialize() Combine layers list must have size = nLayers");
71 ATH_MSG_FATAL(
"initialize() Bin scale list must have size = nLayers");
72 else if (std::any_of(
m_binScale.begin(),
m_binScale.end(), [&](
unsigned i){ return m_imageSize_y % i != 0; }))
73 ATH_MSG_FATAL(
"initialize() The imagesize is not divisible by scale");
76 if (!ok)
return StatusCode::FAILURE;
81 ATH_MSG_WARNING(
"initialize() localMaxWindowSize requires tracing hits, turning on automatically");
88 ATH_MSG_WARNING(
"initialize() idealGeoRoads conflicts with useSectors, switching off FPGATrackSim sector matching");
93 ATH_MSG_WARNING(
"initialize() idealGeoRoads requires tracing hits, turning on automatically");
145 m_HT_sel =
static_cast<int>(DBinPhi0_clean / DBinQApt_clean);
147 return StatusCode::SUCCESS;
163 std::vector<std::pair<unsigned, unsigned>> roadListXY;
164 std::unordered_set<std::shared_ptr<const FPGATrackSimHit>> merged_image;
174 roadListXY.push_back({
x,
y});
182 for (
const auto & roadXY : roadListXY){
183 ATH_MSG_DEBUG(
"roadList x : " << roadXY.first <<
", y : " << roadXY.second);
188 ATH_MSG_DEBUG(
"Trace hits is disabled. Road merge doesn't work with this state at the moment.");
192 for (
const auto & roadXY : roadListXY){
193 unsigned x = roadXY.first;
194 unsigned y = roadXY.second;
198 if(!roadListXY.empty()){
199 std::vector<std::pair<unsigned, unsigned>> mergedRoadListXY;
200 size_t roadCounter = 0;
212 while(!roadListXY.empty()){
213 mergedRoadListXY.push_back(roadListXY[0]);
214 roadListXY.erase(roadListXY.begin());
215 for(
size_t i = 0;
i < mergedRoadListXY.size();
i++){
216 for(
int j = 0; j < std::ssize(roadListXY); j++){
217 if(std::abs(
static_cast<int>(roadListXY[j].
first) -
static_cast<int>(mergedRoadListXY[
i].
first)) < 2 &&
218 std::abs(
static_cast<int>(roadListXY[j].
second) -
static_cast<int>(mergedRoadListXY[
i].
second)) < 2){
219 mergedRoadListXY.push_back(roadListXY[j]);
220 roadListXY.erase(roadListXY.begin() + j);
226 if (!mergedRoadListXY.empty()){
227 ATH_MSG_DEBUG( mergedRoadListXY.size() -1 <<
" roads are merged to road(" << mergedRoadListXY[0].first <<
"," << mergedRoadListXY[0].second <<
")");
229 for (
const auto & roadXY : mergedRoadListXY){
230 unsigned x = roadXY.first;
231 unsigned y = roadXY.second;
235 mergedRoadListXY.clear();
236 merged_image.clear();
239 ATH_MSG_DEBUG(
"There is/are " << roadCounter <<
" roads after road merge");
246 return StatusCode::SUCCESS;
253 for (
const auto& hit :
hits) {
259 switch (hit->getLayer()){
280 ATH_MSG_FATAL(
"Debug: something wrong! layer: " << hit->getLayer());
284 ATH_MSG_ERROR(
"FPGATrackSimHoughTransformTool: array index is negative");
287 int phi_L_bin =
m_h_rfix.at(LUT_layer)->FindBin(hit->getGPhi());
289 if(MSB >= 64) MSB = 63;
290 for(
const auto& LUT_i:
m_LUT.at(LUT_layer).at(MSB)){
291 if(LUT_i.input_begin <= phi_L_bin && phi_L_bin <= LUT_i.input_end){
292 for(
const auto&
pos: LUT_i.output){
299 double r = hit->getR();
301 int input_bins_vector_size = 0;
310 std::vector<std::vector<int>> hitLineQAPtPhi0(input_bins_vector_size);
313 for (
const std::vector<int>& it_base : hitLineQAPtPhi0) {
314 if (it_base[0] != -1 && it_base[1] != -1) {
315 image(it_base[0] , it_base[1]).first++;
317 image(it_base[0] , it_base[1]).second.insert(hit);
325 for (
unsigned y_ = 0; y_ < new_size_y; y_++) {
326 unsigned y_bin_min =
scale * y_;
327 unsigned y_bin_max =
scale * (y_ + 1);
333 for (
unsigned y = y_bin_min;
y < y_bin_max;
y++){
355 if (layerImage(
y,
x).
first > 0)
399 if (
i == 0 && j == 0)
continue;
418 static inline int quant(
double min,
double max,
unsigned nSteps,
double val)
424 static inline double unquant(
double min,
double max,
unsigned nSteps,
int step)
429 template <
typename T>
430 static inline std::string
to_string(
const std::vector<T> &
v)
432 std::ostringstream oss;
436 std::copy(
v.begin(),
v.end()-1, std::ostream_iterator<T>(oss,
", "));
449 double r = hit->
getR();
450 double phi_hit = hit->
getGPhi();
452 x = asin(
r * fpgatracksim::A *
y -
d0 /
r) + phi_hit;
458 ATH_MSG_ERROR(
"yToX() not defined for the current m_par selection");
484 std::vector<int64_t> bins_x_new;
485 std::vector<int64_t> bins_y_new;
498 std::vector<int64_t> qApt_ht_array;
499 std::vector<int64_t> phi0_ht_array;
500 std::vector<int>
zeros = {-1, -1};
502 std::vector<std::vector<int>> hitLineQAPtPhi0;
507 for (
unsigned i = 0;
i < bins_along_phi0_sector;
i++) {
508 qApt_ht_array.push_back(0);
510 for (
unsigned i = 0;
i < bins_along_qApt_sector;
i++) {
511 phi0_ht_array.push_back(0);
514 int64_t phi_conv_offseted = 0;
515 int64_t Delta_qApt_after_first_sector = 0;
516 int64_t qApt_ht = -1;
517 int64_t phi0_ht = -1;
518 int64_t qApt_ht_pre = -1;
519 int64_t phi0_ht_pre = -1;
520 int64_t qApt_ht_sector = 0;
521 int64_t phi0_ht_sector = 0;
522 int64_t Delta_phi0_after_first_sector = 0;
524 std::vector<double> d_qApt_ht_array;
525 std::vector<double> d_phi0_ht_array;
527 for (
unsigned i = 0;
i < bins_along_phi0_sector;
i++) {
528 d_qApt_ht_array.push_back(0.0);
530 for (
unsigned i = 0;
i < bins_along_qApt_sector;
i++) {
531 d_phi0_ht_array.push_back(0.0);
534 double d_phi_conv_offseted = 0.0;
535 double d_qApt_ht_sector = 0.0;
536 double d_phi0_ht_sector = 0.0;
537 double d_one_over_r_by_const = 0.0;
538 double d_phi_conv_over_DBinQApt = 0.0;
539 double d_Delta_phi0_after_first_sector = 0.0;
542 double r = hit->
getR();
543 double phi_hit = hit->
getGPhi();
548 hitLineQAPtPhi0.push_back(
zeros);
559 for (
unsigned j = start_phi0_bins_first_sector; j < end_phi0_bins_first_sector; j++) {
561 int64_t
c2 = (
static_cast<int64_t
>(d_c1) -
static_cast<int64_t
>(d_phi_conv_over_DBinQApt)) *
static_cast<int64_t
>(d_one_over_r_by_const) /
m_one_r_const_twoexp_post_conv;
562 int64_t
c3 = (
static_cast<int64_t
>(d_c1) -
static_cast<int64_t
>(d_phi_conv_over_DBinQApt)) *
static_cast<int64_t
>(d_one_over_r_by_const);
565 double d_c3 = (d_c1 - d_phi_conv_over_DBinQApt) * d_one_over_r_by_const;
568 if (
static_cast<int64_t
>(d_c1) <=
static_cast<int64_t
>(d_phi_conv_over_DBinQApt)) {
574 d_qApt_ht_sector = d_c3;
575 qApt_ht_array[j - start_phi0_bins_first_sector] = qApt_ht_sector;
576 d_qApt_ht_array[j - start_phi0_bins_first_sector] = d_qApt_ht_sector;
578 if (
static_cast<int64_t
>(d_c1) <=
static_cast<int64_t
>(d_phi_conv_over_DBinQApt) && (d_c2 -
m_qApt_min_post_conv - d_c4) >= 0) {
579 qApt_ht = qApt_ht_pre;
581 if (
static_cast<int64_t
>(d_c1) >
static_cast<int64_t
>(d_phi_conv_over_DBinQApt) && 1 + (d_c2 -
m_qApt_min_post_conv - d_c4) >= 0) {
582 qApt_ht = qApt_ht_pre;
586 hitLineQAPtPhi0[j] = {
static_cast<int>(qApt_ht),
static_cast<int>(j)};
595 Delta_phi0_after_first_sector = d_c6 *
static_cast<int64_t
>(d_one_over_r_by_const);
597 d_Delta_phi0_after_first_sector =
static_cast<double>(d_c6) * d_one_over_r_by_const;
598 for (
unsigned psi = start_phi0_bins_first_sector; psi < end_phi0_bins_first_sector; psi++) {
603 if (-qApt_ht_array[psi - start_phi0_bins_first_sector] < Delta_phi0_after_first_sector) {
611 if (qApt_ht >= 0 && qApt_ht <=
m_imageSize_y -1 && (1 + d_c7 - d_c9) >= 0 && (0 + d_c7 - d_c9) >= 0 ) {
622 hitLineQAPtPhi0.push_back(
zeros);
632 for (
unsigned n = start_qApt_bins_first_sector;
n < end_qApt_bins_first_sector;
n++) {
637 phi0_ht_sector = phi_conv_offseted + r_bit *
d2;
638 phi0_ht_array[
n - start_qApt_bins_first_sector] = phi0_ht_sector;
640 d_phi0_ht_sector = d_phi_conv_offseted + r_bit *
static_cast<double>(
d2);
641 d_phi0_ht_array[
n - start_qApt_bins_first_sector] = d_phi0_ht_sector;
643 if (phi0_ht_sector >= 0) {
650 phi0_ht = phi0_ht_pre;
657 if (phi0_ht <= m_imageSize_x - 1 && phi0_ht >= 0){
658 hitLineQAPtPhi0[
n] = {
static_cast<int>(
n),
static_cast<int>(phi0_ht)};
665 Delta_qApt_after_first_sector = r_bit * d3;
666 for (
unsigned tsi = start_qApt_bins_first_sector; tsi < end_qApt_bins_first_sector; tsi++) {
667 phi0_ht = (phi0_ht_array[tsi - start_qApt_bins_first_sector] + Delta_qApt_after_first_sector) /
m_bitwise_qApt_conv;
668 double d_phi0_ht = (d_phi0_ht_array[tsi - start_qApt_bins_first_sector] + Delta_qApt_after_first_sector) /
m_bitwise_qApt_conv;
680 ATH_MSG_ERROR(
"lineGenLay() not defined for the current m_par selection");
684 return hitLineQAPtPhi0;
694 if (x_min > x_max)
std::swap(x_min, x_max);
708 if (x_bin_min < 0) x_bin_min = 0;
711 return { x_bin_min, x_bin_max };
736 r.setHits( std::vector<std::vector<std::shared_ptr<const FPGATrackSimHit>>>(
hits));
746 r.setHitLayers(hitLayers);
755 for (
auto const & hit :
hits)
756 hitLayers |= 1 << hit->getLayer();
768 std::vector<std::shared_ptr<const FPGATrackSimHit>> road_hits;
770 for (
const auto & hit :
hits)
774 unsigned int y_bin_max = y_bin_min +
m_binScale[hit->getLayer()];
780 road_hits.push_back(hit);
781 hitLayers |= 1 << hit->getLayer();
793 const std::string filepath =
m_requirements.value() +
"requirement-" +
tag +
".root";
794 TFile*
fin = TFile::Open(filepath.c_str());
797 TTree* requirement =
fin->Get<TTree>(
"requirement");
799 for(
int i=0;
i<5;
i++){
800 v_LUT.at(
i).resize(64);
801 TH1D*
h = (TH1D*)
fin->Get(Form(
"in%d",
i));
804 int nrequirement = requirement->GetEntries();
810 requirement->SetBranchAddress(
"input_begin", &in_min);
811 requirement->SetBranchAddress(
"input_end", &in_max);
812 requirement->SetBranchAddress(
"output_x", &xi);
813 requirement->SetBranchAddress(
"output_y", &yi);
814 requirement->SetBranchAddress(
"output_l", &ri);
815 for(
int i=0;
i<nrequirement;
i++){
816 requirement->GetEntry(
i);
819 bool toBeFilled =
true;
820 if(!v_LUT.at(ri).at(MSB).empty()){
821 for(
auto& LUT_i : v_LUT.at(ri).at(MSB)){
822 if(LUT_i.input_begin == in_min && LUT_i.input_end == in_max){
823 LUT_i.output.push_back({xi, yi, ri});
834 LUT_i.
output.push_back({xi, yi, ri});
835 v_LUT.at(ri).at(MSB).push_back(std::move(LUT_i));
842 if(!v_LUT.at(ri).at(MSB).empty()){
843 for(
auto& LUT_i : v_LUT.at(ri).at(MSB)){
844 if(LUT_i.input_begin == in_min && LUT_i.input_end == in_max){
845 LUT_i.output.push_back({xi, yi, ri});
856 LUT_i.
output.push_back({xi, yi, ri});
857 v_LUT.at(ri).at(MSB).push_back(LUT_i);
861 int LUTsize = v_LUT.size();
863 for(
int msb = 0; msb < 64; msb++){