Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4ErrorSurfaceTrajState Class Reference

#include <G4ErrorSurfaceTrajState.hh>

Inheritance diagram for G4ErrorSurfaceTrajState:
Collaboration diagram for G4ErrorSurfaceTrajState:

Public Member Functions

 G4ErrorSurfaceTrajState (const G4String &partType, const G4Point3D &pos, const G4Vector3D &mom, const G4Plane3D &plane, const G4ErrorTrajErr &errmat=G4ErrorTrajErr(5, 0))
 
 G4ErrorSurfaceTrajState (const G4String &partType, const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW, const G4ErrorTrajErr &errmat=G4ErrorTrajErr(5, 0))
 
 G4ErrorSurfaceTrajState (G4ErrorFreeTrajState &tpSC, const G4Plane3D &plane)
 
 G4ErrorSurfaceTrajState (G4ErrorFreeTrajState &tpSC, const G4Vector3D &vecV, const G4Vector3D &vecW, G4ErrorMatrix &transfM)
 
 ~G4ErrorSurfaceTrajState ()
 
G4ErrorMatrix BuildErrorMatrix (G4ErrorFreeTrajState &tpSC, const G4Vector3D &vecV, const G4Vector3D &vecW)
 
virtual void Dump (std::ostream &out=G4cout) const
 
G4ErrorSurfaceTrajParam GetParameters () const
 
void SetParameters (const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW)
 
void SetParameters (const G4Point3D &pos, const G4Vector3D &mom, const G4Plane3D &plane)
 
G4Vector3D GetVectorV () const
 
G4Vector3D GetVectorW () const
 
virtual void SetPosition (const G4Point3D pos)
 
virtual void SetMomentum (const G4Vector3D &mom)
 
- Public Member Functions inherited from G4ErrorTrajState
 G4ErrorTrajState ()
 
 G4ErrorTrajState (const G4String &partType, const G4Point3D &pos, const G4Vector3D &mom, const G4ErrorTrajErr &errmat=G4ErrorTrajErr(5, 0))
 
virtual ~G4ErrorTrajState ()
 
void SetData (const G4String &partType, const G4Point3D &pos, const G4Vector3D &mom)
 
void BuildCharge ()
 
virtual G4int PropagateError (const G4Track *)
 
virtual G4int Update (const G4Track *)
 
void UpdatePosMom (const G4Point3D &pos, const G4Vector3D &mom)
 
void DumpPosMomError (std::ostream &out=G4cout) const
 
const G4StringGetParticleType () const
 
void SetParticleType (const G4String &partType)
 
G4Point3D GetPosition () const
 
G4Vector3D GetMomentum () const
 
G4ErrorTrajErr GetError () const
 
virtual void SetError (G4ErrorTrajErr em)
 
G4TrackGetG4Track () const
 
void SetG4Track (G4Track *trk)
 
G4double GetCharge () const
 
void SetCharge (G4double ch)
 
virtual G4eTSType GetTSType () const
 

Friends

std::ostream & operator<< (std::ostream &, const G4ErrorSurfaceTrajState &ts)
 

Additional Inherited Members

- Protected Attributes inherited from G4ErrorTrajState
G4String fParticleType
 
G4Point3D fPosition
 
G4Vector3D fMomentum
 
G4double fCharge
 
G4ErrorTrajErr fError
 
G4eTSType theTSType
 
G4TracktheG4Track
 
G4int iverbose
 

Detailed Description

Definition at line 56 of file G4ErrorSurfaceTrajState.hh.

Constructor & Destructor Documentation

G4ErrorSurfaceTrajState::G4ErrorSurfaceTrajState ( const G4String partType,
const G4Point3D pos,
const G4Vector3D mom,
const G4Plane3D plane,
const G4ErrorTrajErr errmat = G4ErrorTrajErr(5,0) 
)

Definition at line 60 of file G4ErrorSurfaceTrajState.cc.

63  : G4ErrorTrajState( partType, pos, mom, errmat )
64 {
65  Init();
66  fTrajParam = G4ErrorSurfaceTrajParam( pos, mom, plane );
67 
68 }
G4ErrorSurfaceTrajState::G4ErrorSurfaceTrajState ( const G4String partType,
const G4Point3D pos,
const G4Vector3D mom,
const G4Vector3D vecV,
const G4Vector3D vecW,
const G4ErrorTrajErr errmat = G4ErrorTrajErr(5,0) 
)

