BOSS 7.1.2
BESIII Offline Software System
Loading...
Searching...
No Matches
TofTrack Class Reference

#include <TofTrack.h>

Public Member Functions

 TofTrack (int run, int event)
 
 ~TofTrack ()
 
int trackId () const
 
int tofTrackId () const
 
int id1 () const
 
int id2 () const
 
int istrip1 () const
 
int istrip2 () const
 
int dstrip1 () const
 
int dstrip2 () const
 
unsigned int barrel () const
 
ExtTrackCase hitCase () const
 
double p () const
 
double path () const
 
double path1 () const
 
double path2 () const
 
double zrhit1 () const
 
double zrhit2 () const
 
double errzrhit1 () const
 
double errzrhit2 () const
 
double xhit1 () const
 
double yhit1 () const
 
double xhit2 () const
 
double yhit2 () const
 
int kal (unsigned int i) const
 
double zr1 (unsigned int i) const
 
double zr2 (unsigned int i) const
 
double errzr1 (unsigned int i) const
 
double errzr2 (unsigned int i) const
 
std::vector< TofData * > tofData1 () const
 
std::vector< TofData * > tofData2 () const
 
int size1 () const
 
int size2 () const
 
int size3 () const
 
int size4 () const
 
double theta1 () const
 
double theta2 () const
 
double phi1 () const
 
double phi2 () const
 
unsigned int quality1 () const
 
unsigned int quality2 () const
 
unsigned int quality () const
 
int tofId1 () const
 
int tofId2 () const
 
int strip1 () const
 
int strip2 () const
 
double ph11 () const
 
double ph12 () const
 
double ph21 () const
 
double ph22 () const
 
double ph1 () const
 
double ph2 () const
 
double ph () const
 
double tof11 (unsigned int i) const
 
double tof12 (unsigned int i) const
 
double tof21 (unsigned int i) const
 
double tof22 (unsigned int i) const
 
double tof1 (unsigned int i) const
 
double tof2 (unsigned int i) const
 
double tof (unsigned int i) const
 
double qch1 () const
 
double qch2 () const
 
double qch3 () const
 
double qch4 () const
 
double adc1 () const
 
double adc2 () const
 
double adc3 () const
 
double adc4 () const
 
double tdc1 () const
 
double tdc2 () const
 
double tdc3 () const
 
double tdc4 () const
 
double texpInner (unsigned int i) const
 
double texpOuter (unsigned int i) const
 
double texp (unsigned int i) const
 
double ztdc1 () const
 
double ztdc2 () const
 
double zadc1 () const
 
double zadc2 () const
 
double estime () const
 
double tdiff1 () const
 
double tdiff2 () const
 
int t0Stat () const
 
unsigned int flag () const
 
bool isNoHit () const
 
void setQuality (int qual)
 
void setQuality1 (int qual1)
 
void setQuality2 (int qual2)
 
void setFlag (unsigned int flag)
 
void setExtTrack (RecExtTrack *extTrack, double costheta, double p[5], int kal[5], double t0, int t0Stat)
 
void getMultiHit (TofTrack *&)
 
void setTofData (TofDataMap tofDataMap)
 
void tofDataAnalysis (TofData *tof, unsigned int iflag)
 
void getTofData (TofData *tof, unsigned int iflag)
 
void getTofDataEast (TofData *tof, unsigned int iflag)
 
void getTofDataWest (TofData *tof, unsigned int iflag)
 
void getTofDataNohit (unsigned int iflag)
 
void getEtfData (TofData *tof, unsigned int iflag, unsigned int qual)
 
void match (bool forCalibration, std::vector< int > deadId, std::vector< TofTrack * > *&tofTrackVec)
 
void findTofDataBarrel (std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, double zrhit, unsigned int iflag, std::vector< TofTrack * > *&tofTrackVec)
 
TofDatachooseTofData (std::vector< TofData * > tofDataVec, double zrhit)
 
TofDatacompareTofData (TofData *tofData1, TofData *tofData2, double zrhit)
 
void findTofDataEndcap (std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, double zr1[5])
 
TofDatachooseTofDataEndcap (std::vector< TofData * > tofDataVec, double zr1[5])
 
TofDatacompareTofDataEndcap (TofData *tofData1, TofData *tofData2)
 
void findEtfData (std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, std::vector< TofData * > tofDataVec3, double zrhit, unsigned int iflag)
 
TofDatachooseEtfData1 (std::vector< TofData * > tofDataVec, double zrhit)
 
TofDatachooseEtfData2 (std::vector< TofData * > tofDataVec, double zrhit)
 
void tofDataStudy ()
 
void setCalibration ()
 
void convert2RecTofTrackCol (RecTofTrackCol *recTofTrackCol)
 
void buildRecTofTrack (RecTofTrack *track, int layerorend)
 
void setRecTofTrack (RecTofTrack *track, int layerorend)
 
void convert2RecBTofCalHitColBarrel (int runNumber, int eventNumber, RecBTofCalHitCol *btofCalHitCol, std::string calibData)
 
void convert2RecETofCalHitCol (int runNumber, int eventNumber, RecETofCalHitCol *etofCalHitCol, std::string calibData)
 
void convert2RecBTofCalHitColETF (int runNumber, int eventNumber, RecBTofCalHitCol *btofCalHitCol, std::string calibData)
 
void qualityAnalysis ()
 

Detailed Description

Definition at line 22 of file TofTrack.h.

Constructor & Destructor Documentation

◆ TofTrack()

TofTrack::TofTrack ( int run,
int event )

Definition at line 10 of file TofTrack.cxx.

10 {
11 m_run = run;
12 m_event = event;
13 m_trackId = -1;
14 m_tofTrackId = -1;
15 m_id1 = -9;
16 m_id2 = -9;
17 m_istrip1 = -9;
18 m_istrip2 = -9;
19 m_hitCase = NoHit;
20 m_momentum = -99.0;
21 m_path = 0.0;
22 m_path1 = 0.0;
23 m_path2 = 0.0;
24 m_zrhit1 = 99.0;
25 m_errzr1 = 99.0;
26 m_zrhit2 = 99.0;
27 m_errzr2 = 99.0;
28 m_xhit1 = 99.0;
29 m_yhit1 = 99.0;
30 m_xhit2 = 99.0;
31 m_yhit2 = 99.0;
32 for( unsigned int i=0; i<5; i++ ) {
33 m_kal[i] = -1;
34 m_zr1[i] = 99.0;
35 m_zr2[i] = 99.0;
36 m_ezr1[i] = 99.0;
37 m_ezr2[i] = 99.0;
38 }
39 if( m_tofData1.size()>0 ) {
40 m_tofData1.clear();
41 }
42 if( m_tofData2.size()>0 ) {
43 m_tofData2.clear();
44 }
45 if( m_tofData3.size()>0 ) {
46 m_tofData3.clear();
47 }
48 if( m_tofData4.size()>0 ) {
49 m_tofData4.clear();
50 }
51 if( m_tofData5.size()>0 ) {
52 m_tofData5.clear();
53 }
54 if( m_tofData6.size()>0 ) {
55 m_tofData6.clear();
56 }
57 m_quality1 = 0;
58 m_quality2 = 0;
59 m_quality = 10;
60
61 m_delStrip1 = 20;
62 m_delStrip2 = 20;
63
64 m_tofId1 = -9;
65 m_tofId2 = -9;
66 m_strip1 = -9;
67 m_strip2 = -9;
68
69 m_ph11 = -99.0;
70 m_ph12 = -99.0;
71 m_ph21 = -99.0;
72 m_ph22 = -99.0;
73 m_ph1 = -99.0;
74 m_ph2 = -99.0;
75 m_ph = -99.0;
76
77 for( unsigned int i=0; i<5; i++ ) {
78 m_tof11[i] = 0.0;
79 m_tof12[i] = 0.0;
80 m_tof21[i] = 0.0;
81 m_tof22[i] = 0.0;
82 m_tof1[i] = 0.0;
83 m_tof2[i] = 0.0;
84 m_tof[i] = 0.0;
85 }
86
87 m_sigma11 = -99.0;
88 m_sigma12 = -99.0;
89 m_sigma21 = -99.0;
90 m_sigma22 = -99.0;
91 m_sigma1 = -99.0;
92 m_sigma2 = -99.0;
93 m_sigma = -99.0;
94
95 m_qch1 = -99.0;
96 m_qch2 = -99.0;
97 m_qch3 = -99.0;
98 m_qch4 = -99.0;
99 m_adc1 = -99.0;
100 m_adc2 = -99.0;
101 m_adc3 = -99.0;
102 m_adc4 = -99.0;
103 m_tdc1 = -99.0;
104 m_tdc2 = -99.0;
105 m_tdc3 = -99.0;
106 m_tdc4 = -99.0;
107
108 for( unsigned int i=0; i<5; i++ ) {
109 m_texpInner[i] = -99.0;
110 m_texpOuter[i] = -99.0;
111 m_texp[i] = -99.0;
112 }
113
114 m_ztdc1 = -99.0;
115 m_ztdc2 = -99.0;
116 m_zadc1 = -99.0;
117 m_zadc2 = -99.0;
118
119 m_estime = -99.0;
120 m_tdiff1 = -99.0;
121 m_tdiff2 = -99.0;
122
123 m_flag = 0;
124
125 return;
126}
@ NoHit
Definition TofTrack.h:20

◆ ~TofTrack()

TofTrack::~TofTrack ( )

Definition at line 129 of file TofTrack.cxx.

129 {
130 if( m_tofData1.size()>0 ) {
131 m_tofData1.clear();
132 }
133 if( m_tofData2.size()>0 ) {
134 m_tofData2.clear();
135 }
136 if( m_tofData3.size()>0 ) {
137 m_tofData3.clear();
138 }
139 if( m_tofData4.size()>0 ) {
140 m_tofData4.clear();
141 }
142 if( m_tofData5.size()>0 ) {
143 m_tofData5.clear();
144 }
145 if( m_tofData6.size()>0 ) {
146 m_tofData6.clear();
147 }
148 return;
149}

Member Function Documentation

◆ adc1()

double TofTrack::adc1 ( ) const
inline

◆ adc2()

double TofTrack::adc2 ( ) const
inline

◆ adc3()

double TofTrack::adc3 ( ) const
inline

◆ adc4()

double TofTrack::adc4 ( ) const
inline

◆ barrel()

unsigned int TofTrack::barrel ( ) const
inline

◆ buildRecTofTrack()

void TofTrack::buildRecTofTrack ( RecTofTrack * track,
int layerorend )

Definition at line 2406 of file TofTrack.cxx.

2406 {
2407
2408 track->setTofTrackID( m_tofTrackId );
2409 track->setTrackID( m_trackId );
2410
2411 track->setErrTof( 0.0 );
2412 track->setBeta( 0.0 );
2413
2414 double sigma[6];
2415 for( int i=0; i<6; i++ ) {
2416 sigma[i] = 0.0;
2417 }
2418 track->setSigma( sigma );
2419 track->setQuality( m_quality );
2420 track->setT0( m_estime );
2421 track->setErrT0( 0.0 );
2422 track->setPhi( 9999.0 );
2423 track->setErrPhi( 9999.0 );
2424 track->setEnergy( 9999.0 );
2425 track->setErrEnergy( 9999.0 );
2426
2427 if( ( layerorend == 11 ) || ( layerorend == 12 ) || ( layerorend == 1 ) ) {
2428 if( m_strip1<0 ) {
2429 track->setTofID( m_tofId1 ); // scintillator
2430 }
2431 else {
2432 track->setTofID( 12*m_tofId1+m_strip1 ); // MRPC
2433 }
2434 track->setPath( m_path1 );
2435 track->setZrHit( m_zrhit1 );
2436 track->setErrZ( m_errzr1 );
2437 track->setTexp( m_texpInner );
2438
2439 setRecTofTrack( track, layerorend );
2440 }
2441
2442 if( ( layerorend==21 ) || ( layerorend==22 ) || ( layerorend==2 ) ) {
2443 if( m_strip2<0 ) {
2444 track->setTofID( m_tofId2 ); // scintillator
2445 }
2446 else {
2447 track->setTofID( 12*m_tofId2+m_strip2 ); // MRPC
2448 }
2449 track->setPath( m_path2 );
2450 track->setZrHit( m_zrhit2 );
2451 track->setErrZ( m_errzr2 );
2452 track->setTexp( m_texpOuter );
2453
2454 setRecTofTrack( track, layerorend );
2455 }
2456
2457 if( layerorend==0 ) {
2458 track->setTofID( m_tofId1 );
2459 track->setPath( m_path2 );
2460 track->setZrHit( m_zrhit2 );
2461 track->setErrZ( m_errzr2 );
2462 track->setTexp( m_texp );
2463
2464 setRecTofTrack( track, layerorend );
2465 }
2466
2467 if( layerorend == 3 ) {
2468 if( m_strip1<0 ) {
2469 track->setTofID( m_id1 ); // scintillator
2470 }
2471 else {
2472 track->setTofID( 12*m_id1 + m_strip1 ); // mrpc
2473 }
2474 track->setPath( m_path1 );
2475 track->setZrHit( m_zrhit1 );
2476 track->setErrZ( m_errzr1 );
2477 track->setTexp( m_texpInner );
2478 }
2479
2480 return;
2481}
TTree * sigma
void setSigma(double sigma[6])
void setEnergy(double energy)
void setPath(double path)
Definition DstTofTrack.h:95
void setTofTrackID(int tofTrackID)
Definition DstTofTrack.h:90
void setQuality(int quality)
void setZrHit(double zrhit)
Definition DstTofTrack.h:96
void setPhi(double phi)
void setErrTof(double etof)
Definition DstTofTrack.h:99
void setErrT0(double errt0)
void setTexp(double texp[5])
void setBeta(double beta)
void setErrEnergy(double errenergy)
void setTrackID(int trackID)
Definition DstTofTrack.h:91
void setT0(double t0)
void setErrZ(double errz)
void setErrPhi(double errphi)
void setTofID(int tofID)
Definition DstTofTrack.h:92
void setRecTofTrack(RecTofTrack *track, int layerorend)

Referenced by convert2RecTofTrackCol().

◆ chooseEtfData1()

TofData * TofTrack::chooseEtfData1 ( std::vector< TofData * > tofDataVec,
double zrhit )

Definition at line 1480 of file TofTrack.cxx.

1480 {
1481 if( tofDataVec.size() == 0 ) {
1482 return 0;
1483 }
1484 std::vector<TofData*>::iterator igood = tofDataVec.begin();
1485 if( tofDataVec.size() == 1 ) {
1486 return (*igood);
1487 }
1488 else if( tofDataVec.size() > 1 ) {
1489 double deltaZ = 1000.0;
1490 std::vector<TofData*>::iterator iter = tofDataVec.begin();
1491 for( ; iter != tofDataVec.end(); iter++ ) {
1492 if( ( (*iter)->quality() & 0xf ) == 0xf ) {
1493 if( abs( (*iter)->ztdc() - zrhit ) < deltaZ ) {
1494 deltaZ = abs( (*iter)->ztdc() - zrhit );
1495 igood = iter;
1496 }
1497 }
1498 }
1499 // Max Q
1500 if( deltaZ > 999.0 ) {
1501 double maxQ = -1;
1502 iter = tofDataVec.begin();
1503 for( ; iter != tofDataVec.end(); iter++ ) {
1504 if( ( (*iter)->quality() & 0xc ) == 0xc ) {
1505 if( (*iter)->adc1() > maxQ ) {
1506 maxQ = (*iter)->adc1();
1507 igood = iter;
1508 }
1509 }
1510 else if( ( (*iter)->quality() & 0x3 ) == 0x3 ) {
1511 if( (*iter)->adc2() > maxQ ) {
1512 maxQ = (*iter)->adc2();
1513 igood = iter;
1514 }
1515 }
1516 }
1517 }
1518 }
1519
1520 return (*igood);
1521}
EvtStreamInputIterator< typename Generator::result_type > iter(Generator gen, int N=0)

Referenced by findEtfData().

◆ chooseEtfData2()

TofData * TofTrack::chooseEtfData2 ( std::vector< TofData * > tofDataVec,
double zrhit )

Definition at line 1526 of file TofTrack.cxx.

1526 {
1527 if( tofDataVec.size() == 0 ) {
1528 return 0;
1529 }
1530 std::vector<TofData*>::iterator igood = tofDataVec.begin();
1531 if( tofDataVec.size() == 1 ) {
1532 return (*igood);
1533 }
1534 else if( tofDataVec.size() > 1 ) {
1535 double maxQ = -1;
1536 std::vector<TofData*>::iterator iter = tofDataVec.begin();
1537 for( ; iter != tofDataVec.end(); iter++ ) {
1538 if( ( (*iter)->quality() & 0xc ) == 0xc ) {
1539 if( (*iter)->adc1() > maxQ ) {
1540 maxQ = (*iter)->adc1();
1541 igood = iter;
1542 }
1543 }
1544 else if( ( (*iter)->quality() & 0x3 ) == 0x3 ) {
1545 if( (*iter)->adc2() > maxQ ) {
1546 maxQ = (*iter)->adc2();
1547 igood = iter;
1548 }
1549 }
1550 }
1551 }
1552
1553 return (*igood);
1554}

Referenced by findEtfData().

◆ chooseTofData()

TofData * TofTrack::chooseTofData ( std::vector< TofData * > tofDataVec,
double zrhit )

Definition at line 1192 of file TofTrack.cxx.

1192 {
1193 if( tofDataVec.size() == 0 ) {
1194 cout << "TofRec::TofTrack::ChooseTofData: Size of TofData Vector is Zero!" << endl;
1195 return 0;
1196 }
1197 std::vector<TofData*>::iterator igood = tofDataVec.begin();
1198 if( tofDataVec.size() > 1 ) {
1199 double deltaZ = 1000.0;
1200 std::vector<TofData*>::iterator iter = tofDataVec.begin();
1201 // ZTDC compare
1202 for( ; iter != tofDataVec.end(); iter++ ) {
1203 if( ( (*iter)->quality() & 0x5 ) == 0x5 ) {
1204 if( abs( (*iter)->ztdc() - zrhit ) < deltaZ ) {
1205 deltaZ = abs( (*iter)->ztdc() - zrhit );
1206 igood = iter;
1207 }
1208 }
1209 }
1210 // ZADC compare
1211 if( deltaZ > 999.0 ) {
1212 iter = tofDataVec.begin();
1213 for( ; iter != tofDataVec.end(); iter++ ) {
1214 if( ( (*iter)->quality() & 0xa ) == 0xa ) {
1215 if( abs( (*iter)->zadc() - zrhit ) < deltaZ ) {
1216 deltaZ = abs( (*iter)->zadc() - zrhit );
1217 igood = iter;
1218 }
1219 }
1220 }
1221 }
1222 // Max Q
1223 if( deltaZ > 999.0 ) {
1224 unsigned int ibad = 0xf0;
1225 iter = tofDataVec.begin();
1226 for( ; iter != tofDataVec.end(); iter++ ) {
1227 if( ( (*iter)->quality() & 0xf0 ) < ibad ) {
1228 igood = iter;
1229 ibad = ( (*iter)->quality() & 0xf0 );
1230 }
1231 else if( ( (*iter)->quality() & 0xf0 ) == ibad ) {
1232 if( ( (*iter)->adc1() + (*iter)->adc2() ) > ( (*igood)->adc1() + (*igood)->adc2() ) ) {
1233 igood = iter;
1234 ibad = ( (*iter)->quality() & 0xf0 );
1235 }
1236 }
1237 }
1238 }
1239 }
1240
1241 return (*igood);
1242}

Referenced by findTofDataBarrel().

◆ chooseTofDataEndcap()

TofData * TofTrack::chooseTofDataEndcap ( std::vector< TofData * > tofDataVec,
double zr1[5] )

Definition at line 1352 of file TofTrack.cxx.

1352 {
1353 if( tofDataVec.size() == 0 ) {
1354 cout << "TofRec::TofTrack::ChooseTofData: Size of TofData Vector is Zero!" << endl;
1355 return 0;
1356 }
1357 std::vector<TofData*>::iterator igood = tofDataVec.begin();
1358 if( tofDataVec.size() > 1 ) {
1359 bool multihit = false;
1360 std::vector<TofData*>::iterator iter = tofDataVec.begin();
1361 for( ; iter != tofDataVec.end(); iter++ ) {
1362 if( (*iter)->qtimes1()>1 ) { multihit = true; }
1363 }
1364 iter = tofDataVec.begin();
1365 if( multihit ) {
1366 double tcorr = -999.0;
1367 double deltaTMin = 999.0;
1368 for( ; iter != tofDataVec.end(); iter++ ) {
1369 tcorr = tofCaliSvc->ETime( (*iter)->adc(), (*iter)->tdc()-m_estime, zr1[2], (*iter)->tofId() );
1370 for( unsigned int i=0; i<5; i++ ) {
1371 if( abs(tcorr-m_texpInner[i]) < deltaTMin ) {
1372 deltaTMin = abs(tcorr-m_texpInner[i]);
1373 igood = iter;
1374 }
1375 }
1376 }
1377 }
1378 else {
1379 double maxQ = 0.0;
1380 for( ; iter != tofDataVec.end(); iter++ ) {
1381 if( (*iter)->adc() > maxQ ) {
1382 maxQ = (*iter)->adc();
1383 igood = iter;
1384 }
1385 }
1386 }
1387 }
1388 return (*igood);
1389}
ITofCaliSvc * tofCaliSvc
virtual const double ETime(double ADC, double TDC, double rHit, unsigned id)=0
double zr1(unsigned int i) const
Definition TofTrack.h:50

Referenced by findTofDataEndcap().

◆ compareTofData()

TofData * TofTrack::compareTofData ( TofData * tofData1,
TofData * tofData2,
double zrhit )

Definition at line 1249 of file TofTrack.cxx.

1249 {
1250 TofData* tof = tofData1;
1251 // ZTDC compare
1252 if( abs(tofData1->ztdc() - zrhit ) > abs(tofData2->ztdc() - zrhit ) ) {
1253 // SingleEnd/NoT/NoQ compare
1254 if( ( tofData1->quality() & 0xf0 ) >= ( tofData1->quality() & 0xf0 ) ) {
1255 // QDC compare
1256 // if( ( tofData1->adc1() + tofData1->adc2() ) < ( tofData2->adc1() + tofData2->adc2() ) ) {
1257 tof = tofData2;
1258 // }
1259 }
1260 // }
1261 }
1262
1263 return tof;
1264}
std::vector< TofData * > tofData1() const
Definition TofTrack.h:54
double tof(unsigned int i) const
Definition TofTrack.h:87
std::vector< TofData * > tofData2() const
Definition TofTrack.h:55

Referenced by findTofDataBarrel().

◆ compareTofDataEndcap()

TofData * TofTrack::compareTofDataEndcap ( TofData * tofData1,
TofData * tofData2 )

Definition at line 1395 of file TofTrack.cxx.

1395 {
1396 TofData* tof = tofData1;
1397 if( tof->adc() < tofData2->adc() ) {
1398 tof = tofData2;
1399 }
1400 return tof;
1401}

Referenced by findTofDataEndcap().

◆ convert2RecBTofCalHitColBarrel()

void TofTrack::convert2RecBTofCalHitColBarrel ( int runNumber,
int eventNumber,
RecBTofCalHitCol * btofCalHitCol,
std::string calibData )

Definition at line 2568 of file TofTrack.cxx.

2568 {
2569
2570 if( ( m_quality1 & 0x800 ) == 0x800 ) {
2571
2572 RecBTofCalHit* ahit = new RecBTofCalHit;
2573 ahit->setRun( runNumber );
2574 ahit->setEvent( eventNumber );
2575 ahit->setMod( m_tofId1 );
2576 ahit->setQual( m_hitCase );
2577 ahit->setdZHit( 1 );
2578
2579 for( int i=0; i<5; i++ ) {
2580 ahit->setTpred( i, m_texpInner[i] );
2581 }
2582 if( calibData == "Dimu" ) {
2583 ahit->setTpred( m_texpInner[1] );
2584 ahit->setZHit( m_zr1[1] );
2585 // ahit->setdZHit( m_ezr1[1] );
2586 }
2587 else {
2588 ahit->setTpred( m_texpInner[0] );
2589 ahit->setZHit( m_zr1[0] );
2590 // ahit->setdZHit( m_ezr1[0] );
2591 }
2592
2593 ahit->setTdc1( m_tdc1-m_estime );
2594 ahit->setTdc2( m_tdc2-m_estime );
2595 ahit->setAdc1( m_adc1 );
2596 ahit->setAdc2( m_adc2 );
2597 // ahit->setZHit( m_zrhit1 );
2598 // ahit->setdZHit( m_errzr1 );
2599 ahit->setDeltaPhi( m_estime );
2600 ahit->setsinTheta( m_theta1 );
2601 ahit->setP( m_momentum );
2602 ahit->setQ( m_ph1 );
2603 ahit->setPath( m_path1 );
2604
2605 btofCalHitCol->push_back( ahit );
2606
2607 if( ( m_quality2 & 0x800 ) == 0x800 ) {
2608
2609 RecBTofCalHit* bhit = new RecBTofCalHit;
2610 bhit->setRun( runNumber );
2611 bhit->setEvent( eventNumber );
2612 bhit->setMod( m_tofId2 );
2613 bhit->setQual( m_hitCase );
2614 bhit->setdZHit( 1 );
2615
2616 for( int i=0; i<5; i++ ) {
2617 bhit->setTpred( i, m_texpOuter[i] );
2618 }
2619 if( calibData == "Dimu" ) {
2620 bhit->setTpred( m_texpOuter[1] );
2621 bhit->setZHit( m_zr2[1] );
2622 // bhit->setdZHit( m_ezr2[1] );
2623 }
2624 else {
2625 bhit->setTpred( m_texpOuter[0] );
2626 bhit->setZHit( m_zr2[0] );
2627 // bhit->setdZHit( m_ezr2[0] );
2628 }
2629
2630 bhit->setTdc1( m_tdc3-m_estime );
2631 bhit->setTdc2( m_tdc4-m_estime );
2632 bhit->setAdc1( m_adc3 );
2633 bhit->setAdc2( m_adc4 );
2634 // bhit->setZHit( m_zrhit2 );
2635 // bhit->setdZHit( m_errzr2 );
2636 bhit->setDeltaPhi( m_estime );
2637 bhit->setsinTheta( m_theta2 );
2638 bhit->setP( m_momentum );
2639 bhit->setQ( m_ph2 );
2640 bhit->setPath( m_path2 );
2641
2642 ahit->setnext(bhit);
2643
2644 btofCalHitCol->push_back( bhit );
2645
2646 }
2647 }
2648
2649 return;
2650}
void setsinTheta(double sint)
void setnext(RecBTofCalHit *n)
void setTdc2(double tdc2)
void setRun(int run)
void setP(double p)
void setQ(double q)
void setdZHit(double dzHit)
void setPath(double l)
void setAdc2(double adc2)
void setTdc1(double tdc1)
void setMod(int imod)
void setZHit(double zHit)
void setQual(int qual)
void setDeltaPhi(double deltaPhi)
void setEvent(int event)
void setTpred(int idx, double t)
void setAdc1(double adc1)

◆ convert2RecBTofCalHitColETF()

void TofTrack::convert2RecBTofCalHitColETF ( int runNumber,
int eventNumber,
RecBTofCalHitCol * btofCalHitCol,
std::string calibData )

Definition at line 2693 of file TofTrack.cxx.

2693 {
2694
2695 if( ( m_quality1 & 0x800 ) == 0x800 ) {
2696
2697 RecBTofCalHit* ahit = new RecBTofCalHit;
2698 ahit->setRun( runNumber );
2699 ahit->setEvent( eventNumber );
2700 ahit->setMod( m_tofId1 );
2701 ahit->setQual( m_hitCase );
2702 ahit->setdZHit( 0 );
2703
2704 for( int i=0; i<5; i++ ) {
2705 ahit->setTpred( i, m_texpInner[i] );
2706 }
2707 if( calibData == "Dimu" ) {
2708 ahit->setTpred( m_texpInner[1] );
2709 ahit->setZHit( m_zr1[1] );
2710 }
2711 else {
2712 ahit->setTpred( m_texpInner[0] );
2713 ahit->setZHit( m_zr1[0] );
2714 }
2715
2716 ahit->setTdc1( m_tdc1-m_estime );
2717 ahit->setTdc2( m_tdc2-m_estime );
2718 ahit->setAdc1( m_adc1 );
2719 ahit->setAdc2( m_adc2 );
2720 ahit->setDeltaPhi( m_estime );
2721 ahit->setsinTheta( m_strip1 );
2722 ahit->setP( m_momentum );
2723 ahit->setQ( m_ph1 );
2724 ahit->setPath( m_path1 );
2725
2726 btofCalHitCol->push_back( ahit );
2727 }
2728
2729 if( ( m_quality2 & 0x800 ) == 0x800 ) {
2730
2731 RecBTofCalHit* bhit = new RecBTofCalHit;
2732 bhit->setRun( runNumber );
2733 bhit->setEvent( eventNumber );
2734 bhit->setMod( m_tofId2 );
2735 bhit->setQual( m_hitCase );
2736 bhit->setdZHit( 0 );
2737
2738 for( int i=0; i<5; i++ ) {
2739 bhit->setTpred( i, m_texpOuter[i] );
2740 }
2741 if( calibData == "Dimu" ) {
2742 bhit->setTpred( m_texpOuter[1] );
2743 bhit->setZHit( m_zr2[1] );
2744 }
2745 else {
2746 bhit->setTpred( m_texpOuter[0] );
2747 bhit->setZHit( m_zr2[0] );
2748 }
2749
2750 bhit->setTdc1( m_tdc3-m_estime );
2751 bhit->setTdc2( m_tdc4-m_estime );
2752 bhit->setAdc1( m_adc3 );
2753 bhit->setAdc2( m_adc4 );
2754 bhit->setDeltaPhi( m_estime );
2755 bhit->setsinTheta( m_strip2 );
2756 bhit->setP( m_momentum );
2757 bhit->setQ( m_ph2 );
2758 bhit->setPath( m_path2 );
2759
2760 btofCalHitCol->push_back( bhit );
2761 }
2762
2763 return;
2764}

◆ convert2RecETofCalHitCol()

void TofTrack::convert2RecETofCalHitCol ( int runNumber,
int eventNumber,
RecETofCalHitCol * etofCalHitCol,
std::string calibData )

Definition at line 2653 of file TofTrack.cxx.

2653 {
2654
2655 if( ( m_quality1 & 0x800 ) != 0x800 ) return;
2656
2657 RecETofCalHit* chit = new RecETofCalHit;
2658 chit->setRun( runNumber );
2659 chit->setEvent( eventNumber );
2660 chit->setMod( m_tofId1 );
2661 chit->setQual( m_hitCase );
2662
2663 for( int i=0; i<5; i++ ) {
2664 chit->setTpred( i, m_texpInner[i] );
2665 }
2666 if( calibData == "Dimu" ) {
2667 chit->setTpred( m_texpInner[1] );
2668 chit->setRHit( m_zr1[1] );
2669 chit->setdRHit( m_ezr1[1] );
2670 }
2671 else {
2672 chit->setTpred( m_texpInner[0] );
2673 chit->setRHit( m_zr1[0] );
2674 chit->setdRHit( m_ezr1[0] );
2675 }
2676
2677 chit->setTdc( m_tdc1-m_estime );
2678 chit->setAdc( m_adc1 );
2679 // chit->setRHit( m_zrhit1 );
2680 // chit->setdRHit( m_errzr1 );
2681 chit->setDeltaPhi( m_estime );
2682 chit->setcosTheta( m_theta1 );
2683 chit->setQ( m_ph1 );
2684 chit->setP( m_momentum );
2685 chit->setPath( m_path1 );
2686
2687 etofCalHitCol->push_back( chit );
2688
2689 return;
2690}
void setEvent(int event)
void setdRHit(double drHit)
void setDeltaPhi(double deltaPhi)
void setTdc(double tdc)
void setRHit(double rHit)
void setRun(int run)
void setPath(double l)
void setMod(int imod)
void setP(double p)
void setAdc(double adc)
void setcosTheta(double cost)
void setQ(double q)
void setQual(int qual)
void setTpred(int idx, double t)

◆ convert2RecTofTrackCol()

void TofTrack::convert2RecTofTrackCol ( RecTofTrackCol * recTofTrackCol)

Definition at line 2091 of file TofTrack.cxx.

