PamVMC/src/PamVMCPrimaryGenerator.cxx~

// $Id: PamVMCPrimaryGenerator.cxx,v 1.0 2006/06/03 
//

#include <TVirtualMC.h>
#include <TVirtualMCStack.h>
#include <TPDGCode.h>
#include <TDatabasePDG.h>
#include <TParticlePDG.h>
#include <TVector3.h>
#include <TMath.h>

#include "PamVMCPrimaryGenerator.h"

ClassImp(PamVMCPrimaryGenerator)

PamVMCPrimaryGenerator::PamVMCPrimaryGenerator(TVirtualMCStack* stack) 
  : TObject(),
    fStack(stack),
    fPdg(kProton),
    fKinEnergy(100.0),//100 GeV
    fDirX(0.),
    fDirY(0.),
    fDirZ(-1.),
    fPolAngle(0.),
    fNofPrimaries(1)
{
// Standard constructor

}

PamVMCPrimaryGenerator::PamVMCPrimaryGenerator()
  : TObject(),
    fStack(0),
    fPdg(0),
    fKinEnergy(0.),
    fDirX(0.),
    fDirY(0.),
    fDirZ(0.),
    fPolAngle(0.),
    fNofPrimaries(0)
{    
// Default constructor
}

PamVMCPrimaryGenerator::~PamVMCPrimaryGenerator() 
{
// Destructor  
}

// private methods

#include <Riostream.h>

void PamVMCPrimaryGenerator::GeneratePrimary()
{    
// Add one primary particle to the user stack (derived from TVirtualMCStack).
  
  // Track ID (filled by stack)
  Int_t ntr;
 
  // Option: to be tracked
  Int_t toBeDone = 1; 
 
  // Particle type
  Int_t pdg  = fPdg;
  TParticlePDG* particlePDG = TDatabasePDG::Instance()->GetParticle(fPdg);
 
  // Position
  Double_t vx  = 0.; 
  Double_t vy  = 0.; 
  Double_t vz =  80.;
  Double_t tof = 0.;

  // Energy (in GeV)
  Double_t kinEnergy = fKinEnergy; 
  Double_t mass = particlePDG->Mass();      
  Double_t e  = mass + kinEnergy;
 
  // Particle momentum
  Double_t p0, px, py, pz;
  p0 = sqrt(e*e - mass*mass); 
  px = p0 * fDirX;
  py = p0 * fDirY; 
  pz = p0 * fDirZ; 
 
  // Polarization
  TVector3 polar;
//  if ( fPdg == 50000050 ) {
//    TVector3 normal (1., 0., 0.);
//    TVector3 kphoton = TVector3(fDirX, fDirY, fDirZ);
//    TVector3 product = normal.Cross(kphoton); 
//    Double_t modul2  = product*product;
 
//    TVector3 e_perpend (0., 0., 1.);
//    if (modul2 > 0.) e_perpend = (1./sqrt(modul2))*product; 
//    TVector3 e_paralle = e_perpend.Cross(kphoton);
 
///    polar =   TMath::Cos(fPolAngle*TMath::DegToRad())*e_paralle 
  //           + TMath::Sin(fPolAngle*TMath::DegToRad())*e_perpend;
  //  } 
  //else 
  //  Warning("GeneratePrimary",
  //          "The primary particle is not an opticalphoton");

  // Add particle to stack 
  fStack->PushTrack(toBeDone, -1, pdg, px, py, pz, e, vx, vy, vz, tof, 
                   polar.X(), polar.Y(), polar.Z(), 
                   kPPrimary, ntr, 1., 0);
}

// public methods

void PamVMCPrimaryGenerator::GeneratePrimaries()
{    
// Fill the user stack (derived from TVirtualMCStack) with primary particle

  for (Int_t i=0; i<fNofPrimaries; i++) GeneratePrimary();  
}

void PamVMCPrimaryGenerator::SetDirection(
                              Double_t dirX, Double_t dirY, Double_t dirZ) 
{ 
// Set normalized direction

  Double_t norm = TMath::Sqrt(dirX*dirX + dirY*dirY + dirZ*dirZ);
  
  fDirX = dirX/norm; 
  fDirY = dirY/norm;  
  fDirZ = dirZ/norm;
}                                    