Definition at line 48 of file G4ErrorSurfaceTrajState.cc.

51  : G4ErrorTrajState( partType, pos, mom, errmat )
52 {
53  Init();
54  fTrajParam = G4ErrorSurfaceTrajParam( pos, mom, vecU, vecV );
55 }
G4ErrorSurfaceTrajState::G4ErrorSurfaceTrajState ( G4ErrorFreeTrajState tpSC,
const G4Plane3D plane 
)

Definition at line 73 of file G4ErrorSurfaceTrajState.cc.

75  tpSC.GetMomentum() )
76 {
77  // fParticleType = tpSC.GetParticleType();
78  // fPosition = tpSC.GetPosition();
79  // fMomentum = tpSC.GetMomentum();
80  fTrajParam = G4ErrorSurfaceTrajParam( fPosition, fMomentum, plane );
81  Init();
82 
83  //----- Get the error matrix in SC coordinates
85 }
G4ErrorMatrix BuildErrorMatrix(G4ErrorFreeTrajState &tpSC, const G4Vector3D &vecV, const G4Vector3D &vecW)
G4Point3D GetPosition() const
G4Vector3D GetMomentum() const
const G4String & GetParticleType() const

Here is the call graph for this function:

G4ErrorSurfaceTrajState::G4ErrorSurfaceTrajState ( G4ErrorFreeTrajState tpSC,
const G4Vector3D vecV,
const G4Vector3D vecW,
G4ErrorMatrix transfM 
)

Definition at line 90 of file G4ErrorSurfaceTrajState.cc.

93  tpSC.GetMomentum() )
94 {
95  Init(); // needed to define charge sign
96  fTrajParam = G4ErrorSurfaceTrajParam( fPosition, fMomentum, vecU, vecV );
97  //----- Get the error matrix in SC coordinates
98  transfM= BuildErrorMatrix( tpSC, vecU, vecV );
99 }
G4ErrorMatrix BuildErrorMatrix(G4ErrorFreeTrajState &tpSC, const G4Vector3D &vecV, const G4Vector3D &vecW)
G4Point3D GetPosition() const
G4Vector3D GetMomentum() const
const G4String & GetParticleType() const

Here is the call graph for this function:

G4ErrorSurfaceTrajState::~G4ErrorSurfaceTrajState ( )
inline

Definition at line 81 of file G4ErrorSurfaceTrajState.hh.

81 {}

Member Function Documentation

G4ErrorMatrix G4ErrorSurfaceTrajState::BuildErrorMatrix ( G4ErrorFreeTrajState tpSC,
const G4Vector3D vecV,
const G4Vector3D vecW 
)

Definition at line 104 of file G4ErrorSurfaceTrajState.cc.

