Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DMagneticFieldMapSpoiled.cc
Go to the documentation of this file.
1 // $Id$
2 //
3 // File: DMagneticFieldMapSpoiled.cc
4 // Created: Wed Mar 25 04:04:16 EDT 2009
5 // Creator: davidl (on Darwin Amelia.local 9.6.0 i386)
6 //
7 
10 
11 #include <iostream>
12 #include <cmath>
13 using namespace std;
14 
15 #include <JANA/JApplication.h>
16 #include <JANA/JParameterManager.h>
17 using namespace jana;
18 
19 #include <DVector3.h>
20 
21 //---------------------------------
22 // DMagneticFieldMapSpoiled (Constructor)
23 //---------------------------------
24 DMagneticFieldMapSpoiled::DMagneticFieldMapSpoiled(JApplication *japp, unsigned int run_number, string namepath)
25 {
26  bfield = new DMagneticFieldMapCalibDB(japp, run_number, namepath);
27 
28  initialized = false;
29 }
30 
31 //---------------------------------
32 // DMagneticFieldMapSpoiled (Constructor)
33 //---------------------------------
34 DMagneticFieldMapSpoiled::DMagneticFieldMapSpoiled(JCalibration *jcalib, string namepath)
35 {
36  bfield = new DMagneticFieldMapCalibDB(jcalib, namepath);
37 
38  initialized = false;
39 }
40 
41 //---------------------------------
42 // ~DMagneticFieldMapSpoiled (Destructor)
43 //---------------------------------
45 {
46  if(bfield)delete bfield;
47 }
48 
49 //---------------------------------
50 // Init
51 //---------------------------------
53 {
54  phi_amp = 0.0;
55  phi_omega = M_PI/(2.0*M_PI);
56  r_amp = 0.0;
57  r_omega = M_PI/65.0;
58  z_amp = 0.0;
59  z_omega = M_PI/(20.0);
60 
61  if(!gPARMS){
62  _DBG_<<"gPARMS==NULL!"<<endl;
63  return;
64  }
65 
66  gPARMS->SetDefaultParameter("BFIELD:PHI_AMP", phi_amp);
67  gPARMS->SetDefaultParameter("BFIELD:PHI_OMEGA", phi_omega);
68  gPARMS->SetDefaultParameter("BFIELD:R_AMP", r_amp);
69  gPARMS->SetDefaultParameter("BFIELD:R_OMEGA", r_omega);
70  gPARMS->SetDefaultParameter("BFIELD:Z_AMP", z_amp);
71  gPARMS->SetDefaultParameter("BFIELD:Z_OMEGA", z_omega);
72 
73  initialized=true;
74 }
75 
77  double Bx,By,Bz;
78  GetField(pos.x(),pos.y(),pos.z(),Bx,By,Bz);
79  Bout.SetXYZ(Bx,By,Bz);
80 }
81 
82 //---------------------------------
83 // GetField
84 //---------------------------------
85 void DMagneticFieldMapSpoiled::GetField(double x, double y, double z, double &Bx, double &By, double &Bz, int method) const
86 {
87  if(!initialized){
88  DMagneticFieldMapSpoiled *mythis = const_cast<DMagneticFieldMapSpoiled*>(this);
89  mythis->Init();
90  }
91  if(!initialized)return;
92  bfield->GetField(x, y, z, Bx, By, Bz, method);
93 
94  DVector3 B(Bx, By, Bz);
95  double r = sqrt(x*x + y*y);
96  double phi = atan2(y,x);
97  if(phi<0.0)phi+=2.0*M_PI;
98 
99  // Apply phi dependance
100  double fac = phi_amp*sin(phi_omega*phi) + 1.0;
101 
102  // Apply r dependance
103  fac *= r_amp*sin(r_omega*r) + 1.0;
104 
105  // Apply z dependance
106  fac *= z_amp*sin(z_omega*z) + 1.0;
107 
108  B *= fac;
109  Bx = B.X();
110  By = B.Y();
111  Bz = B.Z();
112 
113 }
114 
115 
116 //---------------------------------
117 // get z-component of magnetic field
118 //---------------------------------
119 double DMagneticFieldMapSpoiled::GetBz(double x, double y, double z) const
120 {
121  if(!initialized){
122  DMagneticFieldMapSpoiled *mythis = const_cast<DMagneticFieldMapSpoiled*>(this);
123  mythis->Init();
124  }
125  if(!initialized)return 0.;
126  double Bx,By,Bz;
127  bfield->GetField(x, y, z, Bx, By, Bz);
128 
129  DVector3 B(Bx, By, Bz);
130  double r = sqrt(x*x + y*y);
131  double phi = atan2(y,x);
132  if(phi<0.0)phi+=2.0*M_PI;
133 
134  // Apply phi dependance
135  double fac = phi_amp*sin(phi_omega*phi) + 1.0;
136 
137  // Apply r dependance
138  fac *= r_amp*sin(r_omega*r) + 1.0;
139 
140  // Apply z dependance
141  fac *= z_amp*sin(z_omega*z) + 1.0;
142 
143  B *= fac;
144  return B.Z();
145 
146 }
147 
148 
149 //---------------------------------
150 // GetFieldGradient
151 //---------------------------------
152 void DMagneticFieldMapSpoiled::GetFieldGradient(double x, double y, double z,
153  double &dBxdx, double &dBxdy,
154  double &dBxdz,
155  double &dBydx, double &dBydy,
156  double &dBydz,
157  double &dBzdx, double &dBzdy,
158  double &dBzdz) const
159 {
160  if(!initialized){
161  DMagneticFieldMapSpoiled *mythis = const_cast<DMagneticFieldMapSpoiled*>(this);
162  mythis->Init();
163  }
164  if(!initialized)return;
165  bfield->GetFieldGradient(x, y, z, dBxdx, dBxdy, dBxdz, dBydx, dBydy, dBydz, dBzdx, dBzdy, dBzdz);
166 }
167 
168 
169 void DMagneticFieldMapSpoiled::GetFieldBicubic(double x,double y,double z,
170  double &Bx,double &By,double &Bz) const{
171  if(!initialized){
172  DMagneticFieldMapSpoiled *mythis = const_cast<DMagneticFieldMapSpoiled*>(this);
173  mythis->Init();
174  }
175  if(!initialized)return;
176  bfield->GetFieldBicubic(x,y,z,Bx,By,Bz);
177 }
178 
179 
181  double &Bx,double &By,
182  double &Bz,
183  double &dBxdx, double &dBxdy,
184  double &dBxdz,
185  double &dBydx, double &dBydy,
186  double &dBydz,
187  double &dBzdx, double &dBzdy,
188  double &dBzdz) const{
189  if(!initialized){
190  DMagneticFieldMapSpoiled *mythis = const_cast<DMagneticFieldMapSpoiled*>(this);
191  mythis->Init();
192  }
193  if(!initialized)return;
194  bfield->GetFieldAndGradient(x,y,z,Bx,By,Bz,dBxdx,dBxdy,dBxdz,
195  dBydx,dBydy,dBydz,dBzdx,dBzdy,dBzdz);
196 }
TVector3 DVector3
Definition: DVector3.h:14
DMagneticFieldMapSpoiled(JApplication *japp, unsigned int run_number=1, string namepath="Magnets/Solenoid/solenoid_1500")
Double_t x[NCHANNELS]
Definition: st_tw_resols.C:39
#define y
void GetFieldAndGradient(double x, double y, double z, double &Bx, double &By, double &Bz, double &dBxdx, double &dBxdy, double &dBxdz, double &dBydx, double &dBydy, double &dBydz, double &dBzdx, double &dBzdy, double &dBzdz) const
JApplication * japp
#define _DBG_
Definition: HDEVIO.h:12
void GetFieldBicubic(double x, double y, double z, double &Bx, double &By, double &Bz) const
double sqrt(double)
double sin(double)
void GetFieldGradient(double x, double y, double z, double &dBxdx, double &dBxdy, double &dBxdz, double &dBydx, double &dBydy, double &dBydz, double &dBzdx, double &dBzdy, double &dBzdz) const
double GetBz(double x, double y, double z) const
void GetField(const DVector3 &pos, DVector3 &Bout) const