QwGeant4
QweakSimPrimaryGeneratorAction.cc
Go to the documentation of this file.
1 //=============================================================================
2 //
3 // ---------------------------
4 // | Doxygen File Information |
5 // ---------------------------
6 //
7 /**
8 
9  \file QweakSimPrimaryGeneratorAction.cc
10 
11  $Revision: 1.4 $
12  $Date: 2006/05/05 21:35:07 $
13 
14  \author Klaus Hans Grimm
15 
16 */
17 //=============================================================================
18 
19 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
20 
22 
23 // geant4 includes
24 #include "G4Event.hh"
25 #include "G4ParticleGun.hh"
26 #include "G4ParticleTable.hh"
27 #include "Randomize.hh"
28 
29 // user includes
31 
32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
33 
34 //G4int QweakSimPrimaryGeneratorAction::kActiveOctantNumber = 1;
35 
36 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
37 
39 {
40  G4cout << "###### Calling QweakSimPrimaryGeneratorAction::QweakSimPrimaryGeneratorAction " << G4endl;
41 
42  fPositionX_min = -2.0*mm;
43  fPositionX_max = 2.0*mm;
44  fPositionY_min = -2.0*mm;
45  fPositionY_max = 2.0*mm;
46 
47  fPolarization = "L";
48 
49  myUserInfo = myUI;
50  myEvent = myEPEvent;
51 
52  // get my messenger
54 
55  // initialize my own event counter
56  // myEventCounter = 0;
57 
58  G4int n_particle = 1;
59  particleGun = new G4ParticleGun(n_particle);
60 
61  SetBeamParticleType("e-");
63 
65 
66  G4cout << "###### Leaving QweakSimPrimaryGeneratorAction::QweakSimPrimaryGeneratorAction " << G4endl;
67 }
68 
69 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
70 
72 {
73 
74  G4cout << "###### Calling/Leaving QweakSimPrimaryGeneratorAction::~QweakSimPrimaryGeneratorAction " << G4endl;
75 
76  if (particleGun) delete particleGun;
77 }
78 
79 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
80 
82 {
83  /** \page target_energy_loss target energy loss
84  * Pre-scattering (external) target energy loss is simulated by throwing one primary particle
85  * starting at 30 cm upstream of the target and letting it propagate to the actual z vertex position
86  * in the target. A second primary particle is then thrown with the angle and transverse position
87  * of the first primary at that point. Event time starts when the first primary is thrown.
88  * \todo need to double-check that position and momentum are properly transferred to second primary
89  */
90 
91  G4int myEventCounter = myUserInfo->GetPrimaryEventNumber();
92 
93  if (myEventCounter%1000==0) G4cout << "*=== Event number = " << myEventCounter << " ===*" << G4endl;
94 
95  G4double myPositionX, myPositionY, myPositionZ, myVertexZ;
96  G4double myNormMomentumX, myNormMomentumY, myNormMomentumZ;
97  G4double E_beam; // Energy of the incoming and outgoing particle
98 
99  if (myEventCounter%2 == 0) { // for even myEventCounter
100  // select position in x & y randomly
101  myPositionX = myUserInfo->GetBeamPositionX() + (G4UniformRand()-0.5)*(fPositionX_max-fPositionX_min)+(fPositionX_max+fPositionX_min)/2.0;
102  myPositionY = myUserInfo->GetBeamPositionY() + (G4UniformRand()-0.5)*(fPositionY_max-fPositionY_min)+(fPositionY_max+fPositionY_min)/2.0;
103  // select the z position 30 cm upstream of the target center (this is basically a constant)
104  myPositionZ = myUserInfo->TargetCenterPositionZ - 30.0*cm - G4UniformRand();
105 
106  myNormMomentumX = tan(myUserInfo->GetNormMomentumX()); // = 0
107  myNormMomentumY = tan(myUserInfo->GetNormMomentumY()); // = 0
108  myNormMomentumZ = sqrt(1.0 - myNormMomentumX * myNormMomentumX - myNormMomentumY * myNormMomentumY); // = 1
109 
110  E_beam = myUserInfo->GetBeamEnergy() - 0.511*MeV;
111 
113  myUserInfo->EvtGenStatus = 0; // checked in QweakSimSteppingAction.cc
114 
116  }
117  else{ // for odd myEventCounter
118 
119  myPositionX = myUserInfo->GetOriginVertexPositionX();
120  myPositionY = myUserInfo->GetOriginVertexPositionY();
121  myPositionZ = myUserInfo->GetOriginVertexPositionZ();
122 
123  myNormMomentumX = myUserInfo->GetOriginVertexMomentumDirectionX();
124  myNormMomentumY = myUserInfo->GetOriginVertexMomentumDirectionY();
125  myNormMomentumZ = myUserInfo->GetOriginVertexMomentumDirectionZ();
126 
128 
129  if (myEvent->GetReactionType() == 7) {
130  myVertexZ = myPositionZ;
133 
134  // Project x & y positions from vertex to 1 cm downstream of the target exit window
135  myPositionX += (myPositionZ-myVertexZ)*myNormMomentumX/myNormMomentumZ;
136  myPositionY += (myPositionZ-myVertexZ)*myNormMomentumY/myNormMomentumZ;
137  }
139  }
140 
141 // // Relocate the beam gun to the Cerenkov bar to test the light distributions
142 // G4double PositionX_min = -100.0*cm;
143 // G4double PositionX_max = 100.0*cm;
144 // myPositionX = (G4UniformRand()-0.5)*(PositionX_max-PositionX_min)+(PositionX_max+PositionX_min)/2.0;
145 //
146 // G4double PositionY_min = (328.-9.0)*cm;
147 // G4double PositionY_max = (328.+9.0)*cm;
148 // myPositionY = (G4UniformRand()-0.5)*(PositionY_max-PositionY_min)+(PositionY_max+PositionY_min)/2.0;
149 //
150 // myPositionZ = 560.0*cm;
151 //
152 // myNormMomentumX = 0.0;
153 // myNormMomentumY = 0.0;
154 // myNormMomentumZ = 1.0;
155 // //
156 
157  particleGun->SetParticlePosition(G4ThreeVector(myPositionX,
158  myPositionY,
159  myPositionZ ));
160 
161  particleGun->SetParticleMomentumDirection(G4ThreeVector(myNormMomentumX,
162  myNormMomentumY,
163  myNormMomentumZ));
164 
165  if (fPolarization == "L") {
166  // longitudinal polarization (after generation)
167  particleGun->SetParticlePolarization((G4ThreeVector(myNormMomentumX,
168  myNormMomentumY,
169  myNormMomentumZ)));
170  }
171  if (fPolarization == "V") {
172  // vertical transverse polarization (after generation)
173  particleGun->SetParticlePolarization(G4ThreeVector(myNormMomentumX,
174  myNormMomentumY,
175  myNormMomentumZ)
176  .cross(G4ThreeVector(1,0,0)));
177  }
178  if (fPolarization == "H") {
179  // horizontal transverse polarization (after generation)
180  particleGun->SetParticlePolarization(G4ThreeVector(myNormMomentumX,
181  myNormMomentumY,
182  myNormMomentumZ)
183  .cross(G4ThreeVector(0,1,0)));
184  }
185 
186  particleGun->SetParticleEnergy(E_beam);
187 
188  // finally : fire !!!
189  particleGun->GeneratePrimaryVertex(anEvent); // takes an event, generates primary vertex, and associates primary particles with the vertex
190  myUserInfo->StorePrimaryEventNumber(myEventCounter+1);
191  //myUserInfo->SetBeamEnergy(myEvent->GetBeamEnergy());
192  // rest of userInfo filled in QweakSimSteppingAction.cc
193 // G4cout << "###### Leaving QweakSimPrimaryGeneratorAction::GeneratePrimaries" << G4endl;
194 
195 }
196 
197 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
198 
201 }
202 
203 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
204 
207 }
208 
209 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4double GetOriginVertexPositionZ() const
G4double GetOriginVertexPositionX() const
void SetParticleType(G4ParticleDefinition *particle)
G4double GetOriginVertexMomentumDirectionX() const
G4double GetOriginVertexMomentumDirectionZ() const
G4double GetOriginVertexMomentumDirectionY() const
void GeneratePrimaries(G4Event *anEvent)
QweakSimPrimaryGeneratorActionMessenger * myMessenger
G4double GetVertexZ()
void SetBeamEnergy(G4double energy=1.160 *GeV)
void StoreOriginVertexPositionZ(G4double vz)
G4double GetOriginVertexPositionY() const
G4double GetOriginVertexKineticEnergy() const
QweakSimPrimaryGeneratorAction(QweakSimUserInformation *myUI, QweakSimEPEvent *myEPEvent)