pythiaEvent.h

Go to the documentation of this file.

#ifndef __CINT__
# include <stdlib.h> // Replace deprecated headers with C++ standard ones.
# include <ostream>
# include <TROOT.h>
# include <TRint.h>
# include <TApplication.h>
# include <TFile.h>
# include <TObject.h>
#endif

# include <TPythia6.h>
# include <TClonesArray.h>
# include <TTree.h>
# include <math.h> // Replace deprecated headers with C++ standard ones.
#include <vector>
using std::vector; // Get rid of 'using' statment in header.

#define NDET 16
#define NTYPES 8

#define NOSTARLIB // What is the purpose of this?

#ifndef NOSTARLIB
class StPythiaAnalyzer;
class StJetAssociator;
class StProtoJet;
#endif

class pythiaParticle : public TObject{
public:
   static Long64_t count; // Added for debug
   pythiaParticle(){++count;};
   pythiaParticle(Pyjets_t* pyj, Int_t i);
   virtual ~pythiaParticle(){--count;};
   void dump() const;
   
   Int_t mstatus;        // status of particle          ( LUJETS K[1] )
   Int_t mpdg_id;        // flavour code                ( LUJETS K[2] )
   Int_t mparent;        // parrent's id                ( LUJETS K[3] )
   Int_t mfirstChild;    // id of first child           ( LUJETS K[4] )
   Int_t mlastChild;     // id of last  child           ( LUJETS K[5] )
   Float_t mpx;          // X momenta [GeV/c]           ( LUJETS P[1] )
   Float_t mpy;          // Y momenta [GeV/c]           ( LUJETS P[2] )
   Float_t mpz;          // Z momenta [GeV/c]           ( LUJETS P[3] )
   Float_t menergy;      // energy    [GeV/c]           ( LUJETS P[4] )
   Float_t mmass;        // Mass      [Gev/c^2]         ( LUJETS P[5] )
   Float_t mvx;          // X vertex  [cm];             ( LUJETS V[1] )
   Float_t mvy;          // Y vertex  [cm];             ( LUJETS V[2] )
   Float_t mvz;          // Z vertex  [cm];             ( LUJETS V[3] )
   Float_t mvt;          // time of procuction [cm/c]   ( LUJETS V[4] )
   Float_t mlifetime;    // proper lifetime [cm/c]      ( LUJETS V[5] )
   Float_t met;          // trasnverse energy
   Float_t mpt;          // trasnverse energy
   Float_t meta;         // psudo rapidity
   Float_t mphi;         // phi angle
   Double_t mCharge;     //charge in units of e
   
   Int_t c2JetIndex;     //index into c2 jets
   Int_t cluJetIndex;    //index into clu jets
   
   Bool_t cluLp;         //is it the leading particle in the clu Jet
   Bool_t cluLcp;        //is it the leading charged particle in the clu Jet
   Bool_t cluTpcLcp;     //is it the tpc detected lcp in the clu Jet
   Bool_t c2Lp;          //is it the leading particle in c2 Jet
   Bool_t c2Lcp;         //is it the leading charged particle in c2 Jet
   Bool_t c2TpcLcp;      //is it the tpc detected lcp in the c2 Jet
   
   // What are all these variables?
   // Which detectors the particle passes through?
   Bool_t dZdcE;
   Bool_t dZdcW;
   Bool_t dFpd;
   Bool_t dFpdNF;
   Bool_t dFpdNN;
   Bool_t dFpdS;
   Bool_t dFpdT;
   Bool_t dFpdB;
   Bool_t dFMS;
   Bool_t dFMS2;
   Bool_t dBbcE1;
   Bool_t dBbcE2;
   Bool_t dBbcE3;
   Bool_t dBbcW1;
   Bool_t dBbcW2;
   Bool_t dBbcW3;
   Bool_t dCtb;
   Bool_t dTpc;
   
   double theta() const { return acos(mpz/p()); };
   
   
   double eta() const {
      double arg = tan(theta()/2.);
      return (arg>0.) ? -log(arg) : -999;
   }
   
   double getEta() const {
      double arg = tan(theta()/2.);
      return (arg>0.) ? -log(arg) : -999;
   }
   double px() const {return mpx;}
   double py() const {return mpy;}
   double pz() const {return mpz;}
   double p() const {return sqrt(mpx*mpx + mpy*mpy + mpz*mpz);}
   double pt() const {return sqrt(mpx*mpx+mpy*mpy);}
   double getPt() const {return pt();}
   
   double e() const {return sqrt(fabs(mpx*mpx + mpy*mpy + mpz*mpz + mmass*mmass));}
   double eT() const {return met;}
   double eZ() const {return met*sinh(eta());}
   
