Skip to content

GitLab

Commit cf28cd7b authored by Jérémie Dudouet's avatar Jérémie Dudouet
Browse files

Merge branch 'HistoPrep' into 'preprod' 

add Gated Rates in preprocessing

See merge request IPNL_GAMMA/narval_emulator!49
parents 1a16614d 2ba4eeee
Hide whitespace changes
Inline Side-by-side
Showing with 107 additions and 99 deletions
consumers/Histogramer/AGAPRO_HistoDefs.h
View file @ cf28cd7b
......@@ -8,8 +8,8 @@
#define N_TH1HISTS_PROD 420
#define N_TH2HISTS_PROD 0
//38*ERaw + 38*Baseline + 38*T0 + 2*timestamps + 1*DeltaT
#define N_TH1HISTS_PREP 117
//38*ERaw + 38*Baseline + 38*T0 + 4*timestamps + 1*DeltaT
#define N_TH1HISTS_PREP 119
//1*EvsId
#define N_TH2HISTS_PREP 1
......
consumers/Histogramer/AGAPRO_Histogramer.h
View file @ cf28cd7b
......@@ -37,18 +37,18 @@ using namespace std;
namespace AGAPRO {
struct binning_TH1{
Int_t NBins;
Float_t Min;
Float_t Max;
Int_t NBins = 0;
Float_t Min = 0.;
Float_t Max = 0.;
};
struct binning_TH2{
Int_t NBinsX;
Float_t XMin;
Float_t XMax;
Int_t NBinsY;
Float_t YMin;
Float_t YMax;
Int_t NBinsX = 0;
Float_t XMin = 0.;
Float_t XMax = 0.;
Int_t NBinsY = 0;
Float_t YMin = 0.;
Float_t YMax = 0.;
};
class Histogramer : public ADF::NarvalConsumer
......
consumers/Histogramer/PreproHist/AGAPRO_PreproHist.cpp
View file @ cf28cd7b
......@@ -68,9 +68,11 @@ PreproHist::PreproHist() :
vector<void* > TimeStampList;
fListOfHists.push_back(TimeStampList);
vector<void* > GatedRatesList;
fListOfHists.push_back(GatedRatesList);
fNextCycleSec = 3600.;
fNextCycleMin = 720.;
fNextCycleH = 48.;
fTraceLengthPSA = defTraceLengthPSA;
......@@ -81,7 +83,6 @@ PreproHist::PreproHist() :
PreproHist::~PreproHist()
{
cServer.Finish();
for(int nn = 0; nn < kNSG; nn++) {
if(fTracesSG[nn])
delete [] fTracesSG[nn];
......@@ -90,7 +91,6 @@ PreproHist::~PreproHist()
if(fTracesCC[nn])
delete [] fTracesCC[nn];
}
// in principle not needed ... just in case reset it has not been called by narval
UInt_t error = 0u; process_reset(&error) ;
}
......@@ -103,29 +103,24 @@ UInt_t PreproHist::Process()
Log.SetPID(GetPID());
vector< void* > Amplilist = fListOfHists[(Int_t)HistLib::Ampli];
if(fPlotAmpli && Amplilist.size() == (kNSG+kNCC))
{
for(int nn = 0; nn < kNSG; nn++)
{
if(fPlotAmpli && Amplilist.size() == (kNSG+kNCC)) {
for(int nn = 0; nn < kNSG; nn++) {
if(SegE[nn]<1)
continue;
auto hist = (Hist1F*)Amplilist[nn];
FillHist1F(hist,SegE[nn]);
if(fPlotAmplivsId && fListOfHists[(Int_t)HistLib::Ampli_vs_Id].size()==1)
{
if(fPlotAmplivsId && fListOfHists[(Int_t)HistLib::Ampli_vs_Id].size()==1) {
auto hist2d = (Hist2F*)fListOfHists[(Int_t)HistLib::Ampli_vs_Id][0];
FillHist2F(hist2d,SegE[nn],nn);
}
}
for(int nn = 0; nn < kNCC; nn++)
{
for(int nn = 0; nn < kNCC; nn++) {
auto hist = (Hist1F*)Amplilist[kNSG+nn];
FillHist1F(hist,CoreE[nn]);
if(fPlotAmplivsId && fListOfHists[(Int_t)HistLib::Ampli_vs_Id].size()==1)
{
if(fPlotAmplivsId && fListOfHists[(Int_t)HistLib::Ampli_vs_Id].size()==1) {
auto hist2d = (Hist2F*)fListOfHists[(Int_t)HistLib::Ampli_vs_Id][0];
FillHist2F(hist2d,CoreE[nn],kNSG+nn);
}
......@@ -133,12 +128,10 @@ UInt_t PreproHist::Process()
}
vector< void* > Baselinelist = fListOfHists[(Int_t)HistLib::Baseline];
if(fPlotBaseline && Baselinelist.size() == (kNSG+kNCC))
{
if(fPlotBaseline && Baselinelist.size() == (kNSG+kNCC)) {
int baseLen = (fTraceLengthPSA+4)/5;
baseLen = max(1, baseLen);
for(int nn = 0; nn < kNSG; nn++)
{
for(int nn = 0; nn < kNSG; nn++) {
Float_t baseVal = 0;
Float_t *pW = fTracesSG[nn];
for(int ii = 0; ii < baseLen; ii++)
......@@ -149,8 +142,7 @@ UInt_t PreproHist::Process()
auto hist = (Hist1F*)Baselinelist[nn];
FillHist1F(hist,SegBase[nn]);
}
for(int nn = 0; nn < kNCC; nn++)
{
for(int nn = 0; nn < kNCC; nn++) {
Float_t baseVal = 0;
Float_t *pW = fTracesCC[nn];
for(int ii = 0; ii < baseLen; ii++)
......@@ -164,14 +156,11 @@ UInt_t PreproHist::Process()
}
vector< void* > T0list = fListOfHists[(Int_t)HistLib::T0];
if(fPlotT0 && T0list.size() == (kNSG+kNCC))
{
if(fPlotT0 && T0list.size() == (kNSG+kNCC)) {
vector< void* > T0list = fListOfHists[(Int_t)HistLib::T0];
int T0 = CalcT0();
if(T0>0)
{
if(T0>0) {
auto hist = (Hist1F*)T0list[netChargeSegs[0]];
FillHist1F(hist,T0); // relative time placed at the center of spectrum
......@@ -184,43 +173,57 @@ UInt_t PreproHist::Process()
}
vector< void* > DeltaTlist = fListOfHists[(Int_t)HistLib::DeltaT];
if(fPlotDeltaT && DeltaTlist.size() == 1 && fLastTS != 0)
{
if(fPlotDeltaT && DeltaTlist.size() == 1 && fLastTS != 0) {
auto hist = (Hist1F*)DeltaTlist[0];
FillHist1F(hist,timestamp-fLastTS);
}
vector< void* > TSlist = fListOfHists[(Int_t)HistLib::TimeStamps];
if(fPlotTimeStamps && TSlist.size() == 2)
{
if(fFirstTS == 0)
{
vector< void* > TSGatedlist = fListOfHists[(Int_t)HistLib::GatedRates];
if( (fPlotTimeStamps && TSlist.size() == 2) || (fPlotGatedRates && TSGatedlist.size() == 2)) {
if(fFirstTS == 0) {
fFirstTS = timestamp;
}
if(timestamp != 0 && timestamp>fLastTS)
{
auto hists = (Hist1F*)TSlist[0];
auto histm = (Hist1F*)TSlist[1];
Hist1F* hists = nullptr;
Hist1F* histm = nullptr;
Hist1F* hists_gated = nullptr;
Hist1F* histm_gated = nullptr;
if(fPlotTimeStamps) {
hists = (Hist1F*)TSlist[0];
histm = (Hist1F*)TSlist[1];
}
if(fPlotGatedRates) {
hists_gated = (Hist1F*)TSGatedlist[0];
histm_gated = (Hist1F*)TSGatedlist[1];
}
if(timestamp != 0 && timestamp>fLastTS) {
double TSinSec = fmod((timestamp-fFirstTS)*10*1e-9,3600);
double TSinMin = fmod((timestamp-fFirstTS)*10*1e-9/60.,720);
if((timestamp-fFirstTS)*10*1e-9 > fNextCycleSec)
{
ResetHist1F(hists);
if((timestamp-fFirstTS)*10*1e-9 > fNextCycleSec) {
if(hists) ResetHist1F(hists);
if(hists_gated) ResetHist1F(hists_gated);
fNextCycleSec += 3600.;
}
if((timestamp-fFirstTS)*10*1e-9/60. > fNextCycleMin)
{
ResetHist1F(histm);
if((timestamp-fFirstTS)*10*1e-9/60. > fNextCycleMin) {
if(histm) ResetHist1F(histm);
if(histm_gated) ResetHist1F(histm_gated);
fNextCycleMin += 720.;
}
FillHist1F(hists,TSinSec,(3600./3600.));
FillHist1F(histm,TSinMin,(3600./720.)/60.);
if(fPlotTimeStamps) {
FillHist1F(hists,TSinSec,(3600./3600.));
FillHist1F(histm,TSinMin,(3600./720.)/60.);
}
if(fPlotGatedRates && abs(CoreE[0]-fGatedRateEnergy)<fGatedRateWidth) {
FillHist1F(hists_gated,TSinSec,(3600./3600.));
FillHist1F(histm_gated,TSinMin,(3600./720.)/60.);
}
}
else
{
else {
Log.SetLevel(LogMessage::kWarning);
Log << "First TS : "<<fFirstTS<<" ; Last TS : "<<fLastTS<<" ; TS : "<<timestamp<< dolog;
}
......@@ -237,14 +240,13 @@ void PreproHist::process_initialise (UInt_t *error_code)
Log.SetPID(GetPID());
*error_code = GetParameters(GetConfPath() + gActualClass + ".conf");
if(*error_code)
{
if(*error_code) {
Log << dolog;
return;
}
// to get the input/output frames
if( fFrameCrystal ){
if( fFrameCrystal ) {
delete fFrameCrystal;
fFrameCrystal = NULL;
}
......@@ -253,8 +255,7 @@ void PreproHist::process_initialise (UInt_t *error_code)
std::cout << std::endl;
fFrameCrystal = fTrigger.Add("Agata", "data:ccrystal"); // discard input data
if( ! fFrameCrystal )
{
if( ! fFrameCrystal ) {
rerr |= 8;
}
......@@ -262,8 +263,7 @@ void PreproHist::process_initialise (UInt_t *error_code)
if( !fFrameIO.Register(&fTrigger) )
rerr |= 1;
if(rerr)
{
if(rerr) {
std::cout << "Trigger definition error " << rerr << " in PreproHist::process_initialise()" << std::endl;
*error_code = 100 + rerr;
Log << dolog;
......@@ -279,9 +279,7 @@ void PreproHist::process_initialise (UInt_t *error_code)
// state set to kIdle so that the data can be treated
fFrameIO.SetStatus(BaseFrameIO::kIdle);
cServer.SetCommandFile(GetConfPath() + gActualClass + ".live");
cServer.Start(gActualClass);
InitSMZ(crystal_id);
......@@ -292,7 +290,6 @@ void PreproHist::process_initialise (UInt_t *error_code)
Int_t PreproHist::SetInput()
{
Frame *frame_in = fFrameCrystal->GetFrame();
CrystalInterface *cdata = GetDataPointer<CrystalInterface>(frame_in);
// loads the values from the frame into the ADFObject attached to the frame
......@@ -301,16 +298,14 @@ Int_t PreproHist::SetInput()
// get the real length from the data (assuming all traces have the same length)
fTraceLengthPSA = cdata->GetCore(0)->GetSignal()->GetLength();
for(int iseg = 0; iseg < kNSG; iseg++)
{
for(int iseg = 0; iseg < kNSG; iseg++) {
seg = cdata->GetSegment(iseg);
SegE[iseg] = (Float_t)seg->GetE();
seg->GetSignal()->Get(fTracesSG[iseg], fTraceLengthPSA);
}
for(UShort_t icc=0; icc<kNCC; icc++)
{
for(UShort_t icc=0; icc<kNCC; icc++) {
core = cdata->GetCore(icc);
CoreE[icc] = (Float_t)core->GetE();
core->GetSignal()->Get(fTracesCC[icc], fTraceLengthPSA);
......@@ -342,6 +337,8 @@ UInt_t PreproHist::GetParameters(const std::string& confFile, Bool_t doList)
, &fDeltaTBinning.NBins, &fDeltaTBinning.Min, &fDeltaTBinning.Max);
conf.Add("PlotTimeStamps", "Plot Timestamp spectra, optionnal: define the ref T0 TimeStamp (1 TSRef)"
, &fPlotTimeStamps, &fFirstTS).SetRequired(1);
conf.Add("PlotGatedRates", "Plot Rates gated on the specified energy and width"
, &fGatedRateEnergy, &fGatedRateWidth);
if(doList) {conf.Show(gActualClass); return 0;}
......@@ -354,6 +351,7 @@ UInt_t PreproHist::GetParameters(const std::string& confFile, Bool_t doList)
if(fBaselineBinning.NBins>0) fPlotBaseline = true;
if(fT0Binning.NBins>0) fPlotT0 = true;
if(fDeltaTBinning.NBins>0) fPlotDeltaT = true;
if(fGatedRateEnergy>0.) fPlotGatedRates = true;
// To randomize the refresh time among the different crystals
fRefreshTime = fRefreshTime - 0.5 + ((double) rand() / (RAND_MAX));
......@@ -382,12 +380,10 @@ void PreproHist::InitSMZ(Int_t id)
map_ptr = static_cast<PrepMap*>(CreateNewSMz(mapname.str().data(),sizeof(PrepMap),PortNumber));
if(fPlotAmpli)
{
if(fPlotAmpli) {
vector<void *>& list = fListOfHists[(Int_t)HistLib::Ampli];
for(auto i=0 ; i<kNSG ; i++)
{
for(auto i=0 ; i<kNSG ; i++) {
Name<<"Prepro_Ampli_Cry" << id << "_Seg" << i;
Title<<"Prepro_Ampli_Cry"<<CrystalName.data()<<"_Seg"<<GetSegmentNameFromID(i).data();
XTitle<<"Energy (keV)";
......@@ -397,8 +393,7 @@ void PreproHist::InitSMZ(Int_t id)
Name.str("");Title.str("");XTitle.str("");GruName.str("");
if(hist1F) list.push_back(hist1F);
}
for(auto i=0 ; i<kNCC ; i++)
{
for(auto i=0 ; i<kNCC ; i++) {
Name<<"Prepro_Ampli_Cry" << id << "_Core" << i;
Title<<"Prepro_Ampli_Cry"<<CrystalName.data()<<"_"<<GetCoreNameFromID(i).data();
XTitle<<"Energy (keV)";
......@@ -430,8 +425,7 @@ void PreproHist::InitSMZ(Int_t id)
{
vector<void *>& list = fListOfHists[(Int_t)HistLib::Baseline];
for(auto i=0 ; i<kNSG ; i++)
{
for(auto i=0 ; i<kNSG ; i++) {
Name<<"Prepro_Baseline_Cry" << id << "_Seg" << i;
Title<<"Prepro_Baseline_Cry"<<CrystalName.data()<<"_Seg"<<GetSegmentNameFromID(i).data();
XTitle<<"Corrected Amplitude";
......@@ -441,8 +435,7 @@ void PreproHist::InitSMZ(Int_t id)
Name.str("");Title.str("");XTitle.str("");GruName.str("");
if(hist1F) list.push_back(hist1F);
}
for(auto i=0 ; i<kNCC ; i++)
{
for(auto i=0 ; i<kNCC ; i++) {
Name<<"Prepro_Baseline_Cry" << id << "_Core" << i;
Title<<"Prepro_Baseline_Cry"<<CrystalName.data()<<"_"<<GetCoreNameFromID(i).data();
XTitle<<"Corrected Amplitude";
......@@ -457,8 +450,7 @@ void PreproHist::InitSMZ(Int_t id)
{
vector<void *>& list = fListOfHists[(Int_t)HistLib::T0];
for(auto i=0 ; i<kNSG ; i++)
{
for(auto i=0 ; i<kNSG ; i++) {
Name<<"Prepro_T0_Cry" << id << "_Seg" << i;
Title<<"Prepro_T0_Cry"<<CrystalName.data()<<"_Seg"<<GetSegmentNameFromID(i).data();
XTitle<<"relative time compared to Core0";
......@@ -471,13 +463,11 @@ void PreproHist::InitSMZ(Int_t id)
for(auto i=0 ; i<kNCC ; i++)
{
Name<<"Prepro_T0_Cry" << id << "_Core" << i;
if(i==0)
{
if(i==0) {
Title<<"Prepro_T0_Cry"<<CrystalName.data()<<"_AllSegs";
XTitle<<"relative time compared to Core0 for all segments";
}
else
{
else {
Title<<"Prepro_T0_Cry"<<CrystalName.data()<<"_"<<GetCoreNameFromID(i).data();
XTitle<<"relative time of Core1 compared to Core0";
}
......@@ -489,21 +479,19 @@ void PreproHist::InitSMZ(Int_t id)
if(hist1F) list.push_back(hist1F);
}
}
if(fPlotDeltaT)
{
if(fPlotDeltaT) {
vector<void *>& list = fListOfHists[(Int_t)HistLib::DeltaT];
Name<<"Prepro_DeltaT_Cry" << id;
Title<<"Prepro_DeltaT_Cry"<<CrystalName.data();
XTitle<<"#DeltaT (in ts units) between two successive events";
XTitle<<"DeltaT (in ts units) between two successive events";
GruName<<"Prepro/DeltaT/Cry"<<CrystalName.data();
Hist1F *hist1F = NewHist1F(Name.str(),Title.str(),XTitle.str(),GruName.str(),fDeltaTBinning.NBins,fDeltaTBinning.Min,fDeltaTBinning.Max,map_ptr,PortNumber);
Name.str("");Title.str("");XTitle.str("");GruName.str("");
if(hist1F) list.push_back(hist1F);
}
if(fPlotTimeStamps)
{
if(fPlotTimeStamps) {
vector<void *>& list = fListOfHists[(Int_t)HistLib::TimeStamps];
Name<<"Prepro_TimeStamp_Cry" << id << "_1h";
......@@ -520,8 +508,28 @@ void PreproHist::InitSMZ(Int_t id)
XTitle<<"Elapsed time on the last 12 hour (min)";
hist1F = NewHist1F(Name.str(),Title.str(),XTitle.str(),GruName.str(),3600,0,720,map_ptr,PortNumber);
Name.str("");Title.str("");XTitle.str("");GruName.str("");
if(hist1F) list.push_back(hist1F);
}
if(fPlotGatedRates) {
vector<void *>& list = fListOfHists[(Int_t)HistLib::GatedRates];
Name<<"Prepro_GatedRates_Cry" << id << "_1h";
Title<<"Prepro_GR_Cry"<<CrystalName.data()<<"_1hour (1bin=1s)";
XTitle<<"Elapsed time on the last hour (sec)";
GruName<<"Prepro/GatedRates/Cry"<<CrystalName.data();
Hist1F *hist1F = NewHist1F(Name.str(),Title.str(),XTitle.str(),GruName.str(),3600,0,3600,map_ptr,PortNumber);
Name.str("");Title.str("");XTitle.str("");
if(hist1F) list.push_back(hist1F);
Name<<"Prepro_GatedRates_Cry" << id << "_12h";
Title<<"Prepro_GR_Cry"<<CrystalName.data()<<"_12hour (1bin=12s)";
XTitle<<"Elapsed time on the last 12 hour (min)";
hist1F = NewHist1F(Name.str(),Title.str(),XTitle.str(),GruName.str(),3600,0,720,map_ptr,PortNumber);
Name.str("");Title.str("");XTitle.str("");GruName.str("");
if(hist1F) list.push_back(hist1F);
}
}
......@@ -598,11 +606,10 @@ int PreproHist::CalcT0()
//// calculate the real fold of segments
nSegFold = 0;
float sumSegs = 0;
for(int ns = 0; ns < kNSG; ns++)
{
for(int ns = 0; ns < kNSG; ns++) {
float ener = SegE[ns];
if(ener > 10.) // but it is better to use the individual thresholds
{
// but it is better to use the individual thresholds
if(ener > 10.) {
sumSegs += ener;
netChargeSegs[nSegFold++] = ns;
}
......@@ -629,11 +636,9 @@ int PreproHist::CalcT0()
float tCCSG = (tRefSG-cRefCC)*tScale*tGain;
int relative_time = int(int(tCCSG + speTTlen/2));
CoreT[0] = relative_time; //sum of all the segments
if(tRefSG > 1.f)
return relative_time; // relative time placed at the center of spectrum
else
......
consumers/Histogramer/PreproHist/AGAPRO_PreproHist.h
View file @ cf28cd7b
......@@ -56,17 +56,20 @@ protected: //Evt properties
Float_t CoreT[kNCC];
Float_t DeltaT;
Float_t fGatedRateEnergy = 0.;
Float_t fGatedRateWidth = 0.;
private:
bool fPlotAmpli = false;
bool fPlotAmplivsId = false;
bool fPlotBaseline = false;
bool fPlotT0 = false;
bool fPlotTimeStamps = false;
bool fPlotGatedRates = false;
bool fPlotDeltaT = false;
Float_t fNextCycleSec;
Float_t fNextCycleMin;
Float_t fNextCycleH;
binning_TH1 fAmpliBinning;
binning_TH1 fBaselineBinning;
......@@ -77,7 +80,7 @@ private:
Float_t fMinSegAmpli; // minimum segment amplitude to detect segment multiplicity
enum class HistLib: int { Ampli, Ampli_vs_Id, Baseline, T0, TimeStamps, DeltaT };
enum class HistLib: int { Ampli, Ampli_vs_Id, Baseline, T0, TimeStamps, GatedRates, DeltaT };
/// For T0 calculation (from PreprocessingFilterPSA.cpp
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment