From e9d869994028852b24c23a2024a8c0ec27a8d43b Mon Sep 17 00:00:00 2001
From: sanghwapark <sanghwapark@gmail.com>
Date: Mon, 29 Jan 2018 09:43:30 -0500
Subject: [PATCH] Add new physics module to check BCM current for each event
---
Makefile | 1 +
src/THcBCMCurrent.cxx | 225 ++++++++++++++++++++++++++++++++++++++++++
src/THcBCMCurrent.h | 63 ++++++++++++
3 files changed, 289 insertions(+)
create mode 100644 src/THcBCMCurrent.cxx
create mode 100644 src/THcBCMCurrent.h
diff --git a/Makefile b/Makefile
index b0f8576..a3f35cb 100644
--- a/Makefile
+++ b/Makefile
@@ -38,6 +38,7 @@ SRC = src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \
src/THcHodoEff.cxx \
src/THcTrigApp.cxx src/THcTrigDet.cxx src/THcTrigRawHit.cxx \
src/THcRawAdcHit.cxx src/THcRawTdcHit.cxx \
+ src/THcBCMCurrent.cxx
src/THcDummySpectrometer.cxx
endif
diff --git a/src/THcBCMCurrent.cxx b/src/THcBCMCurrent.cxx
new file mode 100644
index 0000000..b0558a2
--- /dev/null
+++ b/src/THcBCMCurrent.cxx
@@ -0,0 +1,225 @@
+/*
+
+ This physics module does:
+ - Read average BCM beam current values from scaler parameter file.
+ - Write the values into bcm#.AvgCurrent for each event
+ - Compare the current values with the threshold and
+ set event flags (BCM1 and BCM2 only)
+
+ You can set the threshold using SetCurrentCut
+ instead of gBCM_Current_threshold
+
+ */
+
+#include "THcParmList.h"
+#include "THcGlobals.h"
+#include "THcHitList.h"
+
+#include "THcBCMCurrent.h"
+
+using namespace std;
+
+THcBCMCurrent::THcBCMCurrent(const char* name,
+ const char* description) :
+ THaPhysicsModule(name, description), fThreshold(0)
+{
+
+ fBCM1flag = 0;
+ fBCM2flag = 0;
+
+ fBCM1avg = 0;
+ fBCM2avg = 0;
+ fBCM4aavg = 0;
+ fBCM4bavg = 0;
+ fBCM17avg = 0;
+}
+
+//__________________________________________________
+
+THcBCMCurrent::~THcBCMCurrent()
+{
+
+ DefineVariables (kDelete);
+
+ delete [] fiBCM1; fiBCM1 = NULL;
+ delete [] fiBCM2; fiBCM2 = NULL;
+ delete [] fiBCM4a; fiBCM4a = NULL;
+ delete [] fiBCM4b; fiBCM4b = NULL;
+ delete [] fiBCM17; fiBCM17 = NULL;
+ delete [] fEvtNum; fEvtNum = NULL;
+
+}
+
+//__________________________________________________
+
+THaAnalysisObject::EStatus THcBCMCurrent::Init( const TDatime& date )
+{
+
+
+ if( THaPhysicsModule::Init( date ) != kOK )
+ return fStatus;
+
+ return fStatus = kOK;
+}
+
+//__________________________________________________
+
+Int_t THcBCMCurrent::ReadDatabase( const TDatime& date )
+{
+
+ DBRequest list1[] = {
+ {"num_scal_reads", &fNscaler, kInt},
+ {0}
+ };
+ gHcParms->LoadParmValues((DBRequest*)&list1);
+
+ if(fThreshold < 1.e-5)
+ {
+ DBRequest list2[] = {
+ {"gBCM_Current_threshold", &fThreshold, kDouble},
+ {0}
+ };
+ gHcParms->LoadParmValues((DBRequest*)&list2);
+ }
+
+ fiBCM1 = new Double_t[fNscaler];
+ fiBCM2 = new Double_t[fNscaler];
+ fiBCM4a = new Double_t[fNscaler];
+ fiBCM4b = new Double_t[fNscaler];
+ fiBCM17 = new Double_t[fNscaler];
+
+ fEvtNum = new Int_t[fNscaler];
+
+ DBRequest list3[] = {
+ {"scal_read_bcm1_current", fiBCM1, kDouble, (UInt_t) fNscaler},
+ {"scal_read_bcm2_current", fiBCM2, kDouble, (UInt_t) fNscaler},
+ {"scal_read_bcm4a_current", fiBCM4a, kDouble, (UInt_t) fNscaler},
+ {"scal_read_bcm4b_current", fiBCM4b, kDouble, (UInt_t) fNscaler},
+ {"scal_read_bcm17_current", fiBCM17, kDouble, (UInt_t) fNscaler},
+ {"scal_read_event", fEvtNum, kInt, (UInt_t) fNscaler},
+ {0}
+ };
+
+ gHcParms->LoadParmValues((DBRequest*)&list3);
+
+ BCMInfo binfo;
+ for(int i=0; i<fNscaler; i++)
+ {
+ binfo.bcm1_current = fiBCM1[i];
+ binfo.bcm2_current = fiBCM2[i];
+ binfo.bcm4a_current = fiBCM4a[i];
+ binfo.bcm4b_current = fiBCM4b[i];
+ binfo.bcm17_current = fiBCM17[i];
+
+ BCMInfoMap.insert( std::make_pair(fEvtNum[i], binfo) );
+ }
+
+ return kOK;
+
+}
+
+//__________________________________________________
+
+Int_t THcBCMCurrent::DefineVariables( EMode mode )
+{
+
+ if( mode == kDefine && fIsSetup ) return kOK;
+ fIsSetup = ( mode == kDefine );
+
+ RVarDef vars[] = {
+ {"bcm1.currentflag", "BCM1 current flag for good event", "fBCM1flag"},
+ {"bcm2.currentflag", "BCM2 current flag for good event", "fBCM2flag"},
+ {"bcm1.AvgCurrent", "BCM1 average beam current", "fBCM1avg"},
+ {"bcm2.AvgCurrent", "BCM2 average beam current", "fBCM2avg"},
+ {"bcm4a.AvgCurrent", "BCM4a average beam current", "fBCM4aavg"},
+ {"bcm4b.AvgCurrent", "BCM4b average beam current", "fBCM4bavg"},
+ {"bcm17.AvgCurrent", "BCM17 average beam current", "fBCM17avg"},
+ { 0 }
+ };
+
+ return DefineVarsFromList(vars, mode);
+
+}
+
+//__________________________________________________
+
+Int_t THcBCMCurrent::Process( const THaEvData& evdata )
+{
+
+ if( !IsOK() ) return -1;
+
+ int fEventNum = evdata.GetEvNum();
+
+ BCMInfo binfo;
+ Int_t fGetScaler = GetAvgCurrent( fEventNum, binfo );
+
+ if(fGetScaler != kOK)
+ {
+ fBCM1avg = 0;
+ fBCM2avg = 0;
+ fBCM4aavg = 0;
+ fBCM4bavg = 0;
+ fBCM17avg = 0;
+ }
+ else
+ {
+ fBCM1avg = binfo.bcm1_current;
+ fBCM2avg = binfo.bcm2_current;
+ fBCM4aavg = binfo.bcm4a_current;
+ fBCM4bavg = binfo.bcm4b_current;
+ fBCM17avg = binfo.bcm17_current;
+ }
+
+ if(fBCM1avg < fThreshold)
+ fBCM1flag = 0;
+ else
+ fBCM1flag = 1;
+
+ if(fBCM2avg < fThreshold)
+ fBCM2flag = 0;
+ else
+ fBCM2flag = 1;
+
+ return kOK;
+
+}
+
+//__________________________________________________
+
+Int_t THcBCMCurrent::GetAvgCurrent( Int_t fevn, BCMInfo &bcminfo )
+{
+
+ map<int, BCMInfo>::iterator it, next;
+ it = BCMInfoMap.find(fevn);
+ if( it != BCMInfoMap.end() )
+ {
+ bcminfo.bcm1_current = it->second.bcm1_current;
+ bcminfo.bcm2_current = it->second.bcm2_current;
+ bcminfo.bcm4a_current = it->second.bcm4a_current;
+ bcminfo.bcm4b_current = it->second.bcm4b_current;
+ bcminfo.bcm17_current = it->second.bcm17_current;
+
+ return kOK;
+ }
+ else
+ {
+ next = BCMInfoMap.upper_bound(fevn);
+ if( next != BCMInfoMap.end() )
+ {
+ bcminfo.bcm1_current = next->second.bcm1_current;
+ bcminfo.bcm2_current = next->second.bcm2_current;
+ bcminfo.bcm4a_current = next->second.bcm4a_current;
+ bcminfo.bcm4b_current = next->second.bcm4b_current;
+ bcminfo.bcm17_current = next->second.bcm17_current;
+ return kOK;
+ }
+ return kOK+1;
+ }
+
+ return kOK+1;
+
+}
+
+//__________________________________________________
+
+ClassImp(THcBCMCurrent)
diff --git a/src/THcBCMCurrent.h b/src/THcBCMCurrent.h
new file mode 100644
index 0000000..4357a30
--- /dev/null
+++ b/src/THcBCMCurrent.h
@@ -0,0 +1,63 @@
+#ifndef __THCBCMCURRENT_H__
+#define __THCBCMCURRENT_H__
+
+#include "THaPhysicsModule.h"
+#include "THaEvData.h"
+#include "VarDef.h"
+#include "VarType.h"
+
+#include <iostream>
+#include <map>
+
+class THcBCMCurrent : public THaPhysicsModule {
+
+ public:
+
+ THcBCMCurrent(const char* name, const char* description);
+ virtual ~THcBCMCurrent();
+
+ virtual EStatus Init( const TDatime& date);
+ virtual Int_t Process( const THaEvData& );
+
+ void SetCurrentCut(Double_t _threshold){ fThreshold = _threshold; }
+
+ private:
+
+ Int_t fNscaler;
+ Double_t fThreshold;
+ Double_t* fiBCM1;
+ Double_t* fiBCM2;
+ Double_t* fiBCM4a;
+ Double_t* fiBCM4b;
+ Double_t* fiBCM17;
+ Int_t* fEvtNum;
+
+ Int_t fBCM1flag;
+ Int_t fBCM2flag;
+
+ Double_t fBCM1avg;
+ Double_t fBCM2avg;
+ Double_t fBCM4aavg;
+ Double_t fBCM4bavg;
+ Double_t fBCM17avg;
+
+ struct BCMInfo{
+ Double_t bcm1_current;
+ Double_t bcm2_current;
+ Double_t bcm4a_current;
+ Double_t bcm4b_current;
+ Double_t bcm17_current;
+ };
+
+ Int_t GetAvgCurrent( Int_t fevn, BCMInfo &bcminfo );
+ virtual Int_t ReadDatabase( const TDatime& date);
+ virtual Int_t DefineVariables( EMode mode = kDefine );
+
+
+ std::map<Int_t, BCMInfo> BCMInfoMap;
+
+ ClassDef(THcBCMCurrent, 0)
+
+};
+
+#endif
--
GitLab