diff --git a/examples/PARAM/general.param b/examples/PARAM/general.param index 5383e6c1028e2cd9b610b1f215cabc03f070d692..f6d2f37fee455def56a4335216aebf7bce92d22d 100644 --- a/examples/PARAM/general.param +++ b/examples/PARAM/general.param @@ -23,7 +23,7 @@ raddeg=3.14159265/180 #include "PARAM/sdc.pos" #include "PARAM/shodo.pos" #include "PARAM/scal.pos" -#include "PARAM/hhodo.param" +#include "PARAM/hhodo.param_v1" #include "PARAM/haero.param" #include "PARAM/hdc.param" #include "PARAM/hdriftmap.param" diff --git a/src/THcHodoHit.h b/src/THcHodoHit.h index 39ea3a088ce14221e72c442f4ef194a6f0bead50..b6d6aaee38d4c14ee37dabbebf41ec16d79d1559 100644 --- a/src/THcHodoHit.h +++ b/src/THcHodoHit.h @@ -20,7 +20,8 @@ public: THcHodoHit(Int_t postdc, Int_t negtdc, Double_t posadc, Double_t negadc, Int_t ipad, THcScintillatorPlane* sp) : fPosTDC(postdc), fNegTDC(negtdc), fPosADC_Ped(posadc), fNegADC_Ped(negadc), - fPaddleNumber(ipad), fPlane(sp) {}; + fPaddleNumber(ipad), fHasCorrectedTimes(kFALSE), + fTwoGoodTimes(kFALSE), fPlane(sp) {}; virtual ~THcHodoHit() {} @@ -36,6 +37,8 @@ public: Double_t GetPosTOFCorrectedTime() const { return fPosTOFCorrectedTime;} Double_t GetNegTOFCorrectedTime() const { return fNegTOFCorrectedTime;} Double_t GetScinCorrectedTime() const { return fScinCorrectedTime;} + Bool_t GetTwoGoodTimes() const { return fTwoGoodTimes;} + Bool_t GetHasCorrectedTimes() const { return fHasCorrectedTimes;} Int_t GetPaddleNumber() const { return fPaddleNumber; } void SetCorrectedTimes(Double_t pos, Double_t neg, Double_t) { @@ -47,6 +50,10 @@ public: fPosCorrectedTime = pos; fNegCorrectedTime = neg; fPosTOFCorrectedTime = postof; fNegTOFCorrectedTime = negtof; fScinCorrectedTime = timeave; + fHasCorrectedTimes = kTRUE; + } + void SetTwoGoodTimes(Bool_t flag) { + fTwoGoodTimes = flag; } protected: @@ -64,6 +71,10 @@ protected: // based on ADCs. Double_t fPosTOFCorrectedTime; // Times corrected for z position Double_t fNegTOFCorrectedTime; // using nominal beta + + Bool_t fHasCorrectedTimes; + Bool_t fTwoGoodTimes; + THcScintillatorPlane* fPlane; // Pointer to parent scintillator plane diff --git a/src/THcHodoscope.cxx b/src/THcHodoscope.cxx index 39d1c9c0f2a092ef7ddc31fc75a2c0c547d10d37..370c85b59122b681d0d2cce49fb0734b7a64cff1 100644 --- a/src/THcHodoscope.cxx +++ b/src/THcHodoscope.cxx @@ -561,7 +561,8 @@ Int_t THcHodoscope::Decode( const THaEvData& evdata ) */ ClearEvent(); // Get the Hall C style hitlist (fRawHitList) for this event - Int_t nhits = DecodeToHitList(evdata); + fNHits = DecodeToHitList(evdata); + // // GN: print event number so we can cross-check with engine // if (evdata.GetEvNum()>1000) @@ -616,7 +617,7 @@ Int_t THcHodoscope::Decode( const THaEvData& evdata ) } - return nhits; + return fNHits; } //_____________________________________________________________________________ @@ -638,19 +639,24 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) timehist[i] = 0; } Int_t ihit=0; + Int_t nscinhits=0; // Total # hits with at least one good tdc for(Int_t ip=0;ip<fNPlanes;ip++) { Int_t nphits=fPlanes[ip]->GetNScinHits(); + nscinhits += nphits; TClonesArray* hodoHits = fPlanes[ip]->GetHits(); for(Int_t i=0;i<nphits;i++) { - Double_t postime=((THcHodoHit*) hodoHits->At(i))->GetPosTOFCorrectedTime(); - Double_t negtime=((THcHodoHit*) hodoHits->At(i))->GetNegTOFCorrectedTime(); - for (Int_t k=0;k<200;k++) { - Double_t tmin=0.5*(k+1); - if ((postime> tmin) && (postime < tmin+fTofTolerance)) { - timehist[k]++; - } - if ((negtime> tmin) && (negtime < tmin+fTofTolerance)) { - timehist[k]++; + THcHodoHit *hit = (THcHodoHit*)hodoHits->At(i); + if(hit->GetHasCorrectedTimes()) { + Double_t postime=hit->GetPosTOFCorrectedTime(); + Double_t negtime=hit->GetNegTOFCorrectedTime(); + for (Int_t k=0;k<200;k++) { + Double_t tmin=0.5*(k+1); + if ((postime> tmin) && (postime < tmin+fTofTolerance)) { + timehist[k]++; + } + if ((negtime> tmin) && (negtime < tmin+fTofTolerance)) { + timehist[k]++; + } } } } @@ -672,29 +678,37 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) Int_t Ngood_hits_plane=0; Double_t Plane_fptime_sum=0.0; - fNoTrkPlaneInfo.clear(); - fNoTrkHitInfo.clear(); + // fNoTrkPlaneInfo.clear(); + // fNoTrkHitInfo.clear(); + Bool_t goodplanetime[fNPlanes]; + Bool_t twogoodtimes[nscinhits]; for(Int_t ip=0;ip<fNPlanes;ip++) { - fNoTrkPlaneInfo.push_back(NoTrkPlaneInfo()); - fNoTrkPlaneInfo[ip].goodplanetime = kFALSE; + // fNoTrkPlaneInfo.push_back(NoTrkPlaneInfo()); + // fNoTrkPlaneInfo[ip].goodplanetime = kFALSE; + goodplanetime[ip] = kFALSE; Int_t nphits=fPlanes[ip]->GetNScinHits(); TClonesArray* hodoHits = fPlanes[ip]->GetHits(); Ngood_hits_plane=0; Plane_fptime_sum=0.0; for(Int_t i=0;i<nphits;i++) { - fNoTrkHitInfo.push_back(NoTrkHitInfo()); - fNoTrkHitInfo[ihit].goodtwotimes = kFALSE; - fNoTrkHitInfo[ihit].goodscintime = kFALSE; + THcHodoHit *hit = (THcHodoHit*)hodoHits->At(i); + // fNoTrkHitInfo.push_back(NoTrkHitInfo()); + // fNoTrkHitInfo[ihit].goodtwotimes = kFALSE; + // fNoTrkHitInfo[ihit].goodscintime = kFALSE; + twogoodtimes[ihit] = kFALSE; Double_t tmin = 0.5*binmax; - Double_t postime=((THcHodoHit*) hodoHits->At(i))->GetPosTOFCorrectedTime(); - Double_t negtime=((THcHodoHit*) hodoHits->At(i))->GetNegTOFCorrectedTime(); + Double_t postime=hit->GetPosCorrectedTime(); + Double_t negtime=hit->GetNegCorrectedTime(); if ((postime>tmin) && (postime<tmin+fTofTolerance) && (negtime>tmin) && (negtime<tmin+fTofTolerance)) { - fNoTrkHitInfo[ihit].goodtwotimes = kTRUE; - fNoTrkHitInfo[ihit].goodscintime = kTRUE; + hit->SetTwoGoodTimes(kTRUE); + // fNoTrkHitInfo[ihit].goodtwotimes = kTRUE; + // fNoTrkHitInfo[ihit].goodscintime = kTRUE; + twogoodtimes[ihit] = kTRUE; // Both tubes fired - Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1; - Double_t fptime = ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime() + Int_t index=hit->GetPaddleNumber()-1; + // Need to put this in a multihit histo + Double_t fptime = hit->GetScinCorrectedTime() - (fPlanes[ip]->GetZpos()+(index%2)*fPlanes[ip]->GetDzpos()) / (29.979 * fBetaNominal); if(TMath::Abs(fptime-fStartTimeCenter)<=fStartTimeSlop) { @@ -702,8 +716,11 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) Plane_fptime_sum+=fptime; fpTimeSum += fptime; fNfptimes++; - fNoTrkPlaneInfo[ip].goodplanetime = kTRUE; + goodplanetime[ip] = kTRUE; + //fNoTrkPlaneInfo[ip].goodplanetime = kTRUE; } + } else { + hit->SetTwoGoodTimes(kFALSE); } fPlanes[ip]->SetFpTime(Plane_fptime_sum/float(Ngood_hits_plane)); fPlanes[ip]->SetNGoodHits(Ngood_hits_plane); @@ -720,8 +737,10 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) } - if ( ( fNoTrkPlaneInfo[0].goodplanetime || fNoTrkPlaneInfo[1].goodplanetime ) && - ( fNoTrkPlaneInfo[2].goodplanetime || fNoTrkPlaneInfo[3].goodplanetime ) ){ + // if ( ( fNoTrkPlaneInfo[0].goodplanetime || fNoTrkPlaneInfo[1].goodplanetime ) && + // ( fNoTrkPlaneInfo[2].goodplanetime || fNoTrkPlaneInfo[3].goodplanetime ) ){ + if((goodplanetime[0]||goodplanetime[1]) + &&(goodplanetime[2]||goodplanetime[3])) { Double_t sumW = 0.; Double_t sumT = 0.; @@ -737,7 +756,8 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) for(Int_t i=0;i<nphits;i++) { Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1; - if ( fNoTrkHitInfo[ihhit].goodscintime ) { + // if ( fNoTrkHitInfo[ihhit].goodscintime ) { + if(twogoodtimes[ihhit]){ Double_t sigma = 0.5 * ( TMath::Sqrt( TMath::Power( fHodoPosSigma[GetScinIndex(ip,index)],2) + TMath::Power( fHodoNegSigma[GetScinIndex(ip,index)],2) ) ); @@ -774,7 +794,8 @@ void THcHodoscope::EstimateFocalPlaneTime( void ) for(Int_t i=0;i<nphits;i++) { Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1; - if ( fNoTrkHitInfo[ihhit].goodscintime ) { + // if ( fNoTrkHitInfo[ihhit].goodscintime ) { + if(twogoodtimes[ihhit]) { Double_t zPosition = fPlanes[ip]->GetZpos() + (index%2)*fPlanes[ip]->GetDzpos(); Double_t timeDif = ( ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime() - t0 ); @@ -891,6 +912,9 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) fNScinHits[ip] = fPlanes[ip]->GetNScinHits(); TClonesArray* hodoHits = fPlanes[ip]->GetHits(); + + Double_t zPos = fPlanes[ip]->GetZpos(); + Double_t dzPos = fPlanes[ip]->GetDzpos(); // first loop over hits with in a single plane for (Int_t iphit = 0; iphit < fNScinHits[ip]; iphit++ ){ @@ -908,17 +932,16 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) fTOFPInfo[ihhit].hit = hit; fTOFPInfo[ihhit].planeIndex = ip; fTOFPInfo[ihhit].hitNumInPlane = iphit; + fTOFPInfo[ihhit].onTrack = kFALSE; Int_t paddle = hit->GetPaddleNumber()-1; - + Double_t zposition = zPos + (paddle%2)*dzPos; + Double_t xHitCoord = theTrack->GetX() + theTrack->GetTheta() * - ( fPlanes[ip]->GetZpos() + - ( paddle % 2 ) * fPlanes[ip]->GetDzpos() ); // Line 183 + ( zposition ); // Line 183 Double_t yHitCoord = theTrack->GetY() + theTrack->GetPhi() * - ( fPlanes[ip]->GetZpos() + - ( paddle % 2 ) * fPlanes[ip]->GetDzpos() ); // Line 184 - + ( zposition ); // Line 184 Double_t scinTrnsCoord, scinLongCoord; if ( ( ip == 0 ) || ( ip == 2 ) ){ // !x plane. Line 185 @@ -937,42 +960,55 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) Double_t scinCenter = fPlanes[ip]->GetPosCenter(paddle) + fPlanes[ip]->GetPosOffset(); // Index to access the 2d arrays of paddle/scintillator properties - Int_t fPIndex = fNPlanes * paddle + ip; - + Int_t fPIndex = GetScinIndex(ip,paddle); if ( TMath::Abs( scinCenter - scinTrnsCoord ) < ( fPlanes[ip]->GetSize() * 0.5 + fPlanes[ip]->GetHodoSlop() ) ){ // Line 293 fTOFPInfo[ihhit].onTrack = kTRUE; - Double_t pathp = fPlanes[ip]->GetPosLeft() - scinLongCoord; - Double_t timep = hit->GetPosCorrectedTime() - pathp/fHodoVelLight[fPIndex]; - fTOFPInfo[ihhit].scin_pos_time = timep; - timep -= (fPlanes[ip]->GetZpos() + (paddle%2)*fPlanes[ip]->GetDzpos()) - /(29.979*fBetaP)* - TMath::Sqrt(1. + theTrack->GetTheta()*theTrack->GetTheta() - + theTrack->GetPhi()*theTrack->GetPhi()); - fTOFPInfo[ihhit].time_pos = timep; + 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; + Double_t timep = tdc_pos*fScinTdcToTime + - fHodoPosInvAdcOffset[fPIndex] + - pathp/fHodoPosInvAdcLinear[fPIndex] + - fHodoPosInvAdcAdc[fPIndex] + /TMath::Sqrt(TMath::Max(20.0,adc_pos)); + fTOFPInfo[ihhit].scin_pos_time = timep; + timep -= zposition/(29.979*fBetaP)* + TMath::Sqrt(1. + theTrack->GetTheta()*theTrack->GetTheta() + + theTrack->GetPhi()*theTrack->GetPhi()); + 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] ++; + 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 pathn = scinLongCoord - fPlanes[ip]->GetPosRight(); - Double_t timen = hit->GetNegCorrectedTime() - pathn/fHodoVelLight[fPIndex]; - fTOFPInfo[ihhit].scin_neg_time = timen; - timen -= - (fPlanes[ip]->GetZpos() + (paddle%2)*fPlanes[ip]->GetDzpos()) - /(29.979*fBetaP)* - TMath::Sqrt( 1. + theTrack->GetTheta()*theTrack->GetTheta() - + theTrack->GetPhi()*theTrack->GetPhi()); - fTOFPInfo[ihhit].time_neg = timen; + 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(); + Double_t timen = tdc_neg*fScinTdcToTime + - fHodoNegInvAdcOffset[fPIndex] + - pathn/fHodoNegInvAdcLinear[fPIndex] + - fHodoNegInvAdcAdc[fPIndex] + /TMath::Sqrt(TMath::Max(20.0,adc_neg)); + fTOFPInfo[ihhit].scin_neg_time = timen; + timen -= zposition/(29.979*fBetaP)* + TMath::Sqrt( 1. + theTrack->GetTheta()*theTrack->GetTheta() + + theTrack->GetPhi()*theTrack->GetPhi()); + 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] ++; + 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 ihhit++; @@ -1039,7 +1075,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) // Double_t scinTrnsCoord = fTOFPInfo[ihhit].scinTrnsCoord; // Double_t scinLongCoord = fTOFPInfo[ihhit].scinLongCoord; - Int_t fPIndex = fNPlanes * paddle + ip; + Int_t fPIndex = GetScinIndex(ip,paddle); if (fTOFPInfo[ihhit].onTrack) { if ( fTOFPInfo[ihhit].keep_pos ) { // 301 @@ -1056,12 +1092,15 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) if ( fTOFCalc[ihhit].good_tdc_neg ){ fTOFCalc[ihhit].scin_time = ( fTOFPInfo[ihhit].scin_pos_time + fTOFPInfo[ihhit].scin_neg_time ) / 2.; + fTOFCalc[ihhit].scin_time_fp = ( fTOFPInfo[ihhit].time_pos + + fTOFPInfo[ihhit].time_neg ) / 2.; fTOFCalc[ihhit].scin_sigma = TMath::Sqrt( fHodoPosSigma[fPIndex] * fHodoPosSigma[fPIndex] + fHodoNegSigma[fPIndex] * fHodoNegSigma[fPIndex] )/2.; fTOFCalc[ihhit].good_scin_time = kTRUE; fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE; } else{ fTOFCalc[ihhit].scin_time = fTOFPInfo[ihhit].scin_pos_time; + fTOFCalc[ihhit].scin_time_fp = fTOFPInfo[ihhit].time_pos; fTOFCalc[ihhit].scin_sigma = fHodoPosSigma[fPIndex]; fTOFCalc[ihhit].good_scin_time = kTRUE; fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE; @@ -1069,6 +1108,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) } else { if ( fTOFCalc[ihhit].good_tdc_neg ){ fTOFCalc[ihhit].scin_time = fTOFPInfo[ihhit].scin_neg_time; + fTOFCalc[ihhit].scin_time_fp = fTOFPInfo[ihhit].time_neg; fTOFCalc[ihhit].scin_sigma = fHodoNegSigma[fPIndex]; fTOFCalc[ihhit].good_scin_time = kTRUE; fGoodFlags[itrack][ip][iphit].goodScinTime = kTRUE; @@ -1081,9 +1121,10 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks ) // 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[ihhit].time_pos + - fTOFPInfo[ihhit].time_neg ) / 2.; - + // Double_t scin_time_fp = ( fTOFPInfo[ihhit].time_pos + + // fTOFPInfo[ihhit].time_neg ) / 2.; + Double_t scin_time_fp = fTOFCalc[ihhit].scin_time_fp; + sumFPTime = sumFPTime + scin_time_fp; nFPTime ++; diff --git a/src/THcHodoscope.h b/src/THcHodoscope.h index 817f899fc9c15245e6161293ab222d7c51e040b8..c646453aa76910d70e8aab18b96fdfd9568a273d 100644 --- a/src/THcHodoscope.h +++ b/src/THcHodoscope.h @@ -68,6 +68,13 @@ public: Double_t GetHodoPosTimeOffset(Int_t iii) const {return fHodoPosTimeOffset[iii];} Double_t GetHodoNegTimeOffset(Int_t iii) const {return fHodoNegTimeOffset[iii];} Double_t GetHodoVelLight(Int_t iii) const {return fHodoVelLight[iii];} + Double_t GetHodoPosInvAdcOffset(Int_t iii) const {return fHodoPosInvAdcOffset[iii];} + Double_t GetHodoNegInvAdcOffset(Int_t iii) const {return fHodoNegInvAdcOffset[iii];} + Double_t GetHodoPosInvAdcLinear(Int_t iii) const {return fHodoPosInvAdcLinear[iii];} + Double_t GetHodoNegInvAdcLinear(Int_t iii) const {return fHodoNegInvAdcLinear[iii];} + Double_t GetHodoPosInvAdcAdc(Int_t iii) const {return fHodoPosInvAdcAdc[iii];} + Double_t GetHodoNegInvAdcAdc(Int_t iii) const {return fHodoNegInvAdcAdc[iii];} + Double_t GetStartTimeCenter() const {return fStartTimeCenter;} Double_t GetStartTimeSlop() const {return fStartTimeSlop;} Double_t GetBetaNotrk() const {return fBetaNoTrk;} @@ -114,6 +121,8 @@ protected: Int_t fAnalyzePedestals; + Int_t fNHits; + // Calibration // Per-event data @@ -242,19 +251,19 @@ protected: // } fDataDest[NDEST]; // Lookup table for decoder // Inforamtion for each plane - struct NoTrkPlaneInfo { - Bool_t goodplanetime; - NoTrkPlaneInfo () : goodplanetime(kFALSE) {} - }; - std::vector<NoTrkPlaneInfo> fNoTrkPlaneInfo; + // struct NoTrkPlaneInfo { + // Bool_t goodplanetime; + // NoTrkPlaneInfo () : goodplanetime(kFALSE) {} + // }; + // std::vector<NoTrkPlaneInfo> fNoTrkPlaneInfo; // Inforamtion for each plane - struct NoTrkHitInfo { - Bool_t goodtwotimes; - Bool_t goodscintime; - NoTrkHitInfo () : goodtwotimes(kFALSE) {} - }; - std::vector<NoTrkHitInfo> fNoTrkHitInfo; + // struct NoTrkHitInfo { + // Bool_t goodtwotimes; + // Bool_t goodscintime; + // NoTrkHitInfo () : goodtwotimes(kFALSE) {} + // }; + // std::vector<NoTrkHitInfo> fNoTrkHitInfo; // Used in TOF calculation (FineProcess) to hold information about hits @@ -291,6 +300,7 @@ scin_pos_time(0.0), scin_neg_time(0.0) {} Bool_t good_tdc_pos; Bool_t good_tdc_neg; Double_t scin_time; + Double_t scin_time_fp; Double_t scin_sigma; Double_t dedx; TOFCalc() : good_scin_time(kFALSE), good_tdc_pos(kFALSE), diff --git a/src/THcScintillatorPlane.cxx b/src/THcScintillatorPlane.cxx index 6ced1dd70fa2c5094c77d06d87cb0fc937552c32..4d3270d44f8c16e70e830c409c8e295c8f97c853 100644 --- a/src/THcScintillatorPlane.cxx +++ b/src/THcScintillatorPlane.cxx @@ -66,12 +66,19 @@ THcScintillatorPlane::~THcScintillatorPlane() delete [] fScinZpos; delete [] fPosCenter; - delete [] fHodoPosMinPh; fHodoPosMinPh = NULL; - delete [] fHodoNegMinPh; fHodoNegMinPh = NULL; - delete [] fHodoPosPhcCoeff; fHodoPosPhcCoeff = NULL; - delete [] fHodoNegPhcCoeff; fHodoNegPhcCoeff = NULL; - delete [] fHodoPosTimeOffset; fHodoPosTimeOffset = NULL; - delete [] fHodoNegTimeOffset; fHodoNegTimeOffset = NULL; + // delete [] fHodoPosMinPh; fHodoPosMinPh = NULL; + // delete [] fHodoNegMinPh; fHodoNegMinPh = NULL; + // delete [] fHodoPosPhcCoeff; fHodoPosPhcCoeff = NULL; + // delete [] fHodoNegPhcCoeff; fHodoNegPhcCoeff = NULL; + // delete [] fHodoPosTimeOffset; fHodoPosTimeOffset = NULL; + // delete [] fHodoNegTimeOffset; fHodoNegTimeOffset = NULL; + delete [] fHodoPosInvAdcOffset; fHodoPosInvAdcOffset = NULL; + delete [] fHodoNegInvAdcOffset; fHodoNegInvAdcOffset = NULL; + delete [] fHodoPosInvAdcLinear; fHodoPosInvAdcLinear = NULL; + delete [] fHodoNegInvAdcLinear; fHodoNegInvAdcLinear = NULL; + delete [] fHodoPosInvAdcAdc; fHodoPosInvAdcAdc = NULL; + delete [] fHodoNegInvAdcAdc; fHodoNegInvAdcAdc = NULL; + delete [] fHodoVelLight; fHodoVelLight = NULL; delete [] fHodoSigma; fHodoSigma = NULL; @@ -164,22 +171,34 @@ Int_t THcScintillatorPlane::ReadDatabase( const TDatime& date ) fStartTimeCenter=((THcHodoscope *)GetParent())->GetStartTimeCenter(); fStartTimeSlop=((THcHodoscope *)GetParent())->GetStartTimeSlop(); // Parameters for this plane - fHodoPosMinPh = new Double_t[fNelem]; - fHodoNegMinPh = new Double_t[fNelem]; - fHodoPosPhcCoeff = new Double_t[fNelem]; - fHodoNegPhcCoeff = new Double_t[fNelem]; - fHodoPosTimeOffset = new Double_t[fNelem]; - fHodoNegTimeOffset = new Double_t[fNelem]; + // fHodoPosMinPh = new Double_t[fNelem]; + // fHodoNegMinPh = new Double_t[fNelem]; + // fHodoPosPhcCoeff = new Double_t[fNelem]; + // fHodoNegPhcCoeff = new Double_t[fNelem]; + // fHodoPosTimeOffset = new Double_t[fNelem]; + // fHodoNegTimeOffset = new Double_t[fNelem]; fHodoVelLight = new Double_t[fNelem]; + fHodoPosInvAdcOffset = new Double_t[fNelem]; + fHodoNegInvAdcOffset = new Double_t[fNelem]; + fHodoPosInvAdcLinear = new Double_t[fNelem]; + fHodoNegInvAdcLinear = new Double_t[fNelem]; + fHodoPosInvAdcAdc = new Double_t[fNelem]; + fHodoNegInvAdcAdc = new Double_t[fNelem]; fHodoSigma = new Double_t[fNelem]; for(Int_t j=0;j<(Int_t) fNelem;j++) { Int_t index=((THcHodoscope *)GetParent())->GetScinIndex(fPlaneNum-1,j); - fHodoPosMinPh[j] = ((THcHodoscope *)GetParent())->GetHodoPosMinPh(index); - fHodoNegMinPh[j] = ((THcHodoscope *)GetParent())->GetHodoNegMinPh(index); - fHodoPosPhcCoeff[j] = ((THcHodoscope *)GetParent())->GetHodoPosPhcCoeff(index); - fHodoNegPhcCoeff[j] = ((THcHodoscope *)GetParent())->GetHodoNegPhcCoeff(index); - fHodoPosTimeOffset[j] = ((THcHodoscope *)GetParent())->GetHodoPosTimeOffset(index); - fHodoNegTimeOffset[j] = ((THcHodoscope *)GetParent())->GetHodoNegTimeOffset(index); + // fHodoPosMinPh[j] = ((THcHodoscope *)GetParent())->GetHodoPosMinPh(index); + // fHodoNegMinPh[j] = ((THcHodoscope *)GetParent())->GetHodoNegMinPh(index); + // fHodoPosPhcCoeff[j] = ((THcHodoscope *)GetParent())->GetHodoPosPhcCoeff(index); + // fHodoNegPhcCoeff[j] = ((THcHodoscope *)GetParent())->GetHodoNegPhcCoeff(index); + // fHodoPosTimeOffset[j] = ((THcHodoscope *)GetParent())->GetHodoPosTimeOffset(index); + // fHodoNegTimeOffset[j] = ((THcHodoscope *)GetParent())->GetHodoNegTimeOffset(index); + fHodoPosInvAdcOffset[j] = ((THcHodoscope *)GetParent())->GetHodoPosInvAdcOffset(index); + fHodoNegInvAdcOffset[j] = ((THcHodoscope *)GetParent())->GetHodoNegInvAdcOffset(index); + fHodoPosInvAdcLinear[j] = ((THcHodoscope *)GetParent())->GetHodoPosInvAdcLinear(index); + fHodoNegInvAdcLinear[j] = ((THcHodoscope *)GetParent())->GetHodoNegInvAdcLinear(index); + fHodoPosInvAdcAdc[j] = ((THcHodoscope *)GetParent())->GetHodoPosInvAdcAdc(index); + fHodoNegInvAdcAdc[j] = ((THcHodoscope *)GetParent())->GetHodoNegInvAdcAdc(index); fHodoVelLight[j] = ((THcHodoscope *)GetParent())->GetHodoVelLight(index); Double_t possigma = ((THcHodoscope *)GetParent())->GetHodoPosSigma(index); Double_t negsigma = ((THcHodoscope *)GetParent())->GetHodoNegSigma(index); @@ -232,7 +251,7 @@ Int_t THcScintillatorPlane::DefineVariables( EMode mode ) {"posadcval", "List of Positive ADC Values", "frPosADCHits.THcSignalHit.GetData()"}, {"negadcval", "List of Negative ADC Values", "frNegADCHits.THcSignalHit.GetData()"}, {"fptime", "Time at focal plane", "GetFpTime()"}, - {"nhits", "Number of paddle hits (passed TDC Min and Max cuts for both ends)", "GetNScinHits() "}, + {"nhits", "Number of paddle hits (passed TDC Min and Max cuts for either end)", "GetNScinHits() "}, {"ngoodhits", "Number of paddle hits (passed tof tolerance and used to determine the focal plane time )", "GetNGoodHits() "}, { 0 } }; @@ -340,19 +359,24 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit) Int_t index=padnum-1; // Need to be finding first hit in TDC range, not the first hit overall + Double_t adc_pos = hit->GetADCPos()-fPosPed[index]; + Double_t adc_neg = hit->GetADCNeg()-fNegPed[index]; if (hit->fNRawHits[2] > 0) ((THcSignalHit*) frPosTDCHits->ConstructedAt(nrPosTDCHits++))->Set(padnum, hit->GetTDCPos()+fTdcOffset); if (hit->fNRawHits[3] > 0) ((THcSignalHit*) frNegTDCHits->ConstructedAt(nrNegTDCHits++))->Set(padnum, hit->GetTDCNeg()+fTdcOffset); + // For making hit maps, we use >= 50 cut + // For making raw hists, we don't want the cut + // We can use a flag to turn on and off these without 50 cut if ((hit->GetADCPos()-fPosPed[index]) >= 50) - ((THcSignalHit*) frPosADCHits->ConstructedAt(nrPosADCHits++))->Set(padnum, hit->GetADCPos()-fPosPed[index]); + ((THcSignalHit*) frPosADCHits->ConstructedAt(nrPosADCHits++))->Set(padnum, adc_pos); if ((hit->GetADCNeg()-fNegPed[index]) >= 50) - ((THcSignalHit*) frNegADCHits->ConstructedAt(nrNegADCHits++))->Set(padnum, hit->GetADCNeg()-fNegPed[index]); + ((THcSignalHit*) frNegADCHits->ConstructedAt(nrNegADCHits++))->Set(padnum, adc_neg); Bool_t btdcraw_pos=kFALSE; Bool_t btdcraw_neg=kFALSE; - Int_t tdc_pos=-1; - Int_t tdc_neg=-1; + Int_t tdc_pos=-99; + Int_t tdc_neg=-99; // Find first in range hit from multihit tdc for(UInt_t thit=0; thit<hit->fNRawHits[2]; thit++) { tdc_pos = hit->GetTDCPos(thit)+fTdcOffset; @@ -370,44 +394,45 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit) } // Proceed if there is a valid TDC on either end of the bar if(btdcraw_pos || btdcraw_neg) { - Double_t adc_pos = hit->GetADCPos()-fPosPed[index]; - Double_t adc_neg = hit->GetADCNeg()-fNegPed[index]; new( (*fHodoHits)[fNScinHits]) THcHodoHit(tdc_pos, tdc_neg, adc_pos, adc_neg, hit->fCounter, this); + // Do corrections if valid TDC on both ends of bar + if(btdcraw_pos && btdcraw_neg) { - // Do the pulse height correction to the time. (Position dependent corrections later) - Double_t timec_pos = tdc_pos*fScinTdcToTime - fHodoPosPhcCoeff[index]* - TMath::Sqrt(TMath::Max(0.0, - (adc_pos)/fHodoPosMinPh[index]-1.0)) - - fHodoPosTimeOffset[index]; - Double_t timec_neg = tdc_neg*fScinTdcToTime - fHodoNegPhcCoeff[index]* - TMath::Sqrt(TMath::Max(0.0, - (adc_neg)/fHodoNegMinPh[index]-1.0)) - - fHodoNegTimeOffset[index]; - - // Find hit position using ADCs - // If postime larger, then hit was nearer negative side. - Double_t dist_from_center=0.5*(timec_neg-timec_pos)*fHodoVelLight[index]; - Double_t scint_center=0.5*(fPosLeft+fPosRight); - Double_t hit_position=scint_center+dist_from_center; - hit_position=TMath::Min(hit_position,fPosLeft); - hit_position=TMath::Max(hit_position,fPosRight); - Double_t postime=timec_pos-(fPosLeft-hit_position)/fHodoVelLight[index]; - Double_t negtime=timec_neg-(hit_position-fPosRight)/fHodoVelLight[index]; - Double_t scin_corrected_time = 0.5*(postime+negtime); - postime = postime-(fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal); - negtime = negtime-(fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal); - - ((THcHodoHit*) fHodoHits->At(fNScinHits))->SetCorrectedTimes(timec_pos,timec_neg, - postime, negtime, - scin_corrected_time); - fNScinHits++; - } - else { + // Do the pulse height correction to the time. (Position dependent corrections later) + Double_t timec_pos = tdc_pos*fScinTdcToTime + - fHodoPosInvAdcOffset[index] + - fHodoPosInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_pos)); + Double_t timec_neg = tdc_neg*fScinTdcToTime + - fHodoNegInvAdcOffset[index] + - fHodoNegInvAdcAdc[index]/TMath::Sqrt(TMath::Max(20.0,adc_neg)); + // Find hit position using ADCs + // If postime larger, then hit was nearer negative side. + // Some incarnations use fixed velocity of 15 cm/nsec + Double_t vellight=fHodoVelLight[index]; + Double_t dist_from_center=0.5*(timec_neg-timec_pos)*vellight; + Double_t scint_center=0.5*(fPosLeft+fPosRight); + Double_t hit_position=scint_center+dist_from_center; + hit_position=TMath::Min(hit_position,fPosLeft); + hit_position=TMath::Max(hit_position,fPosRight); + // Position depenent time corrections + timec_pos -= (fPosLeft-hit_position)/ + fHodoPosInvAdcLinear[index]; + timec_neg -= (hit_position-fPosRight)/ + fHodoNegInvAdcLinear[index]; + Double_t scin_corrected_time = 0.5*(timec_pos+timec_neg); + Double_t postime = timec_pos - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal); + Double_t negtime = timec_neg - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal); + + ((THcHodoHit*) fHodoHits->At(fNScinHits))->SetCorrectedTimes(timec_pos,timec_neg, + postime, negtime, + scin_corrected_time); + } + fNScinHits++; // One or more good time counter } - ihit++; + ihit++; // Raw hit counter } // cout << "THcScintillatorPlane: ihit = " << ihit << endl; diff --git a/src/THcScintillatorPlane.h b/src/THcScintillatorPlane.h index e9adbcbb3bd394e33210c5328203876b53a5077a..c83fa437ee61e4c80523ece3a05fe45cf264caa2 100644 --- a/src/THcScintillatorPlane.h +++ b/src/THcScintillatorPlane.h @@ -83,8 +83,8 @@ class THcScintillatorPlane : public THaSubDetector { Double_t fZpos; /* z position */ Double_t fDzpos; Double_t fHodoSlop; /* max allowed slop for this plane */ - Double_t fPosLeft; /* NOTE: "left" = "top" for a Y scintillator */ - Double_t fPosRight; /* NOTE: "right" = "bottom" for a Y scintillator */ + Double_t fPosLeft; /* NOTE: "left" = "bottom" for a Y scintillator */ + Double_t fPosRight; /* NOTE: "right" = "top" for a Y scintillator */ Double_t fPosOffset; Double_t *fPosCenter; /* array with centers for all scintillators in the plane */ Double_t fScinTdcMin; @@ -94,13 +94,19 @@ class THcScintillatorPlane : public THaSubDetector { Double_t fScinTdcToTime; Double_t fTofTolerance; Double_t fBetaNominal; - Double_t *fHodoPosMinPh; // Minimum pulse height per paddle for this plane - Double_t *fHodoNegMinPh; // Minimum pulse height per paddle for this plane - Double_t *fHodoPosPhcCoeff; // Pulse height to time coefficient per paddle for this plane - Double_t *fHodoNegPhcCoeff; // Pulse height to time coefficient per paddlefor this plane - Double_t *fHodoPosTimeOffset; - Double_t *fHodoNegTimeOffset; + // Double_t *fHodoPosMinPh; // Minimum pulse height per paddle for this plane + // Double_t *fHodoNegMinPh; // Minimum pulse height per paddle for this plane + // Double_t *fHodoPosPhcCoeff; // Pulse height to time coefficient per paddle for this plane + // Double_t *fHodoNegPhcCoeff; // Pulse height to time coefficient per paddlefor this plane + // Double_t *fHodoPosTimeOffset; + // Double_t *fHodoNegTimeOffset; Double_t *fHodoVelLight; + Double_t *fHodoPosInvAdcOffset; + Double_t *fHodoNegInvAdcOffset; + Double_t *fHodoPosInvAdcLinear; + Double_t *fHodoNegInvAdcLinear; + Double_t *fHodoPosInvAdcAdc; + Double_t *fHodoNegInvAdcAdc; Double_t *fHodoSigma; Double_t fTolerance; /* need this for Focal Plane Time estimation */