StMcTrack.cc

#include <algorithm>
#include <assert.h>
#ifndef ST_NO_NAMESPACES
using std::find;
#endif

#include "StMcTrack.hh"

#include "StMcContainers.hh"
#include "StParticleTable.hh"
#include "StParticleDefinition.hh"

#include "tables/St_g2t_track_Table.h"
#include "tables/St_particle_Table.h"

static const char rcsid[] = "$Id: StMcTrack.cc,v 2.38 2016/05/16 23:47:09 perev Exp $";

ClassImp(StMcTrack);

StMcTrack::StMcTrack() 
{
    initToZero();
}



StMcTrack::StMcTrack(g2t_track_st* trk) {
    initToZero();
    mFourMomentum.setPx(trk->p[0]);
    mFourMomentum.setPy(trk->p[1]);
    mFourMomentum.setPz(trk->p[2]);
    mFourMomentum.setE(trk->e);
    mIsShower = trk->is_shower;
    mGeantId = trk->ge_pid;
    assert(mGeantId>=0 && mGeantId<=0xFFFF);
    mPdgId   = trk->eg_pid;
    mKey     = trk->id;
    mParticleDefinition = StParticleTable::instance()->findParticleByGeantId(trk->ge_pid);
    mEventGenLabel = trk->eg_label;
    mIsPrimary = kFALSE;
    // The information to fill the collections 
    // is not available directly from the tables.  
    // We need to decode from trk->hit_tpc_p, 
    // trk->hit_svt_p, trk->hit_ftp_p, and trk-next_vtx_trk_p 
    // the actual Collection objects, not just the id's stored in the table.
}

StMcTrack::StMcTrack(particle_st* trk) {
    initToZero();
    mFourMomentum.setPx(trk->phep[0]);
    mFourMomentum.setPy(trk->phep[1]);
    mFourMomentum.setPz(trk->phep[2]);
    mFourMomentum.setE(trk->phep[3]);
    mParticleDefinition = StParticleTable::instance()->findParticle(trk->idhep);
    mPdgId = trk->idhep;
    mIsPrimary = kFALSE;
    // This constructor is used for particles coming from the
    // particle table.
}



StMcTrack::~StMcTrack() {
    // Not owner, so we don't have to delete.
    mIntermediateVertices.clear();
    mTpcHits.clear();
    mSvtHits.clear();
    mSsdHits.clear();
    mFtpcHits.clear();
    mRichHits.clear();
    mCtbHits.clear();
    mBemcHits.clear();
    mBprsHits.clear();
    mBsmdeHits.clear();
    mBsmdpHits.clear();
    mTofHits.clear();
    mMtdHits.clear();
    mEemcHits.clear();
    mEprsHits.clear();
    mEsmduHits.clear();
    mEsmdvHits.clear();
    mPxlHits.clear();
    mIstHits.clear();
    mFgtHits.clear();
    mEtrHits.clear();
}

void StMcTrack::initToZero()
{
    
    mStartVertex  = 0;
    mStopVertex   = 0;
    mParticleDefinition = 0;
    mParent = 0;
    mIsShower        = 0;
    mGeantId         = 0;
    mPdgId           = 0;
    mKey             = 0;
    mEventGenLabel   = 0;
    mIsPrimary       = 0;
}

int StMcTrack::operator==(const StMcTrack& t) const
{
    return (t.mFourMomentum == mFourMomentum &&
            t.mStartVertex  == mStartVertex   &&
            t.mStopVertex   == mStopVertex );
}

int StMcTrack::operator!=(const StMcTrack& t) const
{
    return !(t == *this);
}

ostream&  operator<<(ostream& os, const StMcTrack& t)
{
    if (t.particleDefinition())
        os << "Particle      : " << t.particleDefinition()->name().c_str() << endl;
    os << "Four Momentum : " << t.fourMomentum() << endl; 
    os << "Pt            : " << t.pt() << endl;
    os << "Rapidity      : " << t.rapidity() << endl;
    os << "PseudoRapidity: " << t.pseudoRapidity()   << endl;
    os << "No. Tpc   Hits: " << t.tpcHits().size()   << endl;
    os << "No. Svt   Hits: " << t.svtHits().size()   << endl;
    os << "No. Ssd   Hits: " << t.ssdHits().size()   << endl;
    os << "No. Ftpc  Hits: " << t.ftpcHits().size()  << endl;
    os << "No. Rich  Hits: " << t.richHits().size()  << endl;
    os << "No. Ctb   Hits: " << t.ctbHits().size()   << endl;
    os << "No. Bemc  Hits: " << t.bemcHits().size()  << endl;
    os << "No. Bprs  Hits: " << t.bprsHits().size()  << endl;
    os << "No. Bsmde Hits: " << t.bsmdeHits().size() << endl;
    os << "No. Bsmdp Hits: " << t.bsmdpHits().size() << endl;
    os << "No. Tof   Hits: " << t.tofHits().size()   << endl;
    os << "No. Mtd   Hits: " << t.mtdHits().size()   << endl;
    os << "No. Eemc  Hits: " << t.eemcHits().size()  << endl;
    os << "No. Eprs  Hits: " << t.eprsHits().size()  << endl;
    os << "No. Esmdu Hits: " << t.esmduHits().size() << endl;
    os << "No. Esmdv Hits: " << t.esmdvHits().size() << endl;
    os << "No. Pxl Hits  : " << t.pxlHits().size()   << endl;
    os << "No. Ist Hits  : " << t.istHits().size()   << endl;
    os << "No. Fgt Hits  : " << t.fgtHits().size()   << endl;
    os << "No. Fsc Hits  : " << t.fscHits().size()   << endl;
    os << "No. Etr Hits  : " << t.etrHits().size()   << endl;
    os << "Is Shower     : " << t.isShower() << endl;
    os << "Geant Id      : " << t.geantId()  << endl;
    os << "Pdg Code      : " << t.pdgId()  << endl;
    os << "Event Gen. Lab: " << t.eventGenLabel()  << endl;
    os << "Key from g2t  : " << t.key()  << endl;
    
    return os;
}
void StMcTrack::Print(Option_t *option) const {
  TString opt(option);
  if (opt == "") {cout << *this << endl; return;}
  if (opt.Contains("desc",TString::kIgnoreCase)) {
    //       0        1         2         3         4
    //       1234567890123456789012345678901234567890
    cout << "Particle"
         << " Four Momentum                  " 
         << "Pt      " 
         << "Rapidity" 
         << "  Pseudo" 
         << " Id"  
         << "   Pdg"  
         << "Egl"  
         << "Key"
         << "Hit"
         << "Tpc" 
         << "Svt" 
         << "Ssd" 
         << "Fpc" 
         << "Rch" 
         << "Ctb" 
         << "Bmc" 
         << "Bpr" 
         << "Bse" 
         << "Bsp" 
         << "Tof" 
         << "Mtd" 
         << "Emc" 
         << "Eps" 
         << "Esu" 
         << "Esv" 
         << "Pxl" 
         << "Ist"  
         << "Fgt"  
         << "Etr"  
         << "ISh" << endl;
    return;
  }
  TString Name("        ");
  if (particleDefinition()) Name = particleDefinition()->name().c_str();
  Double_t eta = pseudoRapidity();
  if (TMath::Abs(eta) > 999.999) eta = TMath::Sign(999.999, eta);
  Double_t y = rapidity();
  if (TMath::Abs(y) > 999.999) y = TMath::Sign(999.999, y);
    cout << 
      Form("%8s%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%8.3f%3li%6li%3li%3li%6li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li%3li",
           Name.Data(),
           fourMomentum().x(), fourMomentum().y(), fourMomentum().z(), fourMomentum().t(), 
           pt(),
           y,
           eta,
           (long int)(geantId()),
           (long int)(pdgId()),
           (long int)(eventGenLabel()),
           (long int)(key()),
           (long int)(tpcHits().size()),
           (long int)(svtHits().size()),
           (long int)(ssdHits().size()),
           (long int)(ftpcHits().size()),
           (long int)(richHits().size()),
           (long int)(ctbHits().size()),
           (long int)(bemcHits().size()),
           (long int)(bprsHits().size()),
           (long int)(bsmdeHits().size()),
           (long int)(bsmdpHits().size()),
           (long int)(tofHits().size()),
           (long int)(mtdHits().size()),
           (long int)(eemcHits().size()),
           (long int)(eprsHits().size()),
           (long int)(esmduHits().size()),
           (long int)(esmdvHits().size()),
           (long int)(pxlHits().size()),
           (long int)(istHits().size()),
           (long int)(fgtHits().size()),
           (long int)(etrHits().size()),
           (long int)(isShower()))
         << endl;
}
void StMcTrack::setFourMomentum(const StLorentzVectorF& val) { mFourMomentum = val; }

