Geant4 Cross Reference |
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 // G4UniformElectricField implementation 27 // 28 // Created: V.Grichine, 30.01.1997 29 // ------------------------------------------------------------------- 30 31 #include "G4UniformElectricField.hh" 32 #include "G4PhysicalConstants.hh" 33 34 G4UniformElectricField:: 35 G4UniformElectricField(const G4ThreeVector& FieldVector) 36 { 37 fFieldComponents[0] = 0.0; 38 fFieldComponents[1] = 0.0; 39 fFieldComponents[2] = 0.0; 40 fFieldComponents[3] = FieldVector.x(); 41 fFieldComponents[4] = FieldVector.y(); 42 fFieldComponents[5] = FieldVector.z(); 43 } 44 45 G4UniformElectricField::G4UniformElectricField(G4double vField, 46 G4double vTheta, 47 G4double vPhi) 48 { 49 if ( (vField<0) || (vTheta<0) || (vTheta>pi) || (vPhi<0) || (vPhi>twopi) ) 50 { 51 G4Exception("G4UniformElectricField::G4UniformElectricField()", 52 "GeomField0002", FatalException, "Invalid parameters."); 53 } 54 55 fFieldComponents[0] = 0.0; 56 fFieldComponents[1] = 0.0; 57 fFieldComponents[2] = 0.0; 58 fFieldComponents[3] = vField*std::sin(vTheta)*std::cos(vPhi) ; 59 fFieldComponents[4] = vField*std::sin(vTheta)*std::sin(vPhi) ; 60 fFieldComponents[5] = vField*std::cos(vTheta) ; 61 } 62 63 G4UniformElectricField::~G4UniformElectricField() = default; 64 65 G4UniformElectricField:: 66 G4UniformElectricField (const G4UniformElectricField& p) 67 : G4ElectricField(p) 68 { 69 for (auto i=0; i<6; ++i) 70 { 71 fFieldComponents[i] = p.fFieldComponents[i]; 72 } 73 } 74 75 G4UniformElectricField& 76 G4UniformElectricField::operator = (const G4UniformElectricField& p) 77 { 78 if (&p == this) { return *this; } 79 G4ElectricField::operator=(p); 80 for (auto i=0; i<6; ++i) 81 { 82 fFieldComponents[i] = p.fFieldComponents[i]; 83 } 84 return *this; 85 } 86 87 G4Field* G4UniformElectricField::Clone() const 88 { 89 return new G4UniformElectricField( G4ThreeVector(fFieldComponents[3], 90 fFieldComponents[4], 91 fFieldComponents[5])); 92 } 93 94 // ------------------------------------------------------------------------ 95 96 void G4UniformElectricField::GetFieldValue (const G4double[4], 97 G4double* fieldBandE) const 98 { 99 fieldBandE[0] = 0.0; 100 fieldBandE[1] = 0.0; 101 fieldBandE[2] = 0.0; 102 fieldBandE[3] = fFieldComponents[3]; 103 fieldBandE[4] = fFieldComponents[4]; 104 fieldBandE[5] = fFieldComponents[5]; 105 } 106