2091 {
2092
2093 bool barrel = ( ( m_hitCase == InnerLayer ) || ( m_hitCase == OuterLayer ) || ( m_hitCase == DoubleLayer ) );
2094
2095 bool innerEast = ( ( m_quality1 & 0xc ) == 0xc );
2096 bool innerWest = ( ( m_quality1 & 0x3 ) == 0x3 );
2097 bool outerEast = ( ( m_quality2 & 0xc ) == 0xc );
2098 bool outerWest = ( ( m_quality2 & 0x3 ) == 0x3 );
2099
2100 if( barrel ) {
2101
2102 if( innerEast ) {
2103 RecTofTrack* atrack11 = new RecTofTrack;
2104 buildRecTofTrack( atrack11, 11 ); // innerlayer east readout
2105 TofHitStatus* hitStatus11 = new TofHitStatus;
2106 if( innerWest ) {
2107 hitStatus11->setBarrelReadout( 1, true ); // innerlayer east readout
2108 }
2109 else {
2110 if( m_hitCase == InnerLayer ) {
2111 hitStatus11->setBarrelCluster( 11 ); // innerlayer east cluster
2112 }
2113 else if( m_hitCase == DoubleLayer ) {
2114 if( outerEast && outerWest ) {
2115 hitStatus11->setBarrelCounter( 11 ); // innerlayer east counter
2116 }
2117 else {
2118 hitStatus11->setBarrelCluster( 11 ); // innerlayer east cluster
2119 }
2120 }
2121 else {
2122 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 11- Impossible!" << endl;
2123 }
2124 }
2125 atrack11->setStatus( hitStatus11->value() );
2126 delete hitStatus11;
2127 recTofTrackCol->push_back( atrack11 );
2128 }
2129
2130 if( innerWest ) {
2131 RecTofTrack* atrack12 = new RecTofTrack;
2132 buildRecTofTrack( atrack12, 12 ); // innerlayer west readout
2133 TofHitStatus* hitStatus12 = new TofHitStatus;
2134 if( innerEast ) {
2135 hitStatus12->setBarrelReadout( 1, false ); // innerlayer west
2136 }
2137 else {
2138 if( m_hitCase == InnerLayer ) {
2139 hitStatus12->setBarrelCluster( 12 ); // innerlayer west cluster
2140 }
2141 else if( m_hitCase == DoubleLayer ) {
2142 if( outerEast && outerWest ) {
2143 hitStatus12->setBarrelCounter( 12 ); // innerlayer west counter
2144 }
2145 else {
2146 hitStatus12->setBarrelCluster( 12 ); // innerlayer west cluster
2147 }
2148 }
2149 else {
2150 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 12- Impossible!" << endl;
2151 }
2152 }
2153 atrack12->setStatus( hitStatus12->value() );
2154 delete hitStatus12;
2155 recTofTrackCol->push_back( atrack12 );
2156 }
2157
2158 if( innerEast && innerWest ) {
2159 RecTofTrack* atrack1 = new RecTofTrack;
2160 buildRecTofTrack( atrack1, 1 ); // innerlayer counter
2161 TofHitStatus* hitStatus1 = new TofHitStatus;
2162 if( m_hitCase == InnerLayer ) {
2163 hitStatus1->setBarrelCluster( 1 ); // innerlayer cluster and counter
2164 }
2165 else if( m_hitCase == DoubleLayer ) {
2166 if( outerEast && outerWest ) {
2167 hitStatus1->setBarrelCounter( 1 ); // innerlayer counter
2168 }
2169 else {
2170 hitStatus1->setBarrelCluster( 1 ); // innerlayer cluster and counter
2171 }
2172 }
2173 else {
2174 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 1- Impossible!" << endl;
2175 }
2176 atrack1->setStatus( hitStatus1->value() );
2177 delete hitStatus1;
2178 recTofTrackCol->push_back( atrack1 );
2179 }
2180
2181 if( outerEast ) {
2182 RecTofTrack* atrack21 = new RecTofTrack;
2183 buildRecTofTrack( atrack21, 21 ); // outerlayer east readout
2184 TofHitStatus* hitStatus21 = new TofHitStatus;
2185 if( outerWest ) {
2186 hitStatus21->setBarrelReadout( 2, true ); // outerlayer east readout
2187 }
2188 else {
2189 if( m_hitCase == OuterLayer ) {
2190 hitStatus21->setBarrelCluster( 21 ); // outerlayer east cluster
2191 }
2192 else if( m_hitCase == DoubleLayer ) {
2193 if( innerEast || innerWest ) {
2194 hitStatus21->setBarrelCounter( 21 ); // outerlayer east counter
2195 }
2196 // else {
2197 // hitStatus21->setBarrelCluster( 21 ); // outerlayer east cluster
2198 // }
2199 }
2200 else {
2201 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 21- Impossible!" << endl;
2202 }
2203 }
2204 atrack21->setStatus( hitStatus21->value() );
2205 delete hitStatus21;
2206 recTofTrackCol->push_back( atrack21 );
2207 }
2208
2209 if( outerWest ) {
2210 RecTofTrack* atrack22 = new RecTofTrack;
2211 buildRecTofTrack( atrack22, 22 ); // outerlayer west readout
2212 TofHitStatus* hitStatus22 = new TofHitStatus;
2213 if( outerEast ) {
2214 hitStatus22->setBarrelReadout( 2, false ); // outerlayer west readout
2215 }
2216 else {
2217 if( m_hitCase == OuterLayer ) {
2218 hitStatus22->setBarrelCluster( 22 ); // outerlayer west cluster
2219 }
2220 else if( m_hitCase == DoubleLayer ) {
2221 if( innerEast || innerWest ) {
2222 hitStatus22->setBarrelCounter( 22 ); // outerlayer west counter
2223 }
2224 // else {
2225 // hitStatus22->setBarrelCluster( 22 ); // outerlayer west cluster
2226 // }
2227 }
2228 else {
2229 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 22- Impossible!" << endl;
2230 }
2231 }
2232 atrack22->setStatus( hitStatus22->value() );
2233 delete hitStatus22;
2234 recTofTrackCol->push_back( atrack22 );
2235 }
2236
2237 if( outerEast && outerWest ) {
2238 RecTofTrack* atrack2 = new RecTofTrack;
2239 buildRecTofTrack( atrack2, 2 ); // outerlayer counter
2240 TofHitStatus* hitStatus2 = new TofHitStatus;
2241 if( m_hitCase == OuterLayer ) {
2242 hitStatus2->setBarrelCluster( 2 ); // outerlayer cluster and counter
2243 }
2244 else if( m_hitCase == DoubleLayer ) {
2245 if( innerEast && innerWest ) {
2246 hitStatus2->setBarrelCounter( 2 ); // outerlayer counter
2247 }
2248 else {
2249 hitStatus2->setBarrelCluster( 2 ); // outerlayer cluster and counter
2250 }
2251 }
2252 else {
2253 cout << "TofRec::TofTrack:convert2RecTofTrackCol: 2- Impossible!" << endl;
2254 }
2255 atrack2->setStatus( hitStatus2->value() );
2256 delete hitStatus2;
2257 recTofTrackCol->push_back( atrack2 );
2258 }
2259
2260 if( innerEast && innerWest && outerEast && outerWest ) {
2261 RecTofTrack* atrack = new RecTofTrack;
2262 buildRecTofTrack( atrack, 0 ); // doublelayer cluster
2263 TofHitStatus* hitStatus = new TofHitStatus;
2264 hitStatus->setBarrelCluster( 3 ); // doublelayer cluster
2265 atrack->setStatus( hitStatus->value() );
2266 delete hitStatus;
2267 recTofTrackCol->push_back( atrack );
2268 }
2269
2270 }
2271
2272 if( ( m_hitCase == EastEndcap ) || ( m_hitCase == WestEndcap ) ) {
2273 RecTofTrack* atrack = new RecTofTrack;
2274 buildRecTofTrack( atrack, 11 ); // eastendcap counter
2275 TofHitStatus* hitStatus = new TofHitStatus;
2276 if( m_hitCase == EastEndcap ) {
2277 hitStatus->setEndcapCluster( true ); // east endcap cluster counter readout
2278 }
2279 else if( m_hitCase == WestEndcap ) {
2280 hitStatus->setEndcapCluster( false ); // west endcap cluster counter readout
2281 }
2282 else {
2283 cout << "TofRec::TofTrack:convert2RecTofTrackCol: endcap- Impossible!" << endl;
2284 }
2285 atrack->setStatus( hitStatus->value() );
2286 delete hitStatus;
2287 recTofTrackCol->push_back( atrack );
2288 }
2289
2290 if( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) {
2291
2292 if( innerEast || innerWest ) {
2293
2294 if( innerEast ) {
2295 RecTofTrack* atrack1 = new RecTofTrack;
2296 buildRecTofTrack( atrack1, 11 ); // mrpc east readout
2297 TofHitStatus* hitStatus1 = new TofHitStatus;
2298 if( innerWest ) {
2299 hitStatus1->setMRPCReadout( true ); // mrpc east readout
2300 }
2301 else {
2302 hitStatus1->setMRPCCluster( false, true ); // mrpc east cluster
2303 }
2304 atrack1->setStatus( hitStatus1->value() );
2305 delete hitStatus1;
2306 recTofTrackCol->push_back( atrack1 );
2307 }
2308
2309 if( innerWest ) {
2310 RecTofTrack* atrack2 = new RecTofTrack;
2311 buildRecTofTrack( atrack2, 12 ); // mrpc west readout
2312 TofHitStatus* hitStatus2 = new TofHitStatus;
2313 if( innerEast ) {
2314 hitStatus2->setMRPCReadout( false ); // mrpc west readout
2315 }
2316 else {
2317 hitStatus2->setMRPCCluster( false, false ); // mrpc east cluster
2318 }
2319 atrack2->setStatus( hitStatus2->value() );
2320 delete hitStatus2;
2321 recTofTrackCol->push_back( atrack2 );
2322 }
2323
2324 if( innerEast && innerWest ) {
2325 RecTofTrack* atrack = new RecTofTrack;
2326 buildRecTofTrack( atrack, 1 ); // mrpc east readout
2327 TofHitStatus* hitStatus = new TofHitStatus;
2328 if( m_hitCase == EastEndcapMRPC ) {
2329 hitStatus->setMRPCCluster( true, true ); // mrpc east cluster
2330 }
2331 else {
2332 hitStatus->setMRPCCluster( true, false );// mrpc west cluster
2333 }
2334 atrack->setStatus( hitStatus->value() );
2335 delete hitStatus;
2336 recTofTrackCol->push_back( atrack );
2337 }
2338
2339 }
2340 else {
2341
2342 if( outerEast ) {
2343 RecTofTrack* atrack1 = new RecTofTrack;
2344 buildRecTofTrack( atrack1, 21 ); // mrpc east readout
2345 TofHitStatus* hitStatus1 = new TofHitStatus;
2346 if( outerWest ) {
2347 hitStatus1->setMRPCReadout( true ); // mrpc east readout
2348 }
2349 else {
2350 hitStatus1->setMRPCCluster( false, true ); // mrpc east cluster
2351 }
2352 atrack1->setStatus( hitStatus1->value() );
2353 delete hitStatus1;
2354 recTofTrackCol->push_back( atrack1 );
2355 }
2356
2357 if( outerWest ) {
2358 RecTofTrack* atrack2 = new RecTofTrack;
2359 buildRecTofTrack( atrack2, 22 ); // mrpc west readout
2360 TofHitStatus* hitStatus2 = new TofHitStatus;
2361 if( outerEast ) {
2362 hitStatus2->setMRPCReadout( false ); // mrpc west readout
2363 }
2364 else {
2365 hitStatus2->setMRPCCluster( false, false ); // mrpc east cluster
2366 }
2367 atrack2->setStatus( hitStatus2->value() );
2368 delete hitStatus2;
2369 recTofTrackCol->push_back( atrack2 );
2370 }
2371
2372 if( outerEast && outerWest ) {
2373 RecTofTrack* atrack = new RecTofTrack;
2374 buildRecTofTrack( atrack, 2 ); // mrpc east readout
2375 TofHitStatus* hitStatus = new TofHitStatus;
2376 if( m_hitCase == EastEndcapMRPC ) {
2377 hitStatus->setMRPCCluster( true, true ); // mrpc east cluster
2378 }
2379 else {
2380 hitStatus->setMRPCCluster( true, false );// mrpc west cluster
2381 }
2382 atrack->setStatus( hitStatus->value() );
2383 delete hitStatus;
2384 recTofTrackCol->push_back( atrack );
2385 }
2386
2387 }
2388
2389 }
2390
2391 if( m_hitCase == NoHit ) {
2392 RecTofTrack* atrack = new RecTofTrack;
2393 buildRecTofTrack( atrack, 3 ); // no hit
2394 TofHitStatus* hitStatus = new TofHitStatus;
2395 hitStatus->setNoHit(); // no hit
2396 atrack->setStatus( hitStatus->value() );
2397 delete hitStatus;
2398 recTofTrackCol->push_back( atrack );
2399 }
2400
2401 return;
2402}
@ WestEndcap
Definition TofTrack.h:20
@ OuterLayer
Definition TofTrack.h:20
@ EastEndcapMRPC
Definition TofTrack.h:20
@ InnerLayer
Definition TofTrack.h:20
@ EastEndcap
Definition TofTrack.h:20
@ DoubleLayer
Definition TofTrack.h:20
@ WestEndcapMRPC
Definition TofTrack.h:20
void setStatus(unsigned int status)
Definition DstTofTrack.h:93
void setMRPCReadout(bool east)
unsigned int value() const
void setMRPCCluster(bool cluster, bool east)
void setBarrelReadout(unsigned int layer, bool east)
void setEndcapCluster(bool east)
void setBarrelCounter(unsigned int layer)
void setBarrelCluster(unsigned int layer)
void buildRecTofTrack(RecTofTrack *track, int layerorend)
unsigned int barrel() const
Definition TofTrack.h:35

◆ dstrip1()

int TofTrack::dstrip1 ( ) const
inline

Definition at line 33 of file TofTrack.h.

33{ return m_delStrip1; }

Referenced by TofCheckDigi::Fill_TofTrack().

◆ dstrip2()

int TofTrack::dstrip2 ( ) const
inline

Definition at line 34 of file TofTrack.h.

34{ return m_delStrip2; }

Referenced by TofCheckDigi::Fill_TofTrack().

◆ errzr1()

double TofTrack::errzr1 ( unsigned int i) const
inline

Definition at line 52 of file TofTrack.h.

52{ return m_ezr1[i]; }

◆ errzr2()

double TofTrack::errzr2 ( unsigned int i) const
inline

Definition at line 53 of file TofTrack.h.

53{ return m_ezr2[i]; }

◆ errzrhit1()

double TofTrack::errzrhit1 ( ) const
inline

◆ errzrhit2()

double TofTrack::errzrhit2 ( ) const
inline

◆ estime()

double TofTrack::estime ( ) const
inline

Definition at line 111 of file TofTrack.h.

111{ return m_estime; }

◆ findEtfData()

void TofTrack::findEtfData ( std::vector< TofData * > tofDataVec1,
std::vector< TofData * > tofDataVec2,
std::vector< TofData * > tofDataVec3,
double zrhit,
unsigned int iflag )

Definition at line 1408 of file TofTrack.cxx.

1408 {
1409
1410 TofData *tof1 = 0;
1411 TofData *tof2 = 0;
1412 TofData *tof3 = 0;
1413
1414 bool findSignal = false;
1415
1416 if( tofDataVec1.size()==0 && tofDataVec2.size()==0 && tofDataVec3.size()==0 ) {
1417 if( iflag == 1 ) {
1418 m_quality1 = ( m_quality1 | 0x300 );
1419 }
1420 else if( iflag == 2 ) {
1421 m_quality2 = ( m_quality2 | 0x300 );
1422 }
1423 }
1424 else {
1425 if( tofDataVec1.size()>0 ) {
1426 tof1 = chooseEtfData1( tofDataVec1, zrhit );
1427 if( ( tof1->quality() & 0xf ) == 0xf ) {
1428 getEtfData( tof1, iflag, 1 );
1429 findSignal = true;
1430 }
1431 }
1432 if( !findSignal && tofDataVec2.size() > 0 ) {
1433 tof2 = chooseEtfData1( tofDataVec2, zrhit );
1434 if( ( tof2->quality() & 0xf ) == 0xf ) {
1435 getEtfData( tof2, iflag, 2 );
1436 findSignal = true;
1437 }
1438 }
1439 if( !findSignal && tofDataVec3.size() > 0 ) {
1440 tof3 = chooseEtfData2( tofDataVec3, zrhit );
1441 if( ( tof3->quality() & 0xf ) == 0xf ) {
1442 getEtfData( tof3, iflag, 3 );
1443 findSignal = true;
1444 }
1445 }
1446 if( !findSignal && tofDataVec1.size()>0 ) {
1447 if( ( ( tof1->quality() & 0xf ) == 0xc ) || ( ( tof1->quality() & 0xf ) == 0x3 ) ) {
1448 getEtfData( tof1, iflag, 4 );
1449 findSignal = true;
1450 }
1451 }
1452 if( !findSignal && tofDataVec2.size() > 0 ) {
1453 if( ( ( tof2->quality() & 0xf ) == 0xc ) || ( ( tof2->quality() & 0xf ) == 0x3 ) ) {
1454 getEtfData( tof2, iflag, 5 );
1455 findSignal = true;
1456 }
1457 }
1458 if( !findSignal && tofDataVec3.size() > 0 ) {
1459 if( ( ( tof3->quality() & 0xf ) == 0xc ) || ( ( tof3->quality() & 0xf ) == 0x3 ) ) {
1460 getEtfData( tof3, iflag, 6 );
1461 findSignal = true;
1462 }
1463 }
1464 if( findSignal ) {
1465 if( iflag == 1 ) {
1466 m_quality1 = ( m_quality1 | 0x300 );
1467 }
1468 else if( iflag == 2 ) {
1469 m_quality2 = ( m_quality2 | 0x300 );
1470 }
1471 }
1472 }
1473
1474 return;
1475}
unsigned int quality() const
Definition TofData.h:168
double tof2(unsigned int i) const
Definition TofTrack.h:86
void getEtfData(TofData *tof, unsigned int iflag, unsigned int qual)
TofData * chooseEtfData2(std::vector< TofData * > tofDataVec, double zrhit)
double tof1(unsigned int i) const
Definition TofTrack.h:85
TofData * chooseEtfData1(std::vector< TofData * > tofDataVec, double zrhit)

Referenced by match().

◆ findTofDataBarrel()

void TofTrack::findTofDataBarrel ( std::vector< TofData * > tofDataVec1,
std::vector< TofData * > tofDataVec2,
double zrhit,
unsigned int iflag,
std::vector< TofTrack * > *& tofTrackVec )

Definition at line 1020 of file TofTrack.cxx.