void StMcTrack::setStartVertex(StMcVertex* val) { mStartVertex = val; }

void StMcTrack::setStopVertex(StMcVertex* val) { mStopVertex = val; }

void StMcTrack::setIntermediateVertices(StPtrVecMcVertex& val) { mIntermediateVertices = val; }

void StMcTrack::setTpcHits(StPtrVecMcTpcHit& val) { mTpcHits = val; }

void StMcTrack::setSvtHits(StPtrVecMcSvtHit& val) { mSvtHits = val; }

void StMcTrack::setSsdHits(StPtrVecMcSsdHit& val) { mSsdHits = val; }

void StMcTrack::setFtpcHits(StPtrVecMcFtpcHit& val) { mFtpcHits = val; }

void StMcTrack::setRichHits(StPtrVecMcRichHit& val) { mRichHits = val; }

void StMcTrack::setCtbHits(StPtrVecMcCtbHit& val) { mCtbHits = val; }

void StMcTrack::setBemcHits(StPtrVecMcCalorimeterHit& val) { mBemcHits = val; }

void StMcTrack::setBprsHits(StPtrVecMcCalorimeterHit& val) { mBprsHits = val; }

void StMcTrack::setBsmdeHits(StPtrVecMcCalorimeterHit& val) { mBsmdeHits = val; }

void StMcTrack::setBsmdpHits(StPtrVecMcCalorimeterHit& val) { mBsmdpHits = val; }

void StMcTrack::setTofHits(StPtrVecMcTofHit& val) { mTofHits = val; }

void StMcTrack::setMtdHits(StPtrVecMcMtdHit& val) { mMtdHits = val; }

void StMcTrack::setEemcHits(StPtrVecMcCalorimeterHit& val) { mEemcHits = val; }

void StMcTrack::setEprsHits(StPtrVecMcCalorimeterHit& val) { mEprsHits = val; }

void StMcTrack::setEsmduHits(StPtrVecMcCalorimeterHit& val) { mEsmduHits = val; }

void StMcTrack::setEsmdvHits(StPtrVecMcCalorimeterHit& val) { mEsmdvHits = val; }

void StMcTrack::setFscHits(StPtrVecMcCalorimeterHit& val) { mFscHits = val; }

void StMcTrack::setPxlHits(StPtrVecMcPxlHit& val) { mPxlHits = val; }

void StMcTrack::setIstHits(StPtrVecMcIstHit& val) { mIstHits = val; }

void StMcTrack::setFgtHits(StPtrVecMcFgtHit& val) { mFgtHits = val; }

void StMcTrack::setEtrHits(StPtrVecMcEtrHit& val) { mEtrHits = val; }

void StMcTrack::setShower(char val) { mIsShower = val; }

void StMcTrack::setGeantId(long val) { mGeantId = val; assert(mGeantId>=0 && mGeantId<=0xFFFF); }

void StMcTrack::setPdgId(long val) { mPdgId = val; }

void StMcTrack::setKey(long val) { mKey = val; }

void StMcTrack::setEventGenLabel(long val) { mEventGenLabel = val; }

void StMcTrack::setParent(StMcTrack* val) { mParent = val; }

void StMcTrack::addTpcHit(StMcTpcHit* hit)
{
  mTpcHits.push_back(hit);
}

void StMcTrack::addFtpcHit(StMcFtpcHit* hit)
{
  mFtpcHits.push_back(hit);
}

void StMcTrack::addSvtHit(StMcSvtHit* hit)
{
  mSvtHits.push_back(hit);
}

void StMcTrack::addSsdHit(StMcSsdHit* hit)
{
  mSsdHits.push_back(hit);
}

void StMcTrack::addRichHit(StMcRichHit* hit)
{
  mRichHits.push_back(hit);
}

void StMcTrack::addCtbHit(StMcCtbHit* hit)
{
  mCtbHits.push_back(hit);
}

void StMcTrack::addBemcHit(StMcCalorimeterHit* hit)
{
  mBemcHits.push_back(hit);
}

void StMcTrack::addBprsHit(StMcCalorimeterHit* hit)
{
  mBprsHits.push_back(hit);
}

void StMcTrack::addBsmdeHit(StMcCalorimeterHit* hit)
{
  mBsmdeHits.push_back(hit);
}

void StMcTrack::addBsmdpHit(StMcCalorimeterHit* hit)
{
  mBsmdpHits.push_back(hit);
}

void StMcTrack::addTofHit(StMcTofHit* hit)
{
  mTofHits.push_back(hit);
}

void StMcTrack::addMtdHit(StMcMtdHit* hit)
{
  mMtdHits.push_back(hit);
}

void StMcTrack::addEemcHit(StMcCalorimeterHit* hit)
{
  mEemcHits.push_back(hit);
}

void StMcTrack::addEprsHit(StMcCalorimeterHit* hit)
{
  mEprsHits.push_back(hit);
}

void StMcTrack::addEsmduHit(StMcCalorimeterHit* hit)
{
  mEsmduHits.push_back(hit);
}

void StMcTrack::addEsmdvHit(StMcCalorimeterHit* hit)
{
  mEsmdvHits.push_back(hit);
}

void StMcTrack::addFpdHit(StMcCalorimeterHit* hit)
{
  mFpdHits.push_back(hit);
}

void StMcTrack::addFscHit(StMcCalorimeterHit* hit)
{
  mFscHits.push_back(hit);
}

void StMcTrack::addPxlHit(StMcPxlHit* hit)
{
  mPxlHits.push_back(hit);
}

void StMcTrack::addIstHit(StMcIstHit* hit)
{
  mIstHits.push_back(hit);
}

void StMcTrack::addFgtHit(StMcFgtHit* hit)
{
  mFgtHits.push_back(hit);
}

void StMcTrack::addEtrHit(StMcEtrHit* hit)
{
  mEtrHits.push_back(hit);
}

// Not very elegant.  Maybe should have kept all collections as
// vector<StMcHit*> so that then we could have used the same
// routine for all of them... 
void StMcTrack::removeTpcHit(StMcTpcHit* hit)
{
    StMcTpcHitIterator iter = find (mTpcHits.begin(), mTpcHits.end(), hit);
    if (iter != mTpcHits.end()) {
        mTpcHits.erase(iter);
    }
}

void StMcTrack::removeFtpcHit(StMcFtpcHit* hit)
{
    StMcFtpcHitIterator iter = find (mFtpcHits.begin(), mFtpcHits.end(),hit);
    if (iter != mFtpcHits.end()) {
        mFtpcHits.erase(iter);
    }
}

void StMcTrack::removeSvtHit(StMcSvtHit* hit)
{
    StMcSvtHitIterator iter = find(mSvtHits.begin(), mSvtHits.end(), hit);
    if (iter != mSvtHits.end()) {
        mSvtHits.erase(iter);
    }
}

void StMcTrack::removeSsdHit(StMcSsdHit* hit)
{
    StMcSsdHitIterator iter = find(mSsdHits.begin(), mSsdHits.end(), hit);
    if (iter != mSsdHits.end()) {
        mSsdHits.erase(iter);
    }
}

