Skip to content

Write t3_isFake to the debug output #57

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
bucket420 wants to merge 2 commits into from
Closed

Write t3_isFake to the debug output #57

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
05289f8
Address comments
bucket420 Jul 8, 2024
7ff7999
Fix styling
bucket420 Jul 8, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
150 changes: 102 additions & 48 deletions RecoTracker/LSTCore/standalone/code/core/trkCore.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -274,45 +274,50 @@ float runTrackCandidate(SDL::Event<SDL::Acc> *event, bool no_pls_dupclean, bool
// ---------------------------------- =========================================== ----------------------------------------------

//___________________________________________________________________________________________________________________________________________________________________________________________
std::vector<int> matchedSimTrkIdxs(std::vector<int> hitidxs, std::vector<int> hittypes, bool verbose) {
std::vector<int> matchedSimTrkIdxs(std::vector<int> hitidxs, std::vector<int> hittypes, bool verbose, float matchfrac) {
std::vector<unsigned int> hitidxs_(std::begin(hitidxs), std::end(hitidxs));
std::vector<unsigned int> hittypes_(std::begin(hittypes), std::end(hittypes));
return matchedSimTrkIdxs(hitidxs_, hittypes_, verbose);
return matchedSimTrkIdxs(hitidxs_, hittypes_, verbose, matchfrac);
}


//___________________________________________________________________________________________________________________________________________________________________________________________
std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs, std::vector<unsigned int> hittypes, bool verbose, float matchfrac) {
std::vector<int> matched_idxs;
std::vector<float> matched_idx_fracs;
std::tie(matched_idxs, matched_idx_fracs) = matchedSimTrkIdxsAndFracs(hitidxs, hittypes, verbose, matchfrac);
return matched_idxs;
}

//___________________________________________________________________________________________________________________________________________________________________________________________
std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
std::vector<unsigned int> hittypes,
bool verbose) {
std::tuple<std::vector<int>, std::vector<float>> matchedSimTrkIdxsAndFracs(std::vector<unsigned int> hitidxs, std::vector<unsigned int> hittypes, bool verbose, float matchfrac) {
if (hitidxs.size() != hittypes.size()) {
std::cout << "Error: matched_sim_trk_idxs() hitidxs and hittypes have different lengths" << std::endl;
std::cout << "hitidxs.size(): " << hitidxs.size() << std::endl;
std::cout << "hittypes.size(): " << hittypes.size() << std::endl;
}

std::vector<std::pair<unsigned int, unsigned int>> to_check_duplicate;
for (size_t i = 0; i < hitidxs.size(); ++i) {
auto hitidx = hitidxs[i];
auto hittype = hittypes[i];
auto item = std::make_pair(hitidx, hittype);
for (int i = 0; i < hitidxs.size(); i++) {
auto item = std::make_pair(hitidxs[i], hittypes[i]);
if (std::find(to_check_duplicate.begin(), to_check_duplicate.end(), item) == to_check_duplicate.end()) {
to_check_duplicate.push_back(item);
to_check_duplicate.push_back(item);
}
}

int nhits_input = to_check_duplicate.size();

std::vector<vector<int>> simtrk_idxs;
std::vector<int> unique_idxs; // to aggregate which ones to count and test
std::vector<int> unique_idxs; // to aggregate which ones to count and test

if (verbose) {
std::cout << " '------------------------': "
<< "------------------------" << std::endl;
}

for (size_t ihit = 0; ihit < to_check_duplicate.size(); ++ihit) {
auto ihitdata = to_check_duplicate[ihit];
auto &&[hitidx, hittype] = ihitdata;
for (int i = 0; i < nhits_input; i++) {
auto hitidx = to_check_duplicate[i].first;
auto hittype = to_check_duplicate[i].second;

if (verbose) {
std::cout << " hitidx: " << hitidx << " hittype: " << hittype << std::endl;
Expand All @@ -323,8 +328,8 @@ std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
const std::vector<vector<int>> *simHitIdxs = hittype == 4 ? &trk.ph2_simHitIdx() : &trk.pix_simHitIdx();

if (verbose) {
std::cout << " trk.ph2_simHitIdx().size(): " << trk.ph2_simHitIdx().size() << std::endl;
std::cout << " trk.pix_simHitIdx().size(): " << trk.pix_simHitIdx().size() << std::endl;
std::cout << " trk.ph2_simHitIdx().size(): " << trk.ph2_simHitIdx().size() << std::endl;
std::cout << " trk.pix_simHitIdx().size(): " << trk.pix_simHitIdx().size() << std::endl;
}

if (static_cast<const unsigned int>((*simHitIdxs).size()) <= hitidx) {
Expand All @@ -344,23 +349,19 @@ std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
}
int simtrk_idx = trk.simhit_simTrkIdx().at(simhit_idx);
if (verbose) {
std::cout << " hitidx: " << hitidx << " simhit_idx: " << simhit_idx << " simtrk_idx: " << simtrk_idx
<< std::endl;
std::cout << " hitidx: " << hitidx << " simhit_idx: " << simhit_idx << " simtrk_idx: " << simtrk_idx << std::endl;
}
simtrk_idxs_per_hit.push_back(simtrk_idx);
if (std::find(unique_idxs.begin(), unique_idxs.end(), simtrk_idx) == unique_idxs.end())
unique_idxs.push_back(simtrk_idx);
if (std::find(unique_idxs.begin(), unique_idxs.end(), simtrk_idx) == unique_idxs.end()) unique_idxs.push_back(simtrk_idx);
}

if (simtrk_idxs_per_hit.size() == 0) {
if (verbose) {
std::cout << " hitidx: " << hitidx << " -1: " << -1 << std::endl;
}
simtrk_idxs_per_hit.push_back(-1);
if (std::find(unique_idxs.begin(), unique_idxs.end(), -1) == unique_idxs.end())
unique_idxs.push_back(-1);
if (std::find(unique_idxs.begin(), unique_idxs.end(), -1) == unique_idxs.end()) unique_idxs.push_back(-1);
}

simtrk_idxs.push_back(simtrk_idxs_per_hit);
}

Expand All @@ -384,21 +385,21 @@ std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
}

// Compute all permutations
std::function<void(vector<vector<int>> &, vector<int>, size_t, vector<vector<int>> &)> perm =
[&](vector<vector<int>> &result, vector<int> intermediate, size_t n, vector<vector<int>> &va) {
// std::cout << " 'called': " << "called" << std::endl;
if (va.size() > n) {
for (auto x : va[n]) {
// std::cout << " n: " << n << std::endl;
// std::cout << " intermediate.size(): " << intermediate.size() << std::endl;
std::vector<int> copy_intermediate(intermediate);
copy_intermediate.push_back(x);
perm(result, copy_intermediate, n + 1, va);
}
} else {
result.push_back(intermediate);
}
};
std::function<void(vector<vector<int>> &, vector<int>, size_t, vector<vector<int>> &)> perm = [&](vector<vector<int>> &result, vector<int> intermediate, size_t n, vector<vector<int>> &va) {
// std::cout << " 'called': " << "called" << std::endl;
if (va.size() > n) {
for (auto x : va[n]) {
// std::cout << " n: " << n << std::endl;
// std::cout << " intermediate.size(): " << intermediate.size() << std::endl;
std::vector<int> copy_intermediate(intermediate);
copy_intermediate.push_back(x);
perm(result, copy_intermediate, n + 1, va);
}
}
else {
result.push_back(intermediate);
}
};

vector<vector<int>> allperms;
perm(allperms, vector<int>(), 0, simtrk_idxs);
Expand All @@ -408,12 +409,13 @@ std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
for (unsigned iperm = 0; iperm < allperms.size(); ++iperm) {
std::cout << " allperms[iperm].size(): " << allperms[iperm].size() << std::endl;
for (unsigned ielem = 0; ielem < allperms[iperm].size(); ++ielem) {
std::cout << " allperms[iperm][ielem]: " << allperms[iperm][ielem] << std::endl;
std::cout << " allperms[iperm][ielem]: " << allperms[iperm][ielem] << std::endl;
}
}
}
int maxHitMatchCount = 0; // ultimate maximum of the number of matched hits
int maxHitMatchCount = 0; // ultimate maximum of the number of matched hits
std::vector<int> matched_sim_trk_idxs;
std::vector<float> matched_sim_trk_idxs_frac;
for (auto &trkidx_perm : allperms) {
std::vector<int> counts;
for (auto &unique_idx : unique_idxs) {
Expand All @@ -423,18 +425,70 @@ std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
auto result = std::max_element(counts.begin(), counts.end());
int rawidx = std::distance(counts.begin(), result);
int trkidx = unique_idxs[rawidx];
if (trkidx < 0)
continue;
if (counts[rawidx] > (((float)nhits_input) * 0.75))
if (trkidx < 0) continue;
// std::cout << " counts[rawidx]: " << counts[rawidx] << std::endl;
// std::cout << " nhits_input: " << nhits_input << std::endl;
if (verbose) {
float fr = ((float)counts[rawidx]) / ((float)nhits_input);
std::cout << " fr: " << fr << std::endl;
}
// std::cout << " matchfrac: " << matchfrac << std::endl;
if (counts[rawidx] > (((float)nhits_input) * matchfrac)) {
matched_sim_trk_idxs.push_back(trkidx);
matched_sim_trk_idxs_frac.push_back( ((float)counts[rawidx]) / ((float)nhits_input) );
}
maxHitMatchCount = std::max(maxHitMatchCount, *std::max_element(counts.begin(), counts.end()));
}
set<int> s;
map<int, float> pairs;
unsigned size = matched_sim_trk_idxs.size();
for (unsigned i = 0; i < size; ++i)
s.insert(matched_sim_trk_idxs[i]);
matched_sim_trk_idxs.assign(s.begin(), s.end());
return matched_sim_trk_idxs;
for (unsigned i = 0; i < size; ++i) {
int idx = matched_sim_trk_idxs[i];
float frac = matched_sim_trk_idxs_frac[i];
if (pairs.find(idx) != pairs.end()) {
if (pairs[idx] < frac)
pairs[idx] = frac;
}
else {
pairs[idx] = frac;
}
// if (std::find(result.begin(), result.end(), idx) == result.end())
Copy link

@YonsiG YonsiG Jul 9, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same

// {
// result.push_back(idx);
// result_frac.push_back(frac);
// }
}

vector<int> result;
vector<float> result_frac;

// Loop over the map using range-based for loop
for (const auto& pair : pairs) {
result.push_back(pair.first);
result_frac.push_back(pair.second);
}

// Sort indices based on 'values'
// auto indices = sort_indices(result_frac);

// // Reorder 'vec1' and 'vec2' based on the sorted indices
Copy link

@YonsiG YonsiG Jul 9, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

still, why are those commented lines needed? If those are important, please uncomment them. Otherwise please drop them

Copy link
Author

@bucket420 bucket420 Jul 10, 2024

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These lines sort the sim tracks by match fraction. I think that's not needed for our use case but correct me if I'm wrong

// std::vector<int> sorted_result(result.size());
// std::vector<float> sorted_result_frac(result_frac.size());
// for (size_t i = 0; i < indices.size(); ++i)
// {
// sorted_result[i] = result[indices[i]];
// sorted_result_frac[i] = result_frac[indices[i]];
// }

// for (size_t i = 0; i < sorted_result.size(); ++i)
// {
// int sorted_idx = sorted_result.at(i);
// float sorted_frac = sorted_result_frac.at(i);
// std::cout << " i: " << i << " sorted_idx: " << sorted_idx << " sorted_frac: " << sorted_frac << " hitidxs.size(): " << hitidxs.size() << std::endl;
// }

return std::make_tuple(result, result_frac);
// matched_sim_trk_idxs.assign(s.begin(), s.end());
// return matched_sim_trk_idxs;
}

//__________________________________________________________________________________________
Expand Down
7 changes: 3 additions & 4 deletions RecoTracker/LSTCore/standalone/code/core/trkCore.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -28,10 +28,9 @@ float runpT3(SDL::Event<SDL::Acc> *event);

// --------------------- ======================== ---------------------

std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs,
std::vector<unsigned int> hittypes,
bool verbose = false);
std::vector<int> matchedSimTrkIdxs(std::vector<int> hitidxs, std::vector<int> hittypes, bool verbose = false);
std::vector<int> matchedSimTrkIdxs(std::vector<unsigned int> hitidxs, std::vector<unsigned int> hittypes, bool verbose=false, float matchfrac=0.75);
std::vector<int> matchedSimTrkIdxs(std::vector<int> hitidxs, std::vector<int> hittypes, bool verbose=false, float matchfrac=0.75);
std::tuple<std::vector<int>, std::vector<float>> matchedSimTrkIdxsAndFracs(std::vector<unsigned int> hitidxs, std::vector<unsigned int> hittypes, bool verbose=false, float matchfrac=0.75);
int getDenomSimTrkType(int isimtrk);
int getDenomSimTrkType(std::vector<int> simidxs);

Expand Down
33 changes: 20 additions & 13 deletions RecoTracker/LSTCore/standalone/code/core/write_sdl_ntuple.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -133,7 +133,8 @@ void createOptionalOutputBranches() {
ana.tx->createBranch<vector<float>>("t5_nonAnchorChiSquared");

// T3 branches
ana.tx->createBranch<vector<float>>("T3_residual");
ana.tx->createBranch<vector<float>>("t3_residual");
ana.tx->createBranch<vector<int>>("t3_isFake");

// Occupancy branches
ana.tx->createBranch<vector<int>>("module_layers");
Expand Down Expand Up @@ -520,19 +521,25 @@ void setQuintupletOutputBranches(SDL::Event<SDL::Acc>* event) {

//________________________________________________________________________________________________________________________________
void setTripletOutputBranches(SDL::Event<SDL::Acc>* event) {
SDL::tripletsBuffer<alpaka::DevCpu>& tripletsInGPU = (*event->getTriplets());
SDL::objectRangesBuffer<alpaka::DevCpu>& rangesInGPU = (*event->getRanges());
SDL::modulesBuffer<alpaka::DevCpu>& modulesInGPU = (*event->getModules());
for (unsigned int lowerModuleIdx = 0; lowerModuleIdx < *(modulesInGPU.nLowerModules); ++lowerModuleIdx)
{
unsigned int nTriplets = tripletsInGPU.nTriplets[lowerModuleIdx];
for (unsigned int idx = 0; idx < nTriplets; idx++)
{
unsigned int tripletIndex = rangesInGPU.tripletModuleIndices[lowerModuleIdx] + idx;
const float residual = tripletsInGPU.residual[tripletIndex];
ana.tx->pushbackToBranch<float>("T3_residual", residual);
}
SDL::tripletsBuffer<alpaka::DevCpu>& tripletsInGPU = (*event->getTriplets());
SDL::objectRangesBuffer<alpaka::DevCpu>& rangesInGPU = (*event->getRanges());
SDL::modulesBuffer<alpaka::DevCpu>& modulesInGPU = (*event->getModules());
for (unsigned int lowerModuleIdx = 0; lowerModuleIdx < *(modulesInGPU.nLowerModules); ++lowerModuleIdx) {
unsigned int nTriplets = tripletsInGPU.nTriplets[lowerModuleIdx];
for (unsigned int idx = 0; idx < nTriplets; idx++) {
unsigned int tripletIndex = rangesInGPU.tripletModuleIndices[lowerModuleIdx] + idx;

const float residual = tripletsInGPU.residual[tripletIndex];
ana.tx->pushbackToBranch<float>("t3_residual", residual);

std::tuple<std::vector<unsigned int>, std::vector<unsigned int>> hitIdxsAndTypes = getHitIdxsAndHitTypesFromT3(event, tripletIndex);
std::vector<unsigned int> hit_idxs = std::get<0>(hitIdxsAndTypes);
std::vector<unsigned int> hit_types = std::get<1>(hitIdxsAndTypes);

std::tuple<std::vector<int>, std::vector<float>> simIdxsAndFracs = matchedSimTrkIdxsAndFracs(hit_idxs, hit_types);
ana.tx->pushbackToBranch<int>("t3_isFake", static_cast<int>(std::get<0>(simIdxsAndFracs).size() == 0));
}
}
}

//________________________________________________________________________________________________________________________________
Expand Down