diff --git a/src/THcHelicity.cxx b/src/THcHelicity.cxx index 2bb741882b57a0d6ba034038310f897a4d94d98b..6cfdd36ccedb3aa35733b07b1975b2e2fd7a0ccd 100644 --- a/src/THcHelicity.cxx +++ b/src/THcHelicity.cxx @@ -73,6 +73,9 @@ THaAnalysisObject::EStatus THcHelicity::Init(const TDatime& date) { return fStatus; } + fEvNumCheck = 0; + fDisabled = kFALSE; + fStatus = kOK; return fStatus; } @@ -140,6 +143,9 @@ Int_t THcHelicity::ReadDatabase( const TDatime& date ) // Initialize mode to find quartets and then seed } + cout << "Helicity decoder initialized with frequency of " << fFreq + << " Hz and reporting delay of " << fHelDelay << " cycles." << endl; + return kOK; } @@ -256,11 +262,32 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) fActualHelicity = fIsMPS?kUnknown:fReportedHelicity; return 0; } - Long64_t lastlastmpstime = fLastMPSTime; + + if(fDisabled) { + fActualHelicity = kUnknown; + return 0; + } + + fEvNumCheck++; + Int_t evnum = evdata.GetEvNum(); + if(fEvNumCheck!=evnum) { + cout << "THcHelicity: Missed " << evnum-fEvNumCheck << " events at event " << evnum << endl; + cout << " Disabling helicity decoding for rest of run." << endl; + cout << " Make sure \"RawDecode_master in cuts file accepts all physics events." <<endl; + fDisabled = kTRUE; + fActualHelicity = kUnknown; + return 0; + } + + fActualHelicity = -10.0; if(fFirstEvProcessed) { // Normal processing + // cout << evnum << " " << fNCycle << " " << fIsMPS << " " << fFoundMPS << " " << fTITime << " " + // << fLastMPSTime << " " << fNBits << endl; Int_t missed = 0; // Double_t elapsed_time = (fTITime - fFirstEvTime)/250000000.0; if(fIsMPS) { + fActualHelicity = kUnknown; + fPredictedHelicity = kUnknown; if(fFoundMPS) { missed = TMath::Nint(fTITime/fTIPeriod-fLastMPSTime/fTIPeriod); if(missed < 1) { // was <=1 @@ -292,12 +319,6 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) } int quartetphase = (newNCycle-fFirstCycle)%4; // cout << " " << fNCycle << " " << newNCycle << " " << fFirstCycle << " " << quartets_missed << " " << quartetphase << endl; - fQuartetStartHelicity = (fRingSeed_actual&1)?kPlus:kMinus; - fQuartetStartPredictedHelicity = (fRingSeed_reported&1)?kPlus:kMinus; - fActualHelicity = (quartetphase==0||quartetphase==3)? - fQuartetStartHelicity:-fQuartetStartHelicity; - fPredictedHelicity = (quartetphase==0||quartetphase==3)? - fQuartetStartPredictedHelicity:-fQuartetStartPredictedHelicity; // cout << "Cycles " << fNCycle << " " << newNCycle << " " << fFirstCycle // << " skipped " << quartets_missed << " quartets" << endl; fNCycle = newNCycle; @@ -305,6 +326,13 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) // reported helicity. So we don't fail quartet testing. // But only do this if we are calibrated. if(fNBits >= fMAXBIT) { + fQuartetStartHelicity = (fRingSeed_actual&1)?kPlus:kMinus; + fQuartetStartPredictedHelicity = (fRingSeed_reported&1)?kPlus:kMinus; + fActualHelicity = (quartetphase==0||quartetphase==3)? + fQuartetStartHelicity:-fQuartetStartHelicity; + fPredictedHelicity = (quartetphase==0||quartetphase==3)? + fQuartetStartPredictedHelicity:-fQuartetStartPredictedHelicity; + if (((fNCycle - fFirstCycle)%2)==1) { fQuartet[0] = fReportedHelicity; fQuartet[1] = fQuartet[2] = -fQuartet[0]; @@ -313,6 +341,7 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) fQuartet[2] = -fQuartet[1]; } } else { + fActualHelicity = kUnknown; fQuartet[0] = fReportedHelicity; fQuartet[1] = 0; } @@ -320,6 +349,18 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) fLastMPSTime += missed*fTIPeriod; fIsNewCycle = kTRUE; fLastReportedHelicity = fReportedHelicity; + } else { // No missed periods. Get helicities from rings + if(fNBits>=fMAXBIT) { + int quartetphase = (fNCycle-fFirstCycle)%4; + fQuartetStartHelicity = (fRingSeed_actual&1)?kPlus:kMinus; + fQuartetStartPredictedHelicity = (fRingSeed_reported&1)?kPlus:kMinus; + fActualHelicity = (quartetphase==0||quartetphase==3)? + fQuartetStartHelicity:-fQuartetStartHelicity; + fPredictedHelicity = (quartetphase==0||quartetphase==3)? + fQuartetStartPredictedHelicity:-fQuartetStartPredictedHelicity; + } else { + fActualHelicity = 0; + } } if(fIsNewCycle) { fQuartet[3]=fQuartet[2]; fQuartet[2]=fQuartet[1]; fQuartet[1]=fQuartet[0]; @@ -361,6 +402,8 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) // << fPredictedHelicity << " " << fActualHelicity << endl; } // Ignore until a MPS Is found + } else { // No MPS found yet + fActualHelicity = kUnknown; } } else { cout << "Initializing" << endl; @@ -377,10 +420,17 @@ Int_t THcHelicity::Decode( const THaEvData& evdata ) fIsNewCycle = kFALSE; fNBits = 0; } - // cout << setprecision(9) << "HEL " << fTITime/250000000.0 << " " << fNCycle << "(" << (fNCycle-fFirstCycle)%4 << "): " - // << fMPS << " " << fReportedHelicity << " " - // << fPredictedHelicity << " " << fActualHelicity << " " << lastlastmpstime/250000000.0 << endl; - // bitset<32>(v) + + // Some sanity checks + if(fActualHelicity < -5) { + cout << "Actual Helicity never got defined" << endl; + } + if(fNBits < fMAXBIT) { + if(fActualHelicity == -1 || fActualHelicity == 1) { + cout << "Helicity of " << fActualHelicity << " reported prematurely at cycle " << fNCycle << endl; + } + } + fLastActualHelicity = fActualHelicity; return 0; } diff --git a/src/THcHelicity.h b/src/THcHelicity.h index adc598ecad469c0f9469bc7e2a6944eb81bbc6db..cbfa7f51a12deb7d8c32b834046f3809089ad0fe 100644 --- a/src/THcHelicity.h +++ b/src/THcHelicity.h @@ -94,6 +94,9 @@ protected: Double_t fErrorCode; Int_t fEvtype; // Current CODA event type + Int_t fLastActualHelicity; + Int_t fEvNumCheck; + Bool_t fDisabled; static const Int_t NHIST = 2; TH1F* fHisto[NHIST];