BOSS 7.0.5
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 ()
 
void init_params (std::vector< double > current, std::vector< double > beamEnergy, int runNo)
 
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 ()
 
void init_params (std::vector< double > current, std::vector< double > beamEnergy, int runNo)
 
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
 

Public Attributes

bool m_readOneTime
 
int m_runFrom
 
int m_runTo
 
int runNo
 
IDataProviderSvc * m_eventSvc
 
FieldDBUtil::ConnectionDBm_connect_run
 
std::map< int, std::vector< double > > m_mapMagnetInfo
 
std::map< int, std::vector< double > > m_mapBeamEnergy
 
std::vector< double > beamEnergy
 
std::vector< double > current
 

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 declareProperty( "ReadOneTime", m_readOneTime=false);
65 declareProperty( "RunFrom", m_runFrom=8093);
66 declareProperty("RunTo",m_runTo=9025);
67
68 m_Mucfield = new MucMagneticField();
69 if(!m_Mucfield) cout<<"error: can not get MucMagneticField pointer"<<endl;
70
71 m_zfield = -1.0*tesla;
72
73 if(getenv("MAGNETICFIELDROOT") != NULL) {
74 path = std::string(getenv( "MAGNETICFIELDROOT" ));
75 } else {
76 std::cerr<<"Couldn't find MAGNETICFIELDROOT"<<std::endl;
77 }
78
79}

◆ ~MagneticFieldSvc() [1/2]

MagneticFieldSvc::~MagneticFieldSvc ( )
protectedvirtual

Virtual destructor.

Definition at line 112 of file MagneticFieldSvc.cxx.

113{
114 if(m_Mucfield) delete m_Mucfield;
115}

◆ 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 1197 of file MagneticFieldSvc.cxx.

1199{
1200
1201 if(m_turnOffField == true) {
1202 newb[0] = 0.;
1203 newb[1] = 0.;
1204 newb[2] = 0.;
1205#ifndef BEAN
1206 return StatusCode::SUCCESS;
1207#else
1208 return true;
1209#endif
1210 }
1211
1212 // wll added 2012-08-27
1213 if(m_uniField) {
1214 uniFieldVector(newr,newb);
1215#ifndef BEAN
1216 return StatusCode::SUCCESS;
1217#else
1218 return true;
1219#endif
1220 }
1221
1222
1223 //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
1224 double new_x = -newr.x();
1225 double new_y = newr.y();
1226 double new_z = -newr.z();
1227 HepPoint3D r(new_x,new_y,new_z);
1228 HepVector3D b;
1229 b[0]=0.0*tesla;
1230 b[1]=0.0*tesla;
1231 b[2]=0.0*tesla;
1232 // This routine is now dummy. Was previously converting to/from CLHEP units
1233 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)
1234 {
1235 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)
1236 //if(std::abs(r.z())<1.2*m&&std::abs(r.x())<=0.9*m&&std::abs(r.y())<=0.9*m)
1237 {
1238 this->fieldGrid( r, b );
1239 }
1240 else
1241 {
1242 this->fieldGrid_TE( r, b );
1243 }
1244 }
1245 if((fabs(r.z())<=1970*mm && sqrt(r.x()*r.x()+r.y()*r.y()) >= 1740*mm) || (fabs(r.z())>=2050*mm))
1246 {
1247 HepPoint3D mr;
1248 HepVector3D tb;
1249 int part = 0, layer = 0, mat = 0;
1250 double theta;
1251 bool ifbar = false;
1252 bool ifend = false;
1253 if(r.x()!=0.){
1254 theta = atan2(fabs(r.y()),fabs(r.x()));
1255 if(0<=theta&&theta<pi/8) {
1256 mr[0] = fabs(r.x()); mr[1] = -fabs(r.y()); mr[2] = fabs(r.z());
1257 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm){
1258 part = 1;
1259 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1260 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1261 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1262 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1263 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1264 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1265 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1266 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1267 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1268 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1269 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1270 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1271 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1272 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1273 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1274 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1275 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1276 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1277 b[0] = tb[0];
1278 b[1] = -tb[1];
1279 b[2] = tb[2];
1280 ifbar = true;
1281 }
1282 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1283 part = 0; layer = 0; mat = 0;
1284 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1285 b[0] = tb[0];
1286 b[1] = -tb[1];
1287 b[2] = tb[2];
1288 ifend = true;
1289 }
1290 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1291 part = 0; layer = 0; mat = 1;
1292 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1293 b[0] = tb[0];
1294 b[1] = -tb[1];
1295 b[2] = tb[2];
1296 ifend = true;
1297 }
1298 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1299 part = 0; layer = 1; mat = 0;
1300 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1301 b[0] = tb[0];
1302 b[1] = -tb[1];
1303 b[2] = tb[2];
1304 ifend = true;
1305 }
1306 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1307 part = 0; layer = 1; mat = 1;
1308 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1309 b[0] = tb[0];
1310 b[1] = -tb[1];
1311 b[2] = tb[2];
1312 ifend = true;
1313 }
1314 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1315 part = 0; layer = 2; mat = 0;
1316 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1317 b[0] = tb[0];
1318 b[1] = -tb[1];
1319 b[2] = tb[2];
1320 ifend = true;
1321 }
1322 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1323 part = 0; layer = 2; mat = 1;
1324 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1325 b[0] = tb[0];
1326 b[1] = -tb[1];
1327 b[2] = tb[2];
1328 ifend = true;
1329 }
1330 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1331 part = 0; layer = 3; mat = 0;
1332 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1333 b[0] = tb[0];
1334 b[1] = -tb[1];
1335 b[2] = tb[2];
1336 ifend = true;
1337 }
1338 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1339 part = 0; layer = 3; mat = 1;
1340 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1341 b[0] = tb[0];
1342 b[1] = -tb[1];
1343 b[2] = tb[2];
1344 ifend = true;
1345 }
1346 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1347 part = 0; layer = 4; mat = 0;
1348 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1349 b[0] = tb[0];
1350 b[1] = -tb[1];
1351 b[2] = tb[2];
1352 ifend = true;
1353 }
1354 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1355 part = 0; layer = 4; mat = 1;
1356 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1357 b[0] = tb[0];
1358 b[1] = -tb[1];
1359 b[2] = tb[2];
1360 ifend = true;
1361 }
1362 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1363 part = 0; layer = 5; mat = 0;
1364 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1365 b[0] = tb[0];
1366 b[1] = -tb[1];
1367 b[2] = tb[2];
1368 ifend = true;
1369 }
1370 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1371 part = 0; layer = 5; mat =1;
1372 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1373 b[0] = tb[0];
1374 b[1] = -tb[1];
1375 b[2] = tb[2];
1376 ifend = true;
1377 }
1378 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1379 part = 0; layer = 6; mat = 0;
1380 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1381 b[0] = tb[0];
1382 b[1] = -tb[1];
1383 b[2] = tb[2];
1384 ifend = true;
1385 }
1386 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1387 part = 0; layer = 6; mat = 1;
1388 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1389 b[0] = tb[0];
1390 b[1] = -tb[1];
1391 b[2] = tb[2];
1392 ifend = true;
1393 }
1394 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1395 part = 0; layer = 7; mat = 0;
1396 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1397 b[0] = tb[0];
1398 b[1] = -tb[1];
1399 b[2] = tb[2];
1400 ifend = true;
1401 }
1402 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1403 part = 0; layer = 7; mat = 1;
1404 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1405 b[0] = tb[0];
1406 b[1] = -tb[1];
1407 b[2] = tb[2];
1408 ifend = true;
1409 }
1410 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1411 part = 0; layer = 8; mat = 0;
1412 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1413 b[0] = tb[0];
1414 b[1] = -tb[1];
1415 b[2] = tb[2];
1416 ifend = true;
1417 }
1418 }
1419 if(pi/8<=theta&&theta<pi/4) {
1420 mr[0] = fabs(r.x())*cos(pi/4)+fabs(r.y())*sin(pi/4);
1421 mr[1] = -fabs(r.x())*sin(pi/4)+fabs(r.y())*cos(pi/4);
1422 mr[2] = fabs(r.z());
1423 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1424 part = 1;
1425 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1426 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1427 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1428 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1429 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1430 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1431 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1432 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1433 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1434 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1435 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1436 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1437 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1438 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1439 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1440 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1441 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1442 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1443 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1444 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1445 b[2] = tb[2];
1446 ifbar = true;
1447 }
1448 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1449 part = 0; layer = 0; mat = 0;
1450 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1451 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1452 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1453 b[2] = tb[2];
1454 ifend = true;
1455 }
1456 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1457 part = 0; layer = 0; mat = 1;
1458 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1459 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1460 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1461 b[2] = tb[2];
1462 ifend = true;
1463 }
1464 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1465 part = 0; layer = 1; mat = 0;
1466 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1467 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1468 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1469 b[2] = tb[2];
1470 ifend = true;
1471 }
1472 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1473 part = 0; layer = 1; mat = 1;
1474 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1475 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1476 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1477 b[2] = tb[2];
1478 ifend = true;
1479 }
1480 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1481 part = 0; layer = 2; mat = 0;
1482 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1483 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1484 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1485 b[2] = tb[2];
1486 ifend = true;
1487 }
1488 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1489 part = 0; layer = 2; mat = 1;
1490 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1491 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1492 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1493 b[2] = tb[2];
1494 ifend = true;
1495 }
1496 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1497 part = 0; layer = 3; mat = 0;
1498 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1499 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1500 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1501 b[2] = tb[2];
1502 ifend = true;
1503 }
1504 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1505 part = 0; layer = 3; mat = 1;
1506 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1507 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1508 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1509 b[2] = tb[2];
1510 ifend = true;
1511 }
1512 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1513 part = 0; layer = 4; mat = 0;
1514 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1515 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1516 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1517 b[2] = tb[2];
1518 ifend = true;
1519 }
1520 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1521 part = 0; layer = 4; mat = 1;
1522 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1523 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1524 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1525 b[2] = tb[2];
1526 ifend = true;
1527 }
1528 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1529 part = 0; layer = 5; mat = 0;
1530 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1531 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1532 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1533 b[2] = tb[2];
1534 ifend = true;
1535 }
1536 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1537 part = 0; layer = 5; mat =1;
1538 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1539 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1540 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1541 b[2] = tb[2];
1542 ifend = true;
1543 }
1544 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1545 part = 0; layer = 6; mat = 0;
1546 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1547 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1548 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1549 b[2] = tb[2];
1550 ifend = true;
1551 }
1552 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1553 part = 0; layer = 6; mat = 1;
1554 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1555 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1556 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1557 b[2] = tb[2];
1558 ifend = true;
1559 }
1560 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1561 part = 0; layer = 7; mat = 0;
1562 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1563 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1564 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1565 b[2] = tb[2];
1566 ifend = true;
1567 }
1568 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1569 part = 0; layer = 7; mat = 1;
1570 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1571 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1572 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1573 b[2] = tb[2];
1574 ifend = true;
1575 }
1576 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1577 part = 0; layer = 8; mat = 0;
1578 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1579 b[0] = tb[0]*cos(pi/4)-tb[1]*sin(pi/4);
1580 b[1] = tb[0]*sin(pi/4)+tb[1]*cos(pi/4);
1581 b[2] = tb[2];
1582 ifend = true;
1583 }
1584 }
1585 if(pi/4<=theta&&theta<3*pi/8) {
1586 mr[0] = fabs(r.x())*cos(pi/4)+fabs(r.y())*sin(pi/4);
1587 mr[1] = fabs(r.x())*sin(pi/4)-fabs(r.y())*cos(pi/4);
1588 mr[2] = fabs(r.z());
1589 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1590 part = 1;
1591 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1592 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1593 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1594 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1595 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1596 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1597 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1598 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1599 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1600 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1601 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1602 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1603 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1604 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1605 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1606 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1607 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1608 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1609 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1610 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1611 b[2] = tb[2];
1612 ifbar = true;
1613 }
1614 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1615 part = 0; layer = 0; mat = 0;
1616 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1617 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1618 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1619 b[2] = tb[2];
1620 ifend = true;
1621 }
1622 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1623 part = 0; layer = 0; mat = 1;
1624 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1625 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1626 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1627 b[2] = tb[2];
1628 ifend = true;
1629 }
1630 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1631 part = 0; layer = 1; mat = 0;
1632 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1633 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1634 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1635 b[2] = tb[2];
1636 ifend = true;
1637 }
1638 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1639 part = 0; layer = 1; mat = 1;
1640 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1641 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1642 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1643 b[2] = tb[2];
1644 ifend = true;
1645 }
1646 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1647 part = 0; layer = 2; mat = 0;
1648 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1649 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1650 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1651 b[2] = tb[2];
1652 ifend = true;
1653 }
1654 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1655 part = 0; layer = 2; mat = 1;
1656 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1657 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1658 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1659 b[2] = tb[2];
1660 ifend = true;
1661 }
1662 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1663 part = 0; layer = 3; mat = 0;
1664 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1665 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1666 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1667 b[2] = tb[2];
1668 ifend = true;
1669 }
1670 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1671 part = 0; layer = 3; mat = 1;
1672 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1673 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1674 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1675 b[2] = tb[2];
1676 ifend = true;
1677 }
1678 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1679 part = 0; layer = 4; mat = 0;
1680 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1681 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1682 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1683 b[2] = tb[2];
1684 ifend = true;
1685 }
1686 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1687 part = 0; layer = 4; mat = 1;
1688 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1689 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1690 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1691 b[2] = tb[2];
1692 ifend = true;
1693 }
1694 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1695 part = 0; layer = 5; mat = 0;
1696 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1697 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1698 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1699 b[2] = tb[2];
1700 ifend = true;
1701 }
1702 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1703 part = 0; layer = 5; mat =1;
1704 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1705 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1706 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1707 b[2] = tb[2];
1708 ifend = true;
1709 }
1710 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1711 part = 0; layer = 6; mat = 0;
1712 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1713 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1714 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1715 b[2] = tb[2];
1716 ifend = true;
1717 }
1718 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1719 part = 0; layer = 6; mat = 1;
1720 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1721 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1722 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1723 b[2] = tb[2];
1724 ifend = true;
1725 }
1726 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1727 part = 0; layer = 7; mat = 0;
1728 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1729 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1730 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1731 b[2] = tb[2];
1732 ifend = true;
1733 }
1734 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1735 part = 0; layer = 7; mat = 1;
1736 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1737 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1738 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1739 b[2] = tb[2];
1740 ifend = true;
1741 }
1742 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1743 part = 0; layer = 8; mat = 0;
1744 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1745 b[0] = tb[0]*cos(pi/4)+tb[1]*sin(pi/4);
1746 b[1] = tb[0]*sin(pi/4)-tb[1]*cos(pi/4);
1747 b[2] = tb[2];
1748 ifend = true;
1749 }
1750 }
1751 if(3*pi/8<=theta&&theta<pi/2) {
1752 mr[0] = fabs(r.y()); mr[1] = -fabs(r.x()); mr[2] = fabs(r.z());
1753 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1754 part = 1;
1755 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1756 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1757 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1758 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1759 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1760 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1761 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1762 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1763 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1764 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1765 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1766 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1767 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1768 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1769 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1770 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1771 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
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 ifbar = true;
1777 }
1778 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1779 part = 0; layer = 0; 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(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1787 part = 0; layer = 0; 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(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1795 part = 0; layer = 1; 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(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1803 part = 0; layer = 1; 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(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1811 part = 0; layer = 2; 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(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1819 part = 0; layer = 2; 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(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1827 part = 0; layer = 3; 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(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1835 part = 0; layer = 3; 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(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
1843 part = 0; layer = 4; 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(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1851 part = 0; layer = 4; 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(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
1859 part = 0; layer = 5; 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(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1867 part = 0; layer = 5; 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(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
1875 part = 0; layer = 6; 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(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1883 part = 0; layer = 6; 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(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1891 part = 0; layer = 7; 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 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
1899 part = 0; layer = 7; mat = 1;
1900 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1901 b[0] = -tb[1];
1902 b[1] = tb[0];
1903 b[2] = tb[2];
1904 ifend = true;
1905 }
1906 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
1907 part = 0; layer = 8; mat = 0;
1908 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1909 b[0] = -tb[1];
1910 b[1] = tb[0];
1911 b[2] = tb[2];
1912 ifend = true;
1913 }
1914 }
1915 }
1916 if(r.x()==0.) {
1917 mr[0] = fabs(r.y()); mr[1] = -fabs(r.x()); mr[2] = fabs(r.z());
1918 if(mr[2]<=1970*mm&&1740*mm<=mr[0]&&mr[0]<=2620*mm) {
1919 part = 1;
1920 if(1740*mm<=mr[0]&&mr[0]<1770*mm) { layer = 0; mat = 0; }
1921 if(1770*mm<=mr[0]&&mr[0]<1810*mm) { layer = 0; mat = 1; }
1922 if(1810*mm<=mr[0]&&mr[0]<1840*mm) { layer = 1; mat = 0; }
1923 if(1840*mm<=mr[0]&&mr[0]<1880*mm) { layer = 1; mat = 1; }
1924 if(1880*mm<=mr[0]&&mr[0]<1910*mm) { layer = 2; mat = 0; }
1925 if(1910*mm<=mr[0]&&mr[0]<1950*mm) { layer = 2; mat = 1; }
1926 if(1950*mm<=mr[0]&&mr[0]<1990*mm) { layer = 3; mat = 0; }
1927 if(1990*mm<=mr[0]&&mr[0]<2030*mm) { layer = 3; mat = 1; }
1928 if(2030*mm<=mr[0]&&mr[0]<2070*mm) { layer = 4; mat = 0; }
1929 if(2070*mm<=mr[0]&&mr[0]<2110*mm) { layer = 4; mat = 1; }
1930 if(2110*mm<=mr[0]&&mr[0]<2190*mm) { layer = 5; mat = 0; }
1931 if(2190*mm<=mr[0]&&mr[0]<2230*mm) { layer = 5; mat = 1; }
1932 if(2230*mm<=mr[0]&&mr[0]<2310*mm) { layer = 6; mat = 0; }
1933 if(2310*mm<=mr[0]&&mr[0]<2350*mm) { layer = 6; mat = 1; }
1934 if(2350*mm<=mr[0]&&mr[0]<2430*mm) { layer = 7; mat = 0; }
1935 if(2430*mm<=mr[0]&&mr[0]<2470*mm) { layer = 7; mat = 1; }
1936 if(2470*mm<=mr[0]&&mr[0]<=2620*mm) { layer = 8; mat = 0; }
1937 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1938 b[0] = -tb[1];
1939 b[1] = tb[0];
1940 b[2] = tb[2];
1941 ifbar = true;
1942 }
1943 if(2050*mm<=mr[2]&&mr[2]<2090*mm&&1034*mm<=mr[0]&&mr[0]<=2500*mm){
1944 part = 0; layer = 0; mat = 0;
1945 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1946 b[0] = -tb[1];
1947 b[1] = tb[0];
1948 b[2] = tb[2];
1949 ifend = true;
1950 }
1951 if(2090*mm<=mr[2]&&mr[2]<2130*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1952 part = 0; layer = 0; mat = 1;
1953 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1954 b[0] = -tb[1];
1955 b[1] = tb[0];
1956 b[2] = tb[2];
1957 ifend = true;
1958 }
1959 if(2130*mm<=mr[2]&&mr[2]<2170*mm&&1067*mm<=mr[0]&&mr[0]<=2500*mm) {
1960 part = 0; layer = 1; mat = 0;
1961 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1962 b[0] = -tb[1];
1963 b[1] = tb[0];
1964 b[2] = tb[2];
1965 ifend = true;
1966 }
1967 if(2170*mm<=mr[2]&&mr[2]<2210*mm&&1100*mm<=mr[0]&&mr[0]<=2500*mm) {
1968 part = 0; layer = 1; mat = 1;
1969 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1970 b[0] = -tb[1];
1971 b[1] = tb[0];
1972 b[2] = tb[2];
1973 ifend = true;
1974 }
1975 if(2210*mm<=mr[2]&&mr[2]<2240*mm&&1100*mm<mr[0]&&mr[0]<=2500*mm) {
1976 part = 0; layer = 2; mat = 0;
1977 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1978 b[0] = -tb[1];
1979 b[1] = tb[0];
1980 b[2] = tb[2];
1981 ifend = true;
1982 }
1983 if(2240*mm<=mr[2]&&mr[2]<2280*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1984 part = 0; layer = 2; mat = 1;
1985 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1986 b[0] = -tb[1];
1987 b[1] = tb[0];
1988 b[2] = tb[2];
1989 ifend = true;
1990 }
1991 if(2280*mm<=mr[2]&&mr[2]<2310*mm&&1133*mm<=mr[0]&&mr[0]<=2500*mm) {
1992 part = 0; layer = 3; mat = 0;
1993 m_Mucfield->getMucField(part,layer,mat,mr,tb);
1994 b[0] = -tb[1];
1995 b[1] = tb[0];
1996 b[2] = tb[2];
1997 ifend = true;
1998 }
1999 if(2310*mm<=mr[2]&&mr[2]<2350*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
2000 part = 0; layer = 3; mat = 1;
2001 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2002 b[0] = -tb[1];
2003 b[1] = tb[0];
2004 b[2] = tb[2];
2005 ifend = true;
2006 }
2007 if(2350*mm<=mr[2]&&mr[2]<2380*mm&&1167*mm<=mr[0]&&mr[0]<=2500*mm) {
2008 part = 0; layer = 4; mat = 0;
2009 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2010 b[0] = -tb[1];
2011 b[1] = tb[0];
2012 b[2] = tb[2];
2013 ifend = true;
2014 }
2015 if(2380*mm<=mr[2]&&mr[2]<2420*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
2016 part = 0; layer = 4; mat = 1;
2017 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2018 b[0] = -tb[1];
2019 b[1] = tb[0];
2020 b[2] = tb[2];
2021 ifend = true;
2022 }
2023 if(2420*mm<=mr[2]&&mr[2]<2470*mm&&1203*mm<=mr[0]&&mr[0]<=2500*mm) {
2024 part = 0; layer = 5; mat = 0;
2025 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2026 b[0] = -tb[1];
2027 b[1] = tb[0];
2028 b[2] = tb[2];
2029 ifend = true;
2030 }
2031 if(2470*mm<=mr[2]&&mr[2]<2510*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
2032 part = 0; layer = 5; mat =1;
2033 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2034 b[0] = -tb[1];
2035 b[1] = tb[0];
2036 b[2] = tb[2];
2037 ifend = true;
2038 }
2039 if(2510*mm<=mr[2]&&mr[2]<2590*mm&&1241*mm<=mr[0]&&mr[0]<=2500*mm) {
2040 part = 0; layer = 6; mat = 0;
2041 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2042 b[0] = -tb[1];
2043 b[1] = tb[0];
2044 b[2] = tb[2];
2045 ifend = true;
2046 }
2047 if(2590*mm<=mr[2]&&mr[2]<2630*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
2048 part = 0; layer = 6; mat = 1;
2049 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2050 b[0] = -tb[1];
2051 b[1] = tb[0];
2052 b[2] = tb[2];
2053 ifend = true;
2054 }
2055 if(2630*mm<=mr[2]&&mr[2]<2710*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
2056 part = 0; layer = 7; mat = 0;
2057 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2058 b[0] = -tb[1];
2059 b[1] = tb[0];
2060 b[2] = tb[2];
2061 ifend = true;
2062 }
2063 if(2710*mm<=mr[2]&&mr[2]<2750*mm&&1362*mm<=mr[0]&&mr[0]<=2500*mm) {
2064 part = 0; layer = 7; mat = 1;
2065 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2066 b[0] = -tb[1];
2067 b[1] = tb[0];
2068 b[2] = tb[2];
2069 ifend = true;
2070 }
2071 if(2750*mm<=mr[2]&&mr[2]<=2800*mm&&1302*mm<=mr[0]&&mr[0]<=2500*mm) {
2072 part = 0; layer = 8; mat = 0;
2073 m_Mucfield->getMucField(part,layer,mat,mr,tb);
2074 b[0] = -tb[1];
2075 b[1] = tb[0];
2076 b[2] = tb[2];
2077 ifend = true;
2078 }
2079 }
2080 if(ifbar==true||ifend==true) {
2081 if( r.x() < 0. && r.y() >= 0. && r.z() > 0. ){
2082 b[0] = -b[0];
2083 }
2084 else if( r.x() <= 0. && r.y() < 0. && r.z() > 0. ){
2085 b[0] = -b[0];
2086 b[1] = -b[1];
2087 }
2088 else if( r.x() > 0. && r.y() < 0. && r.z() > 0. ){
2089 b[1] = -b[1];
2090 }
2091 else if( r.x() >= 0. && r.y() > 0. && r.z() < 0. ){
2092 b[0] = -b[0];
2093 b[1] = -b[1];
2094 }
2095 else if( r.x() < 0. && r.y() >= 0. && r.z() < 0. ){
2096 b[1] = -b[1];
2097 }
2098 else if( r.x() <= 0. && r.y() < 0. && r.z() < 0. ){
2099 b[0] = b[0];
2100 b[1] = b[1];
2101 }
2102 else if( r.x() > 0. && r.y() <= 0. && r.z() < 0. ){
2103 b[0] = -b[0];
2104 }
2105 }
2106 }
2107
2108 //reference frames defined by simulation and SCM are different.
2109 newb[0] = -b[0] * m_scale;
2110 newb[1] = b[1] * m_scale;
2111 newb[2] = -b[2] * m_scale;
2112/* 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) {
2113 return StatusCode::SUCCESS;
2114 }
2115 else {
2116 newb[0] = newb[0] - 0.10*tesla;
2117 newb[1] = newb[1] + 0.10*tesla;
2118 newb[2] = newb[2] - 0.10*tesla;
2119 }*/
2120
2121 //newb[0] = b[0];
2122 //newb[1] = b[1];
2123 //newb[2] = b[2];
2124
2125#ifndef BEAN
2126 return StatusCode::SUCCESS;
2127#else
2128 return true;
2129#endif
2130}
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)
const double b
Definition: slope.cxx:9

◆ 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 809 of file MagneticFieldSvc.cxx.

810{
811 MsgStream log( msgSvc(), name() );
812 //if(m_connect_run) delete m_connect_run;
813 StatusCode status = Service::finalize();
814
815 if ( status.isSuccess() )
816 log << MSG::INFO << "Service finalized successfully" << endreq;
817 return status;
818}

◆ finalize() [2/2]

virtual StatusCode MagneticFieldSvc::finalize ( )
virtual

Finalise the service.

◆ getReferField() [1/2]

double MagneticFieldSvc::getReferField ( )
virtual

Implements IMagneticFieldSvc.

Definition at line 2180 of file MagneticFieldSvc.cxx.

2181{
2182 if(m_useDB == false) {
2183 if(m_runmode == 8 || m_runmode == 7) m_zfield = -0.9*tesla;
2184 else m_zfield = -1.0*tesla;
2185 }
2186
2187 if(m_turnOffField == true) {
2188 m_zfield = 0.;
2189 }
2190 return m_zfield * m_scale;
2191}

◆ getReferField() [2/2]

virtual double MagneticFieldSvc::getReferField ( )
virtual

Implements IMagneticFieldSvc.

◆ handle() [1/2]

void MagneticFieldSvc::handle ( const Incident &  inc)

Definition at line 841 of file MagneticFieldSvc.cxx.

841 {
842 MsgStream log( messageService(), name() );
843 log << MSG::DEBUG << "handle: " << inc.type() << endreq;
844 if ( inc.type() != "NewRun" ){
845 return;
846 }
847 log << MSG::DEBUG << "Begin New Runcc" << endreq;
848
849 SmartDataPtr<Event::EventHeader> eventHeader(m_eventSvc,"/Event/EventHeader");
850 int new_run = eventHeader->runNumber();
851 //cout << endl << "New Run is: " << new_run << " MagnetSvc= "<<endl;
852 if(new_run<0) new_run=-new_run;
853
854
855 if((m_useDB == true)&&(m_readOneTime == false)) {
857 ConnectionDB::eRet e = m_connect_run->getReadSC_MagnetInfo(current,std::abs(new_run));
858 e = m_connect_run->getBeamEnergy(beamEnergy, std::abs(new_run));
860 //cout<<"beamEnergy is:"<<beamEnergy[0]<<endl;
861 }
862 else if ((m_useDB == true)&&(m_readOneTime == true))
863 {
864 std::vector<double> beamEnergy;
865 std::vector<double> current(3,0.);
866 //if (m_mapBeamEnergy == NULL ||m_mapBeamEnergy.size() == 0))
867 if (m_mapBeamEnergy.size() == 0)
868 {
869 cout<<"Run:"<<new_run<<" BeamEnergy is NULL, please check!!"<<endl;
870 exit(1);
871 }
872 beamEnergy.push_back((m_mapBeamEnergy[new_run])[0]);
873 beamEnergy.push_back((m_mapBeamEnergy[new_run])[1]);
874 //if (m_mapMagnetInfo.size()<1||m_mapMagnetInfo.isEmpty())
875 if(m_mapMagnetInfo.size() == 0)
876 {
877 cout<<"Run:"<<new_run<<" MagnetInfo is NULL, please check!!"<<endl;
878 exit(1);
879 }
880 current[0]=(m_mapMagnetInfo[new_run])[0];
881 current[1]=(m_mapMagnetInfo[new_run])[1];
882 current[2]=(m_mapMagnetInfo[new_run])[2];
884 }
885}
ConnectionDB::eRet getReadSC_MagnetInfo(std::vector< double > &current, int runNo)
ConnectionDB::eRet getBeamEnergy(std::vector< double > &beamE, int runNo)

◆ handle() [2/2]

void MagneticFieldSvc::handle ( const Incident &  )

◆ ifRealField() [1/2]

bool MagneticFieldSvc::ifRealField ( ) const
virtual

Implements IMagneticFieldSvc.

Definition at line 2193 of file MagneticFieldSvc.cxx.

2193 {
2194 return m_ifRealField;
2195}

◆ ifRealField() [2/2]

virtual bool MagneticFieldSvc::ifRealField ( ) const
virtual

Implements IMagneticFieldSvc.

◆ init_params() [1/4]

void MagneticFieldSvc::init_params ( )

Referenced by handle(), and init_params().

◆ init_params() [2/4]

void MagneticFieldSvc::init_params ( )

◆ init_params() [3/4]

void MagneticFieldSvc::init_params ( std::vector< double >  current,
std::vector< double >  beamEnergy,
int  runNo 
)

Definition at line 302 of file MagneticFieldSvc.cxx.

303{
304 MsgStream log(msgSvc(), name());
305#else
306void init_params(int runNo)
307{
308 if( !init_mucMagneticField() ) {
309 cerr << " STOP! " << endl;
310 exit(1);
311 }
312#endif
313
314 m_Q.clear();
315 m_P.clear();
316 m_Q_1.clear();
317 m_P_1.clear();
318 m_Q_2.clear();
319 m_P_2.clear();
320 m_P_TE.clear();
321 m_Q_TE.clear();
322
323 if(m_gridDistance == 5) {
324 m_Q.reserve(201250);
325 m_P.reserve(201250);
326 m_Q_1.reserve(201250);
327 m_P_1.reserve(201250);
328 m_Q_2.reserve(201250);
329 m_P_2.reserve(201250);
330 m_Q_TE.reserve(176608);
331 m_P_TE.reserve(176608);
332 }
333 if(m_gridDistance == 10) {
334 m_Q.reserve(27082);
335 m_P.reserve(27082);
336 m_Q_1.reserve(27082);
337 m_P_1.reserve(27082);
338 m_Q_2.reserve(27082);
339 m_P_2.reserve(27082);
340 m_Q_TE.reserve(23833);
341 m_P_TE.reserve(23833);
342 }
343
344 char BPR_PRB[200];
345 char BER_PRB[200];
346 sprintf(BPR_PRB,"%f ",beamEnergy[0]);
347 sprintf(BER_PRB,"%f ",beamEnergy[1]);
348 setenv("BEPCII_INFO.BPR_PRB", BPR_PRB, 1);
349 setenv("BEPCII_INFO.BER_PRB", BER_PRB, 1);
350 int ssm_curr = int (current[0]);
351 double scql_curr = current[1];
352 double scqr_curr = current[2];
353// cout<<"curent is:"<<ssm_curr<<"scql_curr is:"<<scql_curr<<"scqr_curr is:"<<scqr_currendl;
354
355
356
357#ifndef BEAN
358 log << MSG::INFO << "SSM current: " << current[0] << " SCQL current: " << current[1] << " SCQR current: " << current[2] << " in Run " << runNo << endreq;
359 log << MSG::INFO << "beamEnergy is:"<< beamEnergy[0] <<" BER_PRB:"<< beamEnergy[1] << " in Run " << runNo << endreq;
360#else
361
362 cout << "SSM current11: " << current[0] << " SCQL current: " << current[1]
363 << " SCQR current: " << current[2] << " in Run " << runNo << endl;
364 log << MSG::INFO << "beamEnergy is00:"<< beamEnergy[0] <<" BER_PRB:"<< beamEnergy[1] << " in Run " << runNo << endreq;
365#endif
366
367 //int ssm_curr = 3369;
368 //double scql_curr = 426.2;
369 //double scqr_curr = 426.2;
370 //for the energy less than 1.89GeV
371 if((ssm_curr >= 3367) && (ssm_curr <= 3370) && ((scql_curr+scqr_curr)/2 < m_Cur_SCQ1_89)) {
372 m_zfield = -1.0*tesla;
373 m_filename = path;
374 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode2.dat");
375#ifndef BEAN
376 log << MSG::INFO << "Select field map22: " << m_filename << endreq;
377#else
378 cout << "Select field map33: " << m_filename << endl;
379#endif
380
381 if((former_m_filename == "First Run") || (former_m_filename != m_filename))
382 {
383 former_m_filename = m_filename;
384 StatusCode status = parseFile();
385#ifndef BEAN
386 if ( !status.isSuccess() ) {
387 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
388 }
389#else
390 if ( !status ) {
391 cout << "Magnetic field parse failled" << endl;
392 }
393#endif
394 }
395 if(former_m_filename == m_filename) {
396 //cout<<"No need to open a new map file!!!!"<<endl;
397 }
398 //StatusCode status = parseFile();
399
400
401/*#ifndef BEAN
402 if ( !status.isSuccess() ) {
403 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
404 }
405#else
406 if ( !status ) {
407 cout << "Magnetic field parse failled" << endl;
408 }
409#endif*/
410 m_Q_1 = m_Q;
411 m_P_1 = m_P;
412
413 m_Q.clear();
414 m_P.clear();
415 if(m_gridDistance == 5) {
416 m_Q.reserve(201250);
417 m_P.reserve(201250);
418 }
419 if(m_gridDistance == 10) {
420 m_Q.reserve(27082);
421 m_P.reserve(27082);
422 }
423
424 m_filename = path;
425 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode3.dat");
426#ifndef BEAN
427 log << MSG::INFO << "Select field map44: " << m_filename << endreq;
428#else
429 cout << "Select field map55: " << m_filename << endl;
430#endif
431
432if((former_m_filename == "First Run") || (former_m_filename != m_filename))
433 {
434 former_m_filename = m_filename;
435 StatusCode status = parseFile();
436#ifndef BEAN
437 if ( !status.isSuccess() ) {
438 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
439 }
440#else
441 if ( !status ) {
442 cout << "Magnetic field parse failled" << endl;
443 }
444#endif
445 }
446 if(former_m_filename == m_filename) {
447 // cout<<"No need to open a new map file!!!!"<<endl;
448 }
449
450
451/* status = parseFile();
452
453#ifndef BEAN
454 if ( !status.isSuccess() ) {
455 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
456 }
457#else
458 if ( !status ) {
459 cout << "Magnetic field parse failled" << endl;
460 }
461#endif*/
462 m_Q_2 = m_Q;
463 m_P_2 = m_P;
464
465 m_Q.clear();
466 m_P.clear();
467 if(m_gridDistance == 5) {
468 m_Q.reserve(201250);
469 m_P.reserve(201250);
470 }
471 if(m_gridDistance == 10) {
472 m_Q.reserve(27082);
473 m_P.reserve(27082);
474 }
475 //check
476 if(m_Q_2.size() != m_Q_1.size()) {
477#ifndef BEAN
478 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
479#else
480 cout << "The two field maps used in this run are wrong!!!" << endl;
481#endif
482 exit(1);
483 }
484
485 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
486 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];
487 m_Q.push_back(fieldvalue);
488 m_P.push_back(m_P_2[iQ]);
489 }
490 }
491 //for the energy greater than 1.89GeV
492 if((ssm_curr >= 3367) && (ssm_curr <= 3370) && ((scql_curr+scqr_curr)/2 >= m_Cur_SCQ1_89)) {
493 m_zfield = -1.0*tesla;
494 m_filename = path;
495 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode3.dat");
496#ifndef BEAN
497 log << MSG::INFO << "Select field map66: " << m_filename << endreq;
498#else
499 cout << "Select field map77: " << m_filename << endl;
500#endif
501 if((former_m_filename == "First Run") || (former_m_filename != m_filename))
502 {
503 former_m_filename = m_filename;
504 StatusCode status = parseFile();
505#ifndef BEAN
506 if ( !status.isSuccess() ) {
507 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
508 }
509#else
510 if ( !status ) {
511 cout << "Magnetic field parse failled" << endl;
512 }
513#endif
514 }
515 if(former_m_filename == m_filename) {
516 //cout<<"No need to open a new map file!!!!"<<endl;
517 }
518
519/* StatusCode status = parseFile();
520
521#ifndef BEAN
522 if ( !status.isSuccess() ) {
523 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
524 }
525#else
526 if ( !status ) {
527 cout << "Magnetic field parse failled" << endl;
528 }
529#endif*/
530 m_Q_1 = m_Q;
531 m_P_1 = m_P;
532
533 m_Q.clear();
534 m_P.clear();
535 if(m_gridDistance == 5) {
536 m_Q.reserve(201250);
537 m_P.reserve(201250);
538 }
539 if(m_gridDistance == 10) {
540 m_Q.reserve(27082);
541 m_P.reserve(27082);
542 }
543
544 m_filename = path;
545 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode4.dat");
546#ifndef BEAN
547 log << MSG::INFO << "Select field map88: " << m_filename << endreq;
548#else
549 cout << "Select field map99: " << m_filename << endl;
550#endif
551
552 if((former_m_filename == "First Run") || (former_m_filename != m_filename))
553 {
554 former_m_filename = m_filename;
555 StatusCode status = parseFile();
556#ifndef BEAN
557 if ( !status.isSuccess() ) {
558 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
559 }
560#else
561 if ( !status ) {
562 cout << "Magnetic field parse failled" << endl;
563 }
564#endif
565 }
566 if(former_m_filename == m_filename) {
567 //cout<<"No need to open a new map file!!!!"<<endl;
568 }
569
570/* status = parseFile();
571
572#ifndef BEAN
573 if ( !status.isSuccess() ) {
574 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
575 }
576#else
577 if ( !status ) {
578 cout << "Magnetic field parse failled" << endl;
579 }
580#endif*/
581 m_Q_2 = m_Q;
582 m_P_2 = m_P;
583
584 m_Q.clear();
585 m_P.clear();
586 if(m_gridDistance == 5) {
587 m_Q.reserve(201250);
588 m_P.reserve(201250);
589 }
590 if(m_gridDistance == 10) {
591 m_Q.reserve(27082);
592 m_P.reserve(27082);
593 }
594 //check
595 if(m_Q_2.size() != m_Q_1.size()) {
596#ifndef BEAN
597 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
598#else
599 cout << "The two field maps used in this run are wrong!!!" << endl;
600#endif
601 exit(1);
602 }
603
604 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
605 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];
606 m_Q.push_back(fieldvalue);
607 m_P.push_back(m_P_2[iQ]);
608 }
609 }
610 if((ssm_curr >= 3033) && (ssm_curr <= 3035)) {
611 m_zfield = -0.9*tesla;
612 m_filename = path;
613 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode7.dat");
614#ifndef BEAN
615 log << MSG::INFO << "Select field map100: " << m_filename << endreq;
616#else
617 cout << "Select field map200: " << m_filename << endl;
618#endif
619
620 if((former_m_filename == "First Run") || (former_m_filename != m_filename))
621 {
622 former_m_filename = m_filename;
623 StatusCode status = parseFile();
624#ifndef BEAN
625 if ( !status.isSuccess() ) {
626 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
627 }
628#else
629 if ( !status ) {
630 cout << "Magnetic field parse failled" << endl;
631 }
632#endif
633 }
634 if(former_m_filename == m_filename) {
635 // cout<<"No need to open a new map file!!!!"<<endl;
636 }
637
638/* StatusCode status = parseFile();
639
640#ifndef BEAN
641 if ( !status.isSuccess() ) {
642 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
643 }
644#else
645 if ( !status ) {
646 cout << "Magnetic field parse failled" << endl;
647 }
648#endif*/
649 m_Q_1 = m_Q;
650 m_P_1 = m_P;
651
652 m_Q.clear();
653 m_P.clear();
654 if(m_gridDistance == 5) {
655 m_Q.reserve(201250);
656 m_P.reserve(201250);
657 }
658 if(m_gridDistance == 10) {
659 m_Q.reserve(27082);
660 m_P.reserve(27082);
661 }
662
663 m_filename = path;
664 m_filename += std::string( "/dat/2008-05-27BESIII-field-Mode8.dat");
665#ifndef BEAN
666 log << MSG::INFO << "Select field map300: " << m_filename << endreq;
667#else
668 cout << "Select field map400: " << m_filename << endl;
669#endif
670
671 if((former_m_filename == "First Run") || (former_m_filename != m_filename))
672 {
673 former_m_filename = m_filename;
674 StatusCode status = parseFile();
675#ifndef BEAN
676 if ( !status.isSuccess() ) {
677 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
678 }
679#else
680 if ( !status ) {
681 cout << "Magnetic field parse failled" << endl;
682 }
683#endif
684 }
685 if(former_m_filename == m_filename) {
686 // cout<<"No need to open a new map file!!!!"<<endl;
687 }
688
689/* status = parseFile();
690
691#ifndef BEAN
692 if ( !status.isSuccess() ) {
693 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
694 }
695#else
696 if ( !status ) {
697 cout << "Magnetic field parse failled" << endl;
698 }
699#endif*/
700 m_Q_2 = m_Q;
701 m_P_2 = m_P;
702
703 m_Q.clear();
704 m_P.clear();
705 if(m_gridDistance == 5) {
706 m_Q.reserve(201250);
707 m_P.reserve(201250);
708 }
709 if(m_gridDistance == 10) {
710 m_Q.reserve(27082);
711 m_P.reserve(27082);
712 }
713 //check
714 if(m_Q_2.size() != m_Q_1.size()) {
715#ifndef BEAN
716 log << MSG::FATAL << "The two field maps used in this run are wrong!!!" << endreq;
717#else
718 cout << "The two field maps used in this run are wrong!!!" << endl;
719#endif
720 exit(1);
721 }
722
723 for(unsigned int iQ = 0; iQ < m_Q_1.size(); iQ++) {
724 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];
725 m_Q.push_back(fieldvalue);
726 m_P.push_back(m_P_2[iQ]);
727 }
728 }
729
730 if((!((ssm_curr >= 3367) && (ssm_curr <= 3370)) && !((ssm_curr >= 3033) && (ssm_curr <= 3035))) || scql_curr == 0 || scqr_curr == 0) {
731#ifndef BEAN
732 log << MSG::FATAL << "The current of run " << runNo << " is abnormal in database, please check it! " << endreq;
733 log << MSG::FATAL << "The current of SSM is " << ssm_curr << ", SCQR is " << scqr_curr << ", SCQL is " << scql_curr << endreq;
734#else
735 cout << "The current of run " << runNo
736 t << " is abnormal in database, please check it! " << endl;
737 cout << "The current of SSM is " << ssm_curr
738 << ", SCQR is " << scqr_curr << ", SCQL is " << scql_curr << endl;
739#endif
740 exit(1);
741 }
742
743 if(m_Q.size() == 0) {
744#ifndef BEAN
745 log << MSG::FATAL << "This size of field map vector is ZERO, please check the current of run " << runNo << " in database!" << endreq;
746#else
747 cout << "This size of field map vector is ZERO,"
748 << " please check the current of run " << runNo
749 << " in database!" << endl;
750#endif
751 exit(1);
752 }
753
754 m_filename_TE = path;
755 if(m_gridDistance == 10) m_filename_TE += std::string( "/dat/TEMap10cm.dat");
756 if(m_gridDistance == 5) m_filename_TE += std::string( "/dat/TEMap5cm.dat");
757#ifndef BEAN
758 log << MSG::INFO << "Select field map11: " << m_filename_TE << endreq;
759#else
760 cout << "Select field map00: " << m_filename_TE << endl;
761#endif
762 if((former_m_filename_TE == "First Run") || (former_m_filename_TE != m_filename_TE))
763 {
764 former_m_filename_TE = m_filename_TE;
765 StatusCode status = parseFile_TE();
766 #ifndef BEAN
767 if ( !status.isSuccess() ) {
768 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
769 }
770#else
771 if ( !status ) {
772 cout << "Magnetic field parse failled" << endl;
773 }
774#endif
775 }
776 if(former_m_filename_TE == m_filename_TE) {
777 //cout<<"No need to open a new map file!!"<<endl;
778 }
779/* StatusCode status = parseFile_TE();
780
781#ifndef BEAN
782 if ( !status.isSuccess() ) {
783 log << MSG::DEBUG << "Magnetic field parse failled" << endreq;
784 }
785#else
786 if ( !status ) {
787 cout << "Magnetic field parse failled" << endl;
788 }
789#endif
790*/
791 for (int iC = 0; iC< 2; ++iC ){
792 m_min_FL[iC] = -90.0*cm;
793 m_max_FL[iC] = m_min_FL[iC]+( m_Nxyz[iC]-1 )* m_Dxyz[iC];
794 //for tof and emc
795 m_min_FL_TE[iC] = 0.0*cm;
796 m_max_FL_TE[iC] = m_min_FL_TE[iC]+( m_Nxyz_TE[iC]-1 )* m_Dxyz_TE[iC];
797 } // iC
798 m_min_FL[2] = -120.0*cm;
799 m_max_FL[2] = m_min_FL[2]+( m_Nxyz[2]-1 )* m_Dxyz[2];
800 //for tof and emc
801 m_min_FL_TE[2] = 0.0*cm;
802 m_max_FL_TE[2] = m_min_FL_TE[2]+( m_Nxyz_TE[2]-1 )* m_Dxyz_TE[2];
803}
TTree * t
Definition: binning.cxx:23
sprintf(cut,"kal_costheta0_em>-0.93&&kal_costheta0_em<0.93&&kal_pxy0_em>=0.05+%d*0.1&&kal_pxy0_em<0.15+%d*0.1&&NGch>=2", j, j)

◆ init_params() [4/4]

void MagneticFieldSvc::init_params ( std::vector< double >  current,
std::vector< double >  beamEnergy,
int  runNo 
)

◆ initialize() [1/2]

StatusCode MagneticFieldSvc::initialize ( )
virtual

Initialise the service (Inherited Service overrides)

Definition at line 139 of file MagneticFieldSvc.cxx.

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

◆ 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 823 of file MagneticFieldSvc.cxx.

825{
826 StatusCode sc = StatusCode::FAILURE;
827 if ( ppvInterface ) {
828 *ppvInterface = 0;
829
830 if ( riid == IID_IMagneticFieldSvc ) {
831 *ppvInterface = static_cast<IMagneticFieldSvc*>(this);
832 sc = StatusCode::SUCCESS;
833 addRef();
834 }
835 else
836 sc = Service::queryInterface( riid, ppvInterface );
837 }
838 return sc;
839}

◆ 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 2132 of file MagneticFieldSvc.cxx.

2134{
2135 if(m_runmode == 8 || m_runmode == 7) {
2136 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)
2137 {
2138 b[0]=0.0;
2139 b[1]=0.0;
2140 b[2]=-0.9*tesla;
2141 }
2142 else
2143 {
2144 b[0]=0.0;
2145 b[1]=0.0;
2146 b[2]=0.0;
2147 }
2148 }
2149 else {
2150 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)
2151 {
2152 b[0]=0.0;
2153 b[1]=0.0;
2154 b[2]=-1.0*tesla;
2155 }
2156 else
2157 {
2158 b[0]=0.0;
2159 b[1]=0.0;
2160 b[2]=0.0;
2161 }
2162 }
2163
2164 if(m_turnOffField == true) {
2165 b[0] = 0.;
2166 b[1] = 0.;
2167 b[2] = 0.;
2168 }
2169 //yzhang add 2012-04-25
2170 b[0] *= m_scale;
2171 b[1] *= m_scale;
2172 b[2] *= m_scale;
2173#ifndef BEAN
2174 return StatusCode::SUCCESS;
2175#else
2176 return true;
2177#endif
2178}

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 132 of file InstallArea/include/MagneticField/MagneticField/MagneticFieldSvc.h.

Member Data Documentation

◆ beamEnergy

std::vector< double > MagneticFieldSvc::beamEnergy

◆ current

std::vector< double > MagneticFieldSvc::current

◆ m_connect_run

FieldDBUtil::ConnectionDB * MagneticFieldSvc::m_connect_run

◆ m_eventSvc

IDataProviderSvc * MagneticFieldSvc::m_eventSvc

◆ m_mapBeamEnergy

std::map< int, std::vector< double > > MagneticFieldSvc::m_mapBeamEnergy

◆ m_mapMagnetInfo

std::map< int, std::vector< double > > MagneticFieldSvc::m_mapMagnetInfo

◆ m_readOneTime

bool MagneticFieldSvc::m_readOneTime

◆ m_runFrom

int MagneticFieldSvc::m_runFrom

◆ m_runTo

int MagneticFieldSvc::m_runTo

◆ runNo

int MagneticFieldSvc::runNo

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