   double getPhi() const {return mphi;}
   double phi() const {return mphi;}
   
   double mass() const {return mmass;}
   
   double rapidity() const {
      double num = e()+mpz;
      double den = e()-mpz;
      
      if (den==0.) {return -999;}
      double arg = num/den;
      if (arg<0.) {return -999;}
      return 0.5 * log(arg);
   };
   
   double charge() const;
   
private:
   ClassDef(pythiaParticle,4)
};


class pythiaEvent: public TObject{
public:
   
   pythiaEvent(Int_t n=200, Int_t filterout=0, Int_t debugon=0);
   virtual ~pythiaEvent();
   void dump();
   void clear();
   void fill(Int_t i, TPythia6& pythia);
   void fillntp(Int_t i, int runnum, float zvert, TPythia6& pythia);
   void fillntp_1part(Int_t i, int runnum, float zvert, int pId, float px, float py, float pz, float energy, float mass);
   void closentp(float msel, float seed, float cs, float nevt);
   void STARDetector(float zvert);
   void jetAnalysis();
   
#ifndef NOSTARLIB 

   
   // Which Mike?

   void setPythiaAnalyzer(StPythiaAnalyzer* a) {mAnalyzer = a;}
   void setJetAssociator(StJetAssociator* a) {mAssoc = a;}
   
   bool diJet() const {return out1 && out2;}
   
   //Here are some jet-wise utilities:
   typedef vector<pythiaParticle*> ParticleVec; 
   
   ParticleVec cluJetDaughters(int jetIndex); 
   
   ParticleVec c2JetDaughters(int jetIndex); 
#endif
   
   // Need to figure out what all of these mean...
   
   Int_t         debug;
   Int_t         filter;             // filter output to GEANT
   Int_t         eventNumber;
   Int_t         subprocessType;
   Int_t         nParticle;
   Int_t         nFinalState;
   Float_t       shat;     
   Float_t       that;     
   Float_t       Q2;       
   Float_t       costheta; 
   Int_t         in1;      
   Int_t         in2;      
   Int_t         out1;     
   Int_t         out2;     
   Float_t       x1;       
   Float_t       x2;       
   
   Float_t       eZdcW;
   Float_t       eZdcE;
   Float_t       eSumFms[2][2];
   Float_t       eFms[40];
   Float_t       eFpd[2][16][16];
   Float_t       eFpdEM[2][16][16];
   Float_t       eFpdMax[3];
   Float_t       eFpdNN;
   Float_t       eFpdNF;
   Float_t       eFpdS;
   Float_t       eFpdT;
   Float_t       eFpdB;
   Int_t         nBbcE1;
   Int_t         nBbcE2;
   Int_t         nBbcE3;
   Int_t         nBbcW1;
   Int_t         nBbcW2;
   Int_t         nBbcW3;
   Int_t         nCtb;
   Int_t         nTpc;
   
   Int_t         Lcp;
   Float_t       LcpPt;
   Float_t       LcpEta;
   Float_t       LcpPhi;
   Float_t       sumPt;
   Float_t       vectorSumPt;
   Float_t       weightedEta;
   Float_t       weightedPhi;
   Int_t         nPJets;
   Int_t         nFJets;
   Int_t         nSJets;
   Int_t         nCJets;
   Int_t         nC2Jets;
   Int_t         nCluJets;
   Float_t       dEt[NTYPES][2];
   Float_t       dEta[NTYPES][2];
   Float_t       dPhi[NTYPES][2];
   
   TClonesArray* particles;
   TClonesArray* pJets;
   TClonesArray* fJets;
   TClonesArray* sJets;
   TClonesArray* cJets;
   TClonesArray* c2Jets;
   TClonesArray* cluJets;
   
private:
#ifndef NOSTARLIB
   void setCluIndex(StProtoJet& pj, int ijet);
   void setC2Index(StProtoJet& pj, int ijet);    
   StPythiaAnalyzer* mAnalyzer; 
   StJetAssociator* mAssoc; 
#endif
   
   ClassDef(pythiaEvent,5);
};


class pythiaRun: public TObject{
public:
   pythiaRun();
   virtual ~pythiaRun(){};
   Int_t nEvent;
   Int_t MSEL;
   Float_t energy;
   Float_t xsection;
   Float_t CKIN3;
   Float_t CKIN4;
   float   getXSection(TPythia6 &);
private:
   ClassDef(pythiaRun,1);
};

//functors
// Only used in the .cxx file, move there.
struct PythiaPartPtLessThan
{
   bool operator() (pythiaParticle* lhs, pythiaParticle* rhs) {
      return lhs->pt()<rhs->pt();
   }
};