Geant4  10.02.p01
Par01PiModel.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: Par01PiModel.cc 90093 2015-05-13 11:59:54Z gcosmo $
28 //
29 #include "Par01PiModel.hh"
30 
31 #include "G4PionMinus.hh"
32 #include "G4PionPlus.hh"
33 #include "G4Gamma.hh"
34 
36  G4VFastSimulationModel("Par01PiModel",anEnvelope)
37 {;}
38 
40 {;}
41 
43 {
44  return
45  &particleType == G4PionMinus::PionMinusDefinition() ||
46  &particleType == G4PionPlus::PionPlusDefinition();
47 }
48 
50  //-------------------------------------------------------------
51  // UserTrigger() method: method which has to decide if
52  // the parameterisation has to be applied.
53  // Here ModelTrigger() asks the user (ie you) a 0/1 answer.
54  //
55  // Note that quantities like the local/global position/direction etc..
56  // are available at this level via the fastTrack parameter (allowing
57  // to check distance from boundaries, see below to allow the decision)
58  //--------------------------------------------------------------
59  G4cout << "\nPar01PiModel::ModelTrigger() called:" << G4endl;
60  G4cout << "--------------------------------" << G4endl;
61  G4cout << "(particle is a " << fastTrack.GetPrimaryTrack()->
62  GetDefinition()->GetParticleName() << " )\n" << G4endl;
63 
64  // -- Examples of available informations:
65 
66  // -- position:
67  G4cout << " Track position: " <<
68  fastTrack.GetPrimaryTrack()->GetPosition() << "(global coord.)" <<
69  fastTrack.GetPrimaryTrackLocalPosition() << "(in envelope coord.)"
70  << G4endl;
71 
72  // -- direction:
73  G4cout << " Track direction:" <<
74  fastTrack.GetPrimaryTrack()->GetMomentum().unit() <<
75  "(global coord.)" <<
76  fastTrack.GetPrimaryTrackLocalDirection() << "(in envelope coord.)" <<
77  G4endl;
78 
79  return true;
80 }
81 
82 void Par01PiModel::DoIt(const G4FastTrack& fastTrack,
83  G4FastStep& fastStep)
84  //--------------------------------------------------
85  //
86  // User method to code the parameterisation properly
87  // said.
88  //
89  //--------------------------------------------------
90 {
91 
92  //------------------------------------------------
93  // The primary track continues along its direction.
94  // One secondary (a photon) is added:
95  //------------------------------------------------
96  G4cout << " Pion `model' applied\n" << G4endl;
97 
98  //------------------------------
99  // Primary:
100  // idem as in "DefaultModel":
101  //
102  //------------------------------
104  G4double distance;
105  distance = fastTrack.GetEnvelopeSolid()->
106  DistanceToOut(fastTrack.GetPrimaryTrackLocalPosition(),
107  fastTrack.GetPrimaryTrackLocalDirection());
108  position = fastTrack.GetPrimaryTrackLocalPosition() +
109  distance*fastTrack.GetPrimaryTrackLocalDirection();
110 
111  // -- set final position:
112  fastStep.ProposePrimaryTrackFinalPosition(position);
113 
114  //---------------------------
115  // Secondary:
116  // Adds one "secondary":
117  //---------------------------
118  // -- First, user has to say how many secondaries will be created:
119  fastStep.SetNumberOfSecondaryTracks(1);
120 
121  //------------------------
122  // -- Build the secondary:
123  //------------------------
124  // -- direction:
125  G4ParticleMomentum direction(fastTrack.GetPrimaryTrackLocalDirection());
126  direction.setZ(direction.z()*0.5);
127  direction.setY(direction.y()+direction.z()*0.1);
128  direction = direction.unit(); // necessary ?
129 
130  // -- dynamics (Note that many constructors exists for G4DynamicParticle
131  // -- see prototype/particle+matter/particles/management/include/G4DynamicParticle.hh)
133  direction,
134  fastTrack.GetPrimaryTrack()->
135  GetKineticEnergy()/2.);
136  // -- position:
137  G4double Dist;
138  Dist = fastTrack.GetEnvelopeSolid()->
139  DistanceToOut(fastTrack.GetPrimaryTrackLocalPosition(),
140  direction);
141  G4ThreeVector posi;
142  posi = fastTrack.GetPrimaryTrackLocalPosition() + Dist*direction;
143 
144  //------------------------------------
145  //-- Creation of the secondary Track:
146  //------------------------------------
147  fastStep.CreateSecondaryTrack(dynamique, posi,
148  fastTrack.GetPrimaryTrack()->GetGlobalTime());
149 
150 }
G4ThreeVector GetPrimaryTrackLocalPosition() const
Definition: G4FastTrack.hh:213
const G4Track * GetPrimaryTrack() const
Definition: G4FastTrack.hh:208
G4ThreeVector GetPrimaryTrackLocalDirection() const
Definition: G4FastTrack.hh:223
CLHEP::Hep3Vector G4ThreeVector
virtual G4bool IsApplicable(const G4ParticleDefinition &)
Definition: Par01PiModel.cc:42
const G4ThreeVector & GetPosition() const
void SetNumberOfSecondaryTracks(G4int)
virtual G4bool ModelTrigger(const G4FastTrack &)
Definition: Par01PiModel.cc:49
G4Track * CreateSecondaryTrack(const G4DynamicParticle &, G4ThreeVector, G4ThreeVector, G4double, G4bool localCoordinates=true)
Definition: G4FastStep.cc:202
virtual void DoIt(const G4FastTrack &, G4FastStep &)
Definition: Par01PiModel.cc:82
void ProposePrimaryTrackFinalPosition(const G4ThreeVector &, G4bool localCoordinates=true)
Definition: G4FastStep.cc:98
#define position
Definition: xmlparse.cc:622
G4GLOB_DLL std::ostream G4cout
static G4PionMinus * PionMinusDefinition()
Definition: G4PionMinus.cc:93
bool G4bool
Definition: G4Types.hh:79
static G4PionPlus * PionPlusDefinition()
Definition: G4PionPlus.cc:93
G4double GetGlobalTime() const
G4VSolid * GetEnvelopeSolid() const
Definition: G4FastTrack.hh:203
G4ThreeVector GetMomentum() const
Par01PiModel(G4Region *anEnvelope)
Definition: Par01PiModel.cc:35
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
G4ThreeVector G4ParticleMomentum
static G4Gamma * GammaDefinition()
Definition: G4Gamma.cc:81