Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-07-07 22:33:33

0001 #include "L1Trigger/TrackFindingTracklet/interface/PurgeDuplicate.h"
0002 #include "L1Trigger/TrackFindingTracklet/interface/Tracklet.h"
0003 #include "L1Trigger/TrackFindingTracklet/interface/Globals.h"
0004 #include "L1Trigger/TrackFindingTracklet/interface/CleanTrackMemory.h"
0005 #include "L1Trigger/TrackFindingTracklet/interface/Settings.h"
0006 #include "L1Trigger/TrackFindingTracklet/interface/Stub.h"
0007 #include "L1Trigger/TrackFindingTracklet/interface/Track.h"
0008 
0009 #ifdef USEHYBRID
0010 #include "DataFormats/L1TrackTrigger/interface/TTStub.h"
0011 #include "L1Trigger/TrackFindingTMTT/interface/L1track3D.h"
0012 #include "L1Trigger/TrackFindingTMTT/interface/KFParamsComb.h"
0013 #include "L1Trigger/TrackFindingTracklet/interface/HybridFit.h"
0014 #endif
0015 
0016 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0017 #include "FWCore/Utilities/interface/Exception.h"
0018 
0019 #include <unordered_set>
0020 #include <algorithm>
0021 
0022 using namespace std;
0023 using namespace trklet;
0024 
0025 PurgeDuplicate::PurgeDuplicate(std::string name, Settings const& settings, Globals* global)
0026     : ProcessBase(name, settings, global) {}
0027 
0028 void PurgeDuplicate::addOutput(MemoryBase* memory, std::string output) {
0029   if (settings_.writetrace()) {
0030     edm::LogVerbatim("Tracklet") << "In " << name_ << " adding output to " << memory->getName() << " to output "
0031                                  << output;
0032   }
0033   unordered_set<string> outputs = {"trackout",
0034                                    "trackout1",
0035                                    "trackout2",
0036                                    "trackout3",
0037                                    "trackout4",
0038                                    "trackout5",
0039                                    "trackout6",
0040                                    "trackout7",
0041                                    "trackout8",
0042                                    "trackout9",
0043                                    "trackout10",
0044                                    "trackout11"};
0045   if (outputs.find(output) != outputs.end()) {
0046     auto* tmp = dynamic_cast<CleanTrackMemory*>(memory);
0047     assert(tmp != nullptr);
0048     outputtracklets_.push_back(tmp);
0049     return;
0050   }
0051   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find output: " << output;
0052 }
0053 
0054 void PurgeDuplicate::addInput(MemoryBase* memory, std::string input) {
0055   if (settings_.writetrace()) {
0056     edm::LogVerbatim("Tracklet") << "In " << name_ << " adding input from " << memory->getName() << " to input "
0057                                  << input;
0058   }
0059   unordered_set<string> inputs = {"trackin",
0060                                   "trackin1",
0061                                   "trackin2",
0062                                   "trackin3",
0063                                   "trackin4",
0064                                   "trackin5",
0065                                   "trackin6",
0066                                   "trackin7",
0067                                   "trackin8",
0068                                   "trackin9",
0069                                   "trackin10",
0070                                   "trackin11",
0071                                   "trackin12"};
0072   if (inputs.find(input) != inputs.end()) {
0073     auto* tmp = dynamic_cast<TrackFitMemory*>(memory);
0074     assert(tmp != nullptr);
0075     inputtrackfits_.push_back(tmp);
0076     return;
0077   }
0078   throw cms::Exception("BadConfig") << __FILE__ << " " << __LINE__ << " could not find input: " << input;
0079 }
0080 
0081 void PurgeDuplicate::execute(std::vector<Track>& outputtracks_, unsigned int iSector) {
0082   inputtracklets_.clear();
0083   inputtracks_.clear();
0084 
0085   inputstubidslists_.clear();
0086   inputstublists_.clear();
0087   mergedstubidslists_.clear();
0088 
0089   if (settings_.removalType() != "merge") {
0090     for (auto& inputtrackfit : inputtrackfits_) {
0091       if (inputtrackfit->nTracks() == 0)
0092         continue;
0093       for (unsigned int j = 0; j < inputtrackfit->nTracks(); j++) {
0094         Track* aTrack = inputtrackfit->getTrack(j)->getTrack();
0095         aTrack->setSector(iSector);
0096         inputtracks_.push_back(aTrack);
0097       }
0098     }
0099     if (inputtracks_.empty())
0100       return;
0101   }
0102 
0103   unsigned int numTrk = inputtracks_.size();
0104 
0105   ////////////////////
0106   // Hybrid Removal //
0107   ////////////////////
0108 #ifdef USEHYBRID
0109 
0110   if (settings_.removalType() == "merge") {
0111     std::vector<std::pair<int, bool>> trackInfo;  // Track seed & duplicate flag
0112     // Vector to store the relative rank of the track candidate for merging, based on seed type
0113     std::vector<int> seedRank;
0114 
0115     // Get vectors from TrackFit and save them
0116     // inputtracklets: Tracklet objects from the FitTrack (not actually fit yet)
0117     // inputstublists: L1Stubs for that track
0118     // inputstubidslists: Stub stubIDs for that 3rack
0119     // mergedstubidslists: the same as inputstubidslists, but will be used during duplicate removal
0120     for (unsigned int i = 0; i < inputtrackfits_.size(); i++) {
0121       if (inputtrackfits_[i]->nStublists() == 0)
0122         continue;
0123       if (inputtrackfits_[i]->nStublists() != inputtrackfits_[i]->nTracks())
0124         throw "Number of stublists and tracks don't match up!";
0125       for (unsigned int j = 0; j < inputtrackfits_[i]->nStublists(); j++) {
0126         Tracklet* aTrack = inputtrackfits_[i]->getTrack(j);
0127         inputtracklets_.push_back(inputtrackfits_[i]->getTrack(j));
0128 
0129         std::vector<const Stub*> stublist = inputtrackfits_[i]->getStublist(j);
0130 
0131         inputstublists_.push_back(stublist);
0132 
0133         std::vector<std::pair<int, int>> stubidslist = inputtrackfits_[i]->getStubidslist(j);
0134         inputstubidslists_.push_back(stubidslist);
0135         mergedstubidslists_.push_back(stubidslist);
0136 
0137         // Encoding: L1L2=0, L2L3=1, L3L4=2, L5L6=3, D1D2=4, D3D4=5, L1D1=6, L2D1=7
0138         // Best Guess:          L1L2 > L1D1 > L2L3 > L2D1 > D1D2 > L3L4 > L5L6 > D3D4
0139         // Best Rank:           L1L2 > L3L4 > D3D4 > D1D2 > L2L3 > L2D1 > L5L6 > L1D1
0140         // Rank-Informed Guess: L1L2 > L3L4 > L1D1 > L2L3 > L2D1 > D1D2 > L5L6 > D3D4
0141         unsigned int curSeed = aTrack->seedIndex();
0142         if (curSeed == 0) {
0143           seedRank.push_back(1);
0144         } else if (curSeed == 2) {
0145           seedRank.push_back(2);
0146         } else if (curSeed == 5) {
0147           seedRank.push_back(3);
0148         } else if (curSeed == 4) {
0149           seedRank.push_back(4);
0150         } else if (curSeed == 1) {
0151           seedRank.push_back(5);
0152         } else if (curSeed == 7) {
0153           seedRank.push_back(6);
0154         } else if (curSeed == 3) {
0155           seedRank.push_back(7);
0156         } else if (curSeed == 6) {
0157           seedRank.push_back(8);
0158         } else if (settings_.extended()) {
0159           seedRank.push_back(9);
0160         } else {
0161           throw cms::Exception("LogError") << __FILE__ << " " << __LINE__ << " Seed " << curSeed
0162                                            << " not found in list, and settings->extended() not set.";
0163         }
0164 
0165         if (stublist.size() != stubidslist.size())
0166           throw "Number of stubs and stubids don't match up!";
0167 
0168         trackInfo.emplace_back(i, false);
0169       }
0170     }
0171 
0172     if (inputtracklets_.empty())
0173       return;
0174     unsigned int numStublists = inputstublists_.size();
0175 
0176     // Initialize all-false 2D array of tracks being duplicates to other tracks
0177     bool dupMap[numStublists][numStublists];  // Ends up symmetric
0178     for (unsigned int itrk = 0; itrk < numStublists; itrk++) {
0179       for (unsigned int jtrk = 0; jtrk < numStublists; jtrk++) {
0180         dupMap[itrk][jtrk] = false;
0181       }
0182     }
0183 
0184     // Find duplicates; Fill dupMap by looping over all pairs of "tracks"
0185     // numStublists-1 since last track has no other to compare to
0186     for (unsigned int itrk = 0; itrk < numStublists - 1; itrk++) {
0187       for (unsigned int jtrk = itrk + 1; jtrk < numStublists; jtrk++) {
0188         // Get primary track stubids
0189         const std::vector<std::pair<int, int>>& stubsTrk1 = inputstubidslists_[itrk];
0190 
0191         // Get and count secondary track stubids
0192         const std::vector<std::pair<int, int>>& stubsTrk2 = inputstubidslists_[jtrk];
0193 
0194         // Count number of Unique Regions (UR) that share stubs, and the number of UR that each track hits
0195         unsigned int nShareUR = 0;
0196         unsigned int nURStubTrk1 = 0;
0197         unsigned int nURStubTrk2 = 0;
0198         if (settings_.mergeComparison() == "CompareAll") {
0199           bool URArray[16];
0200           for (auto& i : URArray) {
0201             i = false;
0202           };
0203           for (const auto& st1 : stubsTrk1) {
0204             for (const auto& st2 : stubsTrk2) {
0205               if (st1.first == st2.first && st1.second == st2.second) {
0206                 // Converts region encoded in st1->first to an index in the Unique Region (UR) array
0207                 int i = st1.first;
0208                 int reg = (i > 0 && i < 10) * (i - 1) + (i > 10) * (i - 5) - (i < 0) * i;
0209                 if (!URArray[reg]) {
0210                   nShareUR++;
0211                   URArray[reg] = true;
0212                 }
0213               }
0214             }
0215           }
0216         } else if (settings_.mergeComparison() == "CompareBest") {
0217           std::vector<const Stub*> fullStubslistsTrk1 = inputstublists_[itrk];
0218           std::vector<const Stub*> fullStubslistsTrk2 = inputstublists_[jtrk];
0219 
0220           // Arrays to store the index of the best stub in each region
0221           int URStubidsTrk1[16];
0222           int URStubidsTrk2[16];
0223           for (int i = 0; i < 16; i++) {
0224             URStubidsTrk1[i] = -1;
0225             URStubidsTrk2[i] = -1;
0226           }
0227           // For each stub on the first track, find the stub with the best residual and store its index in the URStubidsTrk1 array
0228           for (unsigned int stcount = 0; stcount < stubsTrk1.size(); stcount++) {
0229             int i = stubsTrk1[stcount].first;
0230             int reg = (i > 0 && i < 10) * (i - 1) + (i > 10) * (i - 5) - (i < 0) * i;
0231             double nres = getPhiRes(inputtracklets_[itrk], fullStubslistsTrk1[stcount]);
0232             double ores = 0;
0233             if (URStubidsTrk1[reg] != -1)
0234               ores = getPhiRes(inputtracklets_[itrk], fullStubslistsTrk1[URStubidsTrk1[reg]]);
0235             if (URStubidsTrk1[reg] == -1 || nres < ores) {
0236               URStubidsTrk1[reg] = stcount;
0237             }
0238           }
0239           // For each stub on the second track, find the stub with the best residual and store its index in the URStubidsTrk1 array
0240           for (unsigned int stcount = 0; stcount < stubsTrk2.size(); stcount++) {
0241             int i = stubsTrk2[stcount].first;
0242             int reg = (i > 0 && i < 10) * (i - 1) + (i > 10) * (i - 5) - (i < 0) * i;
0243             double nres = getPhiRes(inputtracklets_[jtrk], fullStubslistsTrk2[stcount]);
0244             double ores = 0;
0245             if (URStubidsTrk2[reg] != -1)
0246               ores = getPhiRes(inputtracklets_[jtrk], fullStubslistsTrk2[URStubidsTrk2[reg]]);
0247             if (URStubidsTrk2[reg] == -1 || nres < ores) {
0248               URStubidsTrk2[reg] = stcount;
0249             }
0250           }
0251           // For all 16 regions (6 layers and 10 disks), count the number of regions who's best stub on both tracks are the same
0252           for (int i = 0; i < 16; i++) {
0253             int t1i = URStubidsTrk1[i];
0254             int t2i = URStubidsTrk2[i];
0255             if (t1i != -1 && t2i != -1 && stubsTrk1[t1i].first == stubsTrk2[t2i].first &&
0256                 stubsTrk1[t1i].second == stubsTrk2[t2i].second)
0257               nShareUR++;
0258           }
0259           // Calculate the number of unique regions hit by each track, so that this number can be used in calculating the number of independent
0260           // stubs on a track (not enabled/used by default)
0261           for (int i = 0; i < 16; i++) {
0262             if (URStubidsTrk1[i] != -1)
0263               nURStubTrk1++;
0264             if (URStubidsTrk2[i] != -1)
0265               nURStubTrk2++;
0266           }
0267         }
0268 
0269         // Fill duplicate map
0270         if (nShareUR >= settings_.minIndStubs()) {  // For number of shared stub merge condition
0271           dupMap[itrk][jtrk] = true;
0272           dupMap[jtrk][itrk] = true;
0273         }
0274       }
0275     }
0276 
0277     // Merge duplicate tracks
0278     for (unsigned int itrk = 0; itrk < numStublists - 1; itrk++) {
0279       for (unsigned int jtrk = itrk + 1; jtrk < numStublists; jtrk++) {
0280         // Merge a track with its first duplicate found.
0281         if (dupMap[itrk][jtrk]) {
0282           // Set preferred track based on seed rank
0283           int preftrk;
0284           int rejetrk;
0285           if (seedRank[itrk] < seedRank[jtrk]) {
0286             preftrk = itrk;
0287             rejetrk = jtrk;
0288           } else {
0289             preftrk = jtrk;
0290             rejetrk = itrk;
0291           }
0292 
0293           // Get a merged stub list
0294           std::vector<const Stub*> newStubList;
0295           std::vector<const Stub*> stubsTrk1 = inputstublists_[rejetrk];
0296           std::vector<const Stub*> stubsTrk2 = inputstublists_[preftrk];
0297           newStubList = stubsTrk1;
0298           for (unsigned int stub2it = 0; stub2it < stubsTrk2.size(); stub2it++) {
0299             if (find(stubsTrk1.begin(), stubsTrk1.end(), stubsTrk2[stub2it]) == stubsTrk1.end()) {
0300               newStubList.push_back(stubsTrk2[stub2it]);
0301             }
0302           }
0303           // Overwrite stublist of preferred track with merged list
0304           inputstublists_[preftrk] = newStubList;
0305 
0306           std::vector<std::pair<int, int>> newStubidsList;
0307           std::vector<std::pair<int, int>> stubidsTrk1 = mergedstubidslists_[rejetrk];
0308           std::vector<std::pair<int, int>> stubidsTrk2 = mergedstubidslists_[preftrk];
0309           newStubidsList = stubidsTrk1;
0310           for (unsigned int stub2it = 0; stub2it < stubidsTrk2.size(); stub2it++) {
0311             if (find(stubidsTrk1.begin(), stubidsTrk1.end(), stubidsTrk2[stub2it]) == stubidsTrk1.end()) {
0312               newStubidsList.push_back(stubidsTrk2[stub2it]);
0313             }
0314           }
0315           // Overwrite stubidslist of preferred track with merged list
0316           mergedstubidslists_[preftrk] = newStubidsList;
0317 
0318           // Mark that rejected track has been merged into another track
0319           trackInfo[rejetrk].second = true;
0320         }
0321       }
0322     }
0323 
0324     // Make the final track objects, fit with KF, and send to output
0325     for (unsigned int itrk = 0; itrk < numStublists; itrk++) {
0326       Tracklet* tracklet = inputtracklets_[itrk];
0327       std::vector<const Stub*> trackstublist = inputstublists_[itrk];
0328 
0329       HybridFit hybridFitter(iSector, settings_, globals_);
0330       hybridFitter.Fit(tracklet, trackstublist);
0331 
0332       // If the track was accepted (and thus fit), add to output
0333       if (tracklet->fit()) {
0334         // Add track to output if it wasn't merged into another
0335         Track* outtrack = tracklet->getTrack();
0336         outtrack->setSector(iSector);
0337         if (trackInfo[itrk].second == true)
0338           outtrack->setDuplicate(true);
0339         else
0340           outputtracklets_[trackInfo[itrk].first]->addTrack(tracklet);
0341 
0342         // Add all tracks to standalone root file output
0343         outtrack->setStubIDpremerge(inputstubidslists_[itrk]);
0344         outtrack->setStubIDprefit(mergedstubidslists_[itrk]);
0345         outputtracks_.push_back(*outtrack);
0346       }
0347     }
0348   }
0349 #endif
0350 
0351   //////////////////
0352   // Grid removal //
0353   //////////////////
0354   if (settings_.removalType() == "grid") {
0355     // Sort tracks by ichisq/DoF so that removal will keep the lower ichisq/DoF track
0356     std::sort(inputtracks_.begin(), inputtracks_.end(), [](const Track* lhs, const Track* rhs) {
0357       return lhs->ichisq() / lhs->stubID().size() < rhs->ichisq() / rhs->stubID().size();
0358     });
0359     bool grid[35][40] = {{false}};
0360 
0361     for (unsigned int itrk = 0; itrk < numTrk; itrk++) {
0362       if (inputtracks_[itrk]->duplicate())
0363         edm::LogPrint("Tracklet") << "WARNING: Track already tagged as duplicate!!";
0364 
0365       double phiBin = (inputtracks_[itrk]->phi0(settings_) - 2 * M_PI / 27 * iSector) / (2 * M_PI / 9 / 50) + 9;
0366       phiBin = std::max(phiBin, 0.);
0367       phiBin = std::min(phiBin, 34.);
0368 
0369       double ptBin = 1 / inputtracks_[itrk]->pt(settings_) * 40 + 20;
0370       ptBin = std::max(ptBin, 0.);
0371       ptBin = std::min(ptBin, 39.);
0372 
0373       if (grid[(int)phiBin][(int)ptBin])
0374         inputtracks_[itrk]->setDuplicate(true);
0375       grid[(int)phiBin][(int)ptBin] = true;
0376 
0377       double phiTest = inputtracks_[itrk]->phi0(settings_) - 2 * M_PI / 27 * iSector;
0378       if (phiTest < -2 * M_PI / 27)
0379         edm::LogVerbatim("Tracklet") << "track phi too small!";
0380       if (phiTest > 2 * 2 * M_PI / 27)
0381         edm::LogVerbatim("Tracklet") << "track phi too big!";
0382     }
0383   }  // end grid removal
0384 
0385   //////////////////////////
0386   // ichi + nstub removal //
0387   //////////////////////////
0388   if (settings_.removalType() == "ichi" || settings_.removalType() == "nstub") {
0389     for (unsigned int itrk = 0; itrk < numTrk - 1; itrk++) {  // numTrk-1 since last track has no other to compare to
0390 
0391       // If primary track is a duplicate, it cannot veto any...move on
0392       if (inputtracks_[itrk]->duplicate() == 1)
0393         continue;
0394 
0395       unsigned int nStubP = 0;
0396       vector<unsigned int> nStubS(numTrk);
0397       vector<unsigned int> nShare(numTrk);
0398       // Get and count primary stubs
0399       std::map<int, int> stubsTrk1 = inputtracks_[itrk]->stubID();
0400       nStubP = stubsTrk1.size();
0401 
0402       for (unsigned int jtrk = itrk + 1; jtrk < numTrk; jtrk++) {
0403         // Skip duplicate tracks
0404         if (inputtracks_[jtrk]->duplicate() == 1)
0405           continue;
0406 
0407         // Get and count secondary stubs
0408         std::map<int, int> stubsTrk2 = inputtracks_[jtrk]->stubID();
0409         nStubS[jtrk] = stubsTrk2.size();
0410 
0411         // Count shared stubs
0412         for (auto& st : stubsTrk1) {
0413           if (stubsTrk2.find(st.first) != stubsTrk2.end()) {
0414             if (st.second == stubsTrk2[st.first])
0415               nShare[jtrk]++;
0416           }
0417         }
0418       }
0419 
0420       // Tag duplicates
0421       for (unsigned int jtrk = itrk + 1; jtrk < numTrk; jtrk++) {
0422         // Skip duplicate tracks
0423         if (inputtracks_[jtrk]->duplicate() == 1)
0424           continue;
0425 
0426         // Chi2 duplicate removal
0427         if (settings_.removalType() == "ichi") {
0428           if ((nStubP - nShare[jtrk] < settings_.minIndStubs()) ||
0429               (nStubS[jtrk] - nShare[jtrk] < settings_.minIndStubs())) {
0430             if ((int)inputtracks_[itrk]->ichisq() / (2 * inputtracks_[itrk]->stubID().size() - 4) >
0431                 (int)inputtracks_[jtrk]->ichisq() / (2 * inputtracks_[itrk]->stubID().size() - 4)) {
0432               inputtracks_[itrk]->setDuplicate(true);
0433             } else if ((int)inputtracks_[itrk]->ichisq() / (2 * inputtracks_[itrk]->stubID().size() - 4) <=
0434                        (int)inputtracks_[jtrk]->ichisq() / (2 * inputtracks_[itrk]->stubID().size() - 4)) {
0435               inputtracks_[jtrk]->setDuplicate(true);
0436             } else {
0437               edm::LogVerbatim("Tracklet") << "Error: Didn't tag either track in duplicate pair.";
0438             }
0439           }
0440         }  // end ichi removal
0441 
0442         // nStub duplicate removal
0443         if (settings_.removalType() == "nstub") {
0444           if ((nStubP - nShare[jtrk] < settings_.minIndStubs()) && (nStubP < nStubS[jtrk])) {
0445             inputtracks_[itrk]->setDuplicate(true);
0446           } else if ((nStubS[jtrk] - nShare[jtrk] < settings_.minIndStubs()) && (nStubS[jtrk] <= nStubP)) {
0447             inputtracks_[jtrk]->setDuplicate(true);
0448           } else {
0449             edm::LogVerbatim("Tracklet") << "Error: Didn't tag either track in duplicate pair.";
0450           }
0451         }  // end nstub removal
0452 
0453       }  // end tag duplicates
0454 
0455     }  // end loop over primary track
0456 
0457   }  // end ichi + nstub removal
0458 
0459   //Add tracks to output
0460   if (settings_.removalType() != "merge") {
0461     for (unsigned int i = 0; i < inputtrackfits_.size(); i++) {
0462       for (unsigned int j = 0; j < inputtrackfits_[i]->nTracks(); j++) {
0463         if (inputtrackfits_[i]->getTrack(j)->getTrack()->duplicate() == 0) {
0464           if (settings_.writeMonitorData("Seeds")) {
0465             ofstream fout("seeds.txt", ofstream::app);
0466             fout << __FILE__ << ":" << __LINE__ << " " << name_ << "_" << iSector << " "
0467                  << inputtrackfits_[i]->getTrack(j)->getISeed() << endl;
0468             fout.close();
0469           }
0470           outputtracklets_[i]->addTrack(inputtrackfits_[i]->getTrack(j));
0471         }
0472         //For root file:
0473         outputtracks_.push_back(*inputtrackfits_[i]->getTrack(j)->getTrack());
0474       }
0475     }
0476   }
0477 }
0478 
0479 double PurgeDuplicate::getPhiRes(Tracklet* curTracklet, const Stub* curStub) {
0480   double phiproj;
0481   double stubphi;
0482   double phires;
0483   // Get phi position of stub
0484   stubphi = curStub->l1tstub()->phi();
0485   // Get region that the stub is in (Layer 1->6, Disk 1->5)
0486   int Layer = curStub->layerdisk() + 1;
0487   if (Layer > N_LAYER) {
0488     Layer = 0;
0489   }
0490   int Disk = curStub->layerdisk() - (N_LAYER - 1);
0491   if (Disk < 0) {
0492     Disk = 0;
0493   }
0494   // Get phi projection of tracklet
0495   int seedindex = curTracklet->seedIndex();
0496   // If this stub is a seed stub, set projection=phi, so that res=0
0497   if ((seedindex == 0 && (Layer == 1 || Layer == 2)) || (seedindex == 1 && (Layer == 2 || abs(Disk) == 0)) ||
0498       (seedindex == 2 && (Layer == 3 || Layer == 4)) || (seedindex == 3 && (Layer == 5 || Layer == 6)) ||
0499       (seedindex == 4 && (abs(Disk) == 1 || abs(Disk) == 2)) ||
0500       (seedindex == 5 && (abs(Disk) == 3 || abs(Disk) == 4)) || (seedindex == 6 && (Layer == 1 || abs(Disk) == 1)) ||
0501       (seedindex == 7 && (Layer == 2 || abs(Disk) == 1)) ||
0502       (seedindex == 8 && (Layer == 2 || Layer == 3 || Layer == 4)) ||
0503       (seedindex == 9 && (Layer == 4 || Layer == 5 || Layer == 6)) ||
0504       (seedindex == 10 && (Layer == 2 || Layer == 3 || abs(Disk) == 1)) ||
0505       (seedindex == 11 && (Layer == 2 || abs(Disk) == 1 || abs(Disk) == 2))) {
0506     phiproj = stubphi;
0507     // Otherwise, get projection of tracklet
0508   } else {
0509     phiproj = curTracklet->proj(curStub->layerdisk()).phiproj();
0510   }
0511   // Calculate residual
0512   phires = std::abs(stubphi - phiproj);
0513   return phires;
0514 }