KinematicFitPrimitives/src/KinematicTree.cc

0001 #include "RecoVertex/KinematicFitPrimitives/interface/KinematicTree.h"
0002
0003 KinematicTree::KinematicTree() {
0004   empt = true;
0005   treeWalker = nullptr;
0006 }
0007
0008 KinematicTree::~KinematicTree() { delete treeWalker; }
0009
0010 bool KinematicTree::isEmpty() const { return empt; }
0011
0012 bool KinematicTree::isConsistent() const {
0013   movePointerToTheTop();
0014   bool des = false;
0015   if (!treeWalker->nextSibling())
0016     des = true;
0017   return des;
0018 }
0019
0020 void KinematicTree::addParticle(RefCountedKinematicVertex prodVtx,
0021                                 RefCountedKinematicVertex decVtx,
0022                                 RefCountedKinematicParticle part) {
0023   part->setTreePointer(this);
0024   treeGraph.addEdge(prodVtx, decVtx, part);
0025   empt = false;
0026   movePointerToTheTop();
0027   prodVtx->setTreePointer(this);
0028   decVtx->setTreePointer(this);
0029 }
0030
0031 std::vector<RefCountedKinematicParticle> KinematicTree::finalStateParticles() const {
0032   if (isEmpty() || !(isConsistent())) {
0033     throw VertexException("KinematicTree::finalStateParticles; tree is empty or not consistent");
0034   } else {
0035     RefCountedKinematicParticle initial = currentParticle();
0036     std::vector<RefCountedKinematicParticle> rs;
0037     movePointerToTheTop();
0038     if (!(leftFinalParticle())) {
0039       std::cout << "top particle has no daughters, empty vector returned" << std::endl;
0040     } else {
0041       //now pointer is at the  most left final particle
0042       rs.push_back(currentParticle());
0043       bool next_right = true;
0044       bool up = true;
0045       do {
0046         next_right = movePointerToTheNextChild();
0047         if (next_right) {
0048           //if there's a way to the right,
0049           //we go right and down possible
0050           leftFinalParticle();
0051           rs.push_back(currentParticle());
0052         } else {
0053           //once there's no way to right anymore
0054           //trying to find a way upper
0055           up = movePointerToTheMother();
0056         }
0057         //loop stops when we are at the top:
0058         //no way up, no way to the right
0059       } while (up);
0060     }
0061     //getting the pointer back
0062     bool back = findParticle(initial);
0063     if (!back)
0064       throw VertexException("KinematicTree::FinalStateParticles; error occured while getting back");
0065     return rs;
0066   }
0067 }
0068
0069 bool KinematicTree::leftFinalParticle() const {
0070   bool res = false;
0071   if (movePointerToTheFirstChild()) {
0072     res = true;
0073     bool next = true;
0074     do {
0075       next = movePointerToTheFirstChild();
0076     } while (next);
0077   } else {
0078     res = false;
0079   }
0080   return res;
0081 }
0082
0083 RefCountedKinematicParticle KinematicTree::topParticle() const {
0084   if (isEmpty())
0085     throw VertexException("KinematicTree::topParticle; tree is empty!");
0086   //putting pointer to the top of the tree
0087   movePointerToTheTop();
0088   return treeWalker->current().second;
0089 }
0090
0091 RefCountedKinematicVertex KinematicTree::currentDecayVertex() const {
0092   if (isEmpty())
0093     throw VertexException("KinematicTree::currentDecayVertex; tree is empty!");
0094   return treeWalker->current().first;
0095 }
0096
0097 std::pair<bool, RefCountedKinematicParticle> KinematicTree::motherParticle() const {
0098   if (isEmpty())
0099     throw VertexException("KinematicTree::motherParticle; tree is empty!");
0100   bool top = currentProductionVertex()->vertexIsValid();
0101   RefCountedKinematicParticle cr = treeWalker->current().second;
0102   bool up = treeWalker->parent();
0103   std::pair<bool, RefCountedKinematicParticle> res;
0104   if (up && top) {
0105     RefCountedKinematicParticle pr = treeWalker->current().second;
0106
0107     //now putting the pointer back
0108     bool fc = treeWalker->firstChild();
0109     if (!fc)
0110       throw VertexException("KinematicTree::motherParticle; tree is incorrect!");
0111     if (*(treeWalker->current().second) != *cr) {
0112       do {
0113         bool nx = treeWalker->nextSibling();
0114         if (!nx)
0115           throw VertexException("KinematicTree::motherParticle; tree is incorrect!");
0116       } while (*(treeWalker->current().second) != *cr);
0117     }
0118     res = std::pair<bool, RefCountedKinematicParticle>(true, pr);
0119     return res;
0120   } else {
0121     RefCountedKinematicParticle fk;
0122     return std::pair<bool, RefCountedKinematicParticle>(false, fk);
0123   }
0124 }
0125
0126 std::vector<RefCountedKinematicParticle> KinematicTree::daughterParticles() const {
0127   if (isEmpty())
0128     throw VertexException("KinematicTree::daughterParticles; tree is empty!");
0129   std::vector<RefCountedKinematicParticle> sResult;
0130   RefCountedKinematicParticle initial = currentParticle();
0131   bool down = treeWalker->firstChild();
0132   if (down) {
0133     sResult.push_back(treeWalker->current().second);
0134     bool sibling = true;
0135     do {
0136       sibling = treeWalker->nextSibling();
0137       if (sibling)
0138         sResult.push_back(treeWalker->current().second);
0139     } while (sibling);
0140   }
0141
0142   //getting the pointer back to the mother
0143   bool back = findParticle(initial);
0144   if (!back)
0145     throw VertexException("KinematicTree::daughterParticles; error occured while getting back");
0146   return sResult;
0147 }
0148
0149 void KinematicTree::movePointerToTheTop() const {
0150   if (isEmpty())
0151     throw VertexException("KinematicTree::movePointerToTheTop; tree is empty!");
0152   delete treeWalker;
0153   treeWalker = new math::GraphWalker<RefCountedKinematicVertex, RefCountedKinematicParticle>(treeGraph);
0154   //now pointer is a pair: fake vertex and
0155   //icoming 0 pointer to the particle
0156   //moving it to decayed particle
0157   bool move = treeWalker->firstChild();
0158   if (!move)
0159     throw VertexException("KinematicTree::movePointerToTheTop; non consistent tree?");
0160 }
0161
0162 RefCountedKinematicVertex KinematicTree::currentProductionVertex() const {
0163   if (isEmpty())
0164     throw VertexException("KinematicTree::currentProductionVertex; tree is empty!");
0165   //current particle
0166   RefCountedKinematicParticle initial = currentParticle();
0167
0168   bool up;
0169   bool down;
0170   RefCountedKinematicVertex res;
0171   up = movePointerToTheMother();
0172
0173   if (up) {
0174     res = treeWalker->current().first;
0175
0176     //pointer moved so we going back
0177     down = treeWalker->firstChild();
0178
0179     //_down_ variable is always TRUE here, if
0180     //the tree is valid.
0181     if (down) {
0182       if (initial == treeWalker->current().second) {
0183         return res;
0184       } else {
0185         bool next = true;
0186         do {
0187           next = treeWalker->nextSibling();
0188           if (treeWalker->current().second == initial)
0189             next = false;
0190         } while (next);
0191         return res;
0192       }
0193     } else {
0194       throw VertexException("KinematicTree::Navigation failed, tree invalid?");
0195     }
0196   } else {
0197     //very unprobable case. This efectively means that user is
0198     //already out of the tree. Moving back to the top
0199     delete treeWalker;
0200     treeWalker = new math::GraphWalker<RefCountedKinematicVertex, RefCountedKinematicParticle>(treeGraph);
0201     res = treeWalker->current().first;
0202     //now pointer is a pair: fake vertex and
0203     //icoming 0 pointer to the particle
0204     //moving it to decayed particle
0205     bool move = treeWalker->firstChild();
0206     if (!move)
0207       throw VertexException("KinematicTree::movePointerToTheTop; non consistent tree?");
0208     return res;
0209   }
0210 }
0211
0212 RefCountedKinematicParticle KinematicTree::currentParticle() const {
0213   if (isEmpty())
0214     throw VertexException("KinematicTree::currentParticle; tree is empty!");
0215   return treeWalker->current().second;
0216 }
0217
0218 bool KinematicTree::movePointerToTheMother() const {
0219   if (isEmpty())
0220     throw VertexException("KinematicTree::movePointerToTheMother; tree is empty!");
0221   bool up = treeWalker->parent();
0222   bool cr = treeWalker->current().first->vertexIsValid();
0223   return (up && cr);
0224 }
0225
0226 bool KinematicTree::movePointerToTheFirstChild() const {
0227   if (isEmpty())
0228     throw VertexException("KinematicTree::movePointerToTheFirstChild; tree is empty!");
0229   return treeWalker->firstChild();
0230 }
0231
0232 bool KinematicTree::movePointerToTheNextChild() const {
0233   if (isEmpty())
0234     throw VertexException("KinematicTree::movePointerToTheNextChild; tree is empty!");
0235   bool res = treeWalker->nextSibling();
0236   return res;
0237 }
0238
0239 bool KinematicTree::findParticle(const RefCountedKinematicParticle part) const {
0240   if (isEmpty() || !(isConsistent())) {
0241     throw VertexException("KinematicTree::findParticle; tree is empty or not consistent");
0242   } else {
0243     bool res = false;
0244     movePointerToTheTop();
0245     if (currentParticle() == part) {
0246       res = true;
0247     } else if (leftBranchSearch(part)) {
0248       res = true;
0249     } else {
0250       bool found = false;
0251       bool up = true;
0252       bool next_right = false;
0253       do {
0254         //    if(*(currentParticle()) == *part) found = true;
0255         next_right = movePointerToTheNextChild();
0256         if (next_right) {
0257           found = leftBranchSearch(part);
0258         } else {
0259           up = movePointerToTheMother();
0260           if (currentParticle() == part)
0261             found = true;
0262         }
0263       } while (up && !found);
0264       res = found;
0265     }
0266     return res;
0267   }
0268 }
0269
0270 bool KinematicTree::leftBranchSearch(RefCountedKinematicParticle part) const {
0271   bool found = false;
0272   bool next = true;
0273   if (currentParticle() == part) {
0274     found = true;
0275   } else {
0276     do {
0277       next = movePointerToTheFirstChild();
0278       if (currentParticle() == part) {
0279         found = true;
0280       }
0281     } while (next && !found);
0282   }
0283   return found;
0284 }
0285
0286 bool KinematicTree::findDecayVertex(const RefCountedKinematicVertex vert) const {
0287   if (isEmpty() || !(isConsistent())) {
0288     throw VertexException("KinematicTree::findParticle; tree is empty or not consistent");
0289   } else {
0290     bool res = false;
0291     movePointerToTheTop();
0292     if (currentDecayVertex() == vert) {
0293       res = true;
0294     } else if (leftBranchVertexSearch(vert)) {
0295       res = true;
0296     } else {
0297       bool up = true;
0298       bool fnd = false;
0299       do {
0300         if (movePointerToTheNextChild()) {
0301           fnd = leftBranchVertexSearch(vert);
0302         } else {
0303           up = movePointerToTheMother();
0304           if (currentDecayVertex() == vert)
0305             fnd = true;
0306         }
0307       } while (up && !fnd);
0308       res = fnd;
0309     }
0310     return res;
0311   }
0312 }
0313
0314 bool KinematicTree::findDecayVertex(KinematicVertex *vert) const {
0315   if (isEmpty() || !(isConsistent())) {
0316     throw VertexException("KinematicTree::findParticle; tree is empty or not consistent");
0317   } else {
0318     bool res = false;
0319     movePointerToTheTop();
0320     if (*currentDecayVertex() == vert) {
0321       res = true;
0322     } else if (leftBranchVertexSearch(vert)) {
0323       res = true;
0324     } else {
0325       bool up = true;
0326       bool fnd = false;
0327       do {
0328         if (movePointerToTheNextChild()) {
0329           fnd = leftBranchVertexSearch(vert);
0330         } else {
0331           up = movePointerToTheMother();
0332           if (currentDecayVertex() == vert)
0333             fnd = true;
0334         }
0335       } while (up && !fnd);
0336       res = fnd;
0337     }
0338     return res;
0339   }
0340 }
0341
0342 bool KinematicTree::leftBranchVertexSearch(RefCountedKinematicVertex vtx) const {
0343   bool found = false;
0344   if (currentDecayVertex() == vtx) {
0345     found = true;
0346   } else {
0347     bool next = true;
0348     bool res = false;
0349     do {
0350       next = movePointerToTheFirstChild();
0351       if (currentDecayVertex() == vtx)
0352         res = true;
0353     } while (next && !res);
0354     found = res;
0355   }
0356   return found;
0357 }
0358
0359 void KinematicTree::leftBranchAdd(KinematicTree *otherTree, RefCountedKinematicVertex vtx) {
0360   RefCountedKinematicVertex previous_decay = otherTree->currentDecayVertex();
0361
0362   //children of current particle of the
0363   //other tree: in the end the whole
0364   //left branch should be added.
0365   bool next = true;
0366   do {
0367     next = otherTree->movePointerToTheFirstChild();
0368     if (next) {
0369       RefCountedKinematicParticle par = otherTree->currentParticle();
0370       RefCountedKinematicVertex current_decay = otherTree->currentDecayVertex();
0371       addParticle(previous_decay, current_decay, par);
0372       previous_decay = current_decay;
0373     }
0374   } while (next);
0375 }
0376
0377 void KinematicTree::replaceCurrentParticle(RefCountedKinematicParticle newPart) const {
0378   RefCountedKinematicParticle cDParticle = currentParticle();
0379   bool replace = treeGraph.replaceEdge(cDParticle, newPart);
0380   if (!replace)
0381     throw VertexException("KinematicTree::Particle To Replace not found");
0382 }
0383
0384 void KinematicTree::replaceCurrentVertex(RefCountedKinematicVertex newVert) const {
0385   RefCountedKinematicVertex cDVertex = currentDecayVertex();
0386   bool replace = treeGraph.replace(cDVertex, newVert);
0387   if (!replace)
0388     throw VertexException("KinematicTree::Vertex To Replace not found");
0389 }
0390
0391 void KinematicTree::addTree(RefCountedKinematicVertex vtx, KinematicTree *tr) {
0392   //adding new tree to the existing one:
0393   bool fnd = findDecayVertex(vtx);
0394   if (!fnd)
0395     throw VertexException("KinematicTree::addTree; Current tree does not contain the vertex passed");
0396   tr->movePointerToTheTop();
0397
0398   // adding the root of the tree:
0399   RefCountedKinematicParticle mP = tr->currentParticle();
0400   RefCountedKinematicVertex dec_vertex = tr->currentDecayVertex();
0401   addParticle(vtx, dec_vertex, mP);
0402
0403   // adding the left branch if any
0404   leftBranchAdd(tr, dec_vertex);
0405
0406   //now the pointer is at the left down
0407   //edge of the otherTree.
0408   //current tree pointer is where the
0409   //add operation was stoped last time.
0410   bool right = true;
0411   bool up = true;
0412   do {
0413     right = tr->movePointerToTheNextChild();
0414     if (right) {
0415       //production vertex is already at the current tree
0416       //adding current partilce
0417       RefCountedKinematicVertex prodVertex = tr->currentProductionVertex();
0418       RefCountedKinematicParticle cPart = tr->currentParticle();
0419       RefCountedKinematicVertex decVertex = tr->currentDecayVertex();
0420       addParticle(prodVertex, decVertex, cPart);
0421
0422       //adding the rest of the branch
0423       leftBranchAdd(tr, decVertex);
0424     } else {
0425       up = tr->movePointerToTheMother();
0426     }
0427   } while (up);
0428 }