void StMcTrack::removeRichHit(StMcRichHit* hit)
{
    StMcRichHitIterator iter = find(mRichHits.begin(), mRichHits.end(), hit);
    if (iter != mRichHits.end()){
        mRichHits.erase(iter);
    }
}

void StMcTrack::removeCtbHit(StMcCtbHit* hit)
{
    StMcCtbHitIterator iter = find(mCtbHits.begin(), mCtbHits.end(), hit);
    if (iter != mCtbHits.end()){
        mCtbHits.erase(iter);
    }
}

void StMcTrack::removeCalorimeterHit(StPtrVecMcCalorimeterHit& vch, StMcCalorimeterHit* hit)
{
    StMcCalorimeterHitIterator iter = find(vch.begin(), vch.end(), hit);
    if (iter != vch.end()) {
        vch.erase(iter);
    }
    
}
void StMcTrack::removeBemcHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mBemcHits, hit);
}
void StMcTrack::removeBprsHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mBprsHits, hit);
}
void StMcTrack::removeBsmdeHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mBsmdeHits, hit);
}
void StMcTrack::removeBsmdpHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mBsmdpHits, hit);
}

void StMcTrack::removeTofHit(StMcTofHit* hit)
{
    StMcTofHitIterator iter = find(mTofHits.begin(), mTofHits.end(), hit);
    if (iter != mTofHits.end()){
        mTofHits.erase(iter);
    }
}

void StMcTrack::removeMtdHit(StMcMtdHit* hit)
{
    StMcMtdHitIterator iter = find(mMtdHits.begin(), mMtdHits.end(), hit);
    if (iter != mMtdHits.end()){
        mMtdHits.erase(iter);
    }
}

void StMcTrack::removeEemcHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mEemcHits, hit);
}

void StMcTrack::removeEprsHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mEprsHits, hit);
}

void StMcTrack::removeEsmduHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mEsmduHits, hit);
}

void StMcTrack::removeEsmdvHit(StMcCalorimeterHit* hit)
{
    removeCalorimeterHit(mEsmdvHits, hit);
}

void StMcTrack::removePxlHit(StMcPxlHit* hit)
{
    StMcPxlHitIterator iter = find (mPxlHits.begin(), mPxlHits.end(),hit);
    if (iter != mPxlHits.end()) {
        mPxlHits.erase(iter);
    }    
}

void StMcTrack::removeIstHit(StMcIstHit* hit)
{
    StMcIstHitIterator iter = find (mIstHits.begin(), mIstHits.end(),hit);
    if (iter != mIstHits.end()) {
        mIstHits.erase(iter);
    }
}


void StMcTrack::removeFgtHit(StMcFgtHit* hit)
{
    StMcFgtHitIterator iter = find (mFgtHits.begin(), mFgtHits.end(),hit);
    if (iter != mFgtHits.end()) {
        mFgtHits.erase(iter);
    }    
}

void StMcTrack::removeEtrHit(StMcEtrHit* hit)
{
    StMcEtrHitIterator iter = find (mEtrHits.begin(), mEtrHits.end(),hit);
    if (iter != mEtrHits.end()) {
        mEtrHits.erase(iter);
    }    
}