1020 {
1021
1022 unsigned int qual = 0xf;
1023 TofData* tof = 0;
1024 if( tofDataVec2.size() == 0 ) {
1025 if( tofDataVec1.size() == 0 ) {
1026 qual = 0;
1027 }
1028 else if( tofDataVec1.size() == 1 ) {
1029 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1030 tof = (*iter1);
1031 qual = 0x1;
1032 }
1033 else if( tofDataVec1.size() > 1 ) {
1034 tof= chooseTofData( tofDataVec1, zrhit );
1035 qual = 0x2;
1036 }
1037 else {
1038 cout << "TofRec::TofTrack::findTofDataBarrel: 1- Impossible!" << endl;
1039 }
1040 }
1041 else if( ( tofDataVec2.size() == 1 ) ) {
1042 if( tofDataVec1.size() == 0 ) {
1043 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1044 tof = (*iter2);
1045 qual = 0x4;
1046 }
1047 else if( tofDataVec1.size() == 1 ) {
1048 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1049 if( ((*iter1)->quality()&0x1ff)==0x01f && abs((*iter1)->ztdc()-zrhit)<ztdc_Cut ) {
1050 tof = (*iter1);
1051 }
1052 else {
1053 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1054 tof = compareTofData( (*iter1), (*iter2), zrhit );
1055 }
1056 qual = 0x5;
1057 }
1058 else if( tofDataVec1.size() > 1 ) {
1059 TofData* tofData1 = chooseTofData( tofDataVec1, zrhit );
1060 if( (tofData1->quality()&0x1ff)==0x01f && abs(tofData1->ztdc()-zrhit)<ztdc_Cut ) {
1061 tof = tofData1;
1062 }
1063 else {
1064 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1065 tof = compareTofData( tofData1, (*iter2), zrhit );
1066 }
1067 qual = 0x6;
1068 }
1069 else {
1070 cout << "TofRec::TofTrack::findTofDataBarrel: 2- Impossible!" << endl;
1071 }
1072 }
1073 else if( ( tofDataVec2.size() > 1 ) ) {
1074 if( tofDataVec1.size() == 0 ) {
1075 tof = chooseTofData( tofDataVec2, zrhit );
1076 qual = 0x8;
1077 }
1078 else if( tofDataVec1.size() == 1 ) {
1079 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1080 if( ((*iter1)->quality()&0x1ff)==0x01f && abs((*iter1)->ztdc()-zrhit)<ztdc_Cut ) {
1081 tof = (*iter1);
1082 }
1083 else {
1084 TofData* tofData2 = chooseTofData( tofDataVec2, zrhit );
1085 tof = compareTofData( (*iter1), tofData2, zrhit );
1086 }
1087 qual = 0x9;
1088 }
1089 else if( tofDataVec1.size() > 1 ) {
1090 TofData* tofData1 = chooseTofData( tofDataVec1, zrhit );
1091 if( (tofData1->quality()&0x1ff)==0x01f && abs(tofData1->ztdc()-zrhit)<ztdc_Cut ) {
1092 tof = tofData1;
1093 }
1094 else {
1095 TofData* tofData2 = chooseTofData( tofDataVec2, zrhit );
1096 tof = compareTofData( tofData1, tofData2, zrhit );
1097 }
1098 qual = 0xa;
1099 }
1100 else {
1101 cout << "TofRec::TofTrack::findTofDataBarrel: 3- Impossible!" << endl;
1102 }
1103 }
1104
1105 if( qual != 0 ) {
1106 if( !(tof->used()) ) {
1107 getTofData( tof, iflag );
1108 }
1109 else {
1110 bool z1=false, z2=false;
1111 bool zc1=false, zc2=false;
1112 TofTrack* track=0;
1113 if( iflag==1 ) {
1114 z1 = ( abs( m_zrhit1 - tof->ztdc() ) < ztdc_Cut );
1115 zc1 = ( m_zrhit1 > tof->ztdc() );
1116 std::vector<TofTrack*>::iterator iter = tofTrackVec->begin();
1117 for( ; iter!=tofTrackVec->end(); iter++ ) {
1118 if( (*iter)->hitCase()!=InnerLayer && (*iter)->hitCase()!=DoubleLayer ) continue;
1119 if( tof->tofId()==(*iter)->tofId1() ) {
1120 track = (*iter);
1121 z2 = ( abs( (*iter)->zrhit1() - tof->ztdc() ) < ztdc_Cut );
1122 zc2 = ( (*iter)->zrhit1() > tof->ztdc() );
1123 }
1124 }
1125 }
1126 else if( iflag==2 ) {
1127 z1 = ( abs( m_zrhit2 - tof->ztdc() ) < ztdc_Cut );
1128 zc1 = ( m_zrhit2 > tof->ztdc() );
1129 std::vector<TofTrack*>::iterator iter = tofTrackVec->begin();
1130 for( ; iter!=tofTrackVec->end(); iter++ ) {
1131 if( (*iter)->hitCase()!=OuterLayer && (*iter)->hitCase()!=DoubleLayer ) continue;
1132 if( tof->tofId()==(*iter)->tofId2() ) {
1133 track = (*iter);
1134 z2 = ( abs( (*iter)->zrhit2() - tof->ztdc() ) < ztdc_Cut );
1135 zc2 = ( (*iter)->zrhit2() > tof->ztdc() );
1136 }
1137 }
1138 }
1139
1140 if( ( z1 && z2 )||( (!z1) && (!z2) ) ) {
1141 if( zc1 && !zc2 ) {
1142 getTofDataEast( tof, iflag );
1143 track->getTofDataWest( tof, iflag );
1144 }
1145 else if( !zc1 && zc2 ) {
1146 getTofDataWest( tof, iflag );
1147 track->getTofDataEast( tof, iflag );
1148 }
1149 }
1150 else if( z1 && !z2 ) {
1151 getTofData( tof, iflag );
1152 track->getTofDataNohit( iflag );
1153 }
1154 else if( !z1 && z2 ) {
1155 qual = 0;
1156 }
1157 }
1158 }
1159
1160 if( qual == 0 ) {
1161 if( ( iflag == 1 ) || ( iflag == 3 ) ) {
1162 m_quality1 = ( m_quality1 | 0x300 );
1163 }
1164 else if( iflag == 2 ) {
1165 m_quality2 = ( m_quality2 | 0x300 );
1166 }
1167 else {
1168 cout << "TofRec::TofTrack::findTofDataBarrel: the 1- IFLAG is Out of Range!" << endl;
1169 }
1170 }
1171 else {
1172 qual = ( qual << 12 );
1173 if( ( iflag == 1 ) || ( iflag == 3 ) ) {
1174 m_quality1 = ( m_quality1 | qual );
1175 }
1176 else if( iflag == 2 ) {
1177 m_quality2 = ( m_quality2 | qual );
1178 }
1179 else {
1180 cout << "TofRec::TofTrack::findTofDataBarrel: the 2- IFLAG is Out of Range!" << endl;
1181 }
1182 }
1183
1184 return;
1185}
const double ztdc_Cut
Definition TofTrack.h:12
void getTofData(TofData *tof, unsigned int iflag)
TofData * compareTofData(TofData *tofData1, TofData *tofData2, double zrhit)
void getTofDataNohit(unsigned int iflag)
TofData * chooseTofData(std::vector< TofData * > tofDataVec, double zrhit)
void getTofDataWest(TofData *tof, unsigned int iflag)
void getTofDataEast(TofData *tof, unsigned int iflag)

Referenced by match().

◆ findTofDataEndcap()

void TofTrack::findTofDataEndcap ( std::vector< TofData * > tofDataVec1,
std::vector< TofData * > tofDataVec2,
double zr1[5] )

Definition at line 1271 of file TofTrack.cxx.

1271 {
1272
1273 unsigned int iflag = 3;
1274 unsigned int qual = 0xf;
1275
1276 if( tofDataVec2.size() == 0 ) {
1277 if( tofDataVec1.size() == 0 ) {
1278 qual = 0;
1279 }
1280 else if( tofDataVec1.size() == 1 ) {
1281 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1282 getTofData( (*iter1), iflag );
1283 qual = 0x1;
1284 }
1285 else if( tofDataVec1.size() > 1 ) {
1286 getTofData( chooseTofDataEndcap( tofDataVec1, zr1 ), iflag );
1287 qual = 0x2;
1288 }
1289 else {
1290 cout << "TofRec::TofTrack::findTofDataEndcap: 1- Impossible!" << endl;
1291 }
1292 }
1293 else if( ( tofDataVec2.size() == 1 ) ) {
1294 if( tofDataVec1.size() == 0 ) {
1295 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1296 getTofData( (*iter2), iflag );
1297 qual = 0x4;
1298 }
1299 else if( tofDataVec1.size() == 1 ) {
1300 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1301 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1302 getTofData( compareTofDataEndcap( (*iter1), (*iter2) ), iflag );
1303 qual = 0x5;
1304 }
1305 else if( tofDataVec1.size() > 1 ) {
1306 TofData* tofData1 = chooseTofDataEndcap( tofDataVec1, zr1 );
1307 std::vector<TofData*>::iterator iter2 = tofDataVec2.begin();
1308 getTofData( compareTofDataEndcap( tofData1, (*iter2) ), iflag );
1309 qual = 0x6;
1310 }
1311 else {
1312 cout << "TofRec::TofTrack::findTofDataBarrel: 2- Impossible!" << endl;
1313 }
1314 }
1315 else if( ( tofDataVec2.size() > 1 ) ) {
1316 if( tofDataVec1.size() == 0 ) {
1317 getTofData( chooseTofDataEndcap( tofDataVec2, zr1 ), iflag );
1318 qual = 0x8;
1319 }
1320 else if( tofDataVec1.size() == 1 ) {
1321 std::vector<TofData*>::iterator iter1 = tofDataVec1.begin();
1322 TofData* tofData2 = chooseTofDataEndcap( tofDataVec2, zr1 );
1323 getTofData( compareTofDataEndcap( (*iter1), tofData2 ), iflag );
1324 qual = 0x9;
1325 }
1326 else if( tofDataVec1.size() > 1 ) {
1327 TofData* tofData1 = chooseTofDataEndcap( tofDataVec1, zr1 );
1328 TofData* tofData2 = chooseTofDataEndcap( tofDataVec2, zr1 );
1330 qual = 0xa;
1331 }
1332 else {
1333 cout << "TofRec::TofTrack::findTofDataBarrel: 3- Impossible!" << endl;
1334 }
1335 }
1336
1337 if( qual == 0 ) {
1338 m_quality1 = ( m_quality1 | 0x300 );
1339 }
1340 else {
1341 qual = ( qual << 12 );
1342 m_quality1 = ( m_quality1 | qual );
1343 }
1344
1345 return;
1346}
TofData * chooseTofDataEndcap(std::vector< TofData * > tofDataVec, double zr1[5])
TofData * compareTofDataEndcap(TofData *tofData1, TofData *tofData2)

Referenced by match().

◆ flag()

unsigned int TofTrack::flag ( ) const
inline

Definition at line 117 of file TofTrack.h.

117{ return m_flag; }

Referenced by setFlag().

◆ getEtfData()

void TofTrack::getEtfData ( TofData * tof,
unsigned int iflag,
unsigned int qual )

Definition at line 1805 of file TofTrack.cxx.

1805 {
1806
1807 if( iflag == 1 && tof->tofId() != m_id2 ) {
1808 m_tofId1 = tof->tofId();
1809 m_strip1 = tof->strip();
1810 m_qch1 = tof->adcChannelEast();
1811 m_adc1 = tof->adc1();
1812 m_tdc1 = tof->tdc1();
1813 m_qch2 = tof->adcChannelWest();
1814 m_adc2 = tof->adc2();
1815 m_tdc2 = tof->tdc2();
1816 m_ztdc1 = tof->ztdc();
1817 m_zadc1 = tof->zadc();
1818 m_quality1 = ( m_quality1 | ( 0x1f & tof->quality() ) );
1819 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit1 ) > ztdc_Cut ) ) ) {
1820 m_quality1 = ( m_quality1 | 0x100 );
1821 }
1822 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit1 ) > zadc_Cut ) ) ) {
1823 m_quality1 = ( m_quality1 | 0x200 );
1824 }
1825 m_quality1 = ( m_quality1 | ( qual << 12 ) );
1826 tof->setUsed();
1827 if( abs(tof->tofId()-m_id1)== 1 || ( tof->tofId()==0 && m_id1==35 ) || ( tof->tofId()==35 && m_id1==0 ) || ( tof->tofId()==36 && m_id1==71 ) || ( tof->tofId()==71 && m_id1==36 ) ) {
1828 for( unsigned int i=0; i<5; i++ ) {
1829 m_texpInner[i] = m_texpOuter[i];
1830 }
1831 }
1832 }
1833 else if( iflag == 2 && tof->tofId() != m_id1 ) {
1834 m_tofId2 = tof->tofId();
1835 m_strip2 = tof->strip();
1836 m_qch3 = tof->adcChannelEast();
1837 m_adc3 = tof->adc1();
1838 m_tdc3 = tof->tdc1();
1839 m_qch4 = tof->adcChannelWest();
1840 m_adc4 = tof->adc2();
1841 m_tdc4 = tof->tdc2();
1842 m_ztdc2 = tof->ztdc();
1843 m_zadc2 = tof->zadc();
1844 m_quality2 = ( m_quality2 | ( 0x1f & tof->quality() ) );
1845 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit2 ) > ztdc_Cut ) ) ) {
1846 m_quality2 = ( m_quality2 | 0x100 );
1847 }
1848 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit2 ) > zadc_Cut ) ) ) {
1849 m_quality2 = ( m_quality2 | 0x200 );
1850 }
1851 m_quality2 = ( m_quality2 | ( qual << 12 ) );
1852 tof->setUsed();
1853 if( abs(tof->tofId()-m_id2)== 1 || ( tof->tofId()==0 && m_id2==35 ) || ( tof->tofId()==35 && m_id2==0 ) || ( tof->tofId()==36 && m_id2==71 ) || ( tof->tofId()==71 && m_id2==36 ) ) {
1854 for( unsigned int i=0; i<5; i++ ) {
1855 m_texpOuter[i] = m_texpInner[i];
1856 }
1857 }
1858 }
1859
1860 return;
1861}
const double zadc_Cut
Definition TofTrack.h:15

Referenced by findEtfData().

◆ getMultiHit()

void TofTrack::getMultiHit ( TofTrack *& track)

Definition at line 426 of file TofTrack.cxx.

426 {
427 if( m_hitCase == InnerLayer || m_hitCase == OuterLayer || m_hitCase == DoubleLayer ) {
428
429 if( ( m_hitCase==InnerLayer || m_hitCase==DoubleLayer ) && ( track->hitCase()==InnerLayer || track->hitCase()==DoubleLayer ) ) {
430 if( ( abs(m_id1-track->id1())<=1 ) || ( m_id1==0 && track->id1()==87 ) || ( m_id1==87 && track->id1()==0 ) ) {
431 track->setQuality1( ( track->quality1() | 0x400 ) );
432 m_quality1 = ( m_quality1 | 0x400 );
433 }
434 }
435
436 if( ( m_hitCase==OuterLayer || m_hitCase==DoubleLayer ) && ( track->hitCase()==OuterLayer || track->hitCase()==DoubleLayer ) ) {
437 if( ( abs(m_id2-track->id2())<=1 ) || ( m_id2==88 && track->id2()==175 ) || ( m_id2==175 && track->id2()==88 ) ) {
438 track->setQuality2( ( track->quality2() | 0x400 ) );
439 m_quality2 = ( m_quality2 | 0x400 );
440 }
441 }
442
443 }
444 else if( m_hitCase == EastEndcap ) {
445 if( track->hitCase()==EastEndcap ) {
446 if( ( abs(m_id1-track->id1())<=1 ) || ( m_id1==0 && track->id1()==47 ) || ( m_id1==47 && track->id1()==0 ) ) {
447 track->setQuality1( ( track->quality1() | 0x400 ) );
448 m_quality1 = ( m_quality1 | 0x400 );
449 }
450 }
451 }
452 else if( m_hitCase == WestEndcap ) {
453 if( track->hitCase()==WestEndcap ) {
454 if( ( abs(m_id1-track->id1())<=1 ) || ( m_id1==48 && track->id1()==95 ) || ( m_id1==95 && track->id1()==48 ) ) {
455 track->setQuality1( ( track->quality1() | 0x400 ) );
456 m_quality1 = ( m_quality1 | 0x400 );
457 }
458 }
459 }
460 if( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) {
461 if( ( m_hitCase==EastEndcapMRPC && track->hitCase()==EastEndcapMRPC ) || ( m_hitCase==WestEndcapMRPC && track->hitCase()==WestEndcapMRPC ) ) {
462 if( m_id1>=0 ) {
463 if( ( m_id1==track->id1() ) && abs(m_istrip1-track->strip1())<=1 ) {
464 track->setQuality1( ( track->quality1() | 0x400 ) );
465 m_quality1 = ( m_quality1 | 0x400 );
466 }
467 }
468 if( m_id2>=0 ) {
469 if( ( m_id2==track->id2() ) && abs(m_istrip2-track->strip2())<=1 ) {
470 track->setQuality1( ( track->quality1() | 0x400 ) );
471 m_quality1 = ( m_quality1 | 0x400 );
472 }
473 }
474 }
475 }
476
477 return;
478}
ExtTrackCase hitCase() const
Definition TofTrack.h:36
int strip1() const
Definition TofTrack.h:70
unsigned int quality1() const
Definition TofTrack.h:64
int id2() const
Definition TofTrack.h:30
int strip2() const
Definition TofTrack.h:71
void setQuality1(int qual1)
Definition TofTrack.h:122
void setQuality2(int qual2)
Definition TofTrack.h:123
unsigned int quality2() const
Definition TofTrack.h:65
int id1() const
Definition TofTrack.h:29

Referenced by TofRec::execute().

◆ getTofData()

void TofTrack::getTofData ( TofData * tof,
unsigned int iflag )

Definition at line 1559 of file TofTrack.cxx.

1559 {
1560
1561 if( iflag == 1 ) {
1562 m_tofId1 = tof->tofId();
1563 m_strip1 = tof->strip();
1564 if( tofCaliSvc->QElec() ) {
1565 m_qch1 = tof->qtc1();
1566 }
1567 else {
1568 m_qch1 = tof->adcChannelEast();
1569 }
1570 m_adc1 = tof->adc1();
1571 m_tdc1 = tof->tdc1();
1572 if( tofCaliSvc->QElec() ) {
1573 m_qch2 = tof->qtc2();
1574 }
1575 else {
1576 m_qch2 = tof->adcChannelWest();
1577 }
1578 m_adc2 = tof->adc2();
1579 m_tdc2 = tof->tdc2();
1580 m_ztdc1 = tof->ztdc();
1581 m_zadc1 = tof->zadc();
1582 m_quality1 = ( m_quality1 | ( 0x1f & tof->quality() ) );
1583 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit1 ) > ztdc_Cut ) ) ) {
1584 m_quality1 = ( m_quality1 | 0x100 );
1585 }
1586 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit1 ) > zadc_Cut ) ) ) {
1587 m_quality1 = ( m_quality1 | 0x200 );
1588 }
1589 }
1590 else if( iflag == 2 ) {
1591 m_tofId2 = tof->tofId();
1592 m_strip2 = tof->strip();
1593 if( tofCaliSvc->QElec() ) {
1594 m_qch3 = tof->qtc1();
1595 }
1596 else {
1597 m_qch3 = tof->adcChannelEast();
1598 }
1599 m_adc3 = tof->adc1();
1600 m_tdc3 = tof->tdc1();
1601 if( tofCaliSvc->QElec() ) {
1602 m_qch4 = tof->qtc2();
1603 }
1604 else {
1605 m_qch4 = tof->adcChannelWest();
1606 }
1607 m_adc4 = tof->adc2();
1608 m_tdc4 = tof->tdc2();
1609 m_ztdc2 = tof->ztdc();
1610 m_zadc2 = tof->zadc();
1611 m_quality2 = ( m_quality2 | ( 0x1f & tof->quality() ) );
1612 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit2 ) > ztdc_Cut ) ) ) {
1613 m_quality2 = ( m_quality2 | 0x100 );
1614 }
1615 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit2 ) > zadc_Cut ) ) ) {
1616 m_quality2 = ( m_quality2 | 0x200 );
1617 }
1618 }
1619 else if( iflag == 3 ) {
1620 m_tofId1 = tof->tofId();
1621 if( tofCaliSvc->QElec() ) {
1622 m_qch1 = tof->qtc();
1623 }
1624 else {
1625 m_qch1 = tof->adcChannel();
1626 }
1627 m_adc1 = tof->adc();
1628 m_tdc1 = tof->tdc();
1629 m_quality1 = ( m_quality1 | ( 0x1f & tof->quality() ) );
1630 m_quality1 = ( m_quality1 | 0x300 );
1631 }
1632 else {
1633 cout << "TofRec::TofTrack::getTofData: Flag which sign the Barrel/Endcap or Inner/Outer is wrong! Please check it!" << endl;
1634 }
1635 tof->setUsed();
1636 return;
1637}
virtual const int QElec()=0

Referenced by findTofDataBarrel(), and findTofDataEndcap().

◆ getTofDataEast()

void TofTrack::getTofDataEast ( TofData * tof,
unsigned int iflag )

Definition at line 1643 of file TofTrack.cxx.

