diff --git a/src/THcHodoscope.cxx b/src/THcHodoscope.cxx
index 8d0fc57f785bafafc482bfc4f713eecc389ff001..164e95d8c028d3214069dd2c6fe5bd8d5e7765e8 100644
--- a/src/THcHodoscope.cxx
+++ b/src/THcHodoscope.cxx
@@ -173,7 +173,7 @@ THaAnalysisObject::EStatus THcHodoscope::Init( const TDatime& date )
// memcpy( fDataDest, tmp, NDEST*sizeof(DataDest) );
fNScinHits = new Int_t [fNPlanes];
- fGoodPlaneTime = new Bool_t [fNPlanes];
+ fGoodPlaneTime = new Bool_t [fNPlanes];
fNPlaneTime = new Int_t [fNPlanes];
fSumPlaneTime = new Double_t [fNPlanes];
@@ -681,9 +681,7 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
Bool_t goodplanetime[fNPlanes];
Bool_t twogoodtimes[nscinhits];
- Int_t temp_planes=fNPlanes;
- if (fSHMS) temp_planes=3;
- for(Int_t ip=0;ip<temp_planes;ip++) {
+ for(Int_t ip=0;ip<fNPlanes;ip++) {
goodplanetime[ip] = kFALSE;
Int_t nphits=fPlanes[ip]->GetNScinHits();
TClonesArray* hodoHits = fPlanes[ip]->GetHits();
@@ -747,9 +745,7 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
Double_t sumTZ = 0.;
Int_t ihhit = 0;
- temp_planes=fNPlanes;
- if (fSHMS) temp_planes=3;
- for(Int_t ip=0;ip<temp_planes;ip++) {
+ for(Int_t ip=0;ip<fNPlanes;ip++) {
Int_t nphits=fPlanes[ip]->GetNScinHits();
TClonesArray* hodoHits = fPlanes[ip]->GetHits();
@@ -785,7 +781,7 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
fBetaNoTrkChiSq = 0.;
ihhit = 0;
- for (Int_t ip = 0; ip < temp_planes; ip++ ){ // Loop over planes
+ for (Int_t ip = 0; ip < fNPlanes; ip++ ){ // Loop over planes
Int_t nphits=fPlanes[ip]->GetNScinHits();
TClonesArray* hodoHits = fPlanes[ip]->GetHits();
@@ -845,8 +841,6 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
// -------------------------------------------------
// fDumpOut << " ntrack = " << ntracks << endl;
- Int_t temp_planes=fNPlanes;
- if (fSHMS) temp_planes=3;
if (tracks.GetLast()+1 > 0 ) {
@@ -860,7 +854,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
THaTrack* theTrack = dynamic_cast<THaTrack*>( tracks.At(itrack) );
if (!theTrack) return -1;
- for (Int_t ip = 0; ip < temp_planes; ip++ ){
+ for (Int_t ip = 0; ip < fNPlanes; ip++ ){
fGoodPlaneTime[ip] = kFALSE;
fNScinHits[ip] = 0;
fNPlaneTime[ip] = 0;
@@ -899,7 +893,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
fTOFPInfo.clear(); // SAW - combine these two?
Int_t ihhit = 0; // Hit # overall
- for(Int_t ip = 0; ip < temp_planes; ip++ ) {
+ for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
std::vector<GoodFlags> goodflagstmp2;
fGoodFlags[itrack].push_back(goodflagstmp2);
@@ -980,7 +974,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
timep -= fHodoPosInvAdcOffset[fPIndex]
+ pathp/fHodoPosInvAdcLinear[fPIndex]
+ fHodoPosInvAdcAdc[fPIndex]
- /TMath::Sqrt(TMath::Max(20.0,adc_pos));
+ /TMath::Sqrt(TMath::Max(20.0*.020,adc_pos));
} else {
timep -= fHodoPosPhcCoeff[fPIndex]*
TMath::Sqrt(TMath::Max(0.0,adc_pos/fHodoPosMinPh[fPIndex]-1.0))
@@ -1008,7 +1002,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
timen -= fHodoNegInvAdcOffset[fPIndex]
+ pathn/fHodoNegInvAdcLinear[fPIndex]
+ fHodoNegInvAdcAdc[fPIndex]
- /TMath::Sqrt(TMath::Max(20.0,adc_neg));
+ /TMath::Sqrt(TMath::Max(20.0*.020,adc_neg));
} else {
timen -= fHodoNegPhcCoeff[fPIndex]*
TMath::Sqrt(TMath::Max(0.0,adc_neg/fHodoNegMinPh[fPIndex]-1.0))
@@ -1094,6 +1088,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
Int_t fPIndex = GetScinIndex(ip,paddle);
if (fTOFPInfo[ih].onTrack) {
+ fGoodFlags[itrack][ip][iphit].onTrack = kTRUE;
if ( fTOFPInfo[ih].keep_pos ) { // 301
fTOFCalc[ih].good_tdc_pos = kTRUE;
fGoodFlags[itrack][ip][iphit].goodTdcPos = kTRUE;
@@ -1197,7 +1192,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
} // Second loop over hits of a scintillator plane ends here
theTrack->SetGoodPlane3( fGoodPlaneTime[2] ? 1 : 0 );
- theTrack->SetGoodPlane4( fGoodPlaneTime[3] ? 1 : 0 );
+ if (fNPlanes==4) theTrack->SetGoodPlane4( fGoodPlaneTime[3] ? 1 : 0 );
//
//------------------------------------------------------------------------------
//------------------------------------------------------------------------------
@@ -1287,7 +1282,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
Double_t FPTimeSum=0.0;
Int_t nFPTimeSum=0;
- for (Int_t ip = 0; ip < temp_planes; ip++ ){
+ for (Int_t ip = 0; ip < fNPlanes; ip++ ){
if ( fNPlaneTime[ip] != 0 ){
fFPTime[ip] = ( fSumPlaneTime[ip] / fNPlaneTime[ip] );
FPTimeSum += fSumPlaneTime[ip];
@@ -1332,7 +1327,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
// *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++ ) {
+ for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
if (!fPlanes[ip])
return -1;
@@ -1353,7 +1348,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
// *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++ ) {
+ for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
// Planes ip = 0 = 1X
// Planes ip = 2 = 2X
if (!fPlanes[ip]) return -1;
@@ -1398,7 +1393,7 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
// *look for clusters in y planes... (10 scins) !this assume both y planes have same
// *number of scintillators.
- for (Int_t ip = 1; ip < temp_planes; ip +=2 ) {
+ for (Int_t ip = 1; ip < fNPlanes; ip +=2 ) {
// Planes ip = 1 = 1Y
// Planes ip = 3 = 2Y
if (!fPlanes[ip]) return -1;
@@ -1539,13 +1534,42 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
{
Int_t Ntracks = tracks.GetLast()+1; // Number of reconstructed tracks
+ Double_t hitPos;
+ Double_t hitDistance;
for (Int_t itrk=0; itrk<Ntracks; itrk++) {
THaTrack* theTrack = static_cast<THaTrack*>( tracks[itrk] );
if (theTrack->GetIndex()==0) {
fBeta=theTrack->GetBeta();
+ for (Int_t ip = 0; ip < fNPlanes; ip++ ){
+ Double_t pl_xypos=0;
+ Double_t pl_zpos=0;
+ Int_t num_good_pad=0;
+ TClonesArray* hodoHits = fPlanes[ip]->GetHits();
+ for (Int_t iphit = 0; iphit < fPlanes[ip]->GetNScinHits(); iphit++ ){
+ if (fGoodFlags[itrk][ip][iphit].goodScinTime) {
+ Int_t padind = ((THcHodoHit*)hodoHits->At(iphit))->GetPaddleNumber()-1;
+ pl_xypos+=fPlanes[ip]->GetPosCenter(padind)+ fPlanes[ip]->GetPosOffset();
+ pl_zpos+=fPlanes[ip]->GetZpos()+ (padind%2)*fPlanes[ip]->GetDzpos();
+ num_good_pad++;
+ }
+ // cout << ip << " " << iphit << " " << fGoodFlags[itrk][ip][iphit].goodScinTime << " " << fGoodFlags[itrk][ip][iphit].goodTdcPos << " " << fGoodFlags[itrk][ip][iphit].goodTdcNeg << endl;
+ }
+ hitDistance=kBig;
+ if (num_good_pad !=0 ) {
+ pl_xypos=pl_xypos/num_good_pad;
+ pl_zpos=pl_zpos/num_good_pad;
+ hitPos = theTrack->GetY() + theTrack->GetPhi()*pl_zpos ;
+ if (ip%2 == 0) hitPos = theTrack->GetX() + theTrack->GetTheta()*pl_zpos ;
+ hitDistance = hitPos - pl_xypos;
+ fPlanes[ip]->SetTrackXPosition(theTrack->GetX() + theTrack->GetTheta()*pl_zpos );
+ fPlanes[ip]->SetTrackYPosition(theTrack->GetY() + theTrack->GetPhi()*pl_zpos );
+ }
+ // cout << " ip " << ip << " " << hitPos << " " << pl_xypos << " " << pl_zpos << " " << hitDistance << endl;
+ fPlanes[ip]->SetHitDistance(hitDistance);
}
- } //over tracks
-
+ }
+ } //over tracks
+ //
return 0;
}
//_____________________________________________________________________________
diff --git a/src/THcScintillatorPlane.cxx b/src/THcScintillatorPlane.cxx
index 7f621d2819fa357f2b11144c12e1c6fd382ba8d8..61f3dea71c6a422f1bebd2b0711681c62432090b 100644
--- a/src/THcScintillatorPlane.cxx
+++ b/src/THcScintillatorPlane.cxx
@@ -425,6 +425,9 @@ Int_t THcScintillatorPlane::DefineVariables( EMode mode )
{"GoodNegAdcPulseAmp", "List of Negative ADC peak amp (passed TDC && ADC Min and Max cuts for either end)", "fGoodNegAdcPulseAmp"},
{"GoodPosAdcPulseTime", "List of positive ADC time (passed TDC && ADC Min and Max cuts for either end)", "fGoodPosAdcPulseTime"},
{"GoodNegAdcPulseTime", "List of Negative ADC time (passed TDC && ADC Min and Max cuts for either end)", "fGoodNegAdcPulseTime"},
+ {"DiffDisTrack", "Difference between track and scintillator position (cm)", "fHitDistance"},
+ {"TrackXPos", "Track X position at plane (cm)", "fTrackXPosition"},
+ {"TrackYPos", "Track Y position at plane (cm)", "fTrackYPosition"},
//{"ngoodhits", "Number of paddle hits (passed tof tolerance and used to determine the focal plane time )", "GetNGoodHits() "},
{ 0 }
};
@@ -514,6 +517,9 @@ void THcScintillatorPlane::Clear( Option_t* )
fpTime = -1.e4;
+ fHitDistance = kBig;
+ fTrackXPosition = kBig;
+ fTrackYPosition = kBig;
}
//_____________________________________________________________________________
@@ -861,7 +867,7 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
fGoodPosAdcPed.at(padnum-1) = ((THcSignalHit*) frPosAdcPed->ConstructedAt(good_ielem_posadc))->GetData();
fGoodPosAdcPulseInt.at(padnum-1) = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(good_ielem_posadc))->GetData();
fGoodPosAdcPulseAmp.at(padnum-1) = ((THcSignalHit*) frPosAdcPulseAmp->ConstructedAt(good_ielem_posadc))->GetData();
- fGoodPosAdcPulseTime.at(padnum-1) = ((THcSignalHit*) frPosAdcPulseTimeRaw->ConstructedAt(good_ielem_posadc))->GetData();
+ fGoodPosAdcPulseTime.at(padnum-1) = ((THcSignalHit*) frPosAdcPulseTime->ConstructedAt(good_ielem_posadc))->GetData();
}
@@ -878,7 +884,7 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
fGoodNegAdcPed.at(padnum-1) = ((THcSignalHit*) frNegAdcPed->ConstructedAt(good_ielem_negadc))->GetData();
fGoodNegAdcPulseInt.at(padnum-1) = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(good_ielem_negadc))->GetData();
fGoodNegAdcPulseAmp.at(padnum-1) = ((THcSignalHit*) frNegAdcPulseAmp->ConstructedAt(good_ielem_negadc))->GetData();
- fGoodNegAdcPulseTime.at(padnum-1) = ((THcSignalHit*) frNegAdcPulseTimeRaw->ConstructedAt(good_ielem_negadc))->GetData();
+ fGoodNegAdcPulseTime.at(padnum-1) = ((THcSignalHit*) frNegAdcPulseTime->ConstructedAt(good_ielem_negadc))->GetData();
}
//DEFINE THE "GOOD +TDC Multiplicities and Occupancies"
@@ -912,9 +918,9 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPosADCpeak(adc_peak);
adc_peak=hit->GetRawAdcHitNeg().GetPulseAmp();
((THcHodoHit*) fHodoHits->At(fNScinHits))->SetNegADCpeak(adc_peak);
- Double_t time_peak=hit->GetRawAdcHitPos().GetPulseTimeRaw();
+ Double_t time_peak=hit->GetRawAdcHitPos().GetPulseTime();
((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPosADCtime(time_peak);
- time_peak=hit->GetRawAdcHitNeg().GetPulseTimeRaw();
+ time_peak=hit->GetRawAdcHitNeg().GetPulseTime();
((THcHodoHit*) fHodoHits->At(fNScinHits))->SetNegADCtime(time_peak);
//Define GoodTdcUnCorrTime
@@ -935,10 +941,10 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
if(fTofUsingInvAdc) {
timec_pos = tdc_pos*fScinTdcToTime
- fHodoPosInvAdcOffset[index]
- - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_pos));
+ - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adc_pos));
timec_neg = tdc_neg*fScinTdcToTime
- fHodoNegInvAdcOffset[index]
- - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_neg));
+ - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adc_neg));
} else { // Old style
timec_pos = tdc_pos*fScinTdcToTime - fHodoPosPhcCoeff[index]*
TMath::Sqrt(TMath::Max(0.0,adc_pos/fHodoPosMinPh[index]-1.0))
@@ -1005,7 +1011,7 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
if(fTofUsingInvAdc) {
timec_pos = tdc_pos*fScinTdcToTime
- fHodoPosInvAdcOffset[index]
- - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_pos));
+ - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adc_pos));
} else { // Old style
timec_pos = tdc_pos*fScinTdcToTime - fHodoPosPhcCoeff[index]*
TMath::Sqrt(TMath::Max(0.0,adc_pos/fHodoPosMinPh[index]-1.0))
@@ -1016,7 +1022,7 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
if(fTofUsingInvAdc) {
timec_neg = tdc_neg*fScinTdcToTime
- fHodoNegInvAdcOffset[index]
- - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_neg));
+ - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0*.020,adc_neg));
} else { // Old style
timec_neg = tdc_neg*fScinTdcToTime - fHodoNegPhcCoeff[index]*
TMath::Sqrt(TMath::Max(0.0,adc_neg/fHodoNegMinPh[index]-1.0))
diff --git a/src/THcScintillatorPlane.h b/src/THcScintillatorPlane.h
index bc900ed565a4bd07e16efd2c6f5029e1d2b095a9..866c9f0c5ff7d2aedf97f936f0b1d33c58b0775c 100644
--- a/src/THcScintillatorPlane.h
+++ b/src/THcScintillatorPlane.h
@@ -44,6 +44,9 @@ class THcScintillatorPlane : public THaSubDetector {
Int_t GetNelem() {return fNelem;}; // return number of paddles in this plane
Int_t GetNScinHits() {return fNScinHits;}; // Get # hits in plane (that pass min/max TDC cuts)
Int_t GetNGoodHits() {return fNGoodHits;}; // Get # hits in plane (used in determining focal plane time)
+ Double_t GetHitDistance() {return fHitDistance;}; // Distance between track and hit paddle
+ Double_t GetTrackXPosition() {return fTrackXPosition;}; // Distance track X position at plane
+ Double_t GetTrackYPosition() {return fTrackYPosition;}; // Distance track Y position at plane
Double_t GetSpacing() {return fSpacing;}; // spacing of paddles
Double_t GetSize() {return fSize;}; // paddle size
Double_t GetHodoSlop() {return fHodoSlop;}; // hodo slop
@@ -59,6 +62,9 @@ class THcScintillatorPlane : public THaSubDetector {
void SetFpTime(Double_t f) {fFptime=f;};
void SetNGoodHits(Int_t ng) {fNGoodHits=ng;};
+ void SetHitDistance(Double_t f) {fHitDistance=f;}; // Distance between track and hit paddle
+ void SetTrackXPosition(Double_t f) {fTrackXPosition=f;}; // Distance track X position at plane
+ void SetTrackYPosition(Double_t f) {fTrackYPosition=f;}; // Distance track Y position at plane
TClonesArray* fParentHitList;
@@ -152,6 +158,9 @@ class THcScintillatorPlane : public THaSubDetector {
Int_t fDebugAdc;
+ Double_t fHitDistance;
+ Double_t fTrackXPosition;
+ Double_t fTrackYPosition;
Int_t fCosmicFlag; //
Int_t fPlaneNum; /* Which plane am I 1-4 */
UInt_t fTotPlanes; /* so we can read variables that are not indexed by plane id */
@@ -204,7 +213,6 @@ class THcScintillatorPlane : public THaSubDetector {
Double_t *fHodoPosInvAdcAdc;
Double_t *fHodoNegInvAdcAdc;
Double_t *fHodoSigma;
-
Double_t fTolerance; /* need this for Focal Plane Time estimation */
Double_t fFptime;
/* Pedestal Quantities */