// ESAF : Euso Simulation and Analysis Framework
// $Id: ChipTrackingTrgEngine1.cc,v 1.12 2005/11/11 09:08:06 pesce Exp $
// R. Pesce created
//_____________________________________________________________________________
//
// Chip Tracking Trigger Engine 1
// ==============================
//
// Tracking trigger built in front end chip.
// The algorithm searches for sequences of pixels (over the threshold)
// consequtive and adjacent. If the track length is over a threshold,
// trigger occurs.
// The algorithm searches also between track segments of two nearby chips:
// if two track segments are consequtive in time and the total length is
// over a threshold (different from the previous one), trigger occurs
//
// Parameters:
//
// fThresholdType: specify if the threshold on pixel hits in a gtu is absolute
// or relative
// Valid options:
// - absolute: the value specified in fThreshold is the minimum number of hits
// in a gtu to consider
// - relative: trigger threshold depends by the mean value of the background
// hits (B) in each pixel:
// threshold = (B) + fThreshold*sqrt(B)
// fThreshold: threshold on pixel hits in a gtu; if fThresholdType=relative the
// true threshold used by trigger algorithm depends by the mean hits
// of background in each pixel (see above)
// fMinTrackLenght: min length of the valid tracks
// fMaxTrackLength: max length of the valid tracks
// fMinTriggerTrackLength: min track length to have a trigger (single chip)
// fMinTriggerTwoLength: min total track length in 2 nearby chips
// fMaxTwoLength: max track length allowed in 2 nearby chips
//
// fAcceptHole [bool] : an hole in a single track segment is or not accepted
// fOnlyWithSignal [bool] : processes or not only chips ith at least a signal
//
// Special multiple configuration:
// The alorithm counts the numebr of tracks of the various trigger lengths,
// fill a specific trigger word and dump the numbers of the recognized tracks
//
#include "ChipTrackingTrgEngine1.hh"
#include "MacroCellData.hh"
#include "MacroCell.hh"
#include "ElementaryCell.hh"
#include "FrontEndChip.hh"
#include "BoolAlgebra.hh"
#include "ChipTrackSegment.hh"
#include "EusoElectronics.hh"
#include "Photomultiplier.hh"
#include "EusoDetector.hh"
#include "TMath.h"
#include "ERunParameters.hh"
#include "EChipTrackTriggerDataAdder.hh"
#include "EChipTriggParsFiller.hh"
#include "Config.hh"
using namespace TMath;
// trigger pattern types
const Int_t kTrigPattern[32] = {BIT(0),BIT(1),BIT(2),BIT(3),BIT(4),BIT(5),BIT(6),BIT(7),
BIT(8),BIT(9),BIT(10),BIT(11),BIT(12),BIT(13),BIT(14),BIT(15),
BIT(16),BIT(17),BIT(18),BIT(19),BIT(20),BIT(21),BIT(22),BIT(23),
BIT(24),BIT(25),BIT(26),BIT(27),BIT(28),BIT(29),BIT(30),BIT(31)};
ClassImp(ChipTrackingTrgEngine1)
//______________________________________________________________________________
ChipTrackingTrgEngine1::ChipTrackingTrgEngine1() :
TriggerEngine(string("ChipTrackingTrgEngine1"), kChipTrackingTrigger1 ) {
//
// Constructor
//
fMinTrackLength = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fMinTrackLength");
fMaxTrackLength = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fMaxTrackLength");
fMinTriggerTrackLength = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fMinTriggerTrackLength");
fMinTriggerTwoLength = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fMinTriggerTwoLength");
fMaxTwoLength = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fMaxTwoLength");
if ( Conf()->GetStr("ChipTrackingTrgEngine1.fThresholdType")=="absolute" ) fRelativeThreshold = kFALSE;
else if ( Conf()->GetStr("ChipTrackingTrgEngine1.fThresholdType")=="relative" ) fRelativeThreshold = kTRUE;
else Msg(EsafMsg::Panic) << "Wrong cfg value fThresholdType" << MsgDispatch;
fAcceptHole = Conf()->GetBool("ChipTrackingTrgEngine1.fAcceptHole");
fOnlyWithSignal = Conf()->GetBool("ChipTrackingTrgEngine1.fOnlyWithSignal");
fThreshold = (Int_t)Conf()->GetNum("ChipTrackingTrgEngine1.fThreshold");
// panic if we want to simulate >32 trigger patterns
if ( ((fMaxTrackLength-fMinTriggerTrackLength+1)+(fMaxTwoLength-fMinTriggerTwoLength+1)) > 32 )
Msg(EsafMsg::Panic) << "Request for > 32 patterns of trigger." << MsgDispatch;
ConfigFileParser *pConfig = Config::Get()->GetCF("Electronics","MacroCell");
if ( !pConfig->GetBool("MacroCell.fSaveAllChipGtuData") )
Msg(EsafMsg::Warning) << "MacroCell.fSaveAllChipGtuData is not enabled. Some infos are missing" << MsgDispatch;
FillRunPars();
Clear();
}
//______________________________________________________________________________
ChipTrackingTrgEngine1::~ChipTrackingTrgEngine1() {
//
// Destructor
//
Clear();
}
//______________________________________________________________________________
ChipTrackingTrgEngine1::ChipTrackingTrgEngine1( string s, ETriggerTypeIdentifier id ) :
TriggerEngine(s, id) {
//
// Special constructor for child classes
//
}
//______________________________________________________________________________
void ChipTrackingTrgEngine1::Clear() {
//
// Clear data containers
//
// clear track vectors
map<Int_t, vector<ChipTrackSegment> >::iterator it;
for(it=fTrackSegments.begin(); it!=fTrackSegments.end(); it++)
it->second.clear();
// clear map
fTrackSegments.clear();
}
//______________________________________________________________________________
void ChipTrackingTrgEngine1::Simulate( MacroCellData* pData) {
//
// Simulate chip tracking algorythm
//
fPatternCounter = 0;
Int_t nTriggeredTracks = 0;
// get number of channels in one chip
Int_t matsize = pData->Cell()->GetEC(0)->FrontEnd()->Channels();
// build a nearest neighbor matrix and its square
BoolMatrix A( matsize );
A.SetNeighborsMatrix();
BoolMatrix B( matsize );
B = A;
B *= A;
// get map of ChipGtuData elements
map<Int_t, map<Int_t,ChipGtuData*>* > &m1 = pData->fChipData2;
map<Int_t, map<Int_t,ChipGtuData*>* >::const_iterator it1;
// build 3 vectors of state vectors:
// x vector is the list of state vectors v
// y vector is the list of products A*v
// z vector is the list of products A*(A*v)
vector<BoolVector> x;
vector<BoolVector> y;
vector<BoolVector> z;
vector<Int_t> gtu;
// iterate on chips in this macrocell
for(it1=m1.begin(); it1!=m1.end(); it1++) {
Int_t chip_id = it1->first;
map<Int_t,ChipGtuData*> &m2 = *(it1->second);
map<Int_t,ChipGtuData*>::const_iterator it2;
x.clear();
y.clear();
z.clear();
gtu.clear();
// iterate on GTU for the current chip
for(it2=m2.begin(); it2!=m2.end(); it2++) {
gtu.push_back(it2->first);
}
// if we do not have enough activity in this chip, skip to next one
if ( gtu.size() < (UInt_t)fMinTrackLength ) continue;
if (fOnlyWithSignal) {
Bool_t fIsSignal = kFALSE;
for(UInt_t i(0); i<gtu.size(); i++) {
fIsSignal += (Bool_t)(m2[gtu[i]]->GetTotalPureSignal());
}
if (!fIsSignal) continue;
}
Int_t mingtu = gtu[0];
Int_t maxgtu = gtu[0];
for ( UInt_t i(0); i<gtu.size(); i++) {
if ( gtu[i] < mingtu ) mingtu = gtu[i];
if ( gtu[i] > maxgtu ) maxgtu = gtu[i];
}
gtu.clear();
// compute the true threshold to use in trigger algorithm
// if fThresholdType=absolute, is simply the value specified in cfg file
Int_t fTrueThreshold;
if ( fRelativeThreshold ) {
Double_t gtu_length = Config::Get()->GetCF("Electronics","MacroCell")->GetNum("MacroCell.fGtuTimeLength");
// mean number of background hits in a pixel in a gtu
Double_t back = (m2[mingtu]->FrontEnd()->GetNightGlowRate()) * gtu_length;
//Msg(EsafMsg::Info) << "back " << back << " " << m2[mingtu]->FrontEnd()->GetNightGlowRate() << MsgDispatch;
fTrueThreshold = (Int_t)Ceil( back + fThreshold * Sqrt(back) );
} else {
fTrueThreshold = fThreshold;
}
//Msg(EsafMsg::Info) << "Chip " << chip_id << " thresh " << fTrueThreshold << MsgDispatch;
for(Int_t i(0); i<(maxgtu-mingtu+1); i++) {
Int_t curr_gtu = mingtu + i;
gtu.push_back(curr_gtu);
// fill a boolean vector of pixels in this gtu
x.push_back( BoolVector(matsize) );
BoolVector& v1 = x.back();
v1.Zero();
if ( m2.count(curr_gtu) == 1 ) {
ChipGtuData &cgd = *(m2[curr_gtu]);
for(Int_t j=0; j<matsize; j++)
v1(j) = ( cgd.GetCounter(j) >= fTrueThreshold ); //true if j-th channel has at least fTrueThreshold photoelectrons
} else if (m2.count(curr_gtu)>1) {
Msg(EsafMsg::Warning) << "Duplicate ChipGtuData chip " << chip_id << " gtu " << curr_gtu << MsgDispatch;
ChipGtuData &cgd = *(m2[curr_gtu]);
for(Int_t j=0; j<matsize; j++)
v1(j) = ( cgd.GetCounter(j) >= fTrueThreshold ); //true if j-th channel has at least fTrueThreshold photoelectrons
}
y.push_back( BoolVector(matsize) );
BoolVector& v2 = y.back();
v2 = v1;
v2 *= A;
z.push_back( BoolVector(matsize) );
BoolVector& v3 = z.back();
v3 = v1;
v3 *= B;
}
// iterate on x and y to build tracks of length from minimum to maximum
// the algorythm accept at most one hole in the middle of any track
// (if this is required)
Int_t size = gtu.size();
Int_t last = size - fMinTrackLength + 1;
#ifdef DEBUG
Msg(EsafMsg::Debug) << "Number of valid GTU: "<< size << "last= "<<last<< MsgDispatch;
#endif /* DEBUG */
//loop on valid gtus
for(Int_t i=0; i <= last ; i++) {
// loop on valid track length
for(Int_t k=(fMinTrackLength-1); k<=(fMaxTrackLength-1); k++) {
//Bool_t good_track = kFALSE;
Bool_t good_track2 = kFALSE;
Bool_t hole = kFALSE;
BoolVector trackVector(matsize);
trackVector.Zero();
// check if a track of length k can be found
if ( (i+k)>=size ) break;
// loop on hits of a possible track of length k
for(Int_t j=i; j < (i+k); j++) {
// at the beginning there cannot be a hole
if ( j==i ) {
//good_track = (x[j])%(y[j+1]);
trackVector = x[j];
trackVector &= y[j+1];
if ( !trackVector.OR() ) break;
}
// otherwise admit at most one hole
else {
//Bool_t a = x[j]%y[j+1];
if ( x[j]%y[j+1] ) {
//good_track &= a;
trackVector &= y[j+1];
}
else if ( fAcceptHole && !hole && j<(i+k-1)) { //last point of tracks cannot be a hole
//good_track &= x[j]%z[j+2];
trackVector &= z[j+2];
hole = kTRUE;
}
else {
//good_track &= a;
trackVector.Zero();
break;
}
}
}
// good_track: original algorithm (boolean scalar product)
// dood_track2: improved algorithm (bitwise AND + OR of the final vector)
good_track2 = trackVector.OR();
// if we succedeed to build a valid track, store it
if ( good_track2 ) {
fTrackSegments[chip_id].push_back( ChipTrackSegment() );
ChipTrackSegment &seg = fTrackSegments[chip_id].back();
seg.SetGtuStart( gtu[i] );
seg.SetGtuEnd( gtu[i+k] );
seg.SetHasHole( hole );
seg.SetTrackLength(k+1);
seg.SetChipUid( chip_id );
#ifdef DEBUG
Msg(EsafMsg::Debug)<< "Track: len=" << seg.GetTrackLength() << " gtu:" << seg.GetGtuStart()
<< "," << seg.GetGtuEnd() <<" chip_id: "<< chip_id << MsgDispatch;
#endif /* DEBUG */
// if track generate trigger, job is done
if ( seg.GetTrackLength() >= fMinTriggerTrackLength ) {
seg.SetTriggered( kTRUE );
//SetTrigger( true );
SetGtuTrigger( seg.GetGtuStart() );
// check wich pattern of trigger is valid
fPatternCounter = 0; // in trigger word first trigger in a microcell only
for ( Int_t l=fMinTriggerTrackLength; l<=fMaxTrackLength; l++) {
if (seg.GetTrackLength()==l) SetTrigger( kTrigPattern[fPatternCounter] );
fPatternCounter++;
}
nTriggeredTracks++;
}
}
}
}
}
// here try to attach tracks of different length
// in this version of the algorythm, we trigger if:
// a track of length fMinTriggerTrackLnegth is found
// 2 tracks whose total sum is at least fMinTriggerTrackLnegth and the
// end gtu of the first is connected to the first gtu of the second track
// the second condition searches for tracks that are split in two ECs
for(it1=m1.begin(); it1!=m1.end(); it1++) {
Int_t chip_id = it1->first;
if ( !(fTrackSegments.count(chip_id))) continue;
map<Int_t,ChipGtuData*> &m2 = *(it1->second);
map<Int_t,ChipGtuData*>::const_iterator it2;
it2 = m2.begin();
Int_t uidSt = (*(it2->second)).FrontEnd()->UniqueChanId(0);
Int_t uidEnd = (*(it2->second)).FrontEnd()->UniqueChanId(matsize-1);
pair<Int_t,Int_t> rcSt = pData->Cell()->GetRowColUniqueId(uidSt);
pair<Int_t,Int_t> rcEnd = pData->Cell()->GetRowColUniqueId(uidEnd);
Int_t rmin(0),rmax(0),cmin(0),cmax(0);
if (rcSt.first > rcEnd.first ) {
rmax = rcSt.first;
rmin = rcEnd.first;
} else {
rmax = rcEnd.first;
rmin = rcSt.first;
}
if (rcSt.second > rcEnd.second ) {
cmax = rcSt.second;
cmin = rcEnd.second;
} else {
cmax = rcEnd.second;
cmin = rcSt.second;
}
// search the neighbors chips with at least a valid track segment
vector<Int_t> fChipNeighbors;
for(Int_t kk=0; kk<4; kk++) {
Int_t row0(0), col0(0), sr(0), sc(0);
switch(kk) {
case 0: {row0 = rmin; col0=cmin; sr=-1; sc=-1;}
case 1: {row0 = rmin; col0=cmax; sr=-1; sc=1;}
case 2: {row0 = rmax; col0=cmin; sr=1; sc=-1;}
case 3: {row0 = rmax; col0=cmax; sr=1; sc=1;}
}
for( Int_t r(0); r<2; r++ ) {
for( Int_t c(0); c<2; c++ ) {
if ( r == 0 && c==0 ) continue;
Int_t row = row0 + sr*r;
Int_t col = col0 + sc*c;
Int_t id = 0;
Photomultiplier* pmt = GetEusoDetector()->GetEusoElectronics()->PmtId(pData->Cell()->GetUniqueIdRowCol(row,col));
if (pmt) id = pmt->FrontEnd()->Id(); else continue;
if (!id || id == chip_id || fTrackSegments.count(id) == 0) continue;
Bool_t flag = kTRUE;
for ( Int_t i(0); i<(Int_t)fChipNeighbors.size(); i++ ) {
if ( fChipNeighbors[i] == id ) flag=kFALSE;
}
if (flag) fChipNeighbors.push_back(id);
}
}
}
// merge tracks segments between nearby chips
for(Int_t i(0); i<(Int_t)fTrackSegments[chip_id].size(); i++) {
ChipTrackSegment &cseg = fTrackSegments[chip_id][i];
//if ( cseg.GetTriggered()==kTRUE ) continue;
if ( cseg.GetTrackLength() > fMinTriggerTrackLength ) continue; // only for my trigger studies (RP)
for( Int_t k(0); k<(Int_t)fChipNeighbors.size(); k++ ){
Int_t id = fChipNeighbors[k];
for( Int_t j(0); j<(Int_t)fTrackSegments[id].size(); j++ ) {
ChipTrackSegment seg = fTrackSegments[id][j];
Int_t totlen = seg.GetTrackLength() + cseg.GetTrackLength();
if ( totlen < fMinTriggerTwoLength || totlen > fMaxTwoLength ) continue;
if ( TMath::Abs(cseg.GetGtuEnd()-seg.GetGtuEnd()) < (seg.GetTrackLength()+1) ) {
//cseg.SetTriggered(kTRUE);
//SetTrigger( true );
fPatternCounter = fMaxTrackLength - fMinTriggerTrackLength + 1; // skip patterns for a microcell only
for(Int_t l=fMinTriggerTwoLength; l<fMaxTwoLength; l++) {
if ( totlen == l ) SetTrigger( kTrigPattern[fPatternCounter] );
fPatternCounter++;
}
SetGtuTrigger( (seg.GetGtuStart() > cseg.GetGtuStart()) ? cseg.GetGtuStart() : seg.GetGtuStart() );
if ( cseg.GetTriggered() == kFALSE ) {
cseg.SetTriggered(kTRUE);
nTriggeredTracks++;
fTrackSegments[chip_id][i] = cseg; }
if ( seg.GetTriggered()==kFALSE ) {
seg.SetTriggered(kTRUE);
nTriggeredTracks++;
fTrackSegments[id][j] = seg;
}
//break;
}
}
//if ( cseg.GetTriggered()==kTRUE ) break;
}
}
fChipNeighbors.clear();
}
//if ( nTriggeredTracks )
//Msg(EsafMsg::Info) << nTriggeredTracks << " tracks have triggered in cell " << pData->Cell()->Id() << MsgDispatch;
//FillEEvent(pData->Cell()->Id());
//Dump(pData->Cell()->Id());
Clear();
}
//_____________________________________________________________________________________
void ChipTrackingTrgEngine1::FillEEvent( Int_t cell_id ) {
//
// Fill root event
//
map<Int_t, vector<ChipTrackSegment> >::iterator it;
for( it = fTrackSegments.begin(); it != fTrackSegments.end(); it++ ) {
Int_t chip_id = (*it).first;
for( UInt_t i(0); i<((*it).second).size(); i++ ) {
if ( EEvent::GetCurrent() ) {
EChipTrackTriggerDataAdder a( &(((*it).second)[i]), cell_id, chip_id );
EEvent::GetCurrent()->Fill(a);
}
}
}
}
//_____________________________________________________________________________________
void ChipTrackingTrgEngine1::FillRunPars() {
//
// Fill run parameters
//
if (!(EEvent::GetCurrent() && EEvent::GetCurrent()->GetRunPars()) ) return;
ERunParameters *runpars = EEvent::GetCurrent()->GetRunPars();
if (runpars) {
EChipTriggParsFiller ctpfiller;
ctpfiller.SetName(GetName().c_str());
ctpfiller.SetId((ETriggerTypeIdentifier)GetTriggerId());
ctpfiller.SetRelativeThreshold(fRelativeThreshold);
ctpfiller.SetThreshold(fThreshold);
ctpfiller.SetMinTrackLength(fMinTrackLength);
ctpfiller.SetMaxTrackLength(fMaxTrackLength);
ctpfiller.SetMinTriggerTrackLength(fMinTriggerTrackLength);
ctpfiller.SetMinTriggerTwoLength(fMinTriggerTwoLength);
ctpfiller.SetMaxTwoLength(fMaxTwoLength);
ctpfiller.SetAcceptHole(fAcceptHole);
ctpfiller.SetOnlyWithSignal(fOnlyWithSignal);
runpars->Fill(ctpfiller);
}
}
//_____________________________________________________________________________________
void ChipTrackingTrgEngine1::Dump( Int_t cell_id ) {
//
// Dump the number of tracks (triggered and total) for each chip
//
map<Int_t, vector<ChipTrackSegment> >::iterator it;
Msg(EsafMsg::Info) << "***DUMPING CELL " << cell_id << MsgDispatch;
if (fTrackSegments.size()==0) return;
for( it = fTrackSegments.begin(); it != fTrackSegments.end(); it++ ) {
Int_t chip_id = (*it).first;
Int_t tottracks = ((*it).second).size();
Int_t trgtracks(0);
Int_t gtu_min, gtu_max;
gtu_min = ((*it).second)[0].GetGtuStart();
gtu_max = ((*it).second)[0].GetGtuEnd();
for( UInt_t i(0); i<((*it).second).size(); i++ ) {
if ( ((*it).second)[i].GetTriggered() ) trgtracks++;
if ( ((*it).second)[i].GetGtuStart() < gtu_min ) gtu_min = ((*it).second)[i].GetGtuStart();
if ( ((*it).second)[i].GetGtuEnd() > gtu_max ) gtu_max = ((*it).second)[i].GetGtuEnd();
}
Msg(EsafMsg::Info) << "Chip " << chip_id << " Total = " << tottracks << " Triggered = " << trgtracks << MsgDispatch;
Msg(EsafMsg::Info) << "Gtu min = " << gtu_min << " Gtu max = " << gtu_max << MsgDispatch;
}
}