1	nikolas	1.1	// $Id: PamVMCPrimaryGenerator.cxx,v 1.0 2006/06/03
2			//
3
4			#include <TVirtualMC.h>
5			#include <TVirtualMCStack.h>
6			#include <TPDGCode.h>
7			#include <TDatabasePDG.h>
8			#include <TParticlePDG.h>
9			#include <TVector3.h>
10			#include <TMath.h>
11
12			#include "PamVMCPrimaryGenerator.h"
13
14			ClassImp(PamVMCPrimaryGenerator)
15
16			PamVMCPrimaryGenerator::PamVMCPrimaryGenerator(TVirtualMCStack* stack)
17			: TObject(),
18			fStack(stack),
19			fPdg(kProton),
20			fKinEnergy(100.0),//100 GeV
21			fDirX(0.),
22			fDirY(0.),
23			fDirZ(-1.),
24			fPolAngle(0.),
25			fNofPrimaries(1)
26			{
27			// Standard constructor
28
29			}
30
31			PamVMCPrimaryGenerator::PamVMCPrimaryGenerator()
32			: TObject(),
33			fStack(0),
34			fPdg(0),
35			fKinEnergy(0.),
36			fDirX(0.),
37			fDirY(0.),
38			fDirZ(0.),
39			fPolAngle(0.),
40			fNofPrimaries(0)
41			{
42			// Default constructor
43			}
44
45			PamVMCPrimaryGenerator::~PamVMCPrimaryGenerator()
46			{
47			// Destructor
48			}
49
50			// private methods
51
52			#include <Riostream.h>
53
54			void PamVMCPrimaryGenerator::GeneratePrimary()
55			{
56			// Add one primary particle to the user stack (derived from TVirtualMCStack).
57
58			// Track ID (filled by stack)
59			Int_t ntr;
60
61			// Option: to be tracked
62			Int_t toBeDone = 1;
63
64			// Particle type
65			Int_t pdg = fPdg;
66			TParticlePDG* particlePDG = TDatabasePDG::Instance()->GetParticle(fPdg);
67
68			// Position
69			Double_t vx = 0.;
70			Double_t vy = 0.;
71			Double_t vz = 80.;
72			Double_t tof = 0.;
73
74			// Energy (in GeV)
75			Double_t kinEnergy = fKinEnergy;
76			Double_t mass = particlePDG->Mass();
77			Double_t e = mass + kinEnergy;
78
79			// Particle momentum
80			Double_t p0, px, py, pz;
81			p0 = sqrt(ee - massmass);
82			px = p0 * fDirX;
83			py = p0 * fDirY;
84			pz = p0 * fDirZ;
85
86			// Polarization
87			TVector3 polar;
88			// if ( fPdg == 50000050 ) {
89			// TVector3 normal (1., 0., 0.);
90			// TVector3 kphoton = TVector3(fDirX, fDirY, fDirZ);
91			// TVector3 product = normal.Cross(kphoton);
92			// Double_t modul2 = product*product;
93
94			// TVector3 e_perpend (0., 0., 1.);
95			// if (modul2 > 0.) e_perpend = (1./sqrt(modul2))*product;
96			// TVector3 e_paralle = e_perpend.Cross(kphoton);
97
98			/// polar = TMath::Cos(fPolAngleTMath::DegToRad())e_paralle
99			// + TMath::Sin(fPolAngleTMath::DegToRad())e_perpend;
100			// }
101			//else
102			// Warning("GeneratePrimary",
103			// "The primary particle is not an opticalphoton");
104
105			// Add particle to stack
106			fStack->PushTrack(toBeDone, -1, pdg, px, py, pz, e, vx, vy, vz, tof,
107			polar.X(), polar.Y(), polar.Z(),
108			kPPrimary, ntr, 1., 0);
109			}
110
111			// public methods
112
113			void PamVMCPrimaryGenerator::GeneratePrimaries()
114			{
115			// Fill the user stack (derived from TVirtualMCStack) with primary particle
116
117			for (Int_t i=0; i<fNofPrimaries; i++) GeneratePrimary();
118			}
119
120			void PamVMCPrimaryGenerator::SetDirection(
121			Double_t dirX, Double_t dirY, Double_t dirZ)
122			{
123			// Set normalized direction
124
125			Double_t norm = TMath::Sqrt(dirXdirX + dirYdirY + dirZ*dirZ);
126
127			fDirX = dirX/norm;
128			fDirY = dirY/norm;
129			fDirZ = dirZ/norm;
130			}
131
132