//void StMcTrack::setTopologyMap(StTrackTopologyMap& val) { mTopologyMap = val; }
StParticleDefinition* StMcTrack::particleDefinition() { 
  if (mParticleDefinition) return mParticleDefinition; 
  if (mGeantId > 0) mParticleDefinition = StParticleTable::instance()->findParticleByGeantId(mGeantId);
  if (mParticleDefinition) return mParticleDefinition; 
  if (mPdgId)  mParticleDefinition = StParticleTable::instance()->findParticle(mPdgId);
  return mParticleDefinition;
}
//--------------------------------------------------------------------------------
const StPtrVecMcHit *StMcTrack::Hits(StDetectorId Id) const {
  StPtrVecMcHit *coll = 0;
  switch (Id) {
  case kTpcId:                  coll = (StPtrVecMcHit *) &mTpcHits; break;                      
  case kSvtId:                  coll = (StPtrVecMcHit *) &mSvtHits; break;          
  case kRichId:                 coll = (StPtrVecMcHit *) &mRichHits; break;                 
  case kFtpcWestId:             
  case kFtpcEastId:             coll = (StPtrVecMcHit *) &mFtpcHits; break;             
  case kTofId:                  coll = (StPtrVecMcHit *) &mTofHits; break;                       
  case kMtdId:                  coll = (StPtrVecMcHit *) &mMtdHits; break;                       
  case kCtbId:                  coll = (StPtrVecMcHit *) &mCtbHits; break;                       
  case kSsdId:                  coll = (StPtrVecMcHit *) &mSsdHits; break;               
  case kBarrelEmcTowerId:       break; 
  case kBarrelEmcPreShowerId:   break; 
  case kBarrelSmdEtaStripId:    break;
  case kBarrelSmdPhiStripId:    break;    
  case kEndcapEmcTowerId:       break;     
  case kEndcapEmcPreShowerId:   break;     
  case kEndcapSmdUStripId:      break;    
  case kEndcapSmdVStripId:      break;
  case kZdcWestId:              
  case kZdcEastId:              break;                   
  case kMwpcWestId:             
  case kMwpcEastId:             break;              
  case kPhmdCpvId:              break;               
  case kPhmdId:                 break;                  
  case kPxlId:                  coll = (StPtrVecMcHit *) &mPxlHits; break;       
  case kIstId:                  coll = (StPtrVecMcHit *) &mIstHits; break;                  
  case kFgtId:                  coll = (StPtrVecMcHit *) &mFgtHits; break;   
  case kEtrId:                  coll = (StPtrVecMcHit *) &mEtrHits; break;   
  default:                      break;          
  };
  return coll;
}
//--------------------------------------------------------------------------------
const StPtrVecMcCalorimeterHit *StMcTrack::CalorimeterHits(StDetectorId Id) const {
  StPtrVecMcCalorimeterHit *coll = 0;
  switch (Id) {
  case kTpcId:                  break;                  
  case kSvtId:                  break;      
  case kRichId:                 break;                 
  case kFtpcWestId:             
  case kFtpcEastId:             break;             
  case kTofId:                  break;                   
  case kMtdId:                  break;                   
  case kCtbId:                  break;                   
  case kSsdId:                  break;           
  case kBarrelEmcTowerId:       coll = (StPtrVecMcCalorimeterHit *) &mBemcHits;  break; 
  case kBarrelEmcPreShowerId:   coll = (StPtrVecMcCalorimeterHit *) &mBprsHits;  break; 
  case kBarrelSmdEtaStripId:    coll = (StPtrVecMcCalorimeterHit *) &mBsmdeHits; break;
  case kBarrelSmdPhiStripId:    coll = (StPtrVecMcCalorimeterHit *) &mBsmdpHits; break;    
  case kEndcapEmcTowerId:       coll = (StPtrVecMcCalorimeterHit *) &mEemcHits;  break;     
  case kEndcapEmcPreShowerId:   coll = (StPtrVecMcCalorimeterHit *) &mEprsHits;  break;     
  case kEndcapSmdUStripId:      coll = (StPtrVecMcCalorimeterHit *) &mEsmduHits; break;    
  case kEndcapSmdVStripId:      coll = (StPtrVecMcCalorimeterHit *) &mEsmdvHits; break;
  case kZdcWestId:              
  case kZdcEastId:              break;                   
  case kMwpcWestId:             
  case kMwpcEastId:             break;              
  case kPhmdCpvId:              break;               
  case kPhmdId:                 break;                  
  case kPxlId:                  break;       
  case kIstId:                  break;                  
  case kFgtId:                  break;   
  case kEtrId:                  break;   
  default:                      break;          
  };
  return coll;
}