Skip to content
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

Candidate ranking à la CMSSW, modified to maximize efficiency and minimize fake rate #173

Closed
wants to merge 5 commits into from
Closed

Candidate ranking à la CMSSW, modified to maximize efficiency and minimize fake rate #173

Show file tree
Hide file tree
Changes from 4 commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
d11e4e9
cmsswranking
mmasciov Sep 20, 2018
f4741d6
Fixing indentation
mmasciov Oct 2, 2018
471360e
Avoiding useless multiplications in initializations
mmasciov Oct 2, 2018
0afe9dd
Modified cmssw ranking to maximize efficiency and minimize fake rate
mmasciov Oct 11, 2018
03ab24e
Seed based candidate ranking
mmasciov Oct 31, 2018
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
4 changes: 4 additions & 0 deletions Config.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -344,6 +344,10 @@ namespace Config
constexpr float treg_eta[2] = {0.45,1.5};
constexpr float track_ptlow = 0.9;

// sorting config (bonus,penalty)
constexpr float validHitBonus_ = 2.5;
constexpr float missingHitPenalty_ = 20.0;

// Threading
extern int numThreadsFinder;
extern int numThreadsSimulation;
Expand Down
3 changes: 2 additions & 1 deletion TTreeValidation.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -880,7 +880,8 @@ void TTreeValidation::mapRefTkToRecoTks(const TrackVec& evt_tracks, TrackExtraVe
{
tmpMatches.emplace_back(evt_tracks[label]);
}
std::sort(tmpMatches.begin(), tmpMatches.end(), sortByHitsChi2); // sort the tracks
//std::sort(tmpMatches.begin(), tmpMatches.end(), sortByHitsChi2); // sort the tracks
std::sort(tmpMatches.begin(), tmpMatches.end(), sortByScoreCand); // sort the tracks
for (auto itrack = 0; itrack < (int) tmpMatches.size(); itrack++) // loop over sorted tracks, now set the vector of sorted labels match
{
refTkMatches.second[itrack] = tmpMatches[itrack].label();
Expand Down
67 changes: 67 additions & 0 deletions Track.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,18 @@ inline int calculateCharge(const float hit0_x, const float hit0_y,
return ((hit2_y-hit0_y)*(hit2_x-hit1_x)>(hit2_y-hit1_y)*(hit2_x-hit0_x)?1:-1);
}

struct IdxChi2List
{
public:
int trkIdx; // candidate index
int hitIdx; // hit index
int nhits; // number of hits (used for sorting)
int nholes; // number of holes (used for sorting)
float pt; // pt (used for sorting)
float eta; // eta (used for sorting)
float chi2; // total chi2 (used for sorting)
};

struct ReducedTrack // used for cmssw reco track validation
{
public:
Expand Down Expand Up @@ -467,6 +479,61 @@ inline bool sortByHitsChi2(const Track & cand1, const Track & cand2)
return cand1.nFoundHits()>cand2.nFoundHits();
}

inline bool sortByScoreLoop(const int nfoundhits[2],
const int nmisshits[2],
const float chi2[2],
const float pt[2],
const float eta[2])
{
float score[2] = {0.f,0.f};
for(int c=0; c<2; ++c){
// For high pT central tracks: double valid hit bonus
if((pt[0]+pt[1])/2.0f > 2.0f && (eta[0]+eta[1])/2.0f < 1.5f){
Copy link
Collaborator

Choose a reason for hiding this comment

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

is there a good reason for the score of one candidate to depend on the other candidate?

Copy link
Collaborator

Choose a reason for hiding this comment

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

Despite the nice improvement in the efficiency in favoring this over PR #167 , I agree with Slava. This is a bit strange: why not just define a function for this section, calling it once for cand1 and again for cand2? That way you can evaluate the score bonuses based on pT/eta for each candidate independently?

Copy link
Collaborator

Choose a reason for hiding this comment

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

with this coupled dependency in the score definition, I'm trying to think if it's obvious that for A, B, C

  • if (score_A(A,B) > score_B(A,B)) and (score_B(B,C) > score_C(B,C))
  • that score_A(A,C) will be > score_C(A,C).

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

I agree that this is not obvious.
However, my choice was based on observation, and in order to avoid "threshold" effects that may favor one candidate instead of another.
If I do the simple way (which was my first attempt), I get normally lower efficiencies:
eff_averagevssingle_pt0
ratio_singleoaverage_pt0
eff_averagevssingle_pt0p9
ratio_singleoaverage_pt0p9
eff_averagevssingle_pt2p0
ratio_singleoaverage_pt2p0

Copy link
Collaborator Author

Choose a reason for hiding this comment

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

So, despite it's not obvious as Slava points out that if A wins over B, and B wins over C, then A wins over C, I simply went for the best option in terms of efficiency.

score[c] = (Config::validHitBonus_*2.0f)*nfoundhits[c] - Config::missingHitPenalty_*nmisshits[c] - chi2[c];
if(pt[c]<0.9f) score[c] -= 0.5f*(Config::validHitBonus_*2.0f)*nfoundhits[c];
else if(nfoundhits[c]>8) score[c] += (Config::validHitBonus_*2.0f)*nfoundhits[c];
}
// For low pT endcap tracks: half valid hit bonus & half missing hit penalty
else if((pt[0]+pt[1])/2.0f < 0.9f && (eta[0]+eta[1])/2.0f >= 1.5f){
score[c] = (Config::validHitBonus_*0.5f)*nfoundhits[c] - (Config::missingHitPenalty_*0.5f)*nmisshits[c] - chi2[c];
if(pt[c]<0.9f) score[c] -= 0.5f*(Config::validHitBonus_*0.5f)*nfoundhits[c];
else if(nfoundhits[c]>8) score[c] += (Config::validHitBonus_*0.5f)*nfoundhits[c];
}
// For all other tracks: unchanged cmssw bonus and penalty
else{
score[c] = Config::validHitBonus_*nfoundhits[c] - Config::missingHitPenalty_*nmisshits[c] - chi2[c];
if(pt[c]<0.9f) score[c] -= 0.5f*Config::validHitBonus_*nfoundhits[c];
else if(nfoundhits[c]>8) score[c] += Config::validHitBonus_*nfoundhits[c];
}
}
return score[0]>score[1];
}

inline bool sortByScoreCand(const Track& cand1, const Track& cand2)
{
int nfoundhits[2] = {cand1.nFoundHits(),cand2.nFoundHits()};
int nmisshits[2] = {cand1.nTotalHits()-cand1.nFoundHits(),cand2.nTotalHits()-cand2.nFoundHits()};
float chi2[2] = {cand1.chi2(),cand2.chi2()};
float pt[2] = {cand1.pT(),cand2.pT()};
float eta[2] = {std::fabs(cand1.momEta()),std::fabs(cand2.momEta())};
return sortByScoreLoop(nfoundhits,nmisshits,chi2,pt,eta);
}

inline bool sortByScoreStruct(const IdxChi2List& cand1, const IdxChi2List& cand2)
{
int nfoundhits[2] = {cand1.nhits,cand2.nhits};
int nmisshits[2] = {cand1.nholes,cand2.nholes};
float chi2[2] = {cand1.chi2,cand2.chi2};
float pt[2] = {cand1.pt,cand2.pt};
float eta[2] = {cand1.eta,cand2.eta};
return sortByScoreLoop(nfoundhits,nmisshits,chi2,pt,eta);
}

inline bool sortByScoreCandPair(const std::pair<Track, TrackState>& cand1, const std::pair<Track, TrackState>& cand2)
{
return sortByScoreCand(cand1.first,cand2.first);
}

template <typename Vector>
inline void squashPhiGeneral(Vector& v)
{
Expand Down
12 changes: 8 additions & 4 deletions buildtest.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -50,12 +50,14 @@ void processCandidates(const BinInfoMap & segmentMap, Event& ev, candvec& candid
}
if (tmp_candidates.size()>static_cast<size_t>(Config::maxCandsPerSeed)) {
dprint("huge size=" << tmp_candidates.size() << " keeping best "<< Config::maxCandsPerSeed << " only");
std::partial_sort(tmp_candidates.begin(),tmp_candidates.begin()+(Config::maxCandsPerSeed-1),tmp_candidates.end(),sortByHitsChi2);
//std::partial_sort(tmp_candidates.begin(),tmp_candidates.begin()+(Config::maxCandsPerSeed-1),tmp_candidates.end(),sortByHitsChi2);
std::partial_sort(tmp_candidates.begin(),tmp_candidates.begin()+(Config::maxCandsPerSeed-1),tmp_candidates.end(),sortByScoreCand);
tmp_candidates.resize(Config::maxCandsPerSeed); // thread local, so ok not thread safe
} else if (tmp_candidates.size()==0) {
// save the best candidate from the previous iteration and then swap in
// the empty new candidate list; seed will be skipped on future iterations
auto best = std::min_element(candidates.begin(),candidates.end(),sortByHitsChi2);
//auto best = std::min_element(candidates.begin(),candidates.end(),sortByHitsChi2);
auto best = std::min_element(candidates.begin(),candidates.end(),sortByScoreCand);
#ifdef TBB
std::lock_guard<std::mutex> evtguard(evtlock); // should be rare
#endif
Expand Down Expand Up @@ -113,7 +115,8 @@ void buildTracksBySeeds(const BinInfoMap & segmentMap, Event& ev)
if (cand.size()>0) {
// only save one track candidate per seed, one with lowest chi2
//std::partial_sort(cand.begin(),cand.begin()+1,cand.end(),sortByHitsChi2);
auto best = std::min_element(cand.begin(),cand.end(),sortByHitsChi2);
//auto best = std::min_element(cand.begin(),cand.end(),sortByHitsChi2);
auto best = std::min_element(cand.begin(),cand.end(),sortByScoreCand);
best->setLabel(evt_track_candidates.size());
dprint("Saving track from seed " << iseed << " label " << best->label() << " hits "
<< best->nFoundHits() << " parameters " << best->parameters() << std::endl);
Expand Down Expand Up @@ -166,7 +169,8 @@ void buildTracksByLayers(const BinInfoMap & segmentMap, Event& ev)
if (cand.size()>0) {
// only save one track candidate per seed, one with lowest chi2
//std::partial_sort(cand.begin(),cand.begin()+1,cand.end(),sortByHitsChi2);
auto best = std::min_element(cand.begin(),cand.end(),sortByHitsChi2);
//auto best = std::min_element(cand.begin(),cand.end(),sortByHitsChi2);
auto best = std::min_element(cand.begin(),cand.end(),sortByScoreCand);
best->setLabel(evt_track_candidates.size());
evt_track_candidates.push_back(*best);
ev.candidateTracksExtra_.emplace_back(evt_seeds_extra[iseed].seedID());
Expand Down
17 changes: 12 additions & 5 deletions mkFit/CandCloner.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -7,12 +7,18 @@

namespace
{
inline bool sortCandListByHitsChi2(const mkfit::MkFinder::IdxChi2List& cand1,
const mkfit::MkFinder::IdxChi2List& cand2)
inline bool sortCandListByHitsChi2(const mkfit::IdxChi2List& cand1,
const mkfit::IdxChi2List& cand2)
{
if (cand1.nhits == cand2.nhits) return cand1.chi2 < cand2.chi2;
return cand1.nhits > cand2.nhits;
}

inline bool sortCandListByScore(const mkfit::IdxChi2List& cand1,
const mkfit::IdxChi2List& cand2)
{
return mkfit::sortByScoreStruct(cand1,cand2);
}
}

namespace mkfit {
Expand All @@ -30,7 +36,7 @@ void CandCloner::ProcessSeedRange(int is_beg, int is_end)
//1) sort the candidates
for (int is = is_beg; is < is_end; ++is)
{
std::vector<MkFinder::IdxChi2List>& hitsForSeed = m_hits_to_add[is];
std::vector<IdxChi2List>& hitsForSeed = m_hits_to_add[is];

std::vector<CombCandidate> &cands = mp_event_of_comb_candidates->m_candidates;

Expand All @@ -52,7 +58,8 @@ void CandCloner::ProcessSeedRange(int is_beg, int is_end)
if ( ! hitsForSeed.empty())
{
//sort the new hits
std::sort(hitsForSeed.begin(), hitsForSeed.end(), sortCandListByHitsChi2);
//std::sort(hitsForSeed.begin(), hitsForSeed.end(), sortCandListByHitsChi2);
std::sort(hitsForSeed.begin(), hitsForSeed.end(), sortCandListByScore);

int num_hits = std::min((int) hitsForSeed.size(), Config::maxCandsPerSeed);

Expand All @@ -61,7 +68,7 @@ void CandCloner::ProcessSeedRange(int is_beg, int is_end)

for (int ih = 0; ih < num_hits; ih++)
{
const MkFinder::IdxChi2List &h2a = hitsForSeed[ih];
const IdxChi2List &h2a = hitsForSeed[ih];

cv.push_back( cands[ m_start_seed + is ][ h2a.trkIdx ] );

Expand Down
4 changes: 2 additions & 2 deletions mkFit/CandCloner.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -81,7 +81,7 @@ class CandCloner
// Do nothing, "secondary" state vars updated when work completed/assigned.
}

void add_cand(int idx, const MkFinder::IdxChi2List& cand_info)
void add_cand(int idx, const IdxChi2List& cand_info)
{
m_hits_to_add[idx].push_back(cand_info);

Expand Down Expand Up @@ -166,7 +166,7 @@ class CandCloner
// eventually, protected or private

int m_idx_max, m_idx_max_prev;
std::vector<std::vector<MkFinder::IdxChi2List>> m_hits_to_add;
std::vector<std::vector<IdxChi2List>> m_hits_to_add;

EventOfCombCandidates *mp_event_of_comb_candidates;
std::vector<std::pair<int,int>> *mp_kalman_update_list;
Expand Down
20 changes: 15 additions & 5 deletions mkFit/MkBuilder.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -165,6 +165,11 @@ namespace

return cand1.nFoundHits() > cand2.nFoundHits();
}

bool sortCandByScore(const Track & cand1, const Track & cand2)
{
return mkfit::sortByScoreCand(cand1,cand2);
}
}

//------------------------------------------------------------------------------
Expand Down Expand Up @@ -1747,7 +1752,8 @@ void MkBuilder::FindTracksStandard()
for (int is = 0; is < n_seeds; ++is)
{
dprint("dump seed n " << is << " with input candidates=" << tmp_cands[is].size());
std::sort(tmp_cands[is].begin(), tmp_cands[is].end(), sortCandByHitsChi2);
//std::sort(tmp_cands[is].begin(), tmp_cands[is].end(), sortCandByHitsChi2);
std::sort(tmp_cands[is].begin(), tmp_cands[is].end(), sortCandByScore);

if (tmp_cands[is].size() > static_cast<size_t>(Config::maxCandsPerSeed))
{
Expand Down Expand Up @@ -1799,7 +1805,8 @@ void MkBuilder::FindTracksStandard()
{
std::vector<Track>& finalcands = eoccs[iseed];
if (finalcands.size() == 0) continue;
std::sort(finalcands.begin(), finalcands.end(), sortCandByHitsChi2);
//std::sort(finalcands.begin(), finalcands.end(), sortCandByHitsChi2);
std::sort(finalcands.begin(), finalcands.end(), sortCandByScore);
}
}); // end parallel-for over chunk of seeds within region
}); // end of parallel-for-each over eta regions
Expand Down Expand Up @@ -1976,7 +1983,8 @@ void MkBuilder::find_tracks_in_layers(CandCloner &cloner, MkFinder *mkfndr,
dprintf("Extra layer %d start_seed %d, is %d, N=%d (orig size = %d) -- ",
curr_layer, start_seed, is, (int) extra_cands[is].size(), (int) eoccs.m_candidates[start_seed + is].size());

std::sort(extra_cands[is].begin(), extra_cands[is].end(), sortCandByHitsChi2);
//std::sort(extra_cands[is].begin(), extra_cands[is].end(), sortCandByHitsChi2);
std::sort(extra_cands[is].begin(), extra_cands[is].end(), sortCandByScore);

auto src_i = extra_cands[is].begin();
auto src_e = extra_cands[is].end();
Expand All @@ -1995,7 +2003,8 @@ void MkBuilder::find_tracks_in_layers(CandCloner &cloner, MkFinder *mkfndr,
int cnt = 0;
while (dest_i != dest_e && src_i != src_e)
{
while (dest_i != dest_e && sortCandByHitsChi2(*dest_i, *src_i)) ++dest_i;
//while (dest_i != dest_e && sortCandByHitsChi2(*dest_i, *src_i)) ++dest_i;
while (dest_i != dest_e && sortCandByScore(*dest_i, *src_i)) ++dest_i;

if (dest_i != dest_e)
{
Expand All @@ -2021,7 +2030,8 @@ void MkBuilder::find_tracks_in_layers(CandCloner &cloner, MkFinder *mkfndr,
{
std::vector<Track>& finalcands = eoccs.m_candidates[iseed];
if (finalcands.size() == 0) continue;
std::sort(finalcands.begin(), finalcands.end(), sortCandByHitsChi2);
//std::sort(finalcands.begin(), finalcands.end(), sortCandByHitsChi2);
std::sort(finalcands.begin(), finalcands.end(), sortCandByScore);
}
}

Expand Down
8 changes: 7 additions & 1 deletion mkFit/MkFinder.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -71,7 +71,7 @@ void MkFinder::InputTracksAndHitIdx(const std::vector<CombCandidate> & track
}

void MkFinder::InputTracksAndHitIdx(const std::vector<CombCandidate> & tracks,
const std::vector<std::pair<int,MkFinder::IdxChi2List>>& idxs,
const std::vector<std::pair<int,IdxChi2List>>& idxs,
int beg, int end, bool inputProp)
{
// Assign track parameters to initial state and copy hit values in.
Expand Down Expand Up @@ -901,6 +901,9 @@ void MkFinder::FindCandidatesCloneEngine(const LayerOfHits &layer_of_hits, CandC
tmpList.trkIdx = CandIdx(itrack, 0, 0);
tmpList.hitIdx = XHitArr.At(itrack, hit_cnt, 0);
tmpList.nhits = NFoundHits(itrack,0,0) + 1;
tmpList.nholes = num_invalid_hits(itrack);
tmpList.pt = std::abs(1.0f/Par[iP].At(itrack,3,0));
tmpList.eta = std::fabs( getEta( Par[iP].At(itrack,5,0) ) );
tmpList.chi2 = Chi2(itrack, 0, 0) + chi2;
cloner.add_cand(SeedIdx(itrack, 0, 0) - offset, tmpList);
// hitsToAdd[SeedIdx(itrack, 0, 0)-offset].push_back(tmpList);
Expand Down Expand Up @@ -934,6 +937,9 @@ void MkFinder::FindCandidatesCloneEngine(const LayerOfHits &layer_of_hits, CandC
tmpList.trkIdx = CandIdx(itrack, 0, 0);
tmpList.hitIdx = fake_hit_idx;
tmpList.nhits = NFoundHits(itrack,0,0);
tmpList.nholes = num_invalid_hits(itrack)+1;
tmpList.pt = std::abs(1.0f/Par[iP].At(itrack,3,0));
tmpList.eta = std::fabs( getEta( Par[iP].At(itrack,5,0) ) );
tmpList.chi2 = Chi2(itrack, 0, 0);
cloner.add_cand(SeedIdx(itrack, 0, 0) - offset, tmpList);
dprint("adding invalid hit " << fake_hit_idx);
Expand Down
8 changes: 0 additions & 8 deletions mkFit/MkFinder.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -54,14 +54,6 @@ class MkFinder : public MkBase
using MPlexHitIdx = Matriplex::Matriplex<int, MPlexHitIdxMax, 1, NN>;
using MPlexQHoT = Matriplex::Matriplex<HitOnTrack, 1, 1, NN>;

struct IdxChi2List
{
int trkIdx; // candidate index
int hitIdx; // hit index
int nhits; // number of hits (used for sorting)
float chi2; // total chi2 (used for sorting)
};

//----------------------------------------------------------------------------

MPlexQF Chi2;
Expand Down
6 changes: 6 additions & 0 deletions mkFit/buildtestMPlex.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,6 +30,12 @@ inline bool sortByHitsChi2(const std::pair<Track, TrackState>& cand1,
return cand1.first.nFoundHits() > cand2.first.nFoundHits();
}

inline bool sortByScore(const std::pair<Track, TrackState>& cand1,
const std::pair<Track, TrackState>& cand2)
{
return sortByScoreCandPair(cand1,cand2);
}

inline bool sortByPhi(const Hit& hit1, const Hit& hit2)
{
return std::atan2(hit1.y(),hit1.x()) < std::atan2(hit2.y(),hit2.x());
Expand Down