BOSS 7.0.1
BESIII Offline Software System
Loading...
Searching...
No Matches
MagneticFieldSvc Class Reference

#include <MagneticFieldSvc.h>

+ Inheritance diagram for MagneticFieldSvc:

Public Member Functions

virtual StatusCode initialize ()
 Initialise the service (Inherited Service overrides)
 
void init_params ()
 
virtual StatusCode finalize ()
 Finalise the service.
 
virtual StatusCode queryInterface (const InterfaceID &riid, void **ppvInterface)
 
virtual const InterfaceID & type () const
 Service type.
 
void handle (const Incident &)
 
virtual StatusCode fieldVector (const HepPoint3D &xyz, HepVector3D &fvec) const
 
virtual StatusCode uniFieldVector (const HepPoint3D &xyz, HepVector3D &fvec) const
 
virtual double getReferField ()
 
virtual bool ifRealField () const
 
virtual StatusCode initialize ()
 Initialise the service (Inherited Service overrides)
 
void init_params ()
 
virtual StatusCode finalize ()
 Finalise the service.
 
virtual StatusCode queryInterface (const InterfaceID &riid, void **ppvInterface)
 
virtual const InterfaceID & type () const
 Service type.
 
void handle (const Incident &)
 
virtual StatusCode fieldVector (const HepPoint3D &xyz, HepVector3D &fvec) const
 
virtual StatusCode uniFieldVector (const HepPoint3D &xyz, HepVector3D &fvec) const
 
virtual double getReferField ()
 
virtual bool ifRealField () const
 
virtual StatusCode fieldVector (const HepGeom::Point3D< double > &xyz, HepGeom::Vector3D< double > &fvec) const =0
 
virtual StatusCode uniFieldVector (const HepGeom::Point3D< double > &xyz, HepGeom::Vector3D< double > &fvec) const =0
 
virtual double getReferField ()=0
 
virtual bool ifRealField () const =0
 
virtual StatusCode fieldVector (const HepGeom::Point3D< double > &xyz, HepGeom::Vector3D< double > &fvec) const =0
 
virtual StatusCode uniFieldVector (const HepGeom::Point3D< double > &xyz, HepGeom::Vector3D< double > &fvec) const =0
 
virtual double getReferField ()=0
 
virtual bool ifRealField () const =0
 

Protected Member Functions

 MagneticFieldSvc (const std::string &name, ISvcLocator *svc)
 
virtual ~MagneticFieldSvc ()
 Virtual destructor.
 
 MagneticFieldSvc (const std::string &name, ISvcLocator *svc)
 
virtual ~MagneticFieldSvc ()
 Virtual destructor.
 

Friends

class SvcFactory< MagneticFieldSvc >
 Allow SvcFactory to instantiate the service.
 

Additional Inherited Members

- Static Public Member Functions inherited from IMagneticFieldSvc
static const InterfaceID & interfaceID ()
 Retrieve interface ID.
 
static const InterfaceID & interfaceID ()
 Retrieve interface ID.
 

Detailed Description

A service for finding the magnetic field vector at a given point in space.

Definition at line 41 of file InstallArea/include/MagneticField/MagneticField/MagneticFieldSvc.h.

Constructor & Destructor Documentation

◆ MagneticFieldSvc() [1/2]

MagneticFieldSvc::MagneticFieldSvc ( const std::string &  name,
ISvcLocator *  svc 
)
protected

Standard Constructor.

Parameters
nameString with service name
svcPointer to service locator interface

Definition at line 47 of file MagneticFieldSvc.cxx.

48 : Service( name, svc )
49{
50 declareProperty( "TurnOffField", m_turnOffField = false);
51 declareProperty( "FieldMapFile", m_filename );
52 declareProperty( "GridDistance", m_gridDistance = 5);
53 declareProperty( "RunMode", m_runmode = 2);
54 declareProperty( "IfRealField", m_ifRealField = true);
55 declareProperty( "OutLevel", m_outlevel = 1);
56 declareProperty( "Scale", m_scale = 1.0);
57 declareProperty( "UniField", m_uniField = false);
58
59 declareProperty( "Cur_SCQ1_55", m_Cur_SCQ1_55 = 349.4);
60 declareProperty( "Cur_SCQ1_89", m_Cur_SCQ1_89 = 426.2);
61 declareProperty( "Cur_SCQ2_10", m_Cur_SCQ2_10 = 474.2);
62
63 declareProperty( "UseDBFlag", m_useDB = true);
64
65 m_Mucfield = new MucMagneticField();
66 if(!m_Mucfield) cout<<"error: can not get MucMagneticField pointer"<<endl;
67
68 m_zfield = -1.0*tesla;
69
70 if(getenv("MAGNETICFIELDROOT") != NULL) {
71 path = std::string(getenv( "MAGNETICFIELDROOT" ));
72 } else {
73 std::cerr<<"Couldn't find MAGNETICFIELDROOT"<<std::endl;
74 }
75
76}

◆ ~MagneticFieldSvc() [1/2]

MagneticFieldSvc::~MagneticFieldSvc ( )
protectedvirtual

Virtual destructor.

Definition at line 109 of file MagneticFieldSvc.cxx.

110{
111 if(m_Mucfield) delete m_Mucfield;
112}

◆ MagneticFieldSvc() [2/2]

MagneticFieldSvc::MagneticFieldSvc ( const std::string &  name,
ISvcLocator *  svc 
)
protected

Standard Constructor.

Parameters
nameString with service name
svcPointer to service locator interface

◆ ~MagneticFieldSvc() [2/2]

virtual MagneticFieldSvc::~MagneticFieldSvc ( )
protectedvirtual

Virtual destructor.

Member Function Documentation

◆ fieldVector() [1/2]

StatusCode MagneticFieldSvc::fieldVector ( const HepPoint3D xyz,
HepVector3D fvec 
) const
virtual

IMagneticFieldSvc interface.

Parameters
[in]xyzPoint at which magnetic field vector will be given
[out]fvecMagnectic field vector.
Returns
StatusCode SUCCESS if calculation was performed.

Implements IMagneticFieldSvc.

Definition at line 1016 of file MagneticFieldSvc.cxx.

1018{
1019
1020 if(m_turnOffField == true) {
1021 newb[0] = 0.;
1022 newb[1] = 0.;
1023 newb[2] = 0.;
1024#ifndef BEAN
1025 return StatusCode::SUCCESS;
1026#else
1027 return true;
1028#endif
1029 }
1030
1031 // wll added 2012-08-27
1032 if(m_uniField) {
1033 uniFieldVector(newr,newb);
1034#ifndef BEAN
1035 return StatusCode::SUCCESS;
1036#else
1037 return true;
1038#endif
1039 }
1040
1041
1042 //reference frames defined by simulation and SCM are different. In simulation: x--south to north, y--down to up, but in SCM: x--north to south, y--down to up
1043 double new_x = -newr.x();
1044 double new_y = newr.y();
1045 double new_z = -newr.z();
1046 HepPoint3D r(new_x,new_y,new_z);
1047 HepVector3D b;
1048 b[0]=0.0*tesla;
1049 b[1]=0.0*tesla;
1050 b[2]=0.0*tesla;
1051 // This routine is now dummy. Was previously converting to/from CLHEP units
1052 if(-2.1*m<r.z() && r.z()<2.1*m && -1.8*m<r.x() && r.x()<1.8*m && -1.8*m<r.y() && r.y()<1.8*m)
1053 {
1054 if(-1.2*m<r.z()&&r.z()<1.2*m&&0*m<=std::sqrt(r.x()*r.x()+r.y()*r.y())&&std::sqrt(r.x()*r.x()+r.y()*r.y())<0.9*m)
1055 //if(std::abs(r.z())<1.2*m&&std::abs(r.x())<=0.9*m&&std::abs(r.y())<=0.9*m)
1056 {
1057 this->fieldGrid( r, b );
1058 }
1059 else
1060 {
1061 this->fieldGrid_TE( r, b );
1062 }
1063 }
1064 if((fabs(r.z())<=1970*mm && sqrt(r.x()*r.x()+r.y()*r.y()) >= 1740*mm) || (fabs(r.z())>=2050*mm))
1065 {
1066 HepPoint3D mr;
1067 HepVector3D tb;
1068 int part = 0, layer = 0, mat = 0;
1069 double theta;
1070 bool ifbar = false;
1071 bool ifend = false;
1072 if(r.x()!=0.){
1073 theta = atan2(fabs(r.y()),fabs(r.x()));
1074 if(0<=theta&&theta<pi/8) {
1075 mr[0] = fabs(r.x()); mr[1] = -fabs(r.y()); mr[2] = fabs(r.z());
1076 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm){
1077 part = 1;
1078 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1079 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1080 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1081 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1082 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1083 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1084 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1085 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1086 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1087 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1088 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1089 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1090 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1091 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1092 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1093 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1094 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1095 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1096 b[0] = tb[0];
1097 b[1] = -tb[1];
1098 b[2] = tb[2];
1099 ifbar = true;
1100 }
1101 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1102 part = 0; layer = 0; mat = 0;
1103 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1104 b[0] = tb[0];
1105 b[1] = -tb[1];
1106 b[2] = tb[2];
1107 ifend = true;
1108 }
1109 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1110 part = 0; layer = 0; mat = 1;
1111 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1112 b[0] = tb[0];
1113 b[1] = -tb[1];
1114 b[2] = tb[2];
1115 ifend = true;
1116 }
1117 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1118 part = 0; layer = 1; mat = 0;
1119 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1120 b[0] = tb[0];
1121 b[1] = -tb[1];
1122 b[2] = tb[2];
1123 ifend = true;
1124 }
1125 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1126 part = 0; layer = 1; mat = 1;
1127 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1128 b[0] = tb[0];
1129 b[1] = -tb[1];
1130 b[2] = tb[2];
1131 ifend = true;
1132 }
1133 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1134 part = 0; layer = 2; mat = 0;
1135 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1136 b[0] = tb[0];
1137 b[1] = -tb[1];
1138 b[2] = tb[2];
1139 ifend = true;
1140 }
1141 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1142 part = 0; layer = 2; mat = 1;
1143 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1144 b[0] = tb[0];
1145 b[1] = -tb[1];
1146 b[2] = tb[2];
1147 ifend = true;
1148 }
1149 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1150 part = 0; layer = 3; mat = 0;
1151 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1152 b[0] = tb[0];
1153 b[1] = -tb[1];
1154 b[2] = tb[2];
1155 ifend = true;
1156 }
1157 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1158 part = 0; layer = 3; mat = 1;
1159 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1160 b[0] = tb[0];
1161 b[1] = -tb[1];
1162 b[2] = tb[2];
1163 ifend = true;
1164 }
1165 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1166 part = 0; layer = 4; mat = 0;
1167 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1168 b[0] = tb[0];
1169 b[1] = -tb[1];
1170 b[2] = tb[2];
1171 ifend = true;
1172 }
1173 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1174 part = 0; layer = 4; mat = 1;
1175 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1176 b[0] = tb[0];
1177 b[1] = -tb[1];
1178 b[2] = tb[2];
1179 ifend = true;
1180 }
1181 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1182 part = 0; layer = 5; mat = 0;
1183 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1184 b[0] = tb[0];
1185 b[1] = -tb[1];
1186 b[2] = tb[2];
1187 ifend = true;
1188 }
1189 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1190 part = 0; layer = 5; mat =1;
1191 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1192 b[0] = tb[0];
1193 b[1] = -tb[1];
1194 b[2] = tb[2];
1195 ifend = true;
1196 }
1197 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1198 part = 0; layer = 6; mat = 0;
1199 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1200 b[0] = tb[0];
1201 b[1] = -tb[1];
1202 b[2] = tb[2];
1203 ifend = true;
1204 }
1205 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1206 part = 0; layer = 6; mat = 1;
1207 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1208 b[0] = tb[0];
1209 b[1] = -tb[1];
1210 b[2] = tb[2];
1211 ifend = true;
1212 }
1213 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1214 part = 0; layer = 7; mat = 0;
1215 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1216 b[0] = tb[0];
1217 b[1] = -tb[1];
1218 b[2] = tb[2];
1219 ifend = true;
1220 }
1221 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1222 part = 0; layer = 7; mat = 1;
1223 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1224 b[0] = tb[0];
1225 b[1] = -tb[1];
1226 b[2] = tb[2];
1227 ifend = true;
1228 }
1229 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1230 part = 0; layer = 8; mat = 0;
1231 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1232 b[0] = tb[0];
1233 b[1] = -tb[1];
1234 b[2] = tb[2];
1235 ifend = true;
1236 }
1237 }
1238 if(pi/8<=theta&&theta<pi/4) {
1239 mr[0] = fabs(r.x())*cos(pi/4)+fabs(r.y())*sin(pi/4);
1240 mr[1] = -fabs(r.x())*sin(pi/4)+fabs(r.y())*cos(pi/4);
1241 mr[2] = fabs(r.z());
1242 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1243 part = 1;
1244 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1245 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1246 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1247 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1248 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1249 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1250 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1251 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1252 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1253 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1254 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1255 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1256 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1257 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1258 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1259 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1260 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1261 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1262 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1263 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1264 b[2] = tb[2];
1265 ifbar = true;
1266 }
1267 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1268 part = 0; layer = 0; mat = 0;
1269 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1270 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1271 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1272 b[2] = tb[2];
1273 ifend = true;
1274 }
1275 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1276 part = 0; layer = 0; mat = 1;
1277 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1278 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1279 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1280 b[2] = tb[2];
1281 ifend = true;
1282 }
1283 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1284 part = 0; layer = 1; mat = 0;
1285 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1286 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1287 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1288 b[2] = tb[2];
1289 ifend = true;
1290 }
1291 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1292 part = 0; layer = 1; mat = 1;
1293 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1294 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1295 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1296 b[2] = tb[2];
1297 ifend = true;
1298 }
1299 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1300 part = 0; layer = 2; mat = 0;
1301 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1302 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1303 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1304 b[2] = tb[2];
1305 ifend = true;
1306 }
1307 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1308 part = 0; layer = 2; mat = 1;
1309 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1310 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1311 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1312 b[2] = tb[2];
1313 ifend = true;
1314 }
1315 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1316 part = 0; layer = 3; mat = 0;
1317 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1318 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1319 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1320 b[2] = tb[2];
1321 ifend = true;
1322 }
1323 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1324 part = 0; layer = 3; mat = 1;
1325 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1326 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1327 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1328 b[2] = tb[2];
1329 ifend = true;
1330 }
1331 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1332 part = 0; layer = 4; mat = 0;
1333 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1334 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1335 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1336 b[2] = tb[2];
1337 ifend = true;
1338 }
1339 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1340 part = 0; layer = 4; mat = 1;
1341 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1342 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1343 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1344 b[2] = tb[2];
1345 ifend = true;
1346 }
1347 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1348 part = 0; layer = 5; mat = 0;
1349 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1350 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1351 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1352 b[2] = tb[2];
1353 ifend = true;
1354 }
1355 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1356 part = 0; layer = 5; mat =1;
1357 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1358 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1359 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1360 b[2] = tb[2];
1361 ifend = true;
1362 }
1363 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1364 part = 0; layer = 6; mat = 0;
1365 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1366 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1367 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1368 b[2] = tb[2];
1369 ifend = true;
1370 }
1371 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1372 part = 0; layer = 6; mat = 1;
1373 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1374 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1375 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1376 b[2] = tb[2];
1377 ifend = true;
1378 }
1379 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1380 part = 0; layer = 7; mat = 0;
1381 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1382 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1383 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1384 b[2] = tb[2];
1385 ifend = true;
1386 }
1387 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1388 part = 0; layer = 7; mat = 1;
1389 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1390 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1391 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1392 b[2] = tb[2];
1393 ifend = true;
1394 }
1395 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1396 part = 0; layer = 8; mat = 0;
1397 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1398 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1399 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1400 b[2] = tb[2];
1401 ifend = true;
1402 }
1403 }
1404 if(pi/4<=theta&&theta<3*pi/8) {
1405 mr[0] = fabs(r.x())*cos(pi/4)+fabs(r.y())*sin(pi/4);
1406 mr[1] = fabs(r.x())*sin(pi/4)-fabs(r.y())*cos(pi/4);
1407 mr[2] = fabs(r.z());
1408 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1409 part = 1;
1410 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1411 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1412 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1413 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1414 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1415 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1416 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1417 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1418 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1419 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1420 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1421 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1422 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1423 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1424 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1425 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1426 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1427 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1428 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1429 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1430 b[2] = tb[2];
1431 ifbar = true;
1432 }
1433 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1434 part = 0; layer = 0; mat = 0;
1435 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1436 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1437 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1438 b[2] = tb[2];
1439 ifend = true;
1440 }
1441 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1442 part = 0; layer = 0; mat = 1;
1443 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1444 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1445 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1446 b[2] = tb[2];
1447 ifend = true;
1448 }
1449 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1450 part = 0; layer = 1; mat = 0;
1451 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1452 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1453 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1454 b[2] = tb[2];
1455 ifend = true;
1456 }
1457 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1458 part = 0; layer = 1; mat = 1;
1459 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1460 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1461 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1462 b[2] = tb[2];
1463 ifend = true;
1464 }
1465 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1466 part = 0; layer = 2; mat = 0;
1467 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1468 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1469 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1470 b[2] = tb[2];
1471 ifend = true;
1472 }
1473 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1474 part = 0; layer = 2; mat = 1;
1475 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1476 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1477 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1478 b[2] = tb[2];
1479 ifend = true;
1480 }
1481 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1482 part = 0; layer = 3; mat = 0;
1483 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1484 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1485 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1486 b[2] = tb[2];
1487 ifend = true;
1488 }
1489 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1490 part = 0; layer = 3; mat = 1;
1491 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1492 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1493 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1494 b[2] = tb[2];
1495 ifend = true;
1496 }
1497 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1498 part = 0; layer = 4; mat = 0;
1499 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1500 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1501 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1502 b[2] = tb[2];
1503 ifend = true;
1504 }
1505 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1506 part = 0; layer = 4; mat = 1;
1507 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1508 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1509 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1510 b[2] = tb[2];
1511 ifend = true;
1512 }
1513 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1514 part = 0; layer = 5; mat = 0;
1515 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1516 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1517 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1518 b[2] = tb[2];
1519 ifend = true;
1520 }
1521 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1522 part = 0; layer = 5; mat =1;
1523 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1524 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1525 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1526 b[2] = tb[2];
1527 ifend = true;
1528 }
1529 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1530 part = 0; layer = 6; mat = 0;
1531 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1532 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1533 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1534 b[2] = tb[2];
1535 ifend = true;
1536 }
1537 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1538 part = 0; layer = 6; mat = 1;
1539 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1540 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1541 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1542 b[2] = tb[2];
1543 ifend = true;
1544 }
1545 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1546 part = 0; layer = 7; mat = 0;
1547 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1548 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1549 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1550 b[2] = tb[2];
1551 ifend = true;
1552 }
1553 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1554 part = 0; layer = 7; mat = 1;
1555 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1556 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1557 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1558 b[2] = tb[2];
1559 ifend = true;
1560 }
1561 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1562 part = 0; layer = 8; mat = 0;
1563 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1564 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1565 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1566 b[2] = tb[2];
1567 ifend = true;
1568 }
1569 }
1570 if(3*pi/8<=theta&&theta<pi/2) {
1571 mr[0] = fabs(r.y()); mr[1] = -fabs(r.x()); mr[2] = fabs(r.z());
1572 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1573 part = 1;
1574 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1575 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1576 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1577 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1578 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1579 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1580 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1581 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1582 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1583 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1584 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1585 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1586 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1587 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1588 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1589 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1590 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1591 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1592 b[0] = -tb[1];
1593 b[1] = tb[0];
1594 b[2] = tb[2];
1595 ifbar = true;
1596 }
1597 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1598 part = 0; layer = 0; mat = 0;
1599 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1600 b[0] = -tb[1];
1601 b[1] = tb[0];
1602 b[2] = tb[2];
1603 ifend = true;
1604 }
1605 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1606 part = 0; layer = 0; mat = 1;
1607 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1608 b[0] = -tb[1];
1609 b[1] = tb[0];
1610 b[2] = tb[2];
1611 ifend = true;
1612 }
1613 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1614 part = 0; layer = 1; mat = 0;
1615 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1616 b[0] = -tb[1];
1617 b[1] = tb[0];
1618 b[2] = tb[2];
1619 ifend = true;
1620 }
1621 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1622 part = 0; layer = 1; mat = 1;
1623 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1624 b[0] = -tb[1];
1625 b[1] = tb[0];
1626 b[2] = tb[2];
1627 ifend = true;
1628 }
1629 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1630 part = 0; layer = 2; mat = 0;
1631 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1632 b[0] = -tb[1];
1633 b[1] = tb[0];
1634 b[2] = tb[2];
1635 ifend = true;
1636 }
1637 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1638 part = 0; layer = 2; mat = 1;
1639 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1640 b[0] = -tb[1];
1641 b[1] = tb[0];
1642 b[2] = tb[2];
1643 ifend = true;
1644 }
1645 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1646 part = 0; layer = 3; mat = 0;
1647 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1648 b[0] = -tb[1];
1649 b[1] = tb[0];
1650 b[2] = tb[2];
1651 ifend = true;
1652 }
1653 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1654 part = 0; layer = 3; mat = 1;
1655 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1656 b[0] = -tb[1];
1657 b[1] = tb[0];
1658 b[2] = tb[2];
1659 ifend = true;
1660 }
1661 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1662 part = 0; layer = 4; mat = 0;
1663 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1664 b[0] = -tb[1];
1665 b[1] = tb[0];
1666 b[2] = tb[2];
1667 ifend = true;
1668 }
1669 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1670 part = 0; layer = 4; mat = 1;
1671 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1672 b[0] = -tb[1];
1673 b[1] = tb[0];
1674 b[2] = tb[2];
1675 ifend = true;
1676 }
1677 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1678 part = 0; layer = 5; mat = 0;
1679 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1680 b[0] = -tb[1];
1681 b[1] = tb[0];
1682 b[2] = tb[2];
1683 ifend = true;
1684 }
1685 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1686 part = 0; layer = 5; mat =1;
1687 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1688 b[0] = -tb[1];
1689 b[1] = tb[0];
1690 b[2] = tb[2];
1691 ifend = true;
1692 }
1693 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1694 part = 0; layer = 6; mat = 0;
1695 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1696 b[0] = -tb[1];
1697 b[1] = tb[0];
1698 b[2] = tb[2];
1699 ifend = true;
1700 }
1701 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1702 part = 0; layer = 6; mat = 1;
1703 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1704 b[0] = -tb[1];
1705 b[1] = tb[0];
1706 b[2] = tb[2];
1707 ifend = true;
1708 }
1709 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1710 part = 0; layer = 7; mat = 0;
1711 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1712 b[0] = -tb[1];
1713 b[1] = tb[0];
1714 b[2] = tb[2];
1715 ifend = true;
1716 }
1717 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1718 part = 0; layer = 7; mat = 1;
1719 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1720 b[0] = -tb[1];
1721 b[1] = tb[0];
1722 b[2] = tb[2];
1723 ifend = true;
1724 }
1725 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1726 part = 0; layer = 8; mat = 0;
1727 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1728 b[0] = -tb[1];
1729 b[1] = tb[0];
1730 b[2] = tb[2];
1731 ifend = true;
1732 }
1733 }
1734 }
1735 if(r.x()==0.) {
1736 mr[0] = fabs(r.y()); mr[1] = -fabs(r.x()); mr[2] = fabs(r.z());
1737 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1738 part = 1;
1739 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1740 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1741 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1742 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1743 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1744 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1745 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1746 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1747 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1748 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1749 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1750 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1751 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1752 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1753 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1754 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1755 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1756 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1757 b[0] = -tb[1];
1758 b[1] = tb[0];
1759 b[2] = tb[2];
1760 ifbar = true;
1761 }
1762 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1763 part = 0; layer = 0; mat = 0;
1764 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1765 b[0] = -tb[1];
1766 b[1] = tb[0];
1767 b[2] = tb[2];
1768 ifend = true;
1769 }
1770 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1771 part = 0; layer = 0; mat = 1;
1772 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1773 b[0] = -tb[1];
1774 b[1] = tb[0];
1775 b[2] = tb[2];
1776 ifend = true;
1777 }
1778 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1779 part = 0; layer = 1; mat = 0;
1780 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1781 b[0] = -tb[1];
1782 b[1] = tb[0];
1783 b[2] = tb[2];
1784 ifend = true;
1785 }
1786 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1787 part = 0; layer = 1; mat = 1;
1788 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1789 b[0] = -tb[1];
1790 b[1] = tb[0];
1791 b[2] = tb[2];
1792 ifend = true;
1793 }
1794 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1795 part = 0; layer = 2; mat = 0;
1796 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1797 b[0] = -tb[1];
1798 b[1] = tb[0];
1799 b[2] = tb[2];
1800 ifend = true;
1801 }
1802 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1803 part = 0; layer = 2; mat = 1;
1804 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1805 b[0] = -tb[1];
1806 b[1] = tb[0];
1807 b[2] = tb[2];
1808 ifend = true;
1809 }
1810 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1811 part = 0; layer = 3; mat = 0;
1812 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1813 b[0] = -tb[1];
1814 b[1] = tb[0];
1815 b[2] = tb[2];
1816 ifend = true;
1817 }
1818 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1819 part = 0; layer = 3; mat = 1;
1820 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1821 b[0] = -tb[1];
1822 b[1] = tb[0];
1823 b[2] = tb[2];
1824 ifend = true;
1825 }
1826 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1827 part = 0; layer = 4; mat = 0;
1828 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1829 b[0] = -tb[1];
1830 b[1] = tb[0];
1831 b[2] = tb[2];
1832 ifend = true;
1833 }
1834 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1835 part = 0; layer = 4; mat = 1;
1836 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1837 b[0] = -tb[1];
1838 b[1] = tb[0];
1839 b[2] = tb[2];
1840 ifend = true;
1841 }
1842 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1843 part = 0; layer = 5; mat = 0;
1844 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1845 b[0] = -tb[1];
1846 b[1] = tb[0];
1847 b[2] = tb[2];
1848 ifend = true;
1849 }
1850 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1851 part = 0; layer = 5; mat =1;
1852 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1853 b[0] = -tb[1];
1854 b[1] = tb[0];
1855 b[2] = tb[2];
1856 ifend = true;
1857 }
1858 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1859 part = 0; layer = 6; mat = 0;
1860 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1861 b[0] = -tb[1];
1862 b[1] = tb[0];
1863 b[2] = tb[2];
1864 ifend = true;
1865 }
1866 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1867 part = 0; layer = 6; mat = 1;
1868 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1869 b[0] = -tb[1];
1870 b[1] = tb[0];
1871 b[2] = tb[2];
1872 ifend = true;
1873 }
1874 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1875 part = 0; layer = 7; mat = 0;
1876 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1877 b[0] = -tb[1];
1878 b[1] = tb[0];
1879 b[2] = tb[2];
1880 ifend = true;
1881 }
1882 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1883 part = 0; layer = 7; mat = 1;
1884 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1885 b[0] = -tb[1];
1886 b[1] = tb[0];
1887 b[2] = tb[2];
1888 ifend = true;
1889 }
1890 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1891 part = 0; layer = 8; mat = 0;
1892 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1893 b[0] = -tb[1];
1894 b[1] = tb[0];
1895 b[2] = tb[2];
1896 ifend = true;
1897 }
1898 }
1899 if(ifbar==true||ifend==true) {
1900 if( r.x() < 0. && r.y() >= 0. && r.z() > 0. ){
1901 b[0] = -b[0];
1902 }
1903 else if( r.x() <= 0. && r.y() < 0. && r.z() > 0. ){
1904 b[0] = -b[0];
1905 b[1] = -b[1];
1906 }
1907 else if( r.x() > 0. && r.y() < 0. && r.z() > 0. ){
1908 b[1] = -b[1];
1909 }
1910 else if( r.x() >= 0. && r.y() > 0. && r.z() < 0. ){
1911 b[0] = -b[0];
1912 b[1] = -b[1];
1913 }
1914 else if( r.x() < 0. && r.y() >= 0. && r.z() < 0. ){
1915 b[1] = -b[1];
1916 }
1917 else if( r.x() <= 0. && r.y() < 0. && r.z() < 0. ){
1918 b[0] = b[0];
1919 b[1] = b[1];
1920 }
1921 else if( r.x() > 0. && r.y() <= 0. && r.z() < 0. ){
1922 b[0] = -b[0];
1923 }
1924 }
1925 }
1926
1927 //reference frames defined by simulation and SCM are different.
1928 newb[0] = -b[0] * m_scale;
1929 newb[1] = b[1] * m_scale;
1930 newb[2] = -b[2] * m_scale;
1931/* if(-1.2*m<r.z()&&r.z()<1.2*m&&0*m<=std::sqrt(r.x()*r.x()+r.y()*r.y())&&std::sqrt(r.x()*r.x()+r.y()*r.y())<0.9*m) {
1932 return StatusCode::SUCCESS;
1933 }
1934 else {
1935 newb[0] = newb[0] - 0.10*tesla;
1936 newb[1] = newb[1] + 0.10*tesla;
1937 newb[2] = newb[2] - 0.10*tesla;
1938 }*/
1939
1940 //newb[0] = b[0];
1941 //newb[1] = b[1];
1942 //newb[2] = b[2];
1943
1944#ifndef BEAN
1945 return StatusCode::SUCCESS;
1946#else
1947 return true;
1948#endif
1949}
double sin(const BesAngle a)
double cos(const BesAngle a)
virtual StatusCode uniFieldVector(const HepPoint3D &xyz, HepVector3D &fvec) const
void getMucField(int part, int layer, int mat, HepPoint3D &r, HepVector3D &b)

◆ fieldVector() [2/2]

virtual StatusCode MagneticFieldSvc::fieldVector ( const HepPoint3D xyz,
HepVector3D fvec 
) const
virtual

IMagneticFieldSvc interface.

Parameters
[in]xyzPoint at which magnetic field vector will be given
[out]fvecMagnectic field vector.
Returns
StatusCode SUCCESS if calculation was performed.

Implements IMagneticFieldSvc.

◆ finalize() [1/2]

StatusCode MagneticFieldSvc::finalize ( )
virtual

Finalise the service.

Definition at line 662 of file MagneticFieldSvc.cxx.

663{
664 MsgStream log( msgSvc(), name() );
665 //if(m_connect_run) delete m_connect_run;
666 StatusCode status = Service::finalize();
667
668 if ( status.isSuccess() )
669 log << MSG::INFO << "Service finalized successfully" << endreq;
670 return status;
671}

◆ finalize() [2/2]

virtual StatusCode MagneticFieldSvc::finalize ( )
virtual

Finalise the service.

◆ getReferField() [1/2]

double MagneticFieldSvc::getReferField ( )
virtual

Implements IMagneticFieldSvc.

Definition at line 1999 of file MagneticFieldSvc.cxx.

2000{
2001 if(m_useDB == false) {
2002 if(m_runmode == 8 || m_runmode == 7) m_zfield = -0.9*tesla;
2003 else m_zfield = -1.0*tesla;
2004 }
2005
2006 if(m_turnOffField == true) {
2007 m_zfield = 0.;
2008 }
2009 return m_zfield * m_scale;
2010}

◆ getReferField() [2/2]

virtual double MagneticFieldSvc::getReferField ( )
virtual

Implements IMagneticFieldSvc.

◆ handle() [1/2]

void MagneticFieldSvc::handle ( const Incident &  inc)

Definition at line 694 of file MagneticFieldSvc.cxx.

694 {
695 MsgStream log( messageService(), name() );
696 log << MSG::DEBUG << "handle: " << inc.type() << endreq;
697 if ( inc.type() != "NewRun" ){
698 return;
699 }
700 log << MSG::DEBUG << "Begin New Run" << endreq;
701 if(m_useDB == true) {
702 init_params();
703 }
704}

◆ handle() [2/2]

void MagneticFieldSvc::handle ( const Incident &  )

◆ ifRealField() [1/2]

bool MagneticFieldSvc::ifRealField ( ) const
virtual

Implements IMagneticFieldSvc.

Definition at line 2012 of file MagneticFieldSvc.cxx.

2012 {
2013 return m_ifRealField;
2014}

◆ ifRealField() [2/2]

virtual bool MagneticFieldSvc::ifRealField ( ) const
virtual

Implements IMagneticFieldSvc.

◆ init_params() [1/2]

void MagneticFieldSvc::init_params ( )

Definition at line 261 of file MagneticFieldSvc.cxx.

262{
263 MsgStream log(msgSvc(), name());
264#else
266{
267 if( !init_mucMagneticField() ) {
268 cerr << " STOP! " << endl;
269 exit(1);
270 }
271#endif
272
273 m_Q.clear();
274 m_P.clear();
275 m_Q_1.clear();
276 m_P_1.clear();
277 m_Q_2.clear();
278 m_P_2.clear();
279 m_P_TE.clear();
280 m_Q_TE.clear();
281
282 if(m_gridDistance == 5) {
283 m_Q.reserve(201250);
284 m_P.reserve(201250);
285 m_Q_1.reserve(201250);
286 m_P_1.reserve(201250);
287 m_Q_2.reserve(201250);
288 m_P_2.reserve(201250);
289 m_Q_TE.reserve(176608);
290 m_P_TE.reserve(176608);
291 }
292 if(m_gridDistance == 10) {
293 m_Q.reserve(27082);
294 m_P.reserve(27082);
295 m_Q_1.reserve(27082);
296 m_P_1.reserve(27082);
297 m_Q_2.reserve(27082);
298 m_P_2.reserve(27082);
299 m_Q_TE.reserve(23833);
300 m_P_TE.reserve(23833);
301 }
302
303#ifndef BEAN
304 int runNo;
305 SmartDataPtr<Event::EventHeader> evt(m_eventSvc,"/Event/EventHeader");
306 if( evt ){
307 runNo = evt -> runNumber();
308 log << MSG::INFO <<"The runNumber of current event is "<<runNo<<endreq;
309 }
310 else {
311 log << MSG::ERROR <<"Can not get EventHeader!"<<endreq;
312 }
313#else
314 int runNo = run;
315#endif
316
318 std::vector<double> current(3,0.);
319 ConnectionDB::eRet e = m_connect_run->getReadSC_MagnetInfo(current,std::abs(runNo));
320
321 std::vector<double> beamEnergy;
322 beamEnergy.clear();
323 e = m_connect_run->getBeamEnergy(beamEnergy, std::abs(runNo));
324 char BPR_PRB[200];
325 char BER_PRB[200];
326
327 sprintf(BPR_PRB,"%f ",beamEnergy[0]);
328 sprintf(BER_PRB,"%f ",beamEnergy[1]);
329
330 setenv("BEPCII_INFO.BPR_PRB", BPR_PRB, 1);
331 setenv("BEPCII_INFO.BER_PRB", BER_PRB, 1);
332
333 int ssm_curr = int (current[0]);
334 double scql_curr = current[1];
335 double scqr_curr = current[2];
336
337#ifndef BEAN
338 log << MSG::INFO << "SSM current: " << current[0] << " SCQL current: " << current[1] << " SCQR current: " << current[2] << " in Run " << runNo << endreq;
339#else
340 cout << "SSM current: " << current[0] << " SCQL current: " << current[1]
341 << " SCQR current: " << current[2] << " in Run " << runNo << endl;
342#endif
343
344 //int ssm_curr = 3369;
345 //double scql_curr = 426.2;
346 //double scqr_curr = 426.2;
347 //for the energy less than 1.89GeV
348 if((ssm_curr >= 3367) && (ssm_curr <= 3370) && ((scql_curr+scqr_curr)/2 < m_Cur_SCQ1_89)) {
349 m_zfield = -1.0*tesla;
350 m_filename = path;
351 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode2.dat");
352#ifndef BEAN
353 log << MSG::INFO << "Select field map: " << m_filename << endreq;
354#else
355 cout << "Select field map: " << m_filename << endl;
356#endif
357
358 StatusCode status = parseFile();
359
360#ifndef BEAN
361 if ( !status.isSuccess() ) {
362 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
363 }
364#else
365 if ( !status ) {
366 cout << "Magnetic field parse failled" << endl;
367 }
368#endif
369 m_Q_1 = m_Q;
370 m_P_1 = m_P;
371
372 m_Q.clear();
373 m_P.clear();
374 if(m_gridDistance == 5) {
375 m_Q.reserve(201250);
376 m_P.reserve(201250);
377 }
378 if(m_gridDistance == 10) {
379 m_Q.reserve(27082);
380 m_P.reserve(27082);
381 }
382
383 m_filename = path;
384 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode3.dat");
385#ifndef BEAN
386 log << MSG::INFO << "Select field map: " << m_filename << endreq;
387#else
388 cout << "Select field map: " << m_filename << endl;
389#endif
390
391 status = parseFile();
392
393#ifndef BEAN
394 if ( !status.isSuccess() ) {
395 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
396 }
397#else
398 if ( !status ) {
399 cout << "Magnetic field parse failled" << endl;
400 }
401#endif
402 m_Q_2 = m_Q;
403 m_P_2 = m_P;
404
405 m_Q.clear();
406 m_P.clear();
407 if(m_gridDistance == 5) {
408 m_Q.reserve(201250);
409 m_P.reserve(201250);
410 }
411 if(m_gridDistance == 10) {
412 m_Q.reserve(27082);
413 m_P.reserve(27082);
414 }
415 //check
416 if(m_Q_2.size() != m_Q_1.size()) {
417#ifndef BEAN
418 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
419#else
420 cout << "The two field maps used in this run are wrong!!!" << endl;
421#endif
422 exit(1);
423 }
424
425 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
426 double fieldvalue = (m_Q_1[iQ] - m_Q_2[iQ])/(m_Cur_SCQ1_55 - m_Cur_SCQ1_89)*((scql_curr+scqr_curr)/2 - m_Cur_SCQ1_89) + m_Q_2[iQ];
427 m_Q.push_back(fieldvalue);
428 m_P.push_back(m_P_2[iQ]);
429 }
430 }
431 //for the energy greater than 1.89GeV
432 if((ssm_curr >= 3367) && (ssm_curr <= 3370) && ((scql_curr+scqr_curr)/2 >= m_Cur_SCQ1_89)) {
433 m_zfield = -1.0*tesla;
434 m_filename = path;
435 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode3.dat");
436#ifndef BEAN
437 log << MSG::INFO << "Select field map: " << m_filename << endreq;
438#else
439 cout << "Select field map: " << m_filename << endl;
440#endif
441
442 StatusCode status = parseFile();
443
444#ifndef BEAN
445 if ( !status.isSuccess() ) {
446 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
447 }
448#else
449 if ( !status ) {
450 cout << "Magnetic field parse failled" << endl;
451 }
452#endif
453 m_Q_1 = m_Q;
454 m_P_1 = m_P;
455
456 m_Q.clear();
457 m_P.clear();
458 if(m_gridDistance == 5) {
459 m_Q.reserve(201250);
460 m_P.reserve(201250);
461 }
462 if(m_gridDistance == 10) {
463 m_Q.reserve(27082);
464 m_P.reserve(27082);
465 }
466
467 m_filename = path;
468 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode4.dat");
469#ifndef BEAN
470 log << MSG::INFO << "Select field map: " << m_filename << endreq;
471#else
472 cout << "Select field map: " << m_filename << endl;
473#endif
474
475 status = parseFile();
476
477#ifndef BEAN
478 if ( !status.isSuccess() ) {
479 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
480 }
481#else
482 if ( !status ) {
483 cout << "Magnetic field parse failled" << endl;
484 }
485#endif
486 m_Q_2 = m_Q;
487 m_P_2 = m_P;
488
489 m_Q.clear();
490 m_P.clear();
491 if(m_gridDistance == 5) {
492 m_Q.reserve(201250);
493 m_P.reserve(201250);
494 }
495 if(m_gridDistance == 10) {
496 m_Q.reserve(27082);
497 m_P.reserve(27082);
498 }
499 //check
500 if(m_Q_2.size() != m_Q_1.size()) {
501#ifndef BEAN
502 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
503#else
504 cout << "The two field maps used in this run are wrong!!!" << endl;
505#endif
506 exit(1);
507 }
508
509 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
510 double fieldvalue = (m_Q_1[iQ] - m_Q_2[iQ])/(m_Cur_SCQ1_89 - m_Cur_SCQ2_10)*((scql_curr+scqr_curr)/2 - m_Cur_SCQ2_10) + m_Q_2[iQ];
511 m_Q.push_back(fieldvalue);
512 m_P.push_back(m_P_2[iQ]);
513 }
514 }
515 if((ssm_curr >= 3033) && (ssm_curr <= 3035)) {
516 m_zfield = -0.9*tesla;
517 m_filename = path;
518 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode7.dat");
519#ifndef BEAN
520 log << MSG::INFO << "Select field map: " << m_filename << endreq;
521#else
522 cout << "Select field map: " << m_filename << endl;
523#endif
524
525 StatusCode status = parseFile();
526
527#ifndef BEAN
528 if ( !status.isSuccess() ) {
529 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
530 }
531#else
532 if ( !status ) {
533 cout << "Magnetic field parse failled" << endl;
534 }
535#endif
536 m_Q_1 = m_Q;
537 m_P_1 = m_P;
538
539 m_Q.clear();
540 m_P.clear();
541 if(m_gridDistance == 5) {
542 m_Q.reserve(201250);
543 m_P.reserve(201250);
544 }
545 if(m_gridDistance == 10) {
546 m_Q.reserve(27082);
547 m_P.reserve(27082);
548 }
549
550 m_filename = path;
551 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode8.dat");
552#ifndef BEAN
553 log << MSG::INFO << "Select field map: " << m_filename << endreq;
554#else
555 cout << "Select field map: " << m_filename << endl;
556#endif
557
558 status = parseFile();
559
560#ifndef BEAN
561 if ( !status.isSuccess() ) {
562 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
563 }
564#else
565 if ( !status ) {
566 cout << "Magnetic field parse failled" << endl;
567 }
568#endif
569 m_Q_2 = m_Q;
570 m_P_2 = m_P;
571
572 m_Q.clear();
573 m_P.clear();
574 if(m_gridDistance == 5) {
575 m_Q.reserve(201250);
576 m_P.reserve(201250);
577 }
578 if(m_gridDistance == 10) {
579 m_Q.reserve(27082);
580 m_P.reserve(27082);
581 }
582 //check
583 if(m_Q_2.size() != m_Q_1.size()) {
584#ifndef BEAN
585 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
586#else
587 cout << "The two field maps used in this run are wrong!!!" << endl;
588#endif
589 exit(1);
590 }
591
592 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
593 double fieldvalue = (m_Q_1[iQ] - m_Q_2[iQ])/(m_Cur_SCQ1_55 - m_Cur_SCQ1_89)*((scql_curr+scqr_curr)/2 - m_Cur_SCQ1_89) + m_Q_2[iQ];
594 m_Q.push_back(fieldvalue);
595 m_P.push_back(m_P_2[iQ]);
596 }
597 }
598
599 if((!((ssm_curr >= 3367) && (ssm_curr <= 3370)) && !((ssm_curr >= 3033) && (ssm_curr <= 3035))) || scql_curr == 0 || scqr_curr == 0) {
600#ifndef BEAN
601 log << MSG::FATAL << "The current of run " << runNo << " is abnormal in database, please check it! " << endreq;
602 log << MSG::FATAL << "The current of SSM is " << ssm_curr << ", SCQR is " << scqr_curr << ", SCQL is " << scql_curr << endreq;
603#else
604 cout << "The current of run " << runNo
605 << " is abnormal in database, please check it! " << endl;
606 cout << "The current of SSM is " << ssm_curr
607 << ", SCQR is " << scqr_curr << ", SCQL is " << scql_curr << endl;
608#endif
609 exit(1);
610 }
611
612 if(m_Q.size() == 0) {
613#ifndef BEAN
614 log << MSG::FATAL << "This size of field map vector is ZERO, please check the current of run " << runNo << " in database!" << endreq;
615#else
616 cout << "This size of field map vector is ZERO,"
617 << " please check the current of run " << runNo
618 << " in database!" << endl;
619#endif
620 exit(1);
621 }
622
623 m_filename_TE = path;
624 if(m_gridDistance == 10) m_filename_TE += std::string( "/dat/TEMap10cm.dat");
625 if(m_gridDistance == 5) m_filename_TE += std::string( "/dat/TEMap5cm.dat");
626#ifndef BEAN
627 log << MSG::INFO << "Select field map: " << m_filename_TE << endreq;
628#else
629 cout << "Select field map: " << m_filename_TE << endl;
630#endif
631
632 StatusCode status = parseFile_TE();
633
634#ifndef BEAN
635 if ( !status.isSuccess() ) {
636 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
637 }
638#else
639 if ( !status ) {
640 cout << "Magnetic field parse failled" << endl;
641 }
642#endif
643
644 for (int iC = 0; iC< 2; ++iC ){
645 m_min_FL[iC] = -90.0*cm;
646 m_max_FL[iC] = m_min_FL[iC]+( m_Nxyz[iC]-1 )* m_Dxyz[iC];
647 //for tof and emc
648 m_min_FL_TE[iC] = 0.0*cm;
649 m_max_FL_TE[iC] = m_min_FL_TE[iC]+( m_Nxyz_TE[iC]-1 )* m_Dxyz_TE[iC];
650 } // iC
651 m_min_FL[2] = -120.0*cm;
652 m_max_FL[2] = m_min_FL[2]+( m_Nxyz[2]-1 )* m_Dxyz[2];
653 //for tof and emc
654 m_min_FL_TE[2] = 0.0*cm;
655 m_max_FL_TE[2] = m_min_FL_TE[2]+( m_Nxyz_TE[2]-1 )* m_Dxyz_TE[2];
656}
int runNo
Definition: DQA_TO_DB.cxx:12
ConnectionDB::eRet getReadSC_MagnetInfo(std::vector< double > &current, int runNo)
ConnectionDB::eRet getBeamEnergy(std::vector< double > &beamE, int runNo)

Referenced by handle(), and init_params().

◆ init_params() [2/2]

void MagneticFieldSvc::init_params ( )

◆ initialize() [1/2]

StatusCode MagneticFieldSvc::initialize ( )
virtual

Initialise the service (Inherited Service overrides)

Definition at line 136 of file MagneticFieldSvc.cxx.

137{
138#ifndef BEAN
139 MsgStream log(msgSvc(), name());
140 StatusCode status = Service::initialize();
141 // Get the references to the services that are needed by the ApplicationMgr itself
142 IIncidentSvc* incsvc;
143 status = service("IncidentSvc", incsvc);
144 int priority = 100;
145 if( status.isSuccess() ){
146 incsvc -> addListener(this, "NewRun", priority);
147 }
148
149 status = serviceLocator()->service("EventDataSvc", m_eventSvc, true);
150 if (status.isFailure() ) {
151 log << MSG::ERROR << "Unable to find EventDataSvc " << endreq;
152 return status;
153 }
154#else
155 if( !init_mucMagneticField() ) return false;
156 StatusCode status = true;
157#endif
158
159 if(m_useDB == false) {
160 if(m_gridDistance == 5) {
161 m_Q.reserve(201250);
162 m_P.reserve(201250);
163 m_Q_TE.reserve(176608);
164 m_P_TE.reserve(176608);
165 }
166 if(m_gridDistance == 10) {
167 m_Q.reserve(27082);
168 m_P.reserve(27082);
169 m_Q_TE.reserve(23833);
170 m_P_TE.reserve(23833);
171 }
172
173 m_filename = path;
174 m_filename_TE = path;
175 if(m_gridDistance == 10) {
176 m_filename_TE += std::string("/dat/TEMap10cm.dat");
177 if(m_runmode == 2)
178 m_filename += std::string( "/dat/2008-04-03BESIII-field-Mode2.dat");
179 else if(m_runmode == 4)
180 m_filename += std::string( "/dat/2008-04-03BESIII-field-Mode4.dat");
181 else if(m_runmode == 6)
182 m_filename += std::string( "/dat/2008-04-03BESIII-field-Mode6.dat");
183 else if(m_runmode == 7)
184 m_filename += std::string( "/dat/2008-04-03BESIII-field-Mode7.dat");
185 else {
186 m_filename += std::string( "/dat/2007-08-28BESIII-field.dat");
187 std::cout<<"WARNING: YOU ARE USING OLD FIELD MAP!!!!"<<std::endl;
188 }
189 }
190 if(m_gridDistance == 5) {
191 m_filename_TE += std::string("/dat/TEMap5cm.dat");
192 if(m_runmode == 1)
193 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode1.dat");
194 else if(m_runmode == 2)
195 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode2.dat");
196 else if(m_runmode == 3)
197 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode3.dat");
198 else if(m_runmode == 4)
199 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode4.dat");
200 else if(m_runmode == 5)
201 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode5.dat");
202 else if(m_runmode == 6)
203 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode6.dat");
204 else if(m_runmode == 7)
205 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode7.dat");
206 else if(m_runmode == 8)
207 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode8.dat");
208 else {
209 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode2.dat");
210 std::cout<<"WARNING: NO RUN MODE, YOU ARE USING DEFAULT FIELD MAP!!!!"<<std::endl;
211 }
212 }
213 cout<<"*** Read field map: "<<m_filename<<endl;
214 cout<<"*** Read field map: "<<m_filename_TE<<endl;
215
216 status = parseFile();
217 status = parseFile_TE();
218
219#ifndef BEAN
220 if ( status.isSuccess() ) {
221 log << MSG::DEBUG << "Magnetic field parsed successfully" << endreq;
222#else
223 if ( status ) {
224 cout << "Magnetic field parsed successfully" << endl;
225#endif
226
227 for (int iC = 0; iC< 2; ++iC ){
228 m_min_FL[iC] = -90.0*cm;
229 m_max_FL[iC] = m_min_FL[iC]+( m_Nxyz[iC]-1 )* m_Dxyz[iC];
230 //for tof and emc
231 m_min_FL_TE[iC] = 0.0*cm;
232 m_max_FL_TE[iC] = m_min_FL_TE[iC]+( m_Nxyz_TE[iC]-1 )* m_Dxyz_TE[iC];
233 } // iC
234 m_min_FL[2] = -120.0*cm;
235 m_max_FL[2] = m_min_FL[2]+( m_Nxyz[2]-1 )* m_Dxyz[2];
236 //for tof and emc
237 m_min_FL_TE[2] = 0.0*cm;
238 m_max_FL_TE[2] = m_min_FL_TE[2]+( m_Nxyz_TE[2]-1 )* m_Dxyz_TE[2];
239 }
240 else {
241#ifndef BEAN
242 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
243#else
244 cout << "Magnetic field parse failled" << endl;
245#endif
246 }
247 }
248 else {
249 m_connect_run = new FieldDBUtil::ConnectionDB();
250#ifndef BEAN
251 if (!m_connect_run) {
252 log << MSG::ERROR << "Could not open connection to database" << endreq;
253 }
254#endif
255 }
256
257 return status;
258}

◆ initialize() [2/2]

virtual StatusCode MagneticFieldSvc::initialize ( )
virtual

Initialise the service (Inherited Service overrides)

◆ queryInterface() [1/2]

StatusCode MagneticFieldSvc::queryInterface ( const InterfaceID &  riid,
void **  ppvInterface 
)
virtual

Query the available interfaces.

Parameters
riidRequested interface ID
ppvInterfacePointer to requested interface
Returns
StatusCode indicating SUCCESS or FAILURE.

Definition at line 676 of file MagneticFieldSvc.cxx.

678{
679 StatusCode sc = StatusCode::FAILURE;
680 if ( ppvInterface ) {
681 *ppvInterface = 0;
682
683 if ( riid == IID_IMagneticFieldSvc ) {
684 *ppvInterface = static_cast<IMagneticFieldSvc*>(this);
685 sc = StatusCode::SUCCESS;
686 addRef();
687 }
688 else
689 sc = Service::queryInterface( riid, ppvInterface );
690 }
691 return sc;
692}

◆ queryInterface() [2/2]

virtual StatusCode MagneticFieldSvc::queryInterface ( const InterfaceID &  riid,
void **  ppvInterface 
)
virtual

Query the available interfaces.

Parameters
riidRequested interface ID
ppvInterfacePointer to requested interface
Returns
StatusCode indicating SUCCESS or FAILURE.

◆ type() [1/2]

virtual const InterfaceID & MagneticFieldSvc::type ( ) const
inlinevirtual

Service type.

Definition at line 74 of file InstallArea/include/MagneticField/MagneticField/MagneticFieldSvc.h.

static const InterfaceID & interfaceID()
Retrieve interface ID.

◆ type() [2/2]

virtual const InterfaceID & MagneticFieldSvc::type ( ) const
inlinevirtual

◆ uniFieldVector() [1/2]

StatusCode MagneticFieldSvc::uniFieldVector ( const HepPoint3D xyz,
HepVector3D fvec 
) const
virtual

Implements IMagneticFieldSvc.

Definition at line 1951 of file MagneticFieldSvc.cxx.

1953{
1954 if(m_runmode == 8 || m_runmode == 7) {
1955 if(-1.5*m<=r.z()&&r.z()<=1.5*m&&0*m<=std::sqrt(r.x()*r.x()+r.y()*r.y())&&std::sqrt(r.x()*r.x()+r.y()*r.y())<=1.5*m)
1956 {
1957 b[0]=0.0;
1958 b[1]=0.0;
1959 b[2]=-0.9*tesla;
1960 }
1961 else
1962 {
1963 b[0]=0.0;
1964 b[1]=0.0;
1965 b[2]=0.0;
1966 }
1967 }
1968 else {
1969 if(-1.5*m<=r.z()&&r.z()<=1.5*m&&0*m<=std::sqrt(r.x()*r.x()+r.y()*r.y())&&std::sqrt(r.x()*r.x()+r.y()*r.y())<=1.5*m)
1970 {
1971 b[0]=0.0;
1972 b[1]=0.0;
1973 b[2]=-1.0*tesla;
1974 }
1975 else
1976 {
1977 b[0]=0.0;
1978 b[1]=0.0;
1979 b[2]=0.0;
1980 }
1981 }
1982
1983 if(m_turnOffField == true) {
1984 b[0] = 0.;
1985 b[1] = 0.;
1986 b[2] = 0.;
1987 }
1988 //yzhang add 2012-04-25
1989 b[0] *= m_scale;
1990 b[1] *= m_scale;
1991 b[2] *= m_scale;
1992#ifndef BEAN
1993 return StatusCode::SUCCESS;
1994#else
1995 return true;
1996#endif
1997}

Referenced by fieldVector().

◆ uniFieldVector() [2/2]

virtual StatusCode MagneticFieldSvc::uniFieldVector ( const HepPoint3D xyz,
HepVector3D fvec 
) const
virtual

Implements IMagneticFieldSvc.

Friends And Related Function Documentation

◆ SvcFactory< MagneticFieldSvc >

friend class SvcFactory< MagneticFieldSvc >
friend

Allow SvcFactory to instantiate the service.

Definition at line 130 of file InstallArea/include/MagneticField/MagneticField/MagneticFieldSvc.h.


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