Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-05-12 01:51:37

0001 #include <memory>
0002 
0003 #include "FWCore/Framework/interface/ConsumesCollector.h"
0004 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0005 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0006 
0007 #include "TFile.h"
0008 
0009 #include "L1Trigger/DTTriggerPhase2/interface/PseudoBayesGrouping.h"
0010 
0011 using namespace edm;
0012 using namespace std;
0013 using namespace cmsdt;
0014 using namespace dtbayesam;
0015 // ============================================================================
0016 // Constructors and destructor
0017 // ============================================================================
0018 PseudoBayesGrouping::PseudoBayesGrouping(const ParameterSet& pset, edm::ConsumesCollector& iC)
0019     : MotherGrouping(pset, iC),
0020       debug_(pset.getUntrackedParameter<bool>("debug")),
0021       pattern_filename_(pset.getParameter<edm::FileInPath>("pattern_filename").fullPath()),
0022       minNLayerHits_(pset.getParameter<int>("minNLayerHits")),
0023       allowedVariance_(pset.getParameter<int>("allowedVariance")),
0024       allowDuplicates_(pset.getParameter<bool>("allowDuplicates")),
0025       allowUncorrelatedPatterns_(pset.getParameter<bool>("allowUncorrelatedPatterns")),
0026       setLateralities_(pset.getParameter<bool>("setLateralities")),
0027       saveOnPlace_(pset.getParameter<bool>("saveOnPlace")),
0028       minSingleSLHitsMax_(pset.getParameter<int>("minSingleSLHitsMax")),
0029       minSingleSLHitsMin_(pset.getParameter<int>("minSingleSLHitsMin")),
0030       minUncorrelatedHits_(pset.getParameter<int>("minUncorrelatedHits")),
0031       maxPathsPerMatch_(pset.getParameter<int>("maxPathsPerMatch")) {
0032   if (debug_)
0033     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping:: constructor";
0034 }
0035 
0036 PseudoBayesGrouping::~PseudoBayesGrouping() {
0037   if (debug_)
0038     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping:: destructor";
0039 }
0040 
0041 // ============================================================================
0042 // Main methods (initialise, run, finish)
0043 // ============================================================================
0044 void PseudoBayesGrouping::initialise(const edm::EventSetup& iEventSetup) {
0045   if (debug_)
0046     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::initialiase";
0047   if (debug_)
0048     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::initialiase using patterns file " << pattern_filename_;
0049   nPatterns_ = 0;
0050 
0051   TString patterns_folder = "L1Trigger/DTTriggerPhase2/data/";
0052 
0053   // Load all patterns
0054   // MB1
0055   LoadPattern(patterns_folder + "createdPatterns_MB1_left.root", 0, 0);
0056   LoadPattern(patterns_folder + "createdPatterns_MB1_right.root", 0, 2);
0057   // MB2
0058   LoadPattern(patterns_folder + "createdPatterns_MB2_left.root", 1, 0);
0059   LoadPattern(patterns_folder + "createdPatterns_MB2_right.root", 1, 2);
0060   // MB3
0061   LoadPattern(patterns_folder + "createdPatterns_MB3.root", 2, 1);
0062   // MB4
0063   LoadPattern(patterns_folder + "createdPatterns_MB4_left.root", 3, 0);
0064   LoadPattern(patterns_folder + "createdPatterns_MB4.root", 3, 1);
0065   LoadPattern(patterns_folder + "createdPatterns_MB4_right.root", 3, 2);
0066 
0067   prelimMatches_ = std::make_unique<CandidateGroupPtrs>();
0068   allMatches_ = std::make_unique<CandidateGroupPtrs>();
0069   finalMatches_ = std::make_unique<CandidateGroupPtrs>();
0070 }
0071 
0072 void PseudoBayesGrouping::LoadPattern(TString pattern_file_name, int MB_number, int SL_shift) {
0073   TFile* f = TFile::Open(pattern_file_name, "READ");
0074 
0075   std::vector<std::vector<std::vector<int>>>* pattern_reader =
0076       (std::vector<std::vector<std::vector<int>>>*)f->Get("allPatterns");
0077 
0078   for (std::vector<std::vector<std::vector<int>>>::iterator itPattern = (*pattern_reader).begin();
0079        itPattern != (*pattern_reader).end();
0080        ++itPattern) {
0081     if (debug_)
0082       LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::LoadPattern Loading patterns seeded by: "
0083                                       << itPattern->at(0).at(0) << ", " << itPattern->at(0).at(1) << ", "
0084                                       << itPattern->at(0).at(2) << ", ";
0085 
0086     auto p = std::make_shared<DTPattern>();
0087 
0088     bool is_seed = true;
0089     for (const auto& itHits : *itPattern) {
0090       // First entry is the seeding information
0091       if (is_seed) {
0092         p = std::make_shared<DTPattern>(DTPattern(itHits.at(0), itHits.at(1), itHits.at(2)));
0093         is_seed = false;
0094       }
0095       // Other entries are the hits information
0096       else {
0097         if (itHits.begin() == itHits.end())
0098           continue;
0099         // We need to correct the geometry from pattern generation to reconstruction as they use slightly displaced basis
0100         else if (itHits.at(0) % 2 == 0) {
0101           p->addHit(std::make_tuple(itHits.at(0), itHits.at(1), itHits.at(2)));
0102         } else if (itHits.at(0) % 2 == 1) {
0103           p->addHit(std::make_tuple(itHits.at(0), itHits.at(1) - 1, itHits.at(2)));
0104         }
0105       }
0106     }
0107     // Classified by seeding layers for optimized search later
0108     // TODO::This can be vastly improved using std::bitset<8>, for example
0109     if (p->sl1() == 0) {
0110       if (p->sl2() == 7)
0111         L0L7Patterns_[MB_number][SL_shift].push_back(p);
0112       if (p->sl2() == 6)
0113         L0L6Patterns_[MB_number][SL_shift].push_back(p);
0114       if (p->sl2() == 5)
0115         L0L5Patterns_[MB_number][SL_shift].push_back(p);
0116       if (p->sl2() == 4)
0117         L0L4Patterns_[MB_number][SL_shift].push_back(p);
0118       if (p->sl2() == 3)
0119         L0L3Patterns_[MB_number][SL_shift].push_back(p);
0120       if (p->sl2() == 2)
0121         L0L2Patterns_[MB_number][SL_shift].push_back(p);
0122       if (p->sl2() == 1)
0123         L0L1Patterns_[MB_number][SL_shift].push_back(p);
0124     }
0125     if (p->sl1() == 1) {
0126       if (p->sl2() == 7)
0127         L1L7Patterns_[MB_number][SL_shift].push_back(p);
0128       if (p->sl2() == 6)
0129         L1L6Patterns_[MB_number][SL_shift].push_back(p);
0130       if (p->sl2() == 5)
0131         L1L5Patterns_[MB_number][SL_shift].push_back(p);
0132       if (p->sl2() == 4)
0133         L1L4Patterns_[MB_number][SL_shift].push_back(p);
0134       if (p->sl2() == 3)
0135         L1L3Patterns_[MB_number][SL_shift].push_back(p);
0136       if (p->sl2() == 2)
0137         L1L2Patterns_[MB_number][SL_shift].push_back(p);
0138     }
0139     if (p->sl1() == 2) {
0140       if (p->sl2() == 7)
0141         L2L7Patterns_[MB_number][SL_shift].push_back(p);
0142       if (p->sl2() == 6)
0143         L2L6Patterns_[MB_number][SL_shift].push_back(p);
0144       if (p->sl2() == 5)
0145         L2L5Patterns_[MB_number][SL_shift].push_back(p);
0146       if (p->sl2() == 4)
0147         L2L4Patterns_[MB_number][SL_shift].push_back(p);
0148       if (p->sl2() == 3)
0149         L2L3Patterns_[MB_number][SL_shift].push_back(p);
0150     }
0151     if (p->sl1() == 3) {
0152       if (p->sl2() == 7)
0153         L3L7Patterns_[MB_number][SL_shift].push_back(p);
0154       if (p->sl2() == 6)
0155         L3L6Patterns_[MB_number][SL_shift].push_back(p);
0156       if (p->sl2() == 5)
0157         L3L5Patterns_[MB_number][SL_shift].push_back(p);
0158       if (p->sl2() == 4)
0159         L3L4Patterns_[MB_number][SL_shift].push_back(p);
0160     }
0161 
0162     if (p->sl1() == 4) {
0163       if (p->sl2() == 7)
0164         L4L7Patterns_[MB_number][SL_shift].push_back(p);
0165       if (p->sl2() == 6)
0166         L4L6Patterns_[MB_number][SL_shift].push_back(p);
0167       if (p->sl2() == 5)
0168         L4L5Patterns_[MB_number][SL_shift].push_back(p);
0169     }
0170     if (p->sl1() == 5) {
0171       if (p->sl2() == 7)
0172         L5L7Patterns_[MB_number][SL_shift].push_back(p);
0173       if (p->sl2() == 6)
0174         L5L6Patterns_[MB_number][SL_shift].push_back(p);
0175     }
0176     if (p->sl1() == 6) {
0177       if (p->sl2() == 7)
0178         L6L7Patterns_[MB_number][SL_shift].push_back(p);
0179     }
0180 
0181     //Also creating a list of all patterns, needed later for deleting and avoid a memory leak
0182     allPatterns_[MB_number][SL_shift].push_back(p);
0183     nPatterns_++;
0184   }
0185   if (debug_)
0186     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::initialiase Total number of loaded patterns: "
0187                                     << nPatterns_;
0188   f->Close();
0189   delete f;
0190 }
0191 
0192 void PseudoBayesGrouping::run(Event& iEvent,
0193                               const EventSetup& iEventSetup,
0194                               const DTDigiCollection& digis,
0195                               MuonPathPtrs& mpaths) {
0196   //Takes dt digis collection and does the grouping for correlated hits, it is saved in a vector of up to 8 (or 4) correlated hits
0197   if (debug_)
0198     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run";
0199   //Do initial cleaning
0200   CleanDigisByLayer();
0201   //Sort digis by layer
0202   FillDigisByLayer(&digis);
0203 
0204   //Sarch for patterns
0205   DTChamberId chamber_id;
0206   // We just want the chamber ID of the first digi
0207   // as they are all the same --> create a loop and break it
0208   // after the first iteration
0209   for (const auto& detUnitIt : digis) {
0210     const DTLayerId layer_Id = detUnitIt.first;
0211     chamber_id = layer_Id.superlayerId().chamberId();
0212     break;
0213   }
0214 
0215   RecognisePatternsByLayerPairs(chamber_id);
0216 
0217   // Now sort patterns by qualities
0218   std::sort(prelimMatches_->begin(), prelimMatches_->end(), CandPointGreat());
0219   if (debug_ && !prelimMatches_->empty()) {
0220     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run Pattern qualities before cleaning: ";
0221     for (const auto& cand_it : *prelimMatches_) {
0222       LogDebug("PseudoBayesGrouping") << cand_it->nLayerhits() << ", " << cand_it->nisGood() << ", " << cand_it->nhits()
0223                                       << ", " << cand_it->quality() << ", " << cand_it->candId();
0224     }
0225   }
0226   //And ghostbust patterns to retain higher quality ones
0227   ReCleanPatternsAndDigis();
0228   if (debug_)
0229     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run Number of found patterns: " << finalMatches_->size();
0230 
0231   //Last organize candidates information into muonpaths to finalize the grouping
0232   FillMuonPaths(mpaths);
0233   if (debug_)
0234     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run ended run";
0235 }
0236 
0237 void PseudoBayesGrouping::FillMuonPaths(MuonPathPtrs& mpaths) {
0238   //Loop over all selected candidates
0239   for (auto itCand = finalMatches_->begin(); itCand != finalMatches_->end(); itCand++) {
0240     if (debug_)
0241       LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run Create pointers ";
0242 
0243     // Vector of all muon paths we may find
0244     std::vector<DTPrimitivePtrs> ptrPrimitive_vector;
0245 
0246     // We will have at least one muon path
0247     DTPrimitivePtrs ptrPrimitive;
0248     for (int i = 0; i < 8; i++)
0249       ptrPrimitive.push_back(std::make_shared<DTPrimitive>());
0250 
0251     ptrPrimitive_vector.push_back(ptrPrimitive);
0252 
0253     qualitybits qualityDTP;
0254     qualitybits qualityDTP2;
0255     int intHit = 0;
0256 
0257     //And for each candidate loop over all grouped hits
0258     for (auto& itDTP : (*itCand)->candHits()) {
0259       if (debug_)
0260         LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run loop over dt hits to fill pointer";
0261 
0262       int layerHit = (*itDTP).layerId();
0263       //Back to the usual basis for SL
0264       if (layerHit >= 4) {
0265         (*itDTP).setLayerId(layerHit - 4);
0266       }
0267       qualitybits ref8Hit(std::pow(2, layerHit));
0268       //Get the predicted laterality
0269       if (setLateralities_) {
0270         int predLat = (*itCand)->pattern()->latHitIn(layerHit, (*itDTP).channelId(), allowedVariance_);
0271         if (predLat == -10 || predLat == 0) {
0272           (*itDTP).setLaterality(NONE);
0273         } else if (predLat == -1) {
0274           (*itDTP).setLaterality(LEFT);
0275         } else if (predLat == +1) {
0276           (*itDTP).setLaterality(RIGHT);
0277         }
0278       }
0279 
0280       // If one hit is already present in the current layer, for each ptrPrimitive already existing,
0281       // create a new with all its hits. Then, fill it with the new hit and add it to the primitives vector.
0282       // Do not consider more than 2 hits in the same wire or more than maxPathsPerMatch_ total muonpaths per finalMatches_
0283       if (qualityDTP == (qualityDTP | ref8Hit) && qualityDTP2 != (qualityDTP2 | ref8Hit) &&
0284           ptrPrimitive_vector.size() < maxPathsPerMatch_) {
0285         if (debug_)
0286           LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run Creating additional muon paths";
0287 
0288         qualityDTP2 = (qualityDTP2 | ref8Hit);
0289 
0290         int n_prim = ptrPrimitive_vector.size();
0291 
0292         for (int j = 0; j < n_prim; j++) {
0293           DTPrimitivePtrs tmpPrimitive;
0294           for (int i = 0; i < NUM_LAYERS_2SL; i++) {
0295             tmpPrimitive.push_back(ptrPrimitive_vector.at(j).at(i));
0296           }
0297           // Now save the hit in the new path
0298           if (saveOnPlace_) {
0299             //This will save the primitive in a place of the vector equal to its L position
0300             tmpPrimitive.at(layerHit) = std::make_shared<DTPrimitive>((*itDTP));
0301           }
0302           if (!saveOnPlace_) {
0303             //This will save the primitive in order
0304             tmpPrimitive.at(intHit) = std::make_shared<DTPrimitive>((*itDTP));
0305           }
0306           // Now add the new path to the vector of paths
0307           ptrPrimitive_vector.push_back(tmpPrimitive);
0308         }
0309       }
0310 
0311       // If there is not one hit already in the layer, fill the DT primitives pointers
0312       else {
0313         if (debug_)
0314           LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::run Adding hit to muon path";
0315 
0316         qualityDTP = (qualityDTP | ref8Hit);
0317 
0318         // for (all paths --> fill them)
0319         for (auto prim_it = ptrPrimitive_vector.begin(); prim_it != ptrPrimitive_vector.end(); ++prim_it) {
0320           if (saveOnPlace_) {
0321             //This will save the primitive in a place of the vector equal to its L position
0322             prim_it->at(layerHit) = std::make_shared<DTPrimitive>((*itDTP));
0323           }
0324           if (!saveOnPlace_) {
0325             //This will save the primitive in order
0326             intHit++;
0327             prim_it->at(intHit) = std::make_shared<DTPrimitive>((*itDTP));
0328           }
0329         }
0330       }
0331     }
0332 
0333     stringstream ss;
0334 
0335     int n_paths = ptrPrimitive_vector.size();
0336 
0337     for (int n_path = 0; n_path < n_paths; ++n_path) {
0338       mpaths.emplace_back(std::make_shared<MuonPath>(
0339           ptrPrimitive_vector.at(n_path), (short)(*itCand)->nLayerUp(), (short)(*itCand)->nLayerDown()));
0340     }
0341   }
0342 }
0343 
0344 void PseudoBayesGrouping::RecognisePatternsByLayerPairs(DTChamberId chamber_ID) {
0345   // chamber_ID traslated to MB, wheel, sector
0346   int MB = chamber_ID.station() - 1;
0347   int wheel = chamber_ID.wheel();
0348   int sector = chamber_ID.sector();
0349 
0350   // shift of SL3 wrt SL1
0351   int shift = -1;
0352 
0353   // Now define DT geometry depending on its ID
0354 
0355   // MB1
0356   if (MB == 0) {
0357     if (wheel == -1 || wheel == -2)
0358       shift = 2;  // positive (right)
0359     else if (wheel == 1 || wheel == 2)
0360       shift = 0;  // negative (left)
0361     else if (wheel == 0) {
0362       if (sector == 1 || sector == 4 || sector == 5 || sector == 8 || sector == 9 || sector == 12)
0363         shift = 2;  // positive (right)
0364       else
0365         shift = 0;  // negative (left)
0366     }
0367   }
0368   // MB2
0369   else if (MB == 1) {
0370     if (wheel == -1 || wheel == -2)
0371       shift = 0;  // negative (left)
0372     else if (wheel == 1 || wheel == 2)
0373       shift = 2;  // positive (right)
0374     else if (wheel == 0) {
0375       if (sector == 1 || sector == 4 || sector == 5 || sector == 8 || sector == 9 || sector == 12)
0376         shift = 0;  // negative (left)
0377       else
0378         shift = 2;  // positive (right)
0379     }
0380   }
0381   // MB3
0382   else if (MB == 2) {
0383     shift = 1;  // shift is always 0 in MB3
0384   }
0385   // MB4
0386   else if (MB == 3) {
0387     if (wheel == -1 || wheel == -2)
0388       if (sector == 4 || sector == 9 || sector == 11 || sector == 13)
0389         shift = 1;  // no shift
0390       else if (sector == 5 || sector == 6 || sector == 7 || sector == 8 || sector == 14)
0391         shift = 2;  // positive (right)
0392       else
0393         shift = 0;  // negative (left)
0394     else if (wheel == 1 || wheel == 2)
0395       if (sector == 4 || sector == 9 || sector == 11 || sector == 13)
0396         shift = 1;  // no shift
0397       else if (sector == 1 || sector == 2 || sector == 3 || sector == 10 || sector == 12)
0398         shift = 2;  // positive (right)
0399       else
0400         shift = 0;  // negative (left)
0401     else if (wheel == 0)
0402       if (sector == 4 || sector == 9 || sector == 11 || sector == 13)
0403         shift = 1;  // no shift
0404       else if (sector == 2 || sector == 3 || sector == 5 || sector == 8 || sector == 10)
0405         shift = 2;  // positive (right)
0406       else
0407         shift = 0;  // negative (left)
0408     else
0409       return;
0410   }
0411 
0412   //Separated from main run function for clarity. Do all pattern recognition steps
0413   pidx_ = 0;
0414   //Compare L0-L7
0415   RecognisePatterns(digisinL0_, digisinL7_, L0L7Patterns_[MB][shift]);
0416   //Compare L0-L6 and L1-L7
0417   RecognisePatterns(digisinL0_, digisinL6_, L0L6Patterns_[MB][shift]);
0418   RecognisePatterns(digisinL1_, digisinL7_, L1L7Patterns_[MB][shift]);
0419   //Compare L0-L5, L1-L6, L2-L7
0420   RecognisePatterns(digisinL0_, digisinL5_, L0L5Patterns_[MB][shift]);
0421   RecognisePatterns(digisinL1_, digisinL6_, L1L6Patterns_[MB][shift]);
0422   RecognisePatterns(digisinL2_, digisinL7_, L2L7Patterns_[MB][shift]);
0423   //L0-L4, L1-L5, L2-L6, L3-L7
0424   RecognisePatterns(digisinL0_, digisinL4_, L0L4Patterns_[MB][shift]);
0425   RecognisePatterns(digisinL1_, digisinL5_, L1L5Patterns_[MB][shift]);
0426   RecognisePatterns(digisinL2_, digisinL6_, L2L6Patterns_[MB][shift]);
0427   RecognisePatterns(digisinL3_, digisinL7_, L3L7Patterns_[MB][shift]);
0428   //L1-L4, L2-L5, L3-L6
0429   RecognisePatterns(digisinL1_, digisinL4_, L1L4Patterns_[MB][shift]);
0430   RecognisePatterns(digisinL2_, digisinL5_, L2L5Patterns_[MB][shift]);
0431   RecognisePatterns(digisinL3_, digisinL6_, L3L6Patterns_[MB][shift]);
0432   //L2-L4, L3-L5
0433   RecognisePatterns(digisinL2_, digisinL4_, L2L4Patterns_[MB][shift]);
0434   RecognisePatterns(digisinL3_, digisinL5_, L3L5Patterns_[MB][shift]);
0435   //L3-L4
0436   RecognisePatterns(digisinL3_, digisinL4_, L3L4Patterns_[MB][shift]);
0437   //Uncorrelated SL1
0438   RecognisePatterns(digisinL0_, digisinL1_, L0L1Patterns_[MB][shift]);
0439   RecognisePatterns(digisinL0_, digisinL2_, L0L2Patterns_[MB][shift]);
0440   RecognisePatterns(digisinL0_, digisinL3_, L0L3Patterns_[MB][shift]);
0441   RecognisePatterns(digisinL1_, digisinL2_, L1L2Patterns_[MB][shift]);
0442   RecognisePatterns(digisinL1_, digisinL3_, L1L3Patterns_[MB][shift]);
0443   RecognisePatterns(digisinL2_, digisinL3_, L2L3Patterns_[MB][shift]);
0444   //Uncorrelated SL3
0445   RecognisePatterns(digisinL4_, digisinL5_, L4L5Patterns_[MB][shift]);
0446   RecognisePatterns(digisinL4_, digisinL6_, L4L6Patterns_[MB][shift]);
0447   RecognisePatterns(digisinL4_, digisinL7_, L4L7Patterns_[MB][shift]);
0448   RecognisePatterns(digisinL5_, digisinL6_, L5L6Patterns_[MB][shift]);
0449   RecognisePatterns(digisinL5_, digisinL7_, L5L7Patterns_[MB][shift]);
0450   RecognisePatterns(digisinL6_, digisinL7_, L6L7Patterns_[MB][shift]);
0451 }
0452 
0453 void PseudoBayesGrouping::RecognisePatterns(std::vector<DTPrimitive> digisinLDown,
0454                                             std::vector<DTPrimitive> digisinLUp,
0455                                             DTPatternPtrs patterns) {
0456   //Loop over all hits and search for matching patterns (there will be four
0457   // amongst ~60, accounting for possible lateralities)
0458   for (auto dtPD_it = digisinLDown.begin(); dtPD_it != digisinLDown.end(); dtPD_it++) {
0459     int LDown = dtPD_it->layerId();
0460     int wireDown = dtPD_it->channelId();
0461     for (auto dtPU_it = digisinLUp.begin(); dtPU_it != digisinLUp.end(); dtPU_it++) {
0462       int LUp = dtPU_it->layerId();
0463       int wireUp = dtPU_it->channelId();
0464       int diff = wireUp - wireDown;
0465       for (auto pat_it = patterns.begin(); pat_it != patterns.end(); pat_it++) {
0466         //For each pair of hits in the layers search for the seeded patterns
0467         if ((*pat_it)->sl1() != (LDown) || (*pat_it)->sl2() != (LUp) || (*pat_it)->diff() != diff)
0468           continue;
0469         //If we are here a pattern was found and we can start comparing
0470         (*pat_it)->setHitDown(wireDown);
0471         auto cand = std::make_shared<CandidateGroup>(*pat_it);
0472         for (auto dtTest_it = alldigis_.begin(); dtTest_it != alldigis_.end(); dtTest_it++) {
0473           //Find hits matching to the pattern
0474           if (((*pat_it)->latHitIn(dtTest_it->layerId(), dtTest_it->channelId(), allowedVariance_)) != -999) {
0475             if (((*pat_it)->latHitIn(dtTest_it->layerId(), dtTest_it->channelId(), allowedVariance_)) == -10)
0476               cand->addHit((*dtTest_it), dtTest_it->layerId(), false);
0477             else
0478               cand->addHit((*dtTest_it), dtTest_it->layerId(), true);
0479           }
0480         }
0481         if ((cand->nhits() >= minNLayerHits_ &&
0482              (cand->nLayerUp() >= minSingleSLHitsMax_ || cand->nLayerDown() >= minSingleSLHitsMax_) &&
0483              (cand->nLayerUp() >= minSingleSLHitsMin_ && cand->nLayerDown() >= minSingleSLHitsMin_)) ||
0484             (allowUncorrelatedPatterns_ && ((cand->nLayerUp() >= minUncorrelatedHits_ && cand->nLayerDown() == 0) ||
0485                                             (cand->nLayerDown() >= minUncorrelatedHits_ && cand->nLayerUp() == 0)))) {
0486           if (debug_) {
0487             LogDebug("PseudoBayesGrouping")
0488                 << "PseudoBayesGrouping::RecognisePatterns Pattern found for pair in " << LDown << " ," << wireDown
0489                 << " ," << LUp << " ," << wireUp << "\n"
0490                 << "Candidate has " << cand->nhits() << " hits with quality " << cand->quality() << "\n"
0491                 << *(*pat_it);
0492           }
0493           //We currently save everything at this level, might want to be more restrictive
0494           pidx_++;
0495           cand->setCandId(pidx_);
0496           prelimMatches_->push_back(std::move(cand));
0497           allMatches_->push_back(std::move(cand));
0498         }
0499       }
0500     }
0501   }
0502 }
0503 
0504 void PseudoBayesGrouping::FillDigisByLayer(const DTDigiCollection* digis) {
0505   //First we need to have separated lists of digis by layer
0506   if (debug_)
0507     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping::FillDigisByLayer Classifying digis by layer";
0508   //  for (auto dtDigi_It = digis->begin(); dtDigi_It != digis->end(); dtDigi_It++) {
0509   for (const auto& dtDigi_It : *digis) {
0510     const DTLayerId dtLId = dtDigi_It.first;
0511     //Skip digis in SL theta which we are not interested on for the grouping
0512     for (auto digiIt = (dtDigi_It.second).first; digiIt != (dtDigi_It.second).second; digiIt++) {
0513       //Need to change notation slightly here
0514       if (dtLId.superlayer() == 2)
0515         continue;
0516       int layer = dtLId.layer() - 1;
0517       if (dtLId.superlayer() == 3)
0518         layer += 4;
0519       //Use the same format as for InitialGrouping to avoid tons of replicating classes, we will have some not used variables
0520       DTPrimitive dtpAux = DTPrimitive();
0521       dtpAux.setTDCTimeStamp(digiIt->time());
0522       dtpAux.setChannelId(digiIt->wire() - 1);
0523       dtpAux.setLayerId(layer);
0524       dtpAux.setSuperLayerId(dtLId.superlayer());
0525       dtpAux.setCameraId(dtLId.rawId());
0526       if (debug_)
0527         LogDebug("PseudoBayesGrouping") << "Hit in L " << layer << " SL " << dtLId.superlayer() << " WIRE "
0528                                         << digiIt->wire() - 1;
0529       if (layer == 0)
0530         digisinL0_.push_back(dtpAux);
0531       else if (layer == 1)
0532         digisinL1_.push_back(dtpAux);
0533       else if (layer == 2)
0534         digisinL2_.push_back(dtpAux);
0535       else if (layer == 3)
0536         digisinL3_.push_back(dtpAux);
0537       else if (layer == 4)
0538         digisinL4_.push_back(dtpAux);
0539       else if (layer == 5)
0540         digisinL5_.push_back(dtpAux);
0541       else if (layer == 6)
0542         digisinL6_.push_back(dtpAux);
0543       else if (layer == 7)
0544         digisinL7_.push_back(dtpAux);
0545       alldigis_.push_back(dtpAux);
0546     }
0547   }
0548 }
0549 
0550 void PseudoBayesGrouping::ReCleanPatternsAndDigis() {
0551   //GhostbustPatterns that share hits and are of lower quality
0552   if (prelimMatches_->empty()) {
0553     return;
0554   };
0555   while ((prelimMatches_->at(0)->nLayerhits() >= minNLayerHits_ &&
0556           (prelimMatches_->at(0)->nLayerUp() >= minSingleSLHitsMax_ ||
0557            prelimMatches_->at(0)->nLayerDown() >= minSingleSLHitsMax_) &&
0558           (prelimMatches_->at(0)->nLayerUp() >= minSingleSLHitsMin_ &&
0559            prelimMatches_->at(0)->nLayerDown() >= minSingleSLHitsMin_)) ||
0560          (allowUncorrelatedPatterns_ &&
0561           ((prelimMatches_->at(0)->nLayerUp() >= minUncorrelatedHits_ && prelimMatches_->at(0)->nLayerDown() == 0) ||
0562            (prelimMatches_->at(0)->nLayerDown() >= minUncorrelatedHits_ && prelimMatches_->at(0)->nLayerUp() == 0)))) {
0563     finalMatches_->push_back(prelimMatches_->at(0));
0564     auto itSel = finalMatches_->end() - 1;
0565     prelimMatches_->erase(prelimMatches_->begin());
0566     if (prelimMatches_->empty()) {
0567       return;
0568     };
0569     for (auto cand_it = prelimMatches_->begin(); cand_it != prelimMatches_->end(); cand_it++) {
0570       if (*(*cand_it) == *(*itSel) && allowDuplicates_)
0571         continue;
0572       for (const auto& dt_it : (*itSel)->candHits()) {
0573         (*cand_it)->removeHit((*dt_it));
0574       }
0575     }
0576 
0577     std::sort(prelimMatches_->begin(), prelimMatches_->end(), CandPointGreat());
0578     if (debug_) {
0579       LogDebug("PseudoBayesGrouping") << "Pattern qualities: ";
0580       for (const auto& cand_it : *prelimMatches_) {
0581         LogDebug("PseudoBayesGrouping") << cand_it->nLayerhits() << ", " << cand_it->nisGood() << ", "
0582                                         << cand_it->nhits() << ", " << cand_it->quality() << ", " << cand_it->candId()
0583                                         << "\n";
0584       }
0585     }
0586   }
0587 }
0588 
0589 void PseudoBayesGrouping::CleanDigisByLayer() {
0590   digisinL0_.clear();
0591   digisinL1_.clear();
0592   digisinL2_.clear();
0593   digisinL3_.clear();
0594   digisinL4_.clear();
0595   digisinL5_.clear();
0596   digisinL6_.clear();
0597   digisinL7_.clear();
0598   alldigis_.clear();
0599   allMatches_->clear();
0600   prelimMatches_->clear();
0601   finalMatches_->clear();
0602 }
0603 
0604 void PseudoBayesGrouping::finish() {
0605   if (debug_)
0606     LogDebug("PseudoBayesGrouping") << "PseudoBayesGrouping: finish";
0607 };