106 {
107  G4double sclambda = tpSC.GetParameters().GetLambda();
108  G4double scphi = tpSC.GetParameters().GetPhi();
110  sclambda *= -1;
111  scphi += CLHEP::pi;
112  }
113  G4double cosLambda = std::cos( sclambda );
114  G4double sinLambda = std::sin( sclambda );
115  G4double sinPhi = std::sin( scphi );
116  G4double cosPhi = std::cos( scphi );
117 
118 #ifdef G4EVERBOSE
119  if( iverbose >= 4) G4cout << " PM " << fMomentum.mag() << " pLambda " << sclambda << " pPhi " << scphi << G4endl;
120 #endif
121 
122  G4ThreeVector vTN( cosLambda*cosPhi, cosLambda*sinPhi,sinLambda );
123  G4ThreeVector vUN( -sinPhi, cosPhi, 0. );
124  G4ThreeVector vVN( -vTN.z()*vUN.y(),vTN.z()*vUN.x(), cosLambda );
125 
126 #ifdef G4EVERBOSE
127  if( iverbose >= 4) G4cout << " SC2SD: vTN " << vTN << " vUN " << vUN << " vVN " << vVN << G4endl;
128 #endif
129  G4double UJ = vUN*GetVectorV();
130  G4double UK = vUN*GetVectorW();
131  G4double VJ = vVN*GetVectorV();
132  G4double VK = vVN*GetVectorW();
133 
134 
135  //--- Get transformation first
136  G4ErrorMatrix transfM(5, 5, 0 );
137  //--- Get magnetic field
139 
140  G4Vector3D vectorU = GetVectorV().cross( GetVectorW() );
141  G4double T1R = 1. / ( vTN * vectorU );
142 
143 #ifdef G4EVERBOSE
144  if( iverbose >= 4) G4cout << "surf vectors:U,V,W " << vectorU << " " << GetVectorV() << " " << GetVectorW() << " T1R:"<<T1R<<G4endl;
145 #endif
146 
147 
148  if( fCharge != 0 && field ) {
149  G4double pos[3]; pos[0] = fPosition.x()*cm; pos[1] = fPosition.y()*cm; pos[2] = fPosition.z()*cm;
150  G4double Hd[3];
151  field->GetFieldValue( pos, Hd );
152  G4ThreeVector H = G4ThreeVector( Hd[0], Hd[1], Hd[2] ) / tesla *10.; //in kilogauus
153  G4double magH = H.mag();
154  G4double invP = 1./(fMomentum.mag()/GeV);
155  G4double magHM = magH * invP;
156  if( magH != 0. ) {
157  G4double magHM2 = fCharge / magH;
158  G4double Q = -magHM * c_light/ (km/ns);
159 #ifdef G4EVERBOSE
160  if( iverbose >= 4) G4cout << GeV << " Q " << Q << " magHM " << magHM << " c_light/(km/ns) " << c_light/(km/ns) << G4endl;
161 #endif
162 
163  G4double sinz = -H*vUN * magHM2;
164  G4double cosz = H*vVN * magHM2;
165  G4double T3R = Q * std::pow(T1R,3);
166  G4double UI = vUN * vectorU;
167  G4double VI = vVN * vectorU;
168 #ifdef G4EVERBOSE
169  if( iverbose >= 4) {
170  G4cout << " T1R " << T1R << " T3R " << T3R << G4endl;
171  G4cout << " UI " << UI << " VI " << VI << " vectorU " << vectorU << G4endl;
172  G4cout << " UJ " << UJ << " VJ " << VJ << G4endl;
173  G4cout << " UK " << UK << " VK " << VK << G4endl;
174  }
175 #endif
176 
177  transfM[1][3] = -UI*( VK*cosz-UK*sinz)*T3R;
178  transfM[1][4] = -VI*( VK*cosz-UK*sinz)*T3R;
179  transfM[2][3] = UI*( VJ*cosz-UJ*sinz)*T3R;
180  transfM[2][4] = VI*( VJ*cosz-UJ*sinz)*T3R;
181  }
182  }
183 
184  G4double T2R = T1R * T1R;
185  transfM[0][0] = 1.;
186  transfM[1][1] = -UK*T2R;
187  transfM[1][2] = VK*cosLambda*T2R;
188  transfM[2][1] = UJ*T2R;
189  transfM[2][2] = -VJ*cosLambda*T2R;
190  transfM[3][3] = VK*T1R;
191  transfM[3][4] = -UK*T1R;
192  transfM[4][3] = -VJ*T1R;
193  transfM[4][4] = UJ*T1R;
194 
195 #ifdef G4EVERBOSE
196  if( iverbose >= 4) G4cout << " SC2SD transf matrix A= " << transfM << G4endl;
197 #endif
198  fError = G4ErrorTrajErr( tpSC.GetError().similarity( transfM ) );
199 
200 #ifdef G4EVERBOSE
201  if( iverbose >= 1) G4cout << "G4EP: error matrix SC2SD " << fError << G4endl;
202  if( iverbose >= 4) G4cout << "G4ErrorSurfaceTrajState from SC " << *this << G4endl;
203 #endif
204 
205  return transfM; // want to use trasnfM for the reverse transform
206 }
static constexpr double tesla
Definition: G4SIunits.hh:268
G4ErrorSymMatrix G4ErrorTrajErr
static constexpr double km
Definition: G4SIunits.hh:133
G4ErrorSymMatrix similarity(const G4ErrorMatrix &m1) const
CLHEP::Hep3Vector G4ThreeVector
virtual void GetFieldValue(const double Point[4], double *fieldArr) const =0
G4double GetPhi() const
G4double GetLambda() const
G4ErrorTrajErr fError
static double Q[]
G4GLOB_DLL std::ostream G4cout
static constexpr double cm
Definition: G4SIunits.hh:119
static G4TransportationManager * GetTransportationManager()
G4FieldManager * GetFieldManager() const
static constexpr double c_light
static constexpr double GeV
Definition: G4SIunits.hh:217
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76
static G4ErrorPropagatorData * GetErrorPropagatorData()
const G4Field * GetDetectorField() const
double mag() const
#define ns
Definition: xmlparse.cc:614
BasicVector3D< T > cross(const BasicVector3D< T > &v) const
static const G4double pos
static constexpr double pi
Definition: SystemOfUnits.h:54
G4ErrorTrajErr GetError() const
G4ErrorFreeTrajParam GetParameters() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ErrorSurfaceTrajState::Dump ( std::ostream &  out = G4cout) const
virtual

Implements G4ErrorTrajState.

Definition at line 218 of file G4ErrorSurfaceTrajState.cc.

219 {
220  out << *this;
221 }
G4ErrorSurfaceTrajParam G4ErrorSurfaceTrajState::GetParameters ( ) const
inline

Definition at line 93 of file G4ErrorSurfaceTrajState.hh.

94  { return fTrajParam; }

Here is the caller graph for this function:

G4Vector3D G4ErrorSurfaceTrajState::GetVectorV ( ) const
inline

Definition at line 111 of file G4ErrorSurfaceTrajState.hh.

112  { return fTrajParam.GetVectorV(); }

Here is the call graph for this function:

Here is the caller graph for this function:

G4Vector3D G4ErrorSurfaceTrajState::GetVectorW ( ) const
inline

Definition at line 114 of file G4ErrorSurfaceTrajState.hh.

115  { return fTrajParam.GetVectorW(); }

Here is the call graph for this function:

Here is the caller graph for this function:

virtual void G4ErrorSurfaceTrajState::SetMomentum ( const G4Vector3D mom)
inlinevirtual

Reimplemented from G4ErrorTrajState.

Definition at line 120 of file G4ErrorSurfaceTrajState.hh.

121  { SetParameters( fPosition, mom, GetVectorV(), GetVectorW() ); }
void SetParameters(const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW)

Here is the call graph for this function:

void G4ErrorSurfaceTrajState::SetParameters ( const G4Point3D pos,
const G4Vector3D mom,
const G4Vector3D vecV,
const G4Vector3D vecW 
)
inline

Definition at line 95 of file G4ErrorSurfaceTrajState.hh.

97  {
98  fPosition = pos;
99  fMomentum = mom;
100  fTrajParam.SetParameters( pos, mom, vecV, vecW );
101  }
void SetParameters(const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW)
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function:

void G4ErrorSurfaceTrajState::SetParameters ( const G4Point3D pos,
const G4Vector3D mom,
const G4Plane3D plane 
)
inline

Definition at line 103 of file G4ErrorSurfaceTrajState.hh.

105  {
106  fPosition = pos;
107  fMomentum = mom;
108  fTrajParam.SetParameters( pos, mom, plane );
109  }
void SetParameters(const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW)
static const G4double pos

Here is the call graph for this function:

virtual void G4ErrorSurfaceTrajState::SetPosition ( const G4Point3D  pos)
inlinevirtual

Reimplemented from G4ErrorTrajState.

Definition at line 117 of file G4ErrorSurfaceTrajState.hh.

118  { SetParameters( pos, fMomentum, GetVectorV(), GetVectorW() ); }
void SetParameters(const G4Point3D &pos, const G4Vector3D &mom, const G4Vector3D &vecV, const G4Vector3D &vecW)

Here is the call graph for this function:

Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  out,
const G4ErrorSurfaceTrajState ts 
)
friend

Definition at line 225 of file G4ErrorSurfaceTrajState.cc.

226 {
227  std::ios::fmtflags oldFlags = out.flags();
228  out.setf(std::ios::fixed,std::ios::floatfield);
229 
230  ts.DumpPosMomError( out );
231 
232  out << " G4ErrorSurfaceTrajState: Params: " << ts.fTrajParam << G4endl;
233  out.flags(oldFlags);
234  return out;
235 }
void DumpPosMomError(std::ostream &out=G4cout) const
#define G4endl
Definition: G4ios.hh:61

The documentation for this class was generated from the following files: