ParticleUtils.hh

#ifndef RIVET_PARTICLEUTILS_HH
#define RIVET_PARTICLEUTILS_HH
 
#include "Rivet/Particle.hh"
#include "Rivet/Tools/ParticleBaseUtils.hh"
#include "Rivet/Tools/ParticleIdUtils.hh"
 
// Macros to map Rivet::Particle functions to PID:: functions of the same name
#define PARTICLE_TO_PID_BOOLFN(fname) inline bool fname (const Particle& p) { return PID:: fname (p.pid()); }
#define PARTICLE_TO_PID_INTFN(fname) inline int fname (const Particle& p) { return PID:: fname (p.pid()); }
#define PARTICLE_TO_PID_DBLFN(fname) inline double fname (const Particle& p) { return PID:: fname (p.pid()); }
 
namespace Rivet {
 



  inline int pid(const Particle& p) { return p.pid(); }

  inline int abspid(const Particle& p) { return p.abspid(); }
 

  PARTICLE_TO_PID_BOOLFN(isCharged)
 
  
  PARTICLE_TO_PID_BOOLFN(isNeutral)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isNeutrino)
 
  
  PARTICLE_TO_PID_BOOLFN(isChargedLepton)
 
  
  PARTICLE_TO_PID_BOOLFN(isLepton)
 
  
  PARTICLE_TO_PID_BOOLFN(isPhoton)
 
  
  PARTICLE_TO_PID_BOOLFN(isElectron)
 
  
  PARTICLE_TO_PID_BOOLFN(isMuon)
 
  
  PARTICLE_TO_PID_BOOLFN(isTau)
 
  
  PARTICLE_TO_PID_BOOLFN(isHadron)
 
  
  PARTICLE_TO_PID_BOOLFN(isMeson)
 
  
  PARTICLE_TO_PID_BOOLFN(isBaryon)
 
  
  PARTICLE_TO_PID_BOOLFN(isQuark)
 
  
  PARTICLE_TO_PID_BOOLFN(isParton)
 
 
 
  
  PARTICLE_TO_PID_BOOLFN(isWplus)
 
  
  PARTICLE_TO_PID_BOOLFN(isWminus)
 
  
  PARTICLE_TO_PID_BOOLFN(isW)
 
  
  PARTICLE_TO_PID_BOOLFN(isZ)
 
  
  PARTICLE_TO_PID_BOOLFN(isHiggs)
 
  
  PARTICLE_TO_PID_BOOLFN(isStrange)
 
  
  PARTICLE_TO_PID_BOOLFN(isCharm)
 
  
  PARTICLE_TO_PID_BOOLFN(isBottom)
 
  
  PARTICLE_TO_PID_BOOLFN(isTop)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isHeavyFlavour)
 
  
  PARTICLE_TO_PID_BOOLFN(isHeavyParton)
 
  
  PARTICLE_TO_PID_BOOLFN(isLightParton)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isHeavyMeson)
 
  
  PARTICLE_TO_PID_BOOLFN(isHeavyBaryon)
 
  
  PARTICLE_TO_PID_BOOLFN(isHeavyHadron)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isLightMeson)
 
  
  PARTICLE_TO_PID_BOOLFN(isLightBaryon)
 
  
  PARTICLE_TO_PID_BOOLFN(isLightHadron)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isBottomMeson)
 
  
  PARTICLE_TO_PID_BOOLFN(isBottomBaryon)
 
  
  PARTICLE_TO_PID_BOOLFN(isBottomHadron)
 
 
  
  PARTICLE_TO_PID_BOOLFN(isCharmMeson)
 
  
  PARTICLE_TO_PID_BOOLFN(isCharmBaryon)
 
  
  PARTICLE_TO_PID_BOOLFN(isCharmHadron)
 
 
  // /// Determine if the PID is that of a strange meson
  // PARTICLE_TO_PID_BOOLFN(isStrangeMeson)
 
  // /// Determine if the PID is that of a strange baryon
  // PARTICLE_TO_PID_BOOLFN(isStrangeBaryon)
 
  // /// Determine if the PID is that of a strange hadron
  // PARTICLE_TO_PID_BOOLFN(isStrangeHadron)
 
 
 
  
  PARTICLE_TO_PID_BOOLFN(isReggeon)
 
  
  PARTICLE_TO_PID_BOOLFN(isDiquark)
 
  
  PARTICLE_TO_PID_BOOLFN(isPentaquark)
 
  
  PARTICLE_TO_PID_BOOLFN(isSUSY)
 
  
  PARTICLE_TO_PID_BOOLFN(isRhadron)
 
  
  PARTICLE_TO_PID_BOOLFN(isTechnicolor)
 
  
  PARTICLE_TO_PID_BOOLFN(isExcited)
 
  
  PARTICLE_TO_PID_BOOLFN(isKK)
 
  
  PARTICLE_TO_PID_BOOLFN(isGraviton)
 
  
  PARTICLE_TO_PID_BOOLFN(isBSM)
 
 
 
  
  PARTICLE_TO_PID_BOOLFN(isGenSpecific)
 
  
  PARTICLE_TO_PID_BOOLFN(isResonance)
 
  
  PARTICLE_TO_PID_BOOLFN(isTransportable)
 
 
 
  
  PARTICLE_TO_PID_BOOLFN(hasUp)
 
  
  PARTICLE_TO_PID_BOOLFN(hasDown)
 
  
  PARTICLE_TO_PID_BOOLFN(hasStrange)
 
  
  PARTICLE_TO_PID_BOOLFN(hasCharm)
 
  
  PARTICLE_TO_PID_BOOLFN(hasBottom)
 
  
  PARTICLE_TO_PID_BOOLFN(hasTop)
 
 
 
  
  PARTICLE_TO_PID_INTFN(jSpin)
 
  
  PARTICLE_TO_PID_INTFN(sSpin)
 
  
  PARTICLE_TO_PID_INTFN(lSpin)
 
 
  
  PARTICLE_TO_PID_DBLFN(charge)
 
  
  PARTICLE_TO_PID_INTFN(charge3)
 
  
  PARTICLE_TO_PID_DBLFN(abscharge)
 
  
  PARTICLE_TO_PID_INTFN(abscharge3)
 
  
  PARTICLE_TO_PID_INTFN(nuclZ)
 
  
  PARTICLE_TO_PID_INTFN(nuclA)
 
  
  PARTICLE_TO_PID_INTFN(nuclNlambda)
 
  
 
 
  
 
  inline bool isSameSign(const Particle& a, const Particle& b) { return PID::isSameSign(a.pid(), b.pid()); }
  inline bool isOppSign(const Particle& a, const Particle& b) { return PID::isOppSign(a.pid(), b.pid()); }
  inline bool isSameFlav(const Particle& a, const Particle& b) { return PID::isSameFlav(a.pid(), b.pid()); }
  inline bool isOppFlav(const Particle& a, const Particle& b) { return PID::isOppFlav(a.pid(), b.pid()); }
 
  inline bool isOSSF(const Particle& a, const Particle& b) { return PID::isOSSF(a.pid(), b.pid()); }
  inline bool isSSSF(const Particle& a, const Particle& b) { return PID::isSSSF(a.pid(), b.pid()); }
  inline bool isOSOF(const Particle& a, const Particle& b) { return PID::isOSOF(a.pid(), b.pid()); }
  inline bool isSSOF(const Particle& a, const Particle& b) { return PID::isSSOF(a.pid(), b.pid()); }

 


  inline bool oppSign(const Particle& a, const Particle& b) {
    return sign(a.charge3()) == -sign(b.charge3()) && sign(a.charge3()) != ZERO;
  }

  inline bool sameSign(const Particle& a, const Particle& b) {
    return sign(a.charge3()) == sign(b.charge3());
  }

  inline bool oppCharge(const Particle& a, const Particle& b) {
    return a.charge3() == -b.charge3() && a.charge3() != 0;
  }

  inline bool sameCharge(const Particle& a, const Particle& b) {
    return a.charge3() == b.charge3();
  }

  inline bool diffCharge(const Particle& a, const Particle& b) {
    return a.charge3() != b.charge3();
  }

 
 

 
 


  inline bool isFirstWith(const Particle& p, const ParticleSelector& f) {
    return p.isFirstWith(f);
  }

  inline bool isFirstWithout(const Particle& p, const ParticleSelector& f) {
    return p.isFirstWithout(f);
  }
 

  inline bool isLastWith(const Particle& p, const ParticleSelector& f) {
    return p.isLastWith(f);
  }

  inline bool isLastWithout(const Particle& p, const ParticleSelector& f) {
    return p.isLastWithout(f);
  }
 
 

  inline bool hasAncestorWith(const Particle& p, const ParticleSelector& f, bool only_physical=true) {
    return p.hasAncestorWith(f, only_physical);
  }

  inline bool hasAncestorWithout(const Particle& p, const ParticleSelector& f, bool only_physical=true) {
    return p.hasAncestorWithout(f, only_physical);
  }
 

  inline bool hasParentWith(const Particle& p, const ParticleSelector& f) {
    return p.hasParentWith(f);
  }

  inline bool hasParentWithout(const Particle& p, const ParticleSelector& f) {
    return p.hasParentWithout(f);
  }
 

  inline bool hasChildWith(const Particle& p, const ParticleSelector& f) {
    return p.hasChildWith(f);
  }

  inline bool hasChildWithout(const Particle& p, const ParticleSelector& f) {
    return p.hasChildWithout(f);
  }
 

  inline bool hasDescendantWith(const Particle& p, const ParticleSelector& f, bool remove_duplicates=true) {
    return p.hasDescendantWith(f, remove_duplicates);
  }

  inline bool hasDescendantWithout(const Particle& p, const ParticleSelector& f, bool remove_duplicates=true) {
    return p.hasDescendantWithout(f, remove_duplicates);
  }
 

  inline bool hasStableDescendantWith(const Particle& p, const ParticleSelector& f) {
    return p.hasStableDescendantWith(f);
  }

  inline bool hasStableDescendantWithout(const Particle& p, const ParticleSelector& f) {
    return p.hasStableDescendantWithout(f);
  }
 
 

  inline bool isVisible(const Particle& p) { return p.isVisible(); }

  inline bool isDirect(const Particle& p, bool allow_from_direct_tau=false, bool allow_from_direct_mu=false) {
    return p.isDirect(allow_from_direct_tau, allow_from_direct_mu);
  }

  inline bool isPrompt(const Particle& p, bool allow_from_prompt_tau=false, bool allow_from_prompt_mu=false) {
    return p.isPrompt(allow_from_prompt_tau, allow_from_prompt_mu);
  }
 

  inline bool isStable(const Particle& p) { return p.isStable(); }

  inline bool hasHadronicDecay(const Particle& p) {
    if (p.isStable()) return false;
    if (p.hasChildWith(isHadron)) return true;
    return false;
  }

  inline bool hasLeptonicDecay(const Particle& p) {
    if (p.isStable()) return false;
    if (p.hasChildWith(isHadron)) return false;
    return true;
  }
 

  inline bool fromBottom(const Particle& p) { return p.fromBottom(); }

  inline bool fromCharm(const Particle& p) { return p.fromCharm(); }

  inline bool fromHadron(const Particle& p) { return p.fromHadron(); }

  inline bool fromTau(const Particle& p, bool prompt_taus_only=false) {
    return p.fromTau(prompt_taus_only);
  }

  inline bool fromPromptTau(const Particle& p) { return p.fromPromptTau(); }

 


  struct BoolParticleFunctor {
    virtual bool operator()(const Particle& p) const = 0;
    virtual ~BoolParticleFunctor() {}
  };

  struct BoolParticleAND : public BoolParticleFunctor {
    BoolParticleAND(const std::vector<ParticleSelector>& sels) : selectors(sels) {}
    BoolParticleAND(const ParticleSelector& a, const ParticleSelector& b) : selectors({a,b}) {}
    BoolParticleAND(const ParticleSelector& a, const ParticleSelector& b, const ParticleSelector& c) : selectors({a,b,c}) {}
    bool operator()(const Particle& p) const {
      for (const ParticleSelector& sel : selectors) if (!sel(p)) return false;
      return true;
    }
    std::vector<ParticleSelector> selectors;
  };

  inline BoolParticleAND operator && (const ParticleSelector& a, const ParticleSelector& b) {
    return BoolParticleAND(a, b);
  }
 

  struct BoolParticleOR : public BoolParticleFunctor {
    BoolParticleOR(const std::vector<ParticleSelector>& sels) : selectors(sels) {}
    BoolParticleOR(const ParticleSelector& a, const ParticleSelector& b) : selectors({a,b}) {}
    BoolParticleOR(const ParticleSelector& a, const ParticleSelector& b, const ParticleSelector& c) : selectors({a,b,c}) {}
    bool operator()(const Particle& p) const {
      for (const ParticleSelector& sel : selectors) if (sel(p)) return true;
      return false;
    }
    std::vector<ParticleSelector> selectors;
  };

  inline BoolParticleOR operator || (const ParticleSelector& a, const ParticleSelector& b) {
    return BoolParticleOR(a, b);
  }

  struct BoolParticleNOT : public BoolParticleFunctor {
    BoolParticleNOT(const ParticleSelector& sel) : selector(sel) {}
    bool operator()(const Particle& p) const { return !selector(p); }
    ParticleSelector selector;
  };

  inline BoolParticleNOT operator ! (const ParticleSelector& a) {
    return BoolParticleNOT(a);
  }
 

  struct HasPID : public BoolParticleFunctor {
    HasPID(PdgId pid) : targetpids{pid} { }
    HasPID(vector<PdgId> pids) : targetpids{pids} { }
    HasPID(initializer_list<PdgId> pids) : targetpids{pids} { }
    bool operator()(const Particle& p) const { return contains(targetpids, p.pid()); }
    vector<PdgId> targetpids;
  };
  using hasPID = HasPID;

  struct HasAbsPID : public BoolParticleFunctor {
    HasAbsPID(PdgId pid) : targetapids{abs(pid)} { }
    HasAbsPID(vector<PdgId> pids) { for (PdgId pid : pids) targetapids.push_back(abs(pid)); }
    HasAbsPID(initializer_list<PdgId> pids) { for (PdgId pid : pids) targetapids.push_back(abs(pid)); }
    bool operator()(const Particle& p) const { return contains(targetapids, p.abspid()); }
    vector<PdgId> targetapids;
  };
  using hasAbsPID = HasAbsPID;
 

  struct FirstParticleWith : public BoolParticleFunctor {
    FirstParticleWith(const ParticleSelector& f) : fn(f) { }
    FirstParticleWith(const Cut& c);
    bool operator()(const Particle& p) const { return isFirstWith(p, fn); }
    ParticleSelector fn;
  };
  using firstParticleWith = FirstParticleWith;

  struct FirstParticleWithout : public BoolParticleFunctor {
    FirstParticleWithout(const ParticleSelector& f) : fn(f) { }
    FirstParticleWithout(const Cut& c);
    bool operator()(const Particle& p) const { return isFirstWithout(p, fn); }
    ParticleSelector fn;
  };
  using firstParticleWithout = FirstParticleWithout;
 

  struct LastParticleWith : public BoolParticleFunctor {
    template <typename FN>
    LastParticleWith(const FN& f) : fn(f) { }
    LastParticleWith(const Cut& c);
    bool operator()(const Particle& p) const { return isLastWith(p, fn); }
    std::function<bool(const Particle&)> fn;
  };
  using lastParticleWith = LastParticleWith;

  struct LastParticleWithout : public BoolParticleFunctor {
    LastParticleWithout(const ParticleSelector& f) : fn(f) { }
    LastParticleWithout(const Cut& c);
    bool operator()(const Particle& p) const { return isLastWithout(p, fn); }
    ParticleSelector fn;
  };
  using lastParticleWithout = LastParticleWithout;
 

  struct HasParticleAncestorWith : public BoolParticleFunctor {
    HasParticleAncestorWith(const ParticleSelector& f, bool only_physical=true) : fn(f), onlyphysical(only_physical) { }
    HasParticleAncestorWith(const Cut& c, bool only_physical=true);
    bool operator()(const Particle& p) const { return hasAncestorWith(p, fn, onlyphysical); }
    ParticleSelector fn;
    bool onlyphysical;
  };
  using hasParticleAncestorWith = HasParticleAncestorWith;

  struct HasParticleAncestorWithout : public BoolParticleFunctor {
    HasParticleAncestorWithout(const ParticleSelector& f, bool only_physical=true) : fn(f), onlyphysical(only_physical) { }
    HasParticleAncestorWithout(const Cut& c, bool only_physical=true);
    bool operator()(const Particle& p) const { return hasAncestorWithout(p, fn, onlyphysical); }
    ParticleSelector fn;
    bool onlyphysical;
  };
  using hasParticleAncestorWithout = HasParticleAncestorWithout;
 

  struct HasParticleParentWith : public BoolParticleFunctor {
    HasParticleParentWith(const ParticleSelector& f) : fn(f) { }
    HasParticleParentWith(const Cut& c);
    bool operator()(const Particle& p) const { return hasParentWith(p, fn); }
    ParticleSelector fn;
  };
  using hasParticleParentWith = HasParticleParentWith;

  struct HasParticleParentWithout : public BoolParticleFunctor {
    HasParticleParentWithout(const ParticleSelector& f) : fn(f) { }
    HasParticleParentWithout(const Cut& c);
    bool operator()(const Particle& p) const { return hasParentWithout(p, fn); }
    ParticleSelector fn;
  };
  using hasParticleParentWithout = HasParticleParentWithout;
 

  struct HasParticleChildWith : public BoolParticleFunctor {
    HasParticleChildWith(const ParticleSelector& f) : fn(f) { }
    HasParticleChildWith(const Cut& c);
    bool operator()(const Particle& p) const { return hasChildWith(p, fn); }
    ParticleSelector fn;
  };
  using hasParticleChildWith = HasParticleChildWith;

  struct HasParticleChildWithout : public BoolParticleFunctor {
    HasParticleChildWithout(const ParticleSelector& f) : fn(f) { }
    HasParticleChildWithout(const Cut& c);
    bool operator()(const Particle& p) const { return hasChildWithout(p, fn); }
    ParticleSelector fn;
  };
  using hasParticleChildWithout = HasParticleChildWithout;
 

  struct HasParticleDescendantWith : public BoolParticleFunctor {
    HasParticleDescendantWith(const ParticleSelector& f, bool remove_duplicates=true) : fn(f), rmduplicates(remove_duplicates) { }
    HasParticleDescendantWith(const Cut& c, bool remove_duplicates=true);
    bool operator()(const Particle& p) const { return hasDescendantWith(p, fn, rmduplicates); }
    ParticleSelector fn;
    bool rmduplicates;
  };
  using hasParticleDescendantWith = HasParticleDescendantWith;

  struct HasParticleDescendantWithout : public BoolParticleFunctor {
    HasParticleDescendantWithout(const ParticleSelector& f, bool remove_duplicates=true) : fn(f), rmduplicates(remove_duplicates) { }
    HasParticleDescendantWithout(const Cut& c, bool remove_duplicates=true);
    bool operator()(const Particle& p) const { return hasDescendantWithout(p, fn, rmduplicates); }
    ParticleSelector fn;
    bool rmduplicates;
  };
  using hasParticleDescendantWithout = HasParticleDescendantWithout;

 


  Particles& iselect(Particles& particles, const Cut& c);

  inline Particles select(const Particles& particles, const Cut& c) {
    Particles rtn = particles;
    return iselect(rtn, c);
  }

  inline Particles select(const Particles& particles, const Cut& c, Particles& out) {
    out = select(particles, c);
    return out;
  }

   Particles& idiscard(Particles& particles, const Cut& c);

  inline Particles discard(const Particles& particles, const Cut& c) {
    Particles rtn = particles;
    return idiscard(rtn, c);
  }

  inline Particles discard(const Particles& particles, const Cut& c, Particles& out) {
    out = discard(particles, c);
    return out;
  }
 
 
  // inline void ifilterIsolateDeltaR(Particles& particles, const FourMomenta& vecs) {
  //   ifilter_discard(particles,
  // }
 
 
  // inline Particles filterIsolateDeltaR(const Particles& particles, const FourMomenta& vecs) {
  // }

 
 


  inline PdgIdPair pids(const ParticlePair& pp) {
    return make_pair(pp.first.pid(), pp.second.pid());
  }
 
  namespace Kin {
  
    inline pair<double,double> energies(const ParticlePair& pp) {
      return make_pair(pp.first.E(), pp.second.E());
    }
   
    inline pair<FourMomentum,FourMomentum> moms(const ParticlePair& pp) {
      return make_pair(pp.first.mom(), pp.second.mom());
    }

    inline double mass(const ParticlePair& pp) {
      return mass(pp.first, pp.second);
    }

    inline double mass2(const ParticlePair& pp) {
      return mass2(pp.first, pp.second);
    }

    inline double mT(const ParticlePair& pp) {
      return mT(pp.first, pp.second);
    }

    inline double pT(const ParticlePair& pp) {
      return pT(pp.first, pp.second);
    }
 
  }  
  
 

 
  namespace Kin {

    inline double mass(const Particle& p, const FourMomentum& p4) {
      return mass(p.mom(), p4);
    }

    inline double mass(const FourMomentum& p4, const Particle& p) {
      return mass(p4, p.mom());
    }

    inline double mass(const Particle& p1, const Particle& p2) {
      return mass(p1.mom(), p2.mom());
    }

    inline double mass2(const Particle& p, const P4& p4) {
      return mass2(p.mom(), p4);
    }

    inline double mass2(const P4& p4, const Particle& p) {
      return mass2(p4, p.mom());
    }

    inline double mass2(const Particle& p1, const Particle& p2) {
      return mass2(p1.mom(), p2.mom());
    }

    inline double pT(const Particle& p, const P4& p4) {
      return pT(p.mom(), p4);
    }

    inline double pT(const P4& p4, const Particle& p) {
      return pT(p4, p.mom());
    }

    inline double pT(const Particle& p1, const Particle& p2) {
      return pT(p1.mom(), p2.mom());
    }
 
  }

 

  namespace Kin {

    inline double pT(const Particle& p) {
      return p.pT();
    }
 
    inline double sumPt(const Particles& ps) {
      return sum(ps, Kin::pT, 0.0);
    }
 
    inline FourMomentum sumP4(const Particles& ps) {
      return sum(ps, Kin::p4, FourMomentum());
    }
 
    inline Vector3 sumP3(const Particles& ps) {
      return sum(ps, Kin::p3, Vector3());
    }


 
  }
 
 
  // Import Kin namespace into Rivet
  using namespace Kin;
 

  inline bool isSame(const Particle& a, const Particle& b) {
    return a.isSame(b);
  }
  
 

  inline bool containsPID(const Particles& parts, int id, bool abs=false) {
    if (abs) return any(parts, HasAbsPID(id));
    return any(parts, HasPID(id));
  }
 
 

  inline bool isRadiative(const Particle& part) {
    const Particles& parents = part.parents();
    if (parents.size() != 1)
      return false;
 
    const Particle& mother = parents[0];
    return ( part.pid() != mother.pid() ) && ( containsPID(mother.children(), mother.pid()) );
  }
 

  bool cascadeContains( const Particles& parts,
            const vector<int>& pids,
                bool absolute,
            bool ignorephoton);

 
}
 
#endif