Newer
Older
Double_t tdc_pos = hit->GetPosTDC();
if(tdc_pos >=fScinTdcMin && tdc_pos <= fScinTdcMax ) {
Double_t adc_pos = hit->GetPosADC();
Double_t pathp = fPlanes[ip]->GetPosLeft() - scinLongCoord;
fTOFPInfo[ihhit].pathp = pathp;
Double_t timep = tdc_pos*fScinTdcToTime;
if(fTofUsingInvAdc) {
timep -= fHodoPosInvAdcOffset[fPIndex]
+ pathp/fHodoPosInvAdcLinear[fPIndex]
+ fHodoPosInvAdcAdc[fPIndex]
/TMath::Sqrt(TMath::Max(20.0,adc_pos));
} else {
timep -= fHodoPosPhcCoeff[fPIndex]*
TMath::Sqrt(TMath::Max(0.0,adc_pos/fHodoPosMinPh[fPIndex]-1.0))
+ pathp/fHodoVelLight[fPIndex]
+ fHodoPosTimeOffset[fPIndex];
}
fTOFPInfo[ihhit].scin_pos_time = timep;
timep -= zcor;
fTOFPInfo[ihhit].time_pos = timep;
for ( Int_t k = 0; k < 200; k++ ){ // Line 211
Double_t tmin = 0.5 * ( k + 1 ) ;
if ( ( timep > tmin ) && ( timep < ( tmin + fTofTolerance ) ) )
timehist[k] ++;
}
Double_t tdc_neg = hit->GetNegTDC();
if(tdc_neg >=fScinTdcMin && tdc_neg <= fScinTdcMax ) {
Double_t adc_neg = hit->GetNegADC();
Double_t pathn = scinLongCoord - fPlanes[ip]->GetPosRight();
fTOFPInfo[ihhit].pathn = pathn;
Double_t timen = tdc_neg*fScinTdcToTime;
if(fTofUsingInvAdc) {
timen -= fHodoNegInvAdcOffset[fPIndex]
+ pathn/fHodoNegInvAdcLinear[fPIndex]
+ fHodoNegInvAdcAdc[fPIndex]
/TMath::Sqrt(TMath::Max(20.0,adc_neg));
} else {
timen -= fHodoNegPhcCoeff[fPIndex]*
TMath::Sqrt(TMath::Max(0.0,adc_neg/fHodoNegMinPh[fPIndex]-1.0))
+ pathn/fHodoVelLight[fPIndex]
+ fHodoNegTimeOffset[fPIndex];
}
fTOFPInfo[ihhit].scin_neg_time = timen;
fTOFPInfo[ihhit].time_neg = timen;
for ( Int_t k = 0; k < 200; k++ ){ // Line 230
Double_t tmin = 0.5 * ( k + 1 );
if ( ( timen > tmin ) && ( timen < ( tmin + fTofTolerance ) ) )
timehist[k] ++;
}
} // condition for cenetr on a paddle
} // First loop over hits in a plane <---------
//-----------------------------------------------------------------------------------------------
//------------- First large loop over scintillator hits ends here --------------------
//-----------------------------------------------------------------------------------------------
}
Int_t nhits=ihhit;
// Find bin with the most hits
Int_t jmax = 0; // Line 240
Int_t maxhit = 0;
for ( Int_t k = 0; k < 200; k++ ){
if ( timehist[k] > maxhit ){
jmax = k+1;
maxhit = timehist[k];
Double_t tmin = 0.5 * jmax;
for(Int_t ih = 0; ih < nhits; ih++) { // loop over all scintillator hits
if ( ( fTOFPInfo[ih].time_pos > tmin ) && ( fTOFPInfo[ih].time_pos < ( tmin + fTofTolerance ) ) ) {
fTOFPInfo[ih].keep_pos=kTRUE;
if ( ( fTOFPInfo[ih].time_neg > tmin ) && ( fTOFPInfo[ih].time_neg < ( tmin + fTofTolerance ) ) ){
fTOFPInfo[ih].keep_neg=kTRUE;
//---------------------------------------------------------------------------------------------
// ---------------------- Scond loop over scint. hits in a plane ------------------------------
//---------------------------------------------------------------------------------------------
for(Int_t ih=0; ih < nhits; ih++) {
THcHodoHit *hit = fTOFPInfo[ih].hit;
Int_t iphit = fTOFPInfo[ih].hitNumInPlane;
Int_t ip = fTOFPInfo[ih].planeIndex;
// fDumpOut << " looping over hits = " << ih << " plane = " << ip+1 << endl;
GoodFlags flags;
// Flags are used by THcHodoEff
fGoodFlags[itrack][ip].push_back(flags);
fGoodFlags[itrack][ip][iphit].onTrack = kFALSE;
fGoodFlags[itrack][ip][iphit].goodScinTime = kFALSE;
fGoodFlags[itrack][ip][iphit].goodTdcNeg = kFALSE;
fGoodFlags[itrack][ip][iphit].goodTdcPos = kFALSE;
fTOFCalc.push_back(TOFCalc());
// Do we set back to false for each track, or just once per event?
fTOFCalc[ih].good_scin_time = kFALSE;
// These need a track index too to calculate efficiencies
fTOFCalc[ih].good_tdc_pos = kFALSE;
fTOFCalc[ih].good_tdc_neg = kFALSE;
fTOFCalc[ih].pindex = ip;
Int_t paddle = hit->GetPaddleNumber()-1;
fTOFCalc[ih].hit_paddle = paddle;
fTOFCalc[ih].good_raw_pad = paddle;
// Double_t scinCenter = fPlanes[ip]->GetPosCenter(paddle) + fPlanes[ip]->GetPosOffset();
// Double_t scinTrnsCoord = fTOFPInfo[ih].scinTrnsCoord;
// Double_t scinLongCoord = fTOFPInfo[ih].scinLongCoord;
Int_t fPIndex = GetScinIndex(ip,paddle);
if (fTOFPInfo[ih].onTrack) {
if ( fTOFPInfo[ih].keep_pos ) { // 301
fTOFCalc[ih].good_tdc_pos = kTRUE;
fGoodFlags[itrack][ip][iphit].goodTdcPos = kTRUE;
if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
fDumpOut << fixed << setprecision(2);
fDumpOut << showpoint << " 1" << setw(3) << ip+1 << setw(3) << hit->GetPaddleNumber() << setw(10) << hit->GetPosTDC()*fScinTdcToTime << setw(10) << fTOFPInfo[ih].pathp << setw(10) << fTOFPInfo[ih].zcor << setw(10) << fTOFPInfo[ih].time_pos << setw(10) << hit->GetPosADC() << endl;
if ( fTOFPInfo[ih].keep_neg ) { //
fTOFCalc[ih].good_tdc_neg = kTRUE;
fGoodFlags[itrack][ip][iphit].goodTdcNeg = kTRUE;
if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
fDumpOut << fixed << setprecision(2);
fDumpOut << showpoint << " 2" << setw(3) << ip+1 << setw(3) << hit->GetPaddleNumber() << setw(10) << hit->GetNegTDC()*fScinTdcToTime << setw(10) << fTOFPInfo[ih].pathn << setw(10) << fTOFPInfo[ih].zcor << setw(10) << fTOFPInfo[ih].time_neg << setw(10) << hit->GetNegADC() << endl;
// ** Calculate ave time for scin and error.
if ( fTOFCalc[ih].good_tdc_pos ){
if ( fTOFCalc[ih].good_tdc_neg ){
fTOFCalc[ih].scin_time = ( fTOFPInfo[ih].scin_pos_time +
fTOFPInfo[ih].scin_neg_time ) / 2.;
fTOFCalc[ih].scin_time_fp = ( fTOFPInfo[ih].time_pos +
fTOFPInfo[ih].time_neg ) / 2.;
fTOFCalc[ih].scin_sigma = TMath::Sqrt( fHodoPosSigma[fPIndex] * fHodoPosSigma[fPIndex] +
fHodoNegSigma[fPIndex] * fHodoNegSigma[fPIndex] )/2.;
fTOFCalc[ih].good_scin_time = kTRUE;
fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE;
} else{
fTOFCalc[ih].scin_time = fTOFPInfo[ih].scin_pos_time;
fTOFCalc[ih].scin_time_fp = fTOFPInfo[ih].time_pos;
fTOFCalc[ih].scin_sigma = fHodoPosSigma[fPIndex];
fTOFCalc[ih].good_scin_time = kTRUE;
fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE;
if ( fTOFCalc[ih].good_tdc_neg ){
fTOFCalc[ih].scin_time = fTOFPInfo[ih].scin_neg_time;
fTOFCalc[ih].scin_time_fp = fTOFPInfo[ih].time_neg;
fTOFCalc[ih].scin_sigma = fHodoNegSigma[fPIndex];
fTOFCalc[ih].good_scin_time = kTRUE;
fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE;
}
} // In h_tof.f this includes the following if condition for time at focal plane
// // because it is written in FORTRAN code
// c Get time at focal plane
if ( fTOFCalc[ih].good_scin_time ){
// scin_time_fp doesn't need to be an array
// Is this any different than the average of time_pos and time_neg?
// Double_t scin_time_fp = ( fTOFPInfo[ih].time_pos +
// fTOFPInfo[ih].time_neg ) / 2.;
Double_t scin_time_fp = fTOFCalc[ih].scin_time_fp;
sumFPTime = sumFPTime + scin_time_fp;
nFPTime ++;
fSumPlaneTime[ip] = fSumPlaneTime[ip] + scin_time_fp;
fNPlaneTime[ip] ++;
fNScinHit[itrack] ++;
if ( ( fTOFCalc[ih].good_tdc_pos ) && ( fTOFCalc[ih].good_tdc_neg ) ){
nPmtHit[itrack] = nPmtHit[itrack] + 2;
} else {
nPmtHit[itrack] = nPmtHit[itrack] + 1;
}
fdEdX[itrack].push_back(0.0);
// --------------------------------------------------------------------------------------------
if ( fTOFCalc[ih].good_tdc_pos ){
if ( fTOFCalc[ih].good_tdc_neg ){
fdEdX[itrack][fNScinHit[itrack]-1]=
TMath::Sqrt( TMath::Max( 0., hit->GetPosADC() * hit->GetNegADC() ) );
} else{
fdEdX[itrack][fNScinHit[itrack]-1]=
TMath::Max( 0., hit->GetPosADC() );
if ( fTOFCalc[ih].good_tdc_neg ){
fdEdX[itrack][fNScinHit[itrack]-1]=
TMath::Max( 0., hit->GetNegADC() );
} else{
fdEdX[itrack][fNScinHit[itrack]-1]=0.0;
}
// --------------------------------------------------------------------------------------------
} // time at focal plane condition
} // on track condition
// ** See if there are any good time measurements in the plane.
if ( fTOFCalc[ih].good_scin_time ){
fGoodPlaneTime[ip] = kTRUE;
fTOFCalc[ih].dedx = fdEdX[itrack][fNScinHit[itrack]-1];
fTOFCalc[ih].dedx = 0.0;
} // Second loop over hits of a scintillator plane ends here
theTrack->SetGoodPlane3( fGoodPlaneTime[2] ? 1 : 0 );
theTrack->SetGoodPlane4( fGoodPlaneTime[3] ? 1 : 0 );
//
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
// * * Fit beta if there are enough time measurements (one upper, one lower)
// From h_tof_fit
if ( ( ( fGoodPlaneTime[0] ) || ( fGoodPlaneTime[1] ) ) &&
( ( fGoodPlaneTime[2] ) || ( fGoodPlaneTime[3] ) ) ){
Double_t sumW = 0.;
Double_t sumT = 0.;
Double_t sumZ = 0.;
Double_t sumZZ = 0.;
Double_t sumTZ = 0.;
for(Int_t ih=0; ih < nhits; ih++) {
Int_t ip = fTOFPInfo[ih].planeIndex;
if ( fTOFCalc[ih].good_scin_time ) {
Double_t scinWeight = 1 / ( fTOFCalc[ih].scin_sigma * fTOFCalc[ih].scin_sigma );
Double_t zPosition = ( fPlanes[ip]->GetZpos()
+( fTOFCalc[ih].hit_paddle % 2 ) *
fPlanes[ip]->GetDzpos() );
sumW += scinWeight;
sumT += scinWeight * fTOFCalc[ih].scin_time;
sumZ += scinWeight * zPosition;
sumZZ += scinWeight * ( zPosition * zPosition );
sumTZ += scinWeight * zPosition * fTOFCalc[ih].scin_time;
} // condition of good scin time
} // loop over hits
Double_t tmp = sumW * sumZZ - sumZ * sumZ ;
Double_t t0 = ( sumT * sumZZ - sumZ * sumTZ ) / tmp ;
Double_t tmpDenom = sumW * sumTZ - sumZ * sumT;
if ( TMath::Abs( tmpDenom ) > ( 1 / 10000000000.0 ) ) {
betaChiSq = 0.;
for(Int_t ih=0; ih < nhits; ih++) {
Int_t ip = fTOFPInfo[ih].planeIndex;
if ( fTOFCalc[ih].good_scin_time ){
Double_t zPosition = ( fPlanes[ip]->GetZpos() + ( fTOFCalc[ih].hit_paddle % 2 ) *
fPlanes[ip]->GetDzpos() );
Double_t timeDif = ( fTOFCalc[ih].scin_time - t0 );
betaChiSq += ( ( zPosition / beta - timeDif ) *
( zPosition / beta - timeDif ) ) /
( fTOFCalc[ih].scin_sigma * fTOFCalc[ih].scin_sigma );
} // condition for good scin time
} // loop over hits
Double_t pathNorm = TMath::Sqrt( 1. + theTrack->GetTheta() * theTrack->GetTheta() +
theTrack->GetPhi() * theTrack->GetPhi() );
// Take angle into account
beta = beta / 29.979; // velocity / c
} // condition for fTmpDenom
} // else condition for fTmpDenom
timeAtFP[itrack] = ( sumFPTime / nFPTime );
//
// ---------------------------------------------------------------------------
Double_t FPTimeSum=0.0;
Int_t nFPTimeSum=0;
for (Int_t ip = 0; ip < temp_planes; ip++ ){
if ( fNPlaneTime[ip] != 0 ){
fFPTime[ip] = ( fSumPlaneTime[ip] / fNPlaneTime[ip] );
FPTimeSum += fSumPlaneTime[ip];
nFPTimeSum += fNPlaneTime[ip];
fFPTime[ip] = 1000. * ( ip + 1 );
fFPTimeAll = fptime;
Double_t dedx=0.0;
for(UInt_t ih=0;ih<fTOFCalc.size();ih++) {
if(fTOFCalc[ih].good_scin_time) {
dedx = fTOFCalc[ih].dedx;
break;
}
}
theTrack->SetDedx(dedx);
theTrack->SetFPTime(fptime);
theTrack->SetBeta(beta);
theTrack->SetBetaChi2( betaChiSq );
theTrack->SetNPMT(nPmtHit[itrack]);
theTrack->SetFPTime( timeAtFP[itrack]);
} // Main loop over tracks ends here.
} // If condition for at least one track
if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
fDumpOut << "0 " << endl;
}
//-----------------------------------------------------------------------
//
// Trnslation of h_track_tests.f file for tracking efficiency
//
//-----------------------------------------------------------------------
//************************now look at some hodoscope tests
// *second, we move the scintillators. here we use scintillator cuts to see
// *if a track should have been found.
for(Int_t ip = 0; ip < temp_planes; ip++ ) {
if (!fPlanes[ip])
return -1;
TClonesArray* hodoHits = fPlanes[ip]->GetHits();
// TClonesArray* scinPosTDC = fPlanes[ip]->GetPosTDC();
// TClonesArray* scinNegTDC = fPlanes[ip]->GetNegTDC();
fNScinHits[ip] = fPlanes[ip]->GetNScinHits();
for (Int_t iphit = 0; iphit < fNScinHits[ip]; iphit++ ){
Int_t paddle = ((THcHodoHit*)hodoHits->At(iphit))->GetPaddleNumber()-1;
fScinHitPaddle[ip][paddle] = 1;
}
}
// *next, look for clusters of hits in a scin plane. a cluster is a group of
// *adjacent scintillator hits separated by a non-firing scintillator.
// *Wwe count the number of three adjacent scintillators too. (A signle track
// *shouldn't fire three adjacent scintillators.
for(Int_t ip = 0; ip < temp_planes; ip++ ) {
1386
1387
1388
1389
1390
1391
1392
1393
1394
1395
1396
1397
1398
1399
1400
1401
1402
1403
1404
1405
1406
// Planes ip = 0 = 1X
// Planes ip = 2 = 2X
if (!fPlanes[ip]) return -1;
fNClust.push_back(0);
fThreeScin.push_back(0);
}
// *look for clusters in x planes... (16 scins) !this assume both x planes have same
// *number of scintillators.
Int_t icount;
for (Int_t ip = 0; ip < 3; ip +=2 ) {
icount = 0;
if ( fScinHitPaddle[ip][0] == 1 )
icount ++;
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[0] - 1; ipaddle++ ){
// !look for number of clusters of 1 or more hits
if ( ( fScinHitPaddle[ip][ipaddle] == 0 ) &&
( fScinHitPaddle[ip][ipaddle + 1] == 1 ) )
icount ++;
} // Loop over paddles
fNClust[ip] = icount;
icount = 0;
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[0] - 2; ipaddle++ ){
// !look for three or more adjacent hits
if ( ( fScinHitPaddle[ip][ipaddle] == 1 ) &&
( fScinHitPaddle[ip][ipaddle + 1] == 1 ) &&
( fScinHitPaddle[ip][ipaddle + 2] == 1 ) )
icount ++;
} // Second loop over paddles
if ( icount > 0 )
fThreeScin[ip] = 1;
} // Loop over X plane
// *look for clusters in y planes... (10 scins) !this assume both y planes have same
1428
1429
1430
1431
1432
1433
1434
1435
1436
1437
1438
1439
1440
1441
1442
1443
1444
1445
1446
1447
1448
1449
1450
1451
1452
1453
1454
1455
1456
1457
1458
1459
1460
1461
1462
1463
1464
1465
1466
1467
// *number of scintillators.
for (Int_t ip = 1; ip < 4; ip +=2 ) {
// Planes ip = 1 = 1Y
// Planes ip = 3 = 2Y
if (!fPlanes[ip]) return -1;
icount = 0;
if ( fScinHitPaddle[ip][0] == 1 )
icount ++;
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[1] - 1; ipaddle++ ){
// !look for number of clusters of 1 or more hits
if ( ( fScinHitPaddle[ip][ipaddle] == 0 ) &&
( fScinHitPaddle[ip][ipaddle + 1] == 1 ) )
icount ++;
} // Loop over Y paddles
fNClust[ip] = icount;
icount = 0;
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[1] - 2; ipaddle++ ){
// !look for three or more adjacent hits
if ( ( fScinHitPaddle[ip][ipaddle] == 1 ) &&
( fScinHitPaddle[ip][ipaddle + 1] == 1 ) &&
( fScinHitPaddle[ip][ipaddle + 2] == 1 ) )
icount ++;
} // Second loop over Y paddles
if ( icount > 0 )
fThreeScin[ip] = 1;
}// Loop over Y planes
// *now put some "tracking" like cuts on the hslopes, based only on scins...
// *by "slope" here, I mean the difference in the position of scin hits in two
// *like-planes. For example, a track that those great straight through will
// *have a slope of zero. If it moves one scin over from s1x to s2x it has an
// *x-slope of 1... I pick the minimum slope if there are multiple scin hits.
Double_t bestXpScin = 100.0;
Double_t bestYpScin = 100.0;
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[0]; ipaddle++ ){
for (Int_t ipaddle2 = 0; ipaddle2 < (Int_t) fNPaddle[0]; ipaddle2++ ){
if ( ( fScinHitPaddle[0][ipaddle] == 1 ) &&
( fScinHitPaddle[2][ipaddle2] == 1 ) ){
Double_t slope = TMath::Abs(ipaddle - ipaddle2);
if ( slope < bestXpScin ) {
bestXpScin = slope;
} // Second loop over X paddles
} // First loop over X paddles
for (Int_t ipaddle = 0; ipaddle < (Int_t) fNPaddle[1]; ipaddle++ ){
for (Int_t ipaddle2 = 0; ipaddle2 < (Int_t) fNPaddle[1]; ipaddle2++ ){
if ( ( fScinHitPaddle[1][ipaddle] == 1 ) &&
( fScinHitPaddle[3][ipaddle2] == 1 ) ){
Double_t slope = TMath::Abs(ipaddle - ipaddle2);
bestYpScin = slope;
1504
1505
1506
1507
1508
1509
1510
1511
1512
1513
1514
1515
1516
1517
1518
1519
1520
1521
1522
1523
1524
1525
1526
1527
1528
1529
1530
1531
1532
1533
1534
1535
1536
1537
1538
1539
1540
1541
1542
1543
1544
1545
1546
1547
1548
1549
1550
1551
1552
1553
1554
1555
1556
1557
1558
1559
1560
1561
1562
1563
1564
1565
1566
1567
1568
1569
1570
1571
1572
1573
1574
1575
1576
1577
1578
1579
1580
1581
} // Second loop over Y paddles
} // First loop over Y paddles
// *next we mask out the edge scintillators, and look at triggers that happened
// *at the center of the acceptance. To change which scins are in the mask
// *change the values of h*loscin and h*hiscin in htracking.param
// fGoodScinHits = 0;
for (Int_t ifidx = fxLoScin[0]; ifidx < (Int_t) fxHiScin[0]; ifidx ++ ){
fGoodScinHitsX.push_back(0);
}
// *first x plane. first see if there are hits inside the scin region
for (Int_t ifidx = fxLoScin[0]-1; ifidx < fxHiScin[0]; ifidx ++ ){
if ( fScinHitPaddle[0][ifidx] == 1 ){
fHitSweet1X = 1;
fSweet1XScin = ifidx + 1;
}
}
// * next make sure nothing fired outside the good region
for (Int_t ifidx = 0; ifidx < fxLoScin[0]-1; ifidx ++ ){
if ( fScinHitPaddle[0][ifidx] == 1 ){ fHitSweet1X = -1; }
}
for (Int_t ifidx = fxHiScin[0]; ifidx < (Int_t) fNPaddle[0]; ifidx ++ ){
if ( fScinHitPaddle[0][ifidx] == 1 ){ fHitSweet1X = -1; }
}
// *second x plane. first see if there are hits inside the scin region
for (Int_t ifidx = fxLoScin[1]-1; ifidx < fxHiScin[1]; ifidx ++ ){
if ( fScinHitPaddle[2][ifidx] == 1 ){
fHitSweet2X = 1;
fSweet2XScin = ifidx + 1;
}
}
// * next make sure nothing fired outside the good region
for (Int_t ifidx = 0; ifidx < fxLoScin[1]-1; ifidx ++ ){
if ( fScinHitPaddle[2][ifidx] == 1 ){ fHitSweet2X = -1; }
}
for (Int_t ifidx = fxHiScin[1]; ifidx < (Int_t) fNPaddle[2]; ifidx ++ ){
if ( fScinHitPaddle[2][ifidx] == 1 ){ fHitSweet2X = -1; }
}
// *first y plane. first see if there are hits inside the scin region
for (Int_t ifidx = fyLoScin[0]-1; ifidx < fyHiScin[0]; ifidx ++ ){
if ( fScinHitPaddle[1][ifidx] == 1 ){
fHitSweet1Y = 1;
fSweet1YScin = ifidx + 1;
}
}
// * next make sure nothing fired outside the good region
for (Int_t ifidx = 0; ifidx < fyLoScin[0]-1; ifidx ++ ){
if ( fScinHitPaddle[1][ifidx] == 1 ){ fHitSweet1Y = -1; }
}
for (Int_t ifidx = fyHiScin[0]; ifidx < (Int_t) fNPaddle[1]; ifidx ++ ){
if ( fScinHitPaddle[1][ifidx] == 1 ){ fHitSweet1Y = -1; }
}
// *second y plane. first see if there are hits inside the scin region
for (Int_t ifidx = fyLoScin[1]-1; ifidx < fyHiScin[1]; ifidx ++ ){
if ( fScinHitPaddle[3][ifidx] == 1 ){
fHitSweet2Y = 1;
fSweet2YScin = ifidx + 1;
}
}
// * next make sure nothing fired outside the good region
for (Int_t ifidx = 0; ifidx < fyLoScin[1]-1; ifidx ++ ){
if ( fScinHitPaddle[3][ifidx] == 1 ){ fHitSweet2Y = -1; }
}
for (Int_t ifidx = fyHiScin[1]; ifidx < (Int_t) fNPaddle[3]; ifidx ++ ){
if ( fScinHitPaddle[3][ifidx] == 1 ){ fHitSweet2Y = -1; }
}
fTestSum = fHitSweet1X + fHitSweet2X + fHitSweet1Y + fHitSweet2Y;
// * now define a 3/4 or 4/4 trigger of only good scintillators the value
// * is specified in htracking.param...
if ( fTestSum >= fTrackEffTestNScinPlanes ){
fGoodScinHits = 1;
for (Int_t ifidx = fxLoScin[0]; ifidx < fxHiScin[0]; ifidx ++ ){
if ( fSweet1XScin == ifidx )
fGoodScinHitsX[ifidx] = 1;
}
}
// * require front/back hodoscopes be close to each other
if ( ( fGoodScinHits == 1 ) && ( fTrackEffTestNScinPlanes == 4 ) ){
if ( TMath::Abs( fSweet1XScin - fSweet2XScin ) > 3 )
fGoodScinHits = 0;
if ( TMath::Abs( fSweet1YScin - fSweet2YScin ) > 2 )
fGoodScinHits = 0;
}
if ( !fChern || !fShower ) {
return 0;
if ( ( fGoodScinHits == 1 ) && ( fShower->GetNormETot() > fNormETot ) &&
( fChern->GetCerNPE() > fNCerNPE ) )
fScinShould = 1;
if ( ( fGoodScinHits == 1 ) && ( fShower->GetNormETot() > fNormETot ) &&
( fChern->GetCerNPE() > fNCerNPE ) && ( tracks.GetLast() + 1 > 0 ) ) {
fScinDid = 1;
Gabriel Niculescu
committed
//_____________________________________________________________________________
Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
{
Int_t Ntracks = tracks.GetLast()+1; // Number of reconstructed tracks
for (Int_t itrk=0; itrk<Ntracks; itrk++) {
THaTrack* theTrack = static_cast<THaTrack*>( tracks[itrk] );
if (theTrack->GetIndex()==0) {
fBeta=theTrack->GetBeta();
}
} //over tracks
return 0;
}
//_____________________________________________________________________________
Gabriel Niculescu
committed
Int_t THcHodoscope::GetScinIndex( Int_t nPlane, Int_t nPaddle ) {
// GN: Return the index of a scintillator given the plane # and the paddle #
// This assumes that both planes and
// paddles start counting from 0!
// Result also counts from 0.
return fNPlanes*nPaddle+nPlane;
Gabriel Niculescu
committed
}
//_____________________________________________________________________________
Int_t THcHodoscope::GetScinIndex( Int_t nSide, Int_t nPlane, Int_t nPaddle ) {
return nSide*fMaxHodoScin+fNPlanes*nPaddle+nPlane-1;
}
//_____________________________________________________________________________
Double_t THcHodoscope::GetPathLengthCentral() {
return fPathLengthCentral;
}
ClassImp(THcHodoscope)
////////////////////////////////////////////////////////////////////////////////