CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
VolumeBasedMagneticField.cc
Go to the documentation of this file.
4 
6  std::vector<MagBLayer *> theBLayers,
7  std::vector<MagESector *> theESectors,
8  std::vector<MagVolume6Faces*> theBVolumes,
9  std::vector<MagVolume6Faces*> theEVolumes,
10  float rMax, float zMax,
11  const MagneticField* param,
12  bool isParamFieldOwned) :
13  field(new MagGeometry(config,theBLayers,theESectors,theBVolumes,theEVolumes)),
14  maxR(rMax),
15  maxZ(zMax),
16  paramField(param),
17  magGeomOwned(true),
18  paramFieldOwned(isParamFieldOwned)
19 {
21 }
22 
23 
25  field(vbf.field),
26  maxR(vbf.maxR),
27  maxZ(vbf.maxZ),
28  paramField(vbf.paramField),
29  magGeomOwned(false),
30  paramFieldOwned(false),
31  theNominalValue(vbf.theNominalValue) {
32  // std::cout << "VolumeBasedMagneticField::clone() (shallow copy)" << std::endl;
33 }
34 
35 
37  return new VolumeBasedMagneticField(*this);
38 }
39 
40 
42  if(magGeomOwned) delete field;
43  if(paramFieldOwned) delete paramField;
44 }
45 
47 
48  // If parametrization of the inner region is available, use it.
50 
51  // If point is outside magfield map, return 0 field (not an error)
52  if (!isDefined(gp)) return GlobalVector();
53 
54  return field->fieldInTesla(gp);
55 }
56 
58  //same as above, but do not check range
60  return field->fieldInTesla(gp);
61 }
62 
63 
65 {
66  return field->findVolume(gp);
67 }
68 
69 
71  return (fabs(gp.z()) < maxZ && gp.perp() < maxR);
72 }
73 
74 
76  return field->isZSymmetric();
77 }
GlobalVector inTeslaUnchecked(const GlobalPoint &g) const
bool isZSymmetric() const
Definition: MagGeometry.cc:211
virtual GlobalVector inTeslaUnchecked(const GlobalPoint &gp) const
Definition: MagneticField.h:46
T perp() const
Definition: PV3DBase.h:66
virtual bool isDefined(const GlobalPoint &gp) const
True if the point is within the region where the concrete field.
Definition: MagneticField.h:40
virtual int nominalValue() const
The nominal field value for this map in kGauss.
const MagVolume * findVolume(const GlobalPoint &gp) const
virtual MagneticField * clone() const
Returns a shallow copy.
bool isDefined(const GlobalPoint &gp) const
True if the point is within the region where the concrete field.
T z() const
Definition: PV3DBase.h:58
GlobalVector inTesla(const GlobalPoint &g) const
Field value ad specified global point, in Tesla.
const MagneticField * paramField
GlobalVector fieldInTesla(const GlobalPoint &gp) const
Return field vector at the specified global point.
Definition: MagGeometry.cc:76
MagVolume * findVolume(const GlobalPoint &gp, double tolerance=0.) const
Find a volume.
Definition: MagGeometry.cc:147
VolumeBasedMagneticField(const edm::ParameterSet &config, std::vector< MagBLayer * > theBLayers, std::vector< MagESector * > theESectors, std::vector< MagVolume6Faces * > theBVolumes, std::vector< MagVolume6Faces * > theEVolumes, float rMax, float zMax, const MagneticField *param=0, bool isParamFieldOwned=false)
Global3DVector GlobalVector
Definition: GlobalVector.h:10