1643 {
1644
1645 if( iflag == 1 ) {
1646 m_tofId1 = tof->tofId();
1647 m_strip1 = tof->strip();
1648 if( tofCaliSvc->QElec() ) {
1649 m_qch1 = tof->qtc1();
1650 }
1651 else {
1652 m_qch1 = tof->adcChannelEast();
1653 }
1654 m_adc1 = tof->adc1();
1655 m_tdc1 = tof->tdc1();
1656 m_qch2 = -999.0;
1657 m_adc2 = -999.0;
1658 m_tdc2 = -999.0;
1659 m_ztdc1 = tof->ztdc();
1660 m_zadc1 = tof->zadc();
1661 m_quality1 = ( ( m_quality1 & 0xfffffff0 ) | ( 0x1c & tof->quality() ) );
1662 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit1 ) > ztdc_Cut ) ) ) {
1663 m_quality1 = ( m_quality1 | 0x100 );
1664 }
1665 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit1 ) > zadc_Cut ) ) ) {
1666 m_quality1 = ( m_quality1 | 0x200 );
1667 }
1668 }
1669 else if( iflag == 2 ) {
1670 m_tofId2 = tof->tofId();
1671 if( tofCaliSvc->QElec() ) {
1672 m_qch3 = tof->qtc1();
1673 }
1674 else {
1675 m_qch3 = tof->adcChannelEast();
1676 }
1677 m_adc3 = tof->adc1();
1678 m_tdc3 = tof->tdc1();
1679 m_qch4 = -999.0;
1680 m_adc4 = -999.0;
1681 m_tdc4 = -999.0;
1682 m_ztdc2 = tof->ztdc();
1683 m_zadc2 = tof->zadc();
1684 m_quality2 = ( ( m_quality2 & 0xfffffff0 ) | ( 0x1c & tof->quality() ) );
1685 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit2 ) > ztdc_Cut ) ) ) {
1686 m_quality2 = ( m_quality2 | 0x100 );
1687 }
1688 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit2 ) > zadc_Cut ) ) ) {
1689 m_quality2 = ( m_quality2 | 0x200 );
1690 }
1691 }
1692 else {
1693 cout << "TofRec::TofTrack::getTofDataEast: Flag which sign the Barrel/Endcap or Inner/Outer is wrong! Please check it!" << endl;
1694 }
1695 tof->setUsed();
1696 return;
1697}

Referenced by findTofDataBarrel().

◆ getTofDataNohit()

void TofTrack::getTofDataNohit ( unsigned int iflag)

Definition at line 1763 of file TofTrack.cxx.

1763 {
1764
1765 if( iflag == 1 ) {
1766 m_tofId1 = -99;
1767 m_strip1 = -99;
1768 m_qch1 = -999.0;
1769 m_adc1 = -999.0;
1770 m_tdc1 = -999.0;
1771 m_qch2 = -999.0;
1772 m_adc2 = -999.0;
1773 m_tdc2 = -999.0;
1774 m_ztdc1 = -999.0;
1775 m_zadc1 = -999.0;
1776 m_quality1 = ( m_quality1 & 0x700 );
1777 if( m_hitCase == InnerLayer ) { m_hitCase = NoHit; }
1778 else if( m_hitCase == DoubleLayer ) { m_hitCase = OuterLayer; }
1779 }
1780 else if( iflag == 2 ) {
1781 m_tofId2 = -99;
1782 m_qch3 = -999.0;
1783 m_adc3 = -999.0;
1784 m_tdc3 = -999.0;
1785 m_qch4 = -999.0;
1786 m_adc4 = -999.0;
1787 m_tdc4 = -999.0;
1788 m_ztdc2 = -999.0;
1789 m_zadc2 = -999.0;
1790 m_quality2 = ( m_quality2 & 0x700 );
1791 if( m_hitCase == OuterLayer ) { m_hitCase = NoHit; }
1792 else if( m_hitCase == DoubleLayer ) { m_hitCase = InnerLayer; }
1793 }
1794 else {
1795 cout << "TofRec::TofTrack::getTofData: Flag which sign the Barrel/Endcap or Inner/Outer is wrong! Please check it!" << endl;
1796 }
1797
1798 return;
1799}

Referenced by findTofDataBarrel().

◆ getTofDataWest()

void TofTrack::getTofDataWest ( TofData * tof,
unsigned int iflag )

Definition at line 1703 of file TofTrack.cxx.

1703 {
1704
1705 if( iflag == 1 ) {
1706 m_tofId1 = tof->tofId();
1707 m_strip1 = tof->strip();
1708 m_qch1 = -999.0;
1709 m_adc1 = -999.0;
1710 m_tdc1 = -999.0;
1711 if( tofCaliSvc->QElec() ) {
1712 m_qch2 = tof->qtc2();
1713 }
1714 else {
1715 m_qch2 = tof->adcChannelWest();
1716 }
1717 m_adc2 = tof->adc2();
1718 m_tdc2 = tof->tdc2();
1719 m_ztdc1 = tof->ztdc();
1720 m_zadc1 = tof->zadc();
1721 m_quality1 = ( ( m_quality1 & 0xfffffff0 ) | ( 0x13 & tof->quality() ) );
1722 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit1 ) > ztdc_Cut ) ) ) {
1723 m_quality1 = ( m_quality1 | 0x100 );
1724 }
1725 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit1 ) > zadc_Cut ) ) ) {
1726 m_quality1 = ( m_quality1 | 0x200 );
1727 }
1728 }
1729 else if( iflag == 2 ) {
1730 m_tofId2 = tof->tofId();
1731 m_qch3 = -999.0;
1732 m_adc3 = -999.0;
1733 m_tdc3 = -999.0;
1734 if( tofCaliSvc->QElec() ) {
1735 m_qch4 = tof->qtc2();
1736 }
1737 else {
1738 m_qch4 = tof->adcChannelWest();
1739 }
1740 m_adc4 = tof->adc2();
1741 m_tdc4 = tof->tdc2();
1742 m_ztdc2 = tof->ztdc();
1743 m_zadc2 = tof->zadc();
1744 m_quality2 = ( ( m_quality2 & 0xfffffff0 ) | ( 0x13 & tof->quality() ) );
1745 if( ( ( tof->quality() & 0x5 ) != 0x5 ) || ( ( ( tof->quality() & 0x5 ) == 0x5 ) && ( abs( tof->ztdc() - m_zrhit2 ) > ztdc_Cut ) ) ) {
1746 m_quality2 = ( m_quality2 | 0x100 );
1747 }
1748 if( ( ( tof->quality() & 0xa ) != 0xa ) || ( ( ( tof->quality() & 0xa ) == 0xa ) && ( abs( tof->zadc() - m_zrhit2 ) > zadc_Cut ) ) ) {
1749 m_quality2 = ( m_quality2 | 0x200 );
1750 }
1751 }
1752 else {
1753 cout << "TofRec::TofTrack::getTofDataWest: Flag which sign the Barrel/Endcap or Inner/Outer is wrong! Please check it!" << endl;
1754 }
1755 tof->setUsed();
1756 return;
1757}

Referenced by findTofDataBarrel().

◆ hitCase()

◆ id1()

int TofTrack::id1 ( ) const
inline

◆ id2()

int TofTrack::id2 ( ) const
inline

◆ isNoHit()

bool TofTrack::isNoHit ( ) const
inline

Definition at line 120 of file TofTrack.h.

120{ return m_hitCase==NoHit; }

◆ istrip1()

int TofTrack::istrip1 ( ) const
inline

◆ istrip2()

int TofTrack::istrip2 ( ) const
inline

◆ kal()

int TofTrack::kal ( unsigned int i) const
inline

◆ match()

void TofTrack::match ( bool forCalibration,
std::vector< int > deadId,
std::vector< TofTrack * > *& tofTrackVec )

Definition at line 849 of file TofTrack.cxx.

849 {
850
851 if( m_hitCase == NoHit ) return;
852
853 if( m_hitCase == InnerLayer ) {
854 findTofDataBarrel( m_tofData1, m_tofData2, m_zrhit1, 1, tofTrackVec );
855 if( ( m_quality1 & 0x10 ) == 0 ) { m_hitCase = NoHit; }
856 }
857 else if( m_hitCase == OuterLayer ) {
858 findTofDataBarrel( m_tofData3, m_tofData4, m_zrhit2, 2, tofTrackVec );
859 if( ( m_quality2 & 0x10 ) == 0 ) { m_hitCase = NoHit; }
860 }
861 else if( m_hitCase == DoubleLayer ) {
862 findTofDataBarrel( m_tofData1, m_tofData2, m_zrhit1, 1, tofTrackVec );
863 if( ( m_quality1 & 0x10 ) == 0 ) { m_hitCase = OuterLayer; }
864 findTofDataBarrel( m_tofData3, m_tofData4, m_zrhit2, 2, tofTrackVec );
865 if( ( m_quality2 & 0x10 ) == 0 ) {
866 if( m_hitCase == DoubleLayer ) {
867 m_hitCase = InnerLayer;
868 }
869 else if( m_hitCase == OuterLayer ) {
870 m_hitCase = NoHit;
871 }
872 else {
873 cout << "TofRec::TofTrack::match: 2- Impossible!" << endl;
874 }
875 }
876 }
877 else if( ( m_hitCase == EastEndcap ) || ( m_hitCase == WestEndcap ) ) {
878 findTofDataEndcap( m_tofData1, m_tofData2, m_zr1 );
879 }
880 else if( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) {
881 findEtfData( m_tofData1, m_tofData2, m_tofData3, m_zrhit1, 1 );
882 findEtfData( m_tofData4, m_tofData5, m_tofData6, m_zrhit2, 2 );
883 }
884 else {
885 cout << "TofRec::TofTrack::match: 1- Impossible!" << endl;
886 }
887
888 if( forCalibration ) {
889 // set Data Sample for Calibration, double layer, only one hit for counter, T and Q.
890 if( m_hitCase == DoubleLayer ) {
891 if( ( ( m_quality1 & 0xf ) == 0xf ) && ( ( m_quality2 & 0xf ) == 0xf ) ) {
892 m_quality1 = ( m_quality1 | 0x800 ); // Calibration Sample
893 m_quality2 = ( m_quality2 | 0x800 ); // Calibration Sample
894 }
895 else {
896 std::vector<int>::iterator iter = deadId.begin();
897 for( ; iter != deadId.end(); iter++ ) {
898 Identifier iden = Identifier(*iter);
899 int barrel = TofID::barrel_ec(iden);
900 int layer = TofID::layer(iden);
901 int tofId = TofID::phi_module(iden);
902 int east = TofID::end(iden);
903 if( barrel == 1 ) {
904 if( layer==0 ) {
905 if( m_tofId1 == tofId ) {
906 if( ( m_quality2 & 0xf ) == 0xf ) {
907 if( ( ( east == 0 ) && ( ( ( m_quality1 & 0xf ) == 0x3 ) || ( ( m_quality1 & 0xf ) == 0x7 ) || ( ( m_quality1 & 0xf ) == 0xb ) ) ) || ( ( east == 1 ) && ( ( ( m_quality1 & 0xf ) == 0xc ) || ( ( m_quality1 & 0xf ) == 0xd ) || ( ( m_quality1 & 0xf ) == 0xe ) ) ) ) {
908 m_quality1 = ( m_quality1 | 0x800 );
909 m_quality2 = ( m_quality2 | 0x800 );
910
911 }
912 }
913 }
914 }
915 else if( layer == 1 ) {
916 if( m_tofId2 == (tofId+88) ) {
917 if( ( m_quality1 & 0xf ) == 0xf ) {
918 if( ( ( east == 0 ) && ( ( ( m_quality2 & 0xf ) == 0x3 ) || ( ( m_quality2 & 0xf ) == 0x7 ) || ( ( m_quality2 & 0xf ) == 0xb ) ) ) || ( ( east == 1 ) && ( ( ( m_quality2 & 0xf ) == 0xc ) || ( ( m_quality2 & 0xf ) == 0xd ) || ( ( m_quality2 & 0xf ) == 0xe ) ) ) ) {
919 m_quality1 = ( m_quality1 | 0x800 );
920 m_quality2 = ( m_quality2 | 0x800 );
921 }
922 }
923 }
924 }
925 }
926 }
927 }
928 }
929 // set Data Sample for Calibration, only one hit for counter, T and Q.
930 else if( m_hitCase == InnerLayer ) {
931 if( ( m_quality1 & 0xf ) == 0xf ) {
932 m_quality1 = ( m_quality1 | 0x800 ); // Calibration Sample
933 }
934 else {
935 std::vector<int>::iterator iter = deadId.begin();
936 for( ; iter != deadId.end(); iter++ ) {
937 Identifier iden = Identifier(*iter);
938 int barrel = TofID::barrel_ec(iden);
939 int layer = TofID::layer(iden);
940 int tofId = TofID::phi_module(iden);
941 int east = TofID::end(iden);
942 if( barrel == 1 ) {
943 if( layer==0 ) {
944 if( m_tofId1 == tofId ) {
945 if( ( ( east == 0 ) && ( ( ( m_quality1 & 0xf ) == 0x3 ) || ( ( m_quality1 & 0xf ) == 0x7 ) || ( ( m_quality1 & 0xf ) == 0xb ) ) ) || ( ( east == 1 ) && ( ( ( m_quality1 & 0xf ) == 0xc ) || ( ( m_quality1 & 0xf ) == 0xd ) || ( ( m_quality1 & 0xf ) == 0xe ) ) ) ) {
946 m_quality1 = ( m_quality1 | 0x800 );
947 }
948 }
949 }
950 }
951 }
952 }
953 }
954
955 // set Data Sample for Calibration, only one hit for counter, T and Q.
956 if( ( ( m_hitCase == EastEndcap ) || ( m_hitCase == WestEndcap ) ) && ( ( m_quality1 & 0xf ) == 0xc ) ) {
957 m_quality1 = ( m_quality1 | 0x800 ); // Calibration Sample
958 }
959
960 // set Data Sample for Calibration.
961 if( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) {
962 if( ( ( m_quality1 & 0xf000 ) == 0x1000 ) || ( ( m_quality1 & 0xf000 ) == 0x2000 ) ) {
963 if( ( m_quality1 & 0xf ) == 0xf ) {
964 m_quality1 = ( m_quality1 | 0x800 ); // Calibration Sample
965 }
966 }
967 else if( ( ( m_quality1 & 0xf000 ) == 0x4000 ) ) {
968 std::vector<int>::iterator iter = deadId.begin();
969 for( ; iter != deadId.end(); iter++ ) {
970 Identifier iden = Identifier(*iter);
971 int barrel = TofID::barrel_ec( iden );
972 if( barrel == 3 ) {
973 int endcap = TofID::endcap( iden );
974 int tofid = TofID::module( iden );
975 int strip = TofID::strip( iden );
976 int east = TofID::end( iden );
977 if( ( m_tofId1 == ( endcap*36 + tofid ) ) && ( m_strip1 == strip ) ) {
978 if( ( ( east == 0 ) && ( ( ( m_quality1 & 0xf ) == 0x3 ) || ( ( m_quality1 & 0xf ) == 0x7 ) || ( ( m_quality1 & 0xf ) == 0xb ) ) ) || ( ( east == 1 ) && ( ( ( m_quality1 & 0xf ) == 0xc ) || ( ( m_quality1 & 0xf ) == 0xd ) || ( ( m_quality1 & 0xf ) == 0xe ) ) ) ) {
979 m_quality1 = ( m_quality1 | 0x800 );
980 }
981 }
982 }
983 }
984 }
985
986 if( ( ( m_quality2 & 0xf000 ) == 0x1000 ) || ( ( m_quality2 & 0xf000 ) == 0x2000 ) ) {
987 if( ( m_quality2 & 0xf ) == 0xf ) {
988 m_quality2 = ( m_quality2 | 0x800 ); // Calibration Sample
989 }
990 }
991 else if( ( m_quality2 & 0xf000 ) == 0x4000 ) {
992 std::vector<int>::iterator iter = deadId.begin();
993 for( ; iter != deadId.end(); iter++ ) {
994 Identifier iden = Identifier(*iter);
995 int barrel = TofID::barrel_ec( iden );
996 if( barrel == 3 ) {
997 int endcap = TofID::endcap( iden );
998 int tofid = TofID::module( iden );
999 int strip = TofID::strip( iden );
1000 int east = TofID::end( iden );
1001 if( ( m_tofId2 == ( endcap*36 + tofid ) ) && ( m_strip2 == strip ) ) {
1002 if( ( ( east == 0 ) && ( ( ( m_quality2 & 0xf ) == 0x3 ) || ( ( m_quality2 & 0xf ) == 0x7 ) || ( ( m_quality2 & 0xf ) == 0xb ) ) ) || ( ( east == 1 ) && ( ( ( m_quality2 & 0xf ) == 0xc ) || ( ( m_quality2 & 0xf ) == 0xd ) || ( ( m_quality2 & 0xf ) == 0xe ) ) ) ) {
1003 m_quality2 = ( m_quality2 | 0x800 );
1004 }
1005 }
1006 }
1007 }
1008 }
1009 }
1010 }
1011
1012 return;
1013}
static int endcap(const Identifier &id)
Definition TofID.cxx:124
static int strip(const Identifier &id)
Definition TofID.cxx:136
static int end(const Identifier &id)
Definition TofID.cxx:79
static int phi_module(const Identifier &id)
Definition TofID.cxx:73
static int barrel_ec(const Identifier &id)
Values of different levels (failure returns 0)
Definition TofID.cxx:61
static int layer(const Identifier &id)
Definition TofID.cxx:66
void findEtfData(std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, std::vector< TofData * > tofDataVec3, double zrhit, unsigned int iflag)
void findTofDataEndcap(std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, double zr1[5])
void findTofDataBarrel(std::vector< TofData * > tofDataVec1, std::vector< TofData * > tofDataVec2, double zrhit, unsigned int iflag, std::vector< TofTrack * > *&tofTrackVec)

◆ p()

double TofTrack::p ( ) const
inline

◆ path()

double TofTrack::path ( ) const
inline

Definition at line 38 of file TofTrack.h.

38{ return m_path; }

◆ path1()

double TofTrack::path1 ( ) const
inline

◆ path2()

double TofTrack::path2 ( ) const
inline

◆ ph()

double TofTrack::ph ( ) const
inline

◆ ph1()

double TofTrack::ph1 ( ) const
inline

◆ ph11()

double TofTrack::ph11 ( ) const
inline

◆ ph12()

double TofTrack::ph12 ( ) const
inline

◆ ph2()

double TofTrack::ph2 ( ) const
inline

◆ ph21()

double TofTrack::ph21 ( ) const
inline

◆ ph22()

double TofTrack::ph22 ( ) const
inline

◆ phi1()

double TofTrack::phi1 ( ) const
inline

◆ phi2()

double TofTrack::phi2 ( ) const
inline

◆ qch1()

double TofTrack::qch1 ( ) const
inline

◆ qch2()

double TofTrack::qch2 ( ) const
inline

◆ qch3()

double TofTrack::qch3 ( ) const
inline

◆ qch4()

double TofTrack::qch4 ( ) const
inline

◆ quality()

unsigned int TofTrack::quality ( ) const
inline

◆ quality1()

unsigned int TofTrack::quality1 ( ) const
inline

◆ quality2()

unsigned int TofTrack::quality2 ( ) const
inline

◆ qualityAnalysis()

void TofTrack::qualityAnalysis ( )

◆ setCalibration()

void TofTrack::setCalibration ( )

Definition at line 1869 of file TofTrack.cxx.

1869 {
1870
1871 bool barrel = ( ( m_hitCase == InnerLayer ) || ( m_hitCase == OuterLayer ) || ( m_hitCase == DoubleLayer ) );
1872 bool endcap = ( ( m_hitCase == EastEndcap ) || ( m_hitCase == WestEndcap ) );
1873 bool endcapMRPC = ( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) );
1874
1875 bool innerEast = ( ( m_quality1 & 0xc ) == 0xc );
1876 bool innerWest = ( ( m_quality1 & 0x3 ) == 0x3 );
1877 bool outerEast = ( ( m_quality2 & 0xc ) == 0xc );
1878 bool outerWest = ( ( m_quality2 & 0x3 ) == 0x3 );
1879 bool innerLayer = ( ( m_quality1 & 0xf ) == 0xf );
1880 bool outerLayer = ( ( m_quality2 & 0xf ) == 0xf );
1881
1882 bool endcapData = ( ( m_quality1 & 0xc ) == 0xc );
1883
1884 if( m_hitCase == DoubleLayer ) {
1885 for( unsigned int i=0; i<5; i++ ) {
1886 m_texp[i] = tofCaliSvc->BTimeCluster( m_texpInner[i], m_texpOuter[i], m_zr1[i], m_zr2[i], m_tofId1, m_tofId2 );
1887 }
1888 m_path = tofCaliSvc->BTimeCluster( m_path1, m_path2, m_zrhit1, m_zrhit2, m_tofId1, m_tofId2 );
1889 }
1890
1891 if( barrel ) {
1892 if( innerEast ) {
1893 for( unsigned int i=0; i<5; i++ ) {
1894 m_tof11[i] = tofCaliSvc->BTime1( m_adc1, m_tdc1-m_estime, m_zr1[i], m_tofId1, m_estime );
1895 }
1896 m_sigma11 = tofCaliSvc->BSigma1( m_zrhit1, m_tofId1 );
1897 m_ph11 = m_adc1;
1898 }
1899
1900 if( innerWest ) {
1901 for( unsigned int i=0; i<5; i++ ) {
1902 m_tof12[i] = tofCaliSvc->BTime2( m_adc2, m_tdc2-m_estime, m_zr1[i], m_tofId1, m_estime );
1903 }
1904 m_sigma12 = tofCaliSvc->BSigma2( m_zrhit1, m_tofId1 );
1905 m_ph12 = m_adc2;
1906 }
1907
1908 if( innerLayer ) {
1909 for( unsigned int i=0; i<5; i++ ) {
1910 m_tof1[i] = tofCaliSvc->BTimeCounter( m_tof11[i], m_tof12[i], m_zr1[i], m_tofId1 );
1911 }
1912 m_sigma1 = tofCaliSvc->BSigmaCounter( m_zrhit1, m_tofId1 );
1913 m_ph1 = tofCaliSvc->BPulseHeight( m_adc1, m_adc2, m_zrhit1, m_theta1, m_tofId1 );
1914 }
1915
1916 if( outerEast ) {
1917 for( unsigned int i=0; i<5; i++ ) {
1918 m_tof21[i] = tofCaliSvc->BTime1( m_adc3, m_tdc3-m_estime, m_zr2[i], m_tofId2, m_estime );
1919 }
1920 m_sigma21 = tofCaliSvc->BSigma1( m_zrhit2, m_tofId2 );
1921 m_ph21 = m_adc3;
1922 }
1923
1924 if( outerWest ) {
1925 for( unsigned int i=0; i<5; i++ ) {
1926 m_tof22[i] = tofCaliSvc->BTime2( m_adc4, m_tdc4-m_estime, m_zr2[i], m_tofId2, m_estime );
1927 }
1928 m_sigma22 = tofCaliSvc->BSigma2( m_zrhit2, m_tofId2 );
1929 m_ph22 = m_adc4;
1930 }
1931
1932 if( outerLayer ) {
1933 for( unsigned int i=0; i<5; i++ ) {
1934 m_tof2[i] = tofCaliSvc->BTimeCounter( m_tof21[i], m_tof22[i], m_zr2[i], m_tofId2 );
1935 }
1936 m_sigma2 = tofCaliSvc->BSigmaCounter( m_zrhit2, m_tofId2 );
1937 m_ph2 = tofCaliSvc->BPulseHeight( m_adc3, m_adc4, m_zrhit2, m_theta2, m_tofId2 );
1938 }
1939
1940 if( innerLayer && outerLayer ) {
1941 for( unsigned int i=0; i<5; i++ ) {
1942 m_tof[i] = tofCaliSvc->BTimeCluster( m_tof1[i], m_tof2[i], m_zr1[i], m_zr2[i], m_tofId1, m_tofId2 );
1943 }
1944 m_sigma = tofCaliSvc->BSigmaCluster( m_zrhit1, m_zrhit2, m_tofId1, m_tofId2 );
1945 m_ph = tofCaliSvc->BTimeCluster( m_ph1, m_ph2, m_zrhit1, m_zrhit2, m_tofId1, m_tofId2 );
1946 }
1947 }
1948
1949 if( endcap ) {
1950 if( endcapData ) {
1951 for( unsigned int i=0; i<5; i++ ) {
1952 m_tof11[i] = tofCaliSvc->ETime( m_adc1, m_tdc1-m_estime, m_zr1[i], m_tofId1 );
1953 }
1954 m_sigma11 = tofCaliSvc->ESigma( m_zrhit1, m_tofId1 );
1955 m_ph11 = tofCaliSvc->EPulseHeight( m_adc1, m_zrhit1, m_theta1, m_tofId1 );
1956 m_quality = 1;
1957 if( (m_quality1&0xa000)!=0 ) { m_quality = 4; }
1958 }
1959 }
1960
1961 if( endcapMRPC ) {
1962 if( innerEast ) {
1963 if( m_tofId1>-1 ) {
1964 for( unsigned int i=0; i<5; i++ ) {
1965 if( m_run > 0 ) {
1966 m_tof11[i] = tofCaliSvc->EtfTime1( m_adc1, m_tdc1-m_estime, m_zr1[i], m_tofId1, m_strip1, m_estime );
1967 }
1968 else {
1969 m_tof11[i] = tofCaliSvc->EtfTimeMC1( m_adc1, m_tdc1-m_estime, m_zr1[i], m_tofId1, m_strip1, m_estime );
1970 }
1971 }
1972 m_ph11 = m_adc1;
1973 }
1974 }
1975 if( innerWest ) {
1976 if( m_tofId1>-1 ) {
1977 for( unsigned int i=0; i<5; i++ ) {
1978 if( m_run > 0 ) {
1979 m_tof12[i] = tofCaliSvc->EtfTime2( m_adc2, m_tdc2-m_estime, m_zr1[i], m_tofId1, m_strip1, m_estime );
1980 }
1981 else {
1982 m_tof12[i] = tofCaliSvc->EtfTimeMC2( m_adc2, m_tdc2-m_estime, m_zr1[i], m_tofId1, m_strip1, m_estime );
1983 }
1984 }
1985 m_ph12 = m_adc2;
1986 }
1987 }
1988 if( outerEast ) {
1989 if( m_tofId2>-1 ) {
1990 for( unsigned int i=0; i<5; i++ ) {
1991 if( m_run > 0 ) {
1992 m_tof21[i] = tofCaliSvc->EtfTime1( m_adc3, m_tdc3-m_estime, m_zr2[i], m_tofId2, m_strip2, m_estime );
1993 }
1994 else {
1995 m_tof21[i] = tofCaliSvc->EtfTimeMC1( m_adc3, m_tdc3-m_estime, m_zr2[i], m_tofId2, m_strip2, m_estime );
1996 }
1997 }
1998 m_ph21 = m_adc3;
1999 }
2000 }
2001 if( outerWest ) {
2002 if( m_tofId2>-1 ) {
2003 for( unsigned int i=0; i<5; i++ ) {
2004 if( m_run > 0 ) {
2005 m_tof22[i] = tofCaliSvc->EtfTime2( m_adc4, m_tdc4-m_estime, m_zr2[i], m_tofId2, m_strip2, m_estime );
2006 }
2007 else {
2008 m_tof22[i] = tofCaliSvc->EtfTimeMC2( m_adc4, m_tdc4-m_estime, m_zr2[i], m_tofId2, m_strip2, m_estime );
2009 }
2010 }
2011 m_ph22 = m_adc4;
2012 }
2013 }
2014 if( innerLayer ) {
2015 if( m_tofId1>-1 ) {
2016 m_tof1[0] = tofCaliSvc->EtfTime( m_tof11[0], m_tof12[0] );
2017 // if( m_run > 0 ) {
2018 // m_tof1[0] = tofCaliSvc->EtfTime( m_adc1, m_adc2, m_tdc1-m_estime, m_tdc2-m_estime, m_tofId1, m_strip1, m_estime );
2019 // }
2020 // else {
2021 // m_tof1[0] = tofCaliSvc->EtfTimeMC( m_adc1, m_adc2, m_tdc1-m_estime, m_tdc2-m_estime, m_tofId1, m_strip1, m_estime );
2022 // m_tof1[0] = tofCaliSvc->EtfTime( m_tof11[0], m_tof12[0] );
2023 // }
2024 for( unsigned int i=1; i<5; i++ ) {
2025 m_tof1[i] = m_tof1[0];
2026 }
2027 m_ph1 = ( m_adc1 + m_adc2 )/2.0;
2028 }
2029 }
2030 if( outerLayer ) {
2031 if( m_tofId2>-1 ) {
2032 m_tof2[0] = tofCaliSvc->EtfTime( m_tof21[0], m_tof22[0] );
2033 // if( m_run > 0 ) {
2034 // m_tof2[0] = tofCaliSvc->EtfTime( m_adc3, m_adc4, m_tdc3-m_estime, m_tdc4-m_estime, m_tofId2, m_strip2, m_estime );
2035 // }
2036 // else {
2037 // m_tof2[0] = tofCaliSvc->EtfTimeMC( m_adc3, m_adc4, m_tdc3-m_estime, m_tdc4-m_estime, m_tofId2, m_strip2, m_estime );
2038 // }
2039 for( unsigned int i=1; i<5; i++ ) {
2040 m_tof2[i] = m_tof2[0];
2041 }
2042 m_ph2 = ( m_adc3 + m_adc4 )/2.0;
2043 }
2044 }
2045 }
2046
2047 // set Quality of Barrel TOF
2048 if( barrel ) {
2049
2050 // double layer
2051 if( innerLayer && outerLayer ) {
2052 m_quality = 1;
2053 }
2054 else {
2055 // single layer
2056 if( innerLayer || outerLayer ) {
2057 m_quality = 2;
2058 }
2059 else {
2060 // single-end of one layer
2061 if( innerEast || innerWest || outerEast || outerWest ) {
2062 m_quality = 3;
2063 }
2064 }
2065 }
2066
2067 // multi-hit
2068 if( ( (m_quality1&0xa000)!=0 ) || ( (m_quality2&0xa000)!=0 ) ) {
2069 m_quality = m_quality + 3;
2070 }
2071
2072 // ztdc and extrapolated zhit is not matched
2073 if( ( (m_quality1&0x100)==0x100 ) || ( (m_quality2&0x100)==0x100 ) ) {
2074 if( ( m_quality == 1 ) || ( m_quality == 4 ) ) { m_quality = 7; }
2075 else if( ( m_quality == 2 ) || ( m_quality == 5 ) ) { m_quality = 8; }
2076 else if( ( m_quality == 3 ) || ( m_quality == 6 ) ) { m_quality = 9; }
2077 else {
2078 cout << "TofRec::TofTrack::setCalibration: Impossible!" << endl;
2079 }
2080 }
2081
2082 }
2083
2084 return;
2085}
virtual const double BSigmaCounter(double zHit, unsigned id)=0
virtual const double BSigma2(double zHit, unsigned id)=0
virtual const double BTimeCluster(double tlayer1, double tlayer2, double z1, double z2, unsigned id1, unsigned int id2)=0
virtual const double EPulseHeight(double ADC, double rHit, double cost, unsigned int id)=0
virtual const double BTime2(double ADC, double TDC, double zHit, unsigned id, double t0)=0
virtual const double BTime1(double ADC, double TDC, double zHit, unsigned id, double t0)=0
virtual const double EtfTime1(double ADC, double TDC, double zHit, unsigned id, unsigned strip, double t0)=0
virtual const double EtfTime2(double ADC, double TDC, double zHit, unsigned id, unsigned strip, double t0)=0
virtual const double BTimeCounter(double tleft, double tright, double z, unsigned id)=0
virtual const double EtfTimeMC2(double ADC, double TDC, double zHit, unsigned id, unsigned strip, double t0)=0
virtual const double BSigmaCluster(double zHit1, double zHit2, unsigned id1, unsigned id2)=0
virtual const double EtfTime(double ADC1, double ADC2, double TDC1, double TDC2, unsigned int id, unsigned int strip, double t0)=0
virtual const double BPulseHeight(double ADC1, double ADC2, double zHit, double sint, unsigned int id)=0
virtual const double BSigma1(double zHit, unsigned id)=0
virtual const double ESigma(double rHit, unsigned id)=0
virtual const double EtfTimeMC1(double ADC, double TDC, double zHit, unsigned id, unsigned strip, double t0)=0

◆ setExtTrack()

void TofTrack::setExtTrack ( RecExtTrack * extTrack,
double costheta,
double p[5],
int kal[5],
double t0,
int t0Stat )

Definition at line 157 of file TofTrack.cxx.

157 {
158
159 m_estime = t0;
160 m_t0Stat = t0Stat;
161
162 int tofId1 = extTrack->tof1VolumeNumber();
163 int tofId2 = extTrack->tof2VolumeNumber();
164
165 int iExist1 = -1;
166 int iExist2 = -1;
167 if( tofId1<0 ) {
168 if( extTrack->tof1VolumeNumber(3)>-1 ) {
169 iExist1 = 3;
170 }
171 else if( extTrack->tof1VolumeNumber(1)>-1 ) {
172 iExist1 = 1;
173 }
174 else if( extTrack->tof1VolumeNumber(4)>-1 ) {
175 iExist1 = 4;
176 }
177 else if( extTrack->tof1VolumeNumber(0)>-1 ) {
178 iExist1 = 0;
179 }
180 if( iExist1!=-1 ) {
181 tofId1 = extTrack->tof1VolumeNumber(iExist1);
182 }
183 }
184 if( tofId2<0 ) {
185 if( extTrack->tof2VolumeNumber(3)>-1 ) {
186 iExist2 = 3;
187 }
188 else if( extTrack->tof2VolumeNumber(1)>-1 ) {
189 iExist2 = 1;
190 }
191 else if( extTrack->tof2VolumeNumber(4)>-1 ) {
192 iExist2 = 4;
193 }
194 else if( extTrack->tof2VolumeNumber(0)>-1 ) {
195 iExist2 = 0;
196 }
197 if( iExist2!=-1 ) {
198 tofId2 = extTrack->tof2VolumeNumber(iExist2);
199 }
200 }
201 if( iExist1 == -1 ) { iExist1 = 2; }
202 if( iExist2 == -1 ) { iExist2 = 2; }
203
204 m_barrel = 3;
205 if( tofId1>=0 && tofId1<=87 ) {
206 m_id1 = tofId1;
207 m_barrel = 1;
208 m_hitCase = InnerLayer;
209 }
210 else if( tofId1>=176 && tofId1<=223 ) {
211 m_id1 = tofId1 - 176 + 48;
212 m_barrel = 2;
213 m_hitCase = WestEndcap;
214 if( costheta>0.0 ) {
215 m_id1 = -9;
216 m_barrel = 3;
217 m_hitCase = NoHit;
218 }
219 }
220 else if( tofId1>=224 && tofId1<=271 ) {
221 m_id1 = tofId1 - 176 - 48;
222 m_barrel = 0;
223 m_hitCase = EastEndcap;
224 if( costheta<0.0 ) {
225 m_id1 = -9;
226 m_barrel = 3;
227 m_hitCase = NoHit;
228 }
229 }
230 else if( tofId1>=272 && tofId1<=1135 ) {
231 m_id1 = tofId1 - 176 - 96;
232 m_istrip1 = m_id1%12;
233 m_id1 = m_id1/12;
234 if( tofId1>=272 && tofId1<=703 ) {
235 if( costheta>0.0 ) {
236 m_barrel = 4;
237 m_hitCase = EastEndcapMRPC;
238 }
239 else {
240 m_id1 = -9;
241 m_istrip1 = -9;
242 }
243 }
244 else if( tofId1>=704 && tofId1<=1135 ) {
245 if( costheta<0.0 ) {
246 m_barrel = 5;
247 m_hitCase = WestEndcapMRPC;
248 }
249 else {
250 m_id1 = -9;
251 m_istrip1 = -9;
252 }
253 }
254 }
255 else {
256 m_barrel = 3;
257 m_hitCase = NoHit;
258 }
259
260 if( tofId2>=88 && tofId2<=175 ) {
261 m_id2 = tofId2;
262 m_barrel = 1;
263 if( m_hitCase == InnerLayer ) {
264 m_hitCase = DoubleLayer;
265 }
266 else if( m_hitCase==NoHit ) {
267 m_hitCase = OuterLayer;
268 }
269 }
270 else if( tofId2>=272 && tofId2<=1135 ) {
271 m_id2 = tofId2 - 176 - 96;
272 m_istrip2 = m_id2%12;
273 m_id2 = m_id2/12;
274 if( m_hitCase==EastEndcapMRPC || m_hitCase==WestEndcapMRPC || m_hitCase==NoHit ) {
275 if( tofId2>=272 && tofId2<=703 ) {
276 if( costheta>0.0 ) {
277 m_barrel = 4;
278 m_hitCase = EastEndcapMRPC;
279 }
280 else {
281 m_id2 = -9;
282 m_istrip2 = -9;
283 }
284 }
285 else if( tofId2>=704 && tofId2<=1135 ) {
286 if( costheta<0.0 ) {
287 m_barrel = 5;
288 m_hitCase = WestEndcapMRPC;
289 }
290 else {
291 m_id2 = -9;
292 m_istrip2 = -9;
293 }
294 }
295 }
296 }
297
298 m_trackId = extTrack->trackId();
299
300 m_momentum = extTrack->tof1Momentum().r();
301 if( !( tofId1>=272 && tofId1<=1135 ) && ( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) ){
302 m_momentum = extTrack->tof2Momentum().r();
303 }
304
305 if( m_hitCase == InnerLayer || m_hitCase == DoubleLayer || m_hitCase == EastEndcap || m_hitCase == WestEndcap || ( ( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) && ( tofId1>=272 && tofId1<=1135 ) ) ) {
306 m_path1 = extTrack->tof1Path(iExist1);
307 m_theta1 = extTrack->tof1Momentum(iExist1).rho()/extTrack->tof1Momentum(iExist1).r();
308 m_phi1 = extTrack->tof1Position(iExist1).phi();
309 for( unsigned int i=0; i<5; i++ ) {
310 m_texpInner[i] = extTrack->tof1(i);
311 if( fabs(m_texpInner[i]+99.0)<1.0e-6 ) {
312 double beta = p[iExist1]/sqrt(p[iExist1]*p[iExist1]+mass[iExist1]*mass[iExist1]);
313 double betaNew = p[i]/sqrt(p[i]*p[i]+mass[i]*mass[i]);
314 m_texpInner[i] = beta*extTrack->tof1(iExist1)/betaNew;
315 }
316 }
317 if( m_hitCase == InnerLayer || m_hitCase == DoubleLayer ) {
318 m_xhit1 = extTrack->tof1Position(iExist1).x();
319 m_yhit1 = extTrack->tof1Position(iExist1).y();
320 m_zrhit1 = extTrack->tof1Position(iExist1).z();
321 m_errzr1 = extTrack->tof1PosSigmaAlongZ(iExist1);
322 for( unsigned int i=0; i<5; i++ ) {
323 m_zr1[i] = extTrack->tof1Position(i).z();
324 m_ezr1[i] = extTrack->tof1PosSigmaAlongZ(i);
325 if( fabs(m_zr1[i]+99.0)<1.0e-6 ) {
326 m_zr1[i] = m_zrhit1;
327 m_ezr1[i] = m_errzr1;
328 }
329 }
330 }
331 else if( m_hitCase == EastEndcap || m_hitCase == WestEndcap ) {
332 m_xhit1 = extTrack->tof1Position(iExist1).x();
333 m_yhit1 = extTrack->tof1Position(iExist1).y();
334 m_zrhit1 = extTrack->tof1Position(iExist1).rho();
335 m_errzr1 = sqrt( extTrack->tof1PosSigmaAlongX(iExist1)*extTrack->tof1PosSigmaAlongX(iExist1) + extTrack->tof1PosSigmaAlongY(iExist1)*extTrack->tof1PosSigmaAlongY(iExist1) );
336 for( unsigned int i=0; i<5; i++ ) {
337 m_zr1[i] = extTrack->tof1Position(i).rho();
338 m_ezr1[i] = sqrt( extTrack->tof1PosSigmaAlongX(i)*extTrack->tof1PosSigmaAlongX(i) + extTrack->tof1PosSigmaAlongY(i)*extTrack->tof1PosSigmaAlongY(i) );
339 if( fabs(m_zr1[i]+99.0)<1.0e-6 ) {
340 m_zr1[i] = m_zrhit1;
341 m_ezr1[i] = m_errzr1;
342 }
343 }
344 }
345 else if( ( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) && ( tofId1>=272 && tofId1<=1135 ) ) {
346 m_xhit1 = extTrack->tof1Position(iExist1).x();
347 m_yhit1 = extTrack->tof1Position(iExist1).z();
348 m_zrhit1 = extTrack->tof1Position(iExist1).x();
349 m_errzr1 = extTrack->tof1PosSigmaAlongX(iExist1);
350 for( unsigned int i=0; i<5; i++ ) {
351 m_zr1[i] = extTrack->tof1Position(i).x();
352 m_ezr1[i] = extTrack->tof1PosSigmaAlongX(i);
353 if( fabs(m_zr1[i]+99.0)<1.0e-6 ) {
354 m_zr1[i] = m_zrhit1;
355 m_ezr1[i] = m_errzr1;
356 }
357 }
358 if( !( tofId2>=272 && tofId2<=1135 ) ) {
359 for( unsigned int i=0; i<5; i++ ) {
360 m_texpOuter[i] = m_texpInner[i] * 136.573/133.673;
361 }
362 }
363 }
364 }
365
366 if( m_hitCase == OuterLayer || m_hitCase == DoubleLayer || ( ( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) && ( tofId2>=272 && tofId2<=1135 ) ) ) {
367 m_path2 = extTrack->tof2Path(iExist2);
368 m_theta2 = extTrack->tof2Momentum(iExist2).rho()/extTrack->tof2Momentum(iExist2).r();
369 m_phi2 = extTrack->tof2Position(iExist2).phi();
370 for( unsigned int i=0; i<5; i++ ) {
371 m_texpOuter[i] = extTrack->tof2(i);
372 if( fabs(m_texpOuter[i]+99.0)<1.0e-6 ) {
373 double beta = p[iExist2]/sqrt(p[iExist2]*p[iExist2]+mass[iExist2]*mass[iExist2]);
374 double betaNew = p[i]/sqrt(p[i]*p[i]+mass[i]*mass[i]);
375 m_texpOuter[i] = beta*extTrack->tof2(iExist2)/betaNew;
376 }
377 }
378 if( m_hitCase == OuterLayer || m_hitCase == DoubleLayer ) {
379 m_xhit2 = extTrack->tof2Position(iExist2).x();
380 m_yhit2 = extTrack->tof2Position(iExist2).y();
381 m_zrhit2 = extTrack->tof2Position(iExist2).z();
382 m_errzr2 = extTrack->tof2PosSigmaAlongZ(iExist2);
383 for( unsigned int i=0; i<5; i++ ) {
384 m_zr2[i] = extTrack->tof2Position(i).z();
385 m_ezr2[i] = extTrack->tof2PosSigmaAlongZ(i);
386 if( fabs(m_zr2[i]+99.0)<1.0e-6 ) {
387 m_zr2[i] = m_zrhit2;
388 m_ezr2[i] = m_errzr2;
389 }
390 }
391 }
392 else if( ( m_hitCase == EastEndcapMRPC || m_hitCase == WestEndcapMRPC ) && ( tofId2>=272 && tofId2<=1135 ) ) {
393 m_xhit2 = extTrack->tof2Position(iExist2).x();
394 m_yhit2 = extTrack->tof2Position(iExist2).z();
395 m_zrhit2 = extTrack->tof2Position(iExist2).x();
396 m_errzr2 = extTrack->tof2PosSigmaAlongX(iExist2);
397 for( unsigned int i=0; i<5; i++ ) {
398 m_zr2[i] = extTrack->tof2Position(i).x();
399 m_ezr2[i] = extTrack->tof2PosSigmaAlongX(i);
400 if( fabs(m_zr2[i]+99.0)<1.0e-6 ) {
401 m_zr2[i] = m_zrhit2;
402 m_ezr2[i] = m_errzr2;
403 }
404 }
405 if( !( tofId1>=272 && tofId1<=1135 ) ) {
406 for( unsigned int i=0; i<5; i++ ) {
407 m_texpInner[i] = m_texpOuter[i] * 133.673/136.573;
408 }
409 }
410 }
411 }
412
413 if( m_hitCase == NoHit ) { m_quality = 11; }
414
415 for( unsigned int i=0; i<5; i++ ) {
416 m_kal[i] = kal[i];
417 }
418
419 return;
420}
const double mass[5]
Definition TofTrack.h:17
const double tof1Path() const
Definition DstExtTrack.h:68
const Hep3Vector tof1Position() const
Definition DstExtTrack.h:58
const int tof1VolumeNumber() const
Definition DstExtTrack.h:64
const double tof1() const
Definition DstExtTrack.h:66
const Hep3Vector tof2Momentum() const
Definition DstExtTrack.h:96
const Hep3Vector tof1Momentum() const
Definition DstExtTrack.h:60
const double tof2() const
const double tof2PosSigmaAlongZ() const
const double tof1PosSigmaAlongX() const
Definition DstExtTrack.h:74
const double tof2Path() const
const int trackId() const
Definition DstExtTrack.h:43
const double tof1PosSigmaAlongY() const
Definition DstExtTrack.h:76
const double tof2PosSigmaAlongX() const
const double tof1PosSigmaAlongZ() const
Definition DstExtTrack.h:70
const int tof2VolumeNumber() const
const Hep3Vector tof2Position() const
Definition DstExtTrack.h:94
int t0Stat() const
Definition TofTrack.h:115
int kal(unsigned int i) const
Definition TofTrack.h:49
int tofId1() const
Definition TofTrack.h:68
double p() const
Definition TofTrack.h:37
int tofId2() const
Definition TofTrack.h:69
float costheta

Referenced by TofRec::execute().

◆ setFlag()

void TofTrack::setFlag ( unsigned int flag)
inline

Definition at line 124 of file TofTrack.h.

124{ m_flag = flag; }
unsigned int flag() const
Definition TofTrack.h:117

◆ setQuality()

void TofTrack::setQuality ( int qual)
inline

Definition at line 121 of file TofTrack.h.

121{ m_quality = qual; }

◆ setQuality1()

void TofTrack::setQuality1 ( int qual1)
inline

Definition at line 122 of file TofTrack.h.

122{ m_quality1 = qual1; }

Referenced by getMultiHit().

◆ setQuality2()

void TofTrack::setQuality2 ( int qual2)
inline

Definition at line 123 of file TofTrack.h.

123{ m_quality2 = qual2; }

Referenced by getMultiHit().

◆ setRecTofTrack()

void TofTrack::setRecTofTrack ( RecTofTrack * track,
int layerorend )

Definition at line 2484 of file TofTrack.cxx.

2484 {
2485
2486 double toffset[6];
2487 for( unsigned int i=0; i<6; i++ ) {
2488 toffset[i] = 0.0;
2489 }
2490
2491 if( layerorend == 0 ) { // cluster or double layer hit
2492 track->setPh( m_ph );
2493 track->setTof( m_tof[0] );
2494 track->setSigmaElectron( m_sigma );
2495 for( unsigned int i=0; i<5; i++ ) {
2496 toffset[i] = m_tof[0] - m_tof[i];
2497 }
2498 track->setToffset( toffset );
2499 track->setBeta( m_path/m_tof[0]/30.0 );
2500 }
2501 else if( layerorend == 1 ) { // inner layer
2502 track->setPh( m_ph1 );
2503 track->setTof( m_tof1[0] );
2504 track->setSigmaElectron( m_sigma1 );
2505 for( unsigned int i=0; i<5; i++ ) {
2506 toffset[i] = m_tof1[0] - m_tof1[i];
2507 }
2508 track->setToffset( toffset );
2509 track->setBeta( m_path1/m_tof1[0]/30.0 );
2510 }
2511 else if( layerorend == 2 ) { // outer layer
2512 track->setPh( m_ph2 );
2513 track->setTof( m_tof2[0] );
2514 track->setSigmaElectron( m_sigma2 );
2515 for( unsigned int i=0; i<5; i++ ) {
2516 toffset[i] = m_tof2[0] - m_tof2[i];
2517 }
2518 track->setToffset( toffset );
2519 track->setBeta( m_path2/m_tof2[0]/30.0 );
2520 }
2521 else if( layerorend == 11 ) { // inner layer east end readout
2522 track->setPh( m_ph11 );
2523 track->setTof( m_tof11[0] );
2524 track->setSigmaElectron( m_sigma11 );
2525 for( unsigned int i=0; i<5; i++ ) {
2526 toffset[i] = m_tof11[0] - m_tof11[i];
2527 }
2528 track->setToffset( toffset );
2529 track->setBeta( m_path1/m_tof11[0]/30.0 );
2530 }
2531 else if( layerorend == 12 ) { // inner layer west end readout
2532 track->setPh( m_ph12 );
2533 track->setTof( m_tof12[0] );
2534 track->setSigmaElectron( m_sigma12 );
2535 for( unsigned int i=0; i<5; i++ ) {
2536 toffset[i] = m_tof12[0] - m_tof12[i];
2537 }
2538 track->setToffset( toffset );
2539 track->setBeta( m_path1/m_tof12[0]/30.0 );
2540 }
2541 else if( layerorend == 21 ) { // outer layer east end readout
2542 track->setPh( m_ph21 );
2543 track->setTof( m_tof21[0] );
2544 track->setSigmaElectron( m_sigma21 );
2545 for( unsigned int i=0; i<5; i++ ) {
2546 toffset[i] = m_tof21[0] - m_tof21[i];
2547 }
2548 track->setToffset( toffset );
2549 track->setBeta( m_path2/m_tof21[0]/30.0 );
2550 }
2551 else if( layerorend == 22 ) { // outer layer west end readout
2552 track->setPh( m_ph22 );
2553 track->setTof( m_tof22[0] );
2554 track->setSigmaElectron( m_sigma22 );
2555 for( unsigned int i=0; i<5; i++ ) {
2556 toffset[i] = m_tof22[0] - m_tof22[i];
2557 }
2558 track->setToffset( toffset );
2559 track->setBeta( m_path2/m_tof22[0]/30.0 );
2560 }
2561 else{
2562 cout << "TofRec TofTrack::SetRecTofTrack layerorend = " << layerorend << endl;
2563 }
2564 return;
2565}
void setToffset(double toffset[6])
void setPh(double ph)
Definition DstTofTrack.h:97
void setTof(double tof)
Definition DstTofTrack.h:98
void setSigmaElectron(double se)

Referenced by buildRecTofTrack().

◆ setTofData()

void TofTrack::setTofData ( TofDataMap tofDataMap)

Definition at line 486 of file TofTrack.cxx.

486 {
487
488 if( m_hitCase == NoHit ) return;
489
490 unsigned int identify[11];
491 unsigned int count[11];
492 for( unsigned int i=0; i<11; i++ ) {
493 identify[i] = 0x0000c000;
494 count[i] = 0;
495 }
496 unsigned int countTot1 = 0;
497 unsigned int countTot2 = 0;
498
499 if( ( ( m_hitCase == InnerLayer ) || ( m_hitCase == DoubleLayer ) ) && ( m_id1 > -1 ) ) {
500 int tofid0 = m_id1;
501 identify[0] = TofID::getIntID( 1, 0, tofid0, 0 );
502 count[0] = tofDataMap.count( identify[0] );
503 int tofid1 = tofid0 - 1;
504 if( tofid1 == -1 ) { tofid1 = 87; }
505 identify[1] = TofID::getIntID( 1, 0, tofid1, 0 );
506 count[1] = tofDataMap.count( identify[1] );
507 int tofid2 = tofid0 + 1;
508 if( tofid2 == 88 ) { tofid2 = 0; }
509 identify[2] = TofID::getIntID( 1, 0, tofid2, 0 );
510 count[2] = tofDataMap.count( identify[2] );
511 }
512
513 if( ( ( m_hitCase == EastEndcap ) || ( m_hitCase == WestEndcap ) )&& ( m_id1 > -1 ) ) {
514 unsigned int whichEndcap = 0;
515 int tofid0 = m_id1;
516 if( m_hitCase == WestEndcap ) {
517 whichEndcap = 2;
518 tofid0 = m_id1 - 48;
519 }
520 identify[0] = TofID::getIntID( whichEndcap, 0, tofid0, 0 );
521 count[0] = tofDataMap.count( identify[0] );
522 int tofid1 = tofid0 - 1;
523 if( tofid1 == -1 ) { tofid1 = 47; }
524 identify[1] = TofID::getIntID( whichEndcap, 0, tofid1, 0 );
525 count[1] = tofDataMap.count( identify[1] );
526 int tofid2 = tofid0 + 1;
527 if( tofid2 == 48 ) { tofid2 = 0; }
528 identify[2] = TofID::getIntID( whichEndcap, 0, tofid2, 0 );
529 count[2] = tofDataMap.count( identify[2] );
530 }
531
532 if( ( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) && ( ( m_id1 > -1 ) && ( m_istrip1 > -1 ) ) ) {
533 IterTofDataMap iter = tofDataMap.begin();
534 for( ; iter != tofDataMap.end(); iter++ ) {
535 Identifier iden = TofID::cell_id( (*iter).first );
536 if( TofID::is_mrpc( iden ) ) {
537 TofData* tof = (*iter).second;
538 if( m_id1 == tof->tofId() && abs( m_istrip1 - tof->strip() )<=abs(m_delStrip1) ) {
539 m_delStrip1 = m_istrip1 - tof->strip();
540 }
541 if( ( abs( m_id1 - tof->tofId() )==1 || ( m_id1==0 && tof->tofId()==35 ) || ( m_id1==35 && tof->tofId()==0 ) || ( m_id1==36 && tof->tofId()==71 ) || ( m_id1==71 && tof->tofId()==36 ) ) && abs( m_istrip1 - tof->strip() )<=abs(m_delStrip1) ) {
542 m_delStrip2 = m_istrip1 - tof->strip();
543 }
544 }
545 }
546
547 unsigned int whichEndcap = 0;
548 int tofid0 = m_id1;
549 if( m_hitCase == WestEndcapMRPC ) {
550 whichEndcap = 1;
551 tofid0 = m_id1 - 36;
552 }
553 int strip0 = m_istrip1;
554 int strip1 = strip0 - 1;
555 int strip2 = strip0 + 1;
556 int strip3 = strip0 - 2;
557 int strip4 = strip0 + 2;
558 int tofid1 = tofid0 - 1;
559 if( tofid1 == -1 ) { tofid1 = 35; }
560 int tofid2 = tofid0 + 1;
561 if( tofid2 == 36 ) { tofid2 = 0; }
562
563 identify[0] = TofID::getIntID( 3, whichEndcap, tofid0, strip0, 0 );
564 count[0] = tofDataMap.count( identify[0] );
565 identify[5] = TofID::getIntID( 3, whichEndcap, tofid1, strip0, 0 );
566 count[5] = tofDataMap.count( identify[5] );
567 identify[6] = TofID::getIntID( 3, whichEndcap, tofid2, strip0, 0 );
568 count[6] = tofDataMap.count( identify[6] );
569
570 if( strip1 == -1 ) {
571 count[1] = 0;
572 count[7] = 0;
573 count[9] = 0;
574 }
575 else {
576 identify[1] = TofID::getIntID( 3, whichEndcap, tofid0, strip1, 0 );
577 count[1] = tofDataMap.count( identify[1] );
578 identify[7] = TofID::getIntID( 3, whichEndcap, tofid1, strip1, 0 );
579 count[7] = tofDataMap.count( identify[7] );
580 identify[9] = TofID::getIntID( 3, whichEndcap, tofid2, strip1, 0 );
581 count[9] = tofDataMap.count( identify[9] );
582 }
583
584 if( strip2 == 12 ) {
585 count[2] = 0;
586 count[8] = 0;
587 count[10] = 0;
588 }
589 else {
590 identify[2] = TofID::getIntID( 3, whichEndcap, tofid0, strip2, 0 );
591 count[2] = tofDataMap.count( identify[2] );
592 identify[8] = TofID::getIntID( 3, whichEndcap, tofid1, strip2, 0 );
593 count[8] = tofDataMap.count( identify[8] );
594 identify[10] = TofID::getIntID( 3, whichEndcap, tofid2, strip2, 0 );
595 count[10] = tofDataMap.count( identify[10] );
596 }
597 if( strip3 == -1 || strip3 == -2 ) { count[3] = 0; }
598 else {
599 identify[3] = TofID::getIntID( 3, whichEndcap, tofid0, strip3, 0 );
600 count[3] = tofDataMap.count( identify[3] );
601 }
602 if( strip4 == 12 || strip4 == 13 ) { count[4] = 0; }
603 else {
604 identify[4] = TofID::getIntID( 3, whichEndcap, tofid0, strip4, 0 );
605 count[4] = tofDataMap.count( identify[4] );
606 }
607 }
608
609 for( unsigned int i=0; i<11; i++ ) {
610 if( count[i] > 0 ) {
611 pair< IterTofDataMap, IterTofDataMap > range = tofDataMap.equal_range( identify[i] );
612 IterTofDataMap iter = range.first;
613 for( unsigned int j=0; j<count[i]; j++,iter++ ) {
614 if( i==0 ) {
615 tofDataAnalysis( (*iter).second, 1 );
616 }
617 else if( i==1 || i==2 ) {
618 tofDataAnalysis( (*iter).second, 2 );
619 }
620 else {
621 tofDataAnalysis( (*iter).second, 3 );
622 }
623 }
624 }
625 countTot1 = countTot1 + count[i];
626 }
627
628 if( countTot1 == 0 ) {
629 if( m_hitCase == DoubleLayer ) {
630 m_hitCase = OuterLayer;
631 }
632 else if( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) {
633 }
634 else {
635 m_hitCase = NoHit;
636 m_quality = 12;
637 }
638 }
639
640
641 for( unsigned int i=0; i<11; i++ ) {
642 identify[i] = 0x0000c000;
643 count[i] = 0;
644 }
645 if( ( ( m_hitCase == OuterLayer ) || ( m_hitCase == DoubleLayer ) ) && ( m_id2 > 87 ) ) {
646 int tofid0 = m_id2 - 88;
647 identify[0] = TofID::getIntID( 1, 1, tofid0, 0 );
648 count[0] = tofDataMap.count( identify[0] );
649 int tofid1 = tofid0 - 1;
650 if( tofid1 == -1 ) { tofid1 = 87; }
651 identify[1] = TofID::getIntID( 1, 1, tofid1, 0 );
652 count[1] = tofDataMap.count( identify[1] );
653 int tofid2 = tofid0 + 1;
654 if( tofid2 == 88 ) { tofid2 = 0; }
655 identify[2] = TofID::getIntID( 1, 1, tofid2, 0 );
656 count[2] = tofDataMap.count( identify[2] );
657
658 for( unsigned int i=0; i<3; i++ ) {
659 if( count[i] > 0 ) {
660 pair< IterTofDataMap, IterTofDataMap > range = tofDataMap.equal_range( identify[i] );
661 IterTofDataMap iter = range.first;
662 for( unsigned int j=0; j<count[i]; j++,iter++ ) {
663 if( i==0 ) {
664 tofDataAnalysis( (*iter).second, 3 );
665 }
666 else {
667 tofDataAnalysis( (*iter).second, 4 );
668 }
669 }
670 }
671 countTot2 = countTot2 + count[i];
672 }
673
674 if( countTot2 == 0 ) {
675 if( m_hitCase != DoubleLayer ) {
676 m_hitCase = NoHit;
677 m_quality = 12;
678 }
679 else {
680 m_hitCase = InnerLayer;
681 }
682 }
683 }
684
685 if( ( ( m_hitCase == EastEndcapMRPC ) || ( m_hitCase == WestEndcapMRPC ) ) && ( ( m_id2 > -1 ) && ( m_istrip2 > -1 ) ) ) {
686 IterTofDataMap iter = tofDataMap.begin();
687 for( ; iter != tofDataMap.end(); iter++ ) {
688 Identifier iden = TofID::cell_id( (*iter).first );
689 if( TofID::is_mrpc( iden ) ) {
690 TofData* tof = (*iter).second;
691 if( m_id2 == tof->tofId() && abs( m_istrip2 - tof->strip() )<=abs(m_delStrip1) ) {
692 m_delStrip1 = m_istrip2 - tof->strip();
693 }
694 if( ( abs( m_id2 - tof->tofId() )==1 || ( m_id2==0 && tof->tofId()==35 ) || ( m_id2==35 && tof->tofId()==0 ) || ( m_id2==36 && tof->tofId()==71 ) || ( m_id2==71 && tof->tofId()==36 ) ) && abs( m_istrip2 - tof->strip() )<=abs(m_delStrip2) ) {
695 m_delStrip2 = m_istrip2 - tof->strip();
696 }
697 }
698 }
699
700 unsigned int whichEndcap = 0;
701 int tofid0 = m_id2;
702 if( m_hitCase == WestEndcapMRPC ) {
703 whichEndcap = 1;
704 tofid0 = m_id2 - 36;
705 }
706 int strip0 = m_istrip2;
707 int strip1 = strip0 - 1;
708 int strip2 = strip0 + 1;
709 int strip3 = strip0 - 2;
710 int strip4 = strip0 + 2;
711 int tofid1 = tofid0 - 1;
712 if( tofid1 == -1 ) { tofid1 = 35; }
713 int tofid2 = tofid0 + 1;
714 if( tofid2 == 36 ) { tofid2 = 0; }
715
716 identify[0] = TofID::getIntID( 3, whichEndcap, tofid0, strip0, 0 );
717 count[0] = tofDataMap.count( identify[0] );
718 identify[5] = TofID::getIntID( 3, whichEndcap, tofid1, strip0, 0 );
719 count[5] = tofDataMap.count( identify[5] );
720 identify[6] = TofID::getIntID( 3, whichEndcap, tofid2, strip0, 0 );
721 count[6] = tofDataMap.count( identify[6] );
722
723 if( strip1 == -1 ) {
724 count[1] = 0;
725 count[7] = 0;
726 count[9] = 0;
727 }
728 else {
729 identify[1] = TofID::getIntID( 3, whichEndcap, tofid0, strip1, 0 );
730 count[1] = tofDataMap.count( identify[1] );
731 identify[7] = TofID::getIntID( 3, whichEndcap, tofid1, strip1, 0 );
732 count[7] = tofDataMap.count( identify[7] );
733 identify[9] = TofID::getIntID( 3, whichEndcap, tofid2, strip1, 0 );
734 count[9] = tofDataMap.count( identify[9] );
735 }
736
737 if( strip2 == 12 ) {
738 count[2] = 0;
739 count[8] = 0;
740 count[10] = 0;
741 }
742 else {
743 identify[2] = TofID::getIntID( 3, whichEndcap, tofid0, strip2, 0 );
744 count[2] = tofDataMap.count( identify[2] );
745 identify[8] = TofID::getIntID( 3, whichEndcap, tofid1, strip2, 0 );
746 count[8] = tofDataMap.count( identify[8] );
747 identify[10] = TofID::getIntID( 3, whichEndcap, tofid2, strip2, 0 );
748 count[10] = tofDataMap.count( identify[10] );
749 }
750 if( strip3 == -1 || strip3 == -2 ) { count[3] = 0; }
751 else {
752 identify[3] = TofID::getIntID( 3, whichEndcap, tofid0, strip3, 0 );
753 count[3] = tofDataMap.count( identify[3] );
754 }
755 if( strip4 == 12 || strip4 == 13 ) { count[4] = 0; }
756 else {
757 identify[4] = TofID::getIntID( 3, whichEndcap, tofid0, strip4, 0 );
758 count[4] = tofDataMap.count( identify[4] );
759 }
760
761 for( unsigned int i=0; i<11; i++ ) {
762 if( count[i] > 0 ) {
763 pair< IterTofDataMap, IterTofDataMap > range = tofDataMap.equal_range( identify[i] );
764 IterTofDataMap iter = range.first;
765 for( unsigned int j=0; j<count[i]; j++,iter++ ) {
766 if( i==0 ) {
767 tofDataAnalysis( (*iter).second, 4 );
768 }
769 else if( i==1 || i==2 ) {
770 tofDataAnalysis( (*iter).second, 5 );
771 }
772 else {
773 tofDataAnalysis( (*iter).second, 6 );
774 }
775 }
776 }
777 countTot2 = countTot2 + count[i];
778 }
779 }
780
781 if( countTot1==0 && countTot2==0 ) {
782 m_hitCase = NoHit;
783 m_quality = 12;
784 }
785
786 return;
787}
DOUBLE_PRECISION count[3]
std::multimap< unsignedint, TofData * >::iterator IterTofDataMap
Definition TofData.h:245
static Identifier cell_id(int barrel_ec, int layer, int phi_module, int end)
For a single crystal.
Definition TofID.cxx:143
static bool is_mrpc(const Identifier &id)
Definition TofID.cxx:113
static value_type getIntID(int barrel_ec, int layer, int phi_module, int end)
Definition TofID.cxx:178
void tofDataAnalysis(TofData *tof, unsigned int iflag)
Definition TofTrack.cxx:793

◆ size1()

int TofTrack::size1 ( ) const
inline

Definition at line 56 of file TofTrack.h.

56{ return m_tofData1.size(); }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ size2()

int TofTrack::size2 ( ) const
inline

Definition at line 57 of file TofTrack.h.

57{ return m_tofData2.size(); }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ size3()

int TofTrack::size3 ( ) const
inline

Definition at line 58 of file TofTrack.h.

58{ return m_tofData3.size(); }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ size4()

int TofTrack::size4 ( ) const
inline

Definition at line 59 of file TofTrack.h.

59{ return m_tofData4.size(); }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ strip1()

int TofTrack::strip1 ( ) const
inline

◆ strip2()

int TofTrack::strip2 ( ) const
inline

◆ t0Stat()

int TofTrack::t0Stat ( ) const
inline

Definition at line 115 of file TofTrack.h.

115{ return m_t0Stat; }

Referenced by setExtTrack().

◆ tdc1()

double TofTrack::tdc1 ( ) const
inline

◆ tdc2()

double TofTrack::tdc2 ( ) const
inline

◆ tdc3()

double TofTrack::tdc3 ( ) const
inline

◆ tdc4()

double TofTrack::tdc4 ( ) const
inline

◆ tdiff1()

double TofTrack::tdiff1 ( ) const
inline

◆ tdiff2()

double TofTrack::tdiff2 ( ) const
inline

◆ texp()

double TofTrack::texp ( unsigned int i) const
inline

◆ texpInner()

double TofTrack::texpInner ( unsigned int i) const
inline

Definition at line 102 of file TofTrack.h.

102{ return m_texpInner[i]; }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ texpOuter()

double TofTrack::texpOuter ( unsigned int i) const
inline

Definition at line 103 of file TofTrack.h.

103{ return m_texpOuter[i]; }

Referenced by TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), TofCheckDigi::Fill_TofTrack(), and TofCheckDigi::Fill_TofTrack().

◆ theta1()

double TofTrack::theta1 ( ) const
inline

◆ theta2()

double TofTrack::theta2 ( ) const
inline

◆ tof()

◆ tof1()

double TofTrack::tof1 ( unsigned int i) const
inline

◆ tof11()

double TofTrack::tof11 ( unsigned int i) const
inline

◆ tof12()

double TofTrack::tof12 ( unsigned int i) const
inline

◆ tof2()

double TofTrack::tof2 ( unsigned int i) const
inline

◆ tof21()

double TofTrack::tof21 ( unsigned int i) const
inline

◆ tof22()

double TofTrack::tof22 ( unsigned int i) const
inline

◆ tofData1()

std::vector< TofData * > TofTrack::tofData1 ( ) const
inline

Definition at line 54 of file TofTrack.h.

54{ return m_tofData1; }

Referenced by compareTofData(), compareTofDataEndcap(), findTofDataBarrel(), and findTofDataEndcap().

◆ tofData2()

std::vector< TofData * > TofTrack::tofData2 ( ) const
inline

Definition at line 55 of file TofTrack.h.

55{ return m_tofData2; }

Referenced by compareTofData(), compareTofDataEndcap(), findTofDataBarrel(), and findTofDataEndcap().

◆ tofDataAnalysis()

void TofTrack::tofDataAnalysis ( TofData * tof,
unsigned int iflag )

Definition at line 793 of file TofTrack.cxx.

793 {
794
795 unsigned int qual = tof->quality();
796
797 if( ( qual & 0x10 ) == 0 ) {
798 qual = ( qual | 0x10 ); // zadc, ztdc unmatched, and track matched
799 if( tof->barrel() || tof->is_mrpc() ) { // Barrel, Endcap has been done
800 if( ( tof->quality() == 0x7 ) || ( tof->quality() == 0xd ) ) {
801 qual = ( qual | 0x20 ); // lost one Q
802 }
803
804 if( ( tof->quality() == 0xb ) || ( tof->quality() == 0xe ) ) {
805 qual = ( qual | 0x40 ); // lost one T
806 }
807
808 if( ( tof->quality() == 0x3 ) || ( tof->quality() == 0xc ) ) {
809 qual = ( qual | 0x80 ); // single end
810 }
811 }
812 if( tof->barrel() ) {
813 if( ( tof->quality() & 0x5 ) == 0x5 ) {
814 double ztdc = tofCaliSvc->ZTDC( tof->tdc1(), tof->tdc2(), tof->tofId() );
815 tof->setZTdc( ztdc );
816 }
817
818 if( ( tof->quality() & 0xa ) == 0xa ) {
819 double zadc = tofCaliSvc->ZADC( tof->adc1(), tof->adc2(), tof->tofId() );
820 tof->setZAdc( zadc );
821 }
822 }
823 if( tof->is_mrpc() ) {
824 if( ( tof->quality() & 0x5 ) == 0x5 ) {
825 double ztdc = tofCaliSvc->EtfZTDC( tof->tdc1(), tof->tdc2(), tof->tofId(), tof->strip() );
826 tof->setZTdc( ztdc );
827 }
828 }
829 tof->setQuality( qual );
830 }
831
832 if( iflag == 1 ) { m_tofData1.push_back( tof ); }
833 else if( iflag == 2 ) { m_tofData2.push_back( tof ); }
834 else if( iflag == 3 ) { m_tofData3.push_back( tof ); }
835 else if( iflag == 4 ) { m_tofData4.push_back( tof ); }
836 else if( iflag == 5 ) { m_tofData5.push_back( tof ); }
837 else if( iflag == 6 ) { m_tofData6.push_back( tof ); }
838 else {
839 cout << "TofRec::TofTrack::TofDataAnalylsis: the Flag should be 1-4, out of the Range!" << endl;
840 }
841
842 return;
843}
virtual const double ZTDC(double tleft, double tright, unsigned id)=0
virtual const double ZADC(double qleft, double qright, unsigned id)=0
virtual const double EtfZTDC(double tleft, double tright, unsigned int id, unsigned int strip)=0

Referenced by setTofData().

◆ tofDataStudy()

void TofTrack::tofDataStudy ( )

◆ tofId1()

int TofTrack::tofId1 ( ) const
inline

◆ tofId2()

int TofTrack::tofId2 ( ) const
inline

◆ tofTrackId()

int TofTrack::tofTrackId ( ) const
inline

◆ trackId()

int TofTrack::trackId ( ) const
inline

◆ xhit1()

double TofTrack::xhit1 ( ) const
inline

◆ xhit2()

double TofTrack::xhit2 ( ) const
inline

◆ yhit1()

double TofTrack::yhit1 ( ) const
inline

◆ yhit2()

double TofTrack::yhit2 ( ) const
inline

◆ zadc1()

double TofTrack::zadc1 ( ) const
inline

◆ zadc2()

double TofTrack::zadc2 ( ) const
inline

◆ zr1()

double TofTrack::zr1 ( unsigned int i) const
inline

◆ zr2()

double TofTrack::zr2 ( unsigned int i) const
inline

◆ zrhit1()

double TofTrack::zrhit1 ( ) const
inline

◆ zrhit2()

double TofTrack::zrhit2 ( ) const
inline

◆ ztdc1()

double TofTrack::ztdc1 ( ) const
inline

◆ ztdc2()

double TofTrack::ztdc2 ( ) const
inline

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