stream  0.10.0
stream analysis framework
Example with second.C and reading of created ROOT file

Here normal first.C

#include "base/ProcMgr.h"
#include "hadaq/HldProcessor.h"
#include "hadaq/TdcProcessor.h"
#include "hadaq/TrbProcessor.h"
void first()
{
// base::ProcMgr::instance()->SetRawAnalysis(true);
// this limits used for liner calibrations when nothing else is available
// default channel numbers and edges mask
// 1 - use only rising edge, falling edge is ignore
// 2 - falling edge enabled and fully independent from rising edge
// 3 - falling edge enabled and uses calibration from rising edge
// 4 - falling edge enabled and common statistic is used for calibration
// About time calibration - there are two possibilities
// 1) automatic calibration after N hits in every enabled channel.
// Just use SetAutoCalibrations method for this
// 2) generate calibration on base of provided data and than use it later statically for analysis
// Than one makes special run with SetWriteCalibrations() enabled.
// Later one reuse such calibrations enabling only LoadCalibrations() call
hadaq::TrbProcessor* trb3_1 = new hadaq::TrbProcessor(0x8000, hld);
trb3_1->SetHistFilling(4);
trb3_1->SetCrossProcess(true);
trb3_1->CreateTDC(0x1000, 0x1001, 0x1002, 0x1003);
// enable automatic calibration, specify required number of hits in each channel
// trb3_1->SetAutoCalibrations(50000);
// calculate and write static calibration at the end of the run
// trb3_1->SetWriteCalibrations("run1");
// load static calibration at the beginning of the run
trb3_1->LoadCalibrations("run1");
//hadaq::TrbProcessor* trb3_2 = new hadaq::TrbProcessor(0x8001, hld);
//trb3_2->SetHistFilling(4);
//trb3_2->SetCrossProcess(true);
//trb3_2->CreateTDC(0x1010, 0x1011, 0x1012, 0x1013);
//trb3_2->SetAutoCalibrations(50000);
//trb3_2->SetWriteCalibrations("run1");
//trb3_2->LoadCalibrations("run1");
// this is array with available TDCs ids
int tdcmap[4] = { 0x1000, 0x1001, 0x1002, 0x1003 };
// TDC subevent header id
for (int cnt=0;cnt<4;cnt++) {
hadaq::TdcProcessor* tdc = hld->FindTDC(tdcmap[cnt]);
if (tdc==0) continue;
// enable storage of first channel, which contains absolute time
// all other channels times are relative to the first channel
tdc->SetCh0Enabled(true);
// specify reference channel
for (int ch=2;ch<=32;ch++)
tdc->SetRefChannel(ch, ch-1, 0xffff, 2000, -10., 90., true);
// specify reference channel from other TDC
//if (cnt > 0) {
// tdc->SetRefChannel(0, 0, 0xc000, 2000, -10., 10., true);
// tdc->SetRefChannel(5, 5, 0xc000, 2000, -10., 10., true);
// tdc->SetRefChannel(7, 7, 0xc000, 2000, -10., 10., true);
// }
// for old FPGA code one should have epoch for each hit, no longer necessary
// tdc->SetEveryEpoch(true);
// When enabled, time of last hit will be used for reference calculations
// tdc->SetUseLastHit(true);
}
}
void SetTriggeredAnalysis(bool on=true)
Enabled/disable triggered analysis is configured.
Definition: ProcMgr.h:192
virtual bool CreateStore(const char *storename)
create data store, for the moment - ROOT tree
Definition: ProcMgr.h:162
static ProcMgr * instance()
Return global instance of processor manager, provided by framework.
Definition: ProcMgr.cxx:46
void SetStoreEnabled(bool on=true)
Enable store - set store kind 1.
Definition: base/Processor.h:239
void SetHistFilling(int lvl=99)
Set histogram filling level.
Definition: base/Processor.h:225
HLD processor.
Definition: HldProcessor.h:93
TdcProcessor * FindTDC(unsigned tdcid) const
Find TDC by id.
Definition: HldProcessor.cxx:136
TDC processor.
Definition: TdcProcessor.h:32
void SetCh0Enabled(bool on=true)
Enable/disable store of channel 0 in output event.
Definition: TdcProcessor.h:708
void SetRefChannel(unsigned ch, unsigned refch, unsigned reftdc=0xffff, int npoints=5000, double left=-10., double right=10., bool twodim=false)
Set reference signal for the TDC channel ch.
Definition: TdcProcessor.cxx:516
TRB processor.
Definition: TrbProcessor.h:46
bool LoadCalibrations(const char *fileprefix)
Load calibrations for all existing TDCs as argument file prefix (without TDC id) should be specified.
Definition: TrbProcessor.cxx:300
static void SetDefaults(unsigned numch=65, unsigned edges=0x1, bool ignore_sync=true)
Set defaults for the next creation of TDC processors.
Definition: TrbProcessor.cxx:35
int CreateTDC(unsigned id1, unsigned id2=0, unsigned id3=0, unsigned id4=0)
Create up to 4 TDC with pre-configured default parameters.
Definition: TrbProcessor.cxx:222
void SetCrossProcess(bool on=true)
Enable Cross-processing - hits correlation between different TDCs.
Definition: TrbProcessor.cxx:522
static void SetFineLimits(unsigned min, unsigned max)
Method set static limits, which are used for simple interpolation of time for fine counter.
Definition: TdcMessage.h:294

Then second.C

#include <cstdio>
#include "TTree.h"
#include "base/EventProc.h"
#include "base/Event.h"
#include "hadaq/TdcSubEvent.h"
class SecondProc : public base::EventProc {
protected:
std::string fTdcId;
double fHits[8];
base::H1handle hNumHits;
public:
SecondProc(const char* procname, const char* _tdcid) :
base::EventProc(procname),
fTdcId(_tdcid),
hNumHits(0)
{
printf("Create %s for %s\n", GetName(), fTdcId.c_str());
hNumHits = MakeH1("NumHits","Number of hits", 32, 0, 32, "number");
// enable storing already in constructor
}
virtual void CreateBranch(TTree* t)
{
// only called when tree is created in first.C
// one can ignore
t->Branch(GetName(), fHits, "hits[8]/D");
}
virtual bool Process(base::Event* ev)
{
for (unsigned n=0;n<8;n++) fHits[n] = 0.;
dynamic_cast<hadaq::TdcSubEvent*> (ev->GetSubEvent(fTdcId));
// printf("%s process sub %p %s\n", GetName(), sub, fTdcId.c_str());
if (sub==0) return false;
double num(0), ch0tm(0);
for (unsigned cnt=0;cnt<sub->Size();cnt++) {
const hadaq::TdcMessageExt& ext = sub->msg(cnt);
unsigned chid = ext.msg().getHitChannel();
if (chid==0) { ch0tm = ext.GetGlobalTime(); continue; }
// full time
double tm = ext.GetGlobalTime() + ch0tm;
// printf(" ch:%3d tm:%f\n", chid, tm);
num+=1;
}
FillH1(hNumHits, num);
return true;
}
};
void second()
{
new SecondProc("Second1", "TDC_1000");
}
Abstract processor of build events.
Definition: EventProc.h:12
virtual bool Process(Event *)
Generic event processing If returns false, processing will be aborted and event will not be stored.
Definition: EventProc.h:26
Event - collection of several subevents.
Definition: Event.h:17
base::SubEvent * GetSubEvent(const std::string &name) const
Return subevent by name.
Definition: Event.h:69
Extended message - any message plus global time stamp.
Definition: base/SubEvent.h:36
const MsgClass & msg() const
message
Definition: base/SubEvent.h:81
double GetGlobalTime() const
global time stamp
Definition: base/SubEvent.h:84
H1handle MakeH1(const char *name, const char *title, int nbins, double left, double right, const char *xtitle=0)
Adds processor prefix to histogram name and calls base::ProcMgr::MakeH1 method.
Definition: base/Processor.cxx:132
virtual void CreateBranch(TTree *)
Create branch.
Definition: base/Processor.h:204
const char * GetName() const
Get processor name.
Definition: base/Processor.h:220
void FillH1(H1handle h1, double x, double weight=1.)
Fill 1-D histogram.
Definition: base/Processor.h:103
Subevent with vector of extended messages.
Definition: base/SubEvent.h:91
unsigned Size() const
Returns number of messages.
Definition: base/SubEvent.h:107
MsgClass & msg(unsigned indx)
Returns message with specified index.
Definition: base/SubEvent.h:116

And script read.C to read data from created ROOT file

// macro reads tree with data, stored by TDC processor
#include "TTree.h"
#include "TFile.h"
void read()
{
TFile* f = TFile::Open("tree.root");
if (f==0) return;
TTree* t = (TTree*) f->Get("T");
if (t==0) return;
std::vector<hadaq::TdcMessageExt> *vect = new std::vector<hadaq::TdcMessageExt>;
t->SetBranchAddress("TDC_1000", &vect);
Long64_t entry(0), total(0);
while (t->GetEntry(entry)>0) {
entry++;
double ch0tm(0);
for (unsigned n=0;n<vect->size();n++) {
hadaq::TdcMessageExt& ext = vect->at(n);
int chid = ext.msg().getHitChannel();
if (chid==0) { ch0tm = ext.GetGlobalTime(); continue; }
double tm = ext.GetGlobalTime() + ch0tm;
// printf(" ch:%3d tm %10.9f\n", chid, tm);
total += 1;
}
}
printf("Read %ld entries %ld messages\n", entry, total);
}