Commit 28ee968a authored by dino's avatar dino
Browse files

Keyword NumGeDets in TrackingFilter to override the default number of...

Keyword NumGeDets in TrackingFilter to override the default number of detectors NUMGEDET defined in commonDef.h

git-svn-id: svn://gal-serv.lnl.infn.it/agata/trunk/narval_emulator@1131 170316e4-aea8-4b27-aad4-0380ec0519c9
parent b2eda57b
......@@ -656,7 +656,7 @@ void EventBuilder::GetParameters(UInt_t *error_code)
forceTailSlash(fOdirPrefix);
KEND;
KPAR("BuilderType", "%s [%s]", "builder type(TimeStamp | EventNumber) window [lower]");
KPAR("BuilderType", "%s %d [%d]", "builder type(TimeStamp | EventNumber) coincidence_window start_from(default = 0)");
stringSplit(data, keyw, value);
if( stringEq(keyw, "TimeStamp") ) {
useEventNumber = false;
......
......@@ -405,7 +405,7 @@ bool adetParams::ReadXtalkCoeffs(std::string cname, bool verbose)
}
cout << "Found all " << nSG*nSG << " expected values" << endl;
doecalF2 = true;
doEcalXT = true;
return true;
}
......
......@@ -41,7 +41,7 @@ static const int NSLI = 6; // number of slices (faster index)
static const int NSEC = 6; // number of sectors (slower index)
public:
adetParams() : vtrig(-1000.f), itrig(0), nCC(NCC), nSG(NSG),
doecalF2(false), bSegEmin(false), vSegEmin(0),
doEcalXT(false), bSegEmin(false), vSegEmin(0),
tntFactor(2097152.f) { // 2^21
for(int nn = 0; nn < nCC; nn++)
pCC[nn].myID = nn;
......@@ -64,10 +64,11 @@ public:
float xTalkProp[NSG][NSG]; // proportional cross-talk
float xTalkDiff[NSG][NSG]; // differential cross-talk
float xTalkAver[NSG]; // average proportional xTalk
bool doecalF2; // flag to enable this correction
bool doEcalXT; // flag to enable this correction
bool bSegEmin;
float vSegEmin;
float tntFactor; // usually 2^21 = 2097152
void SetTriggerLevel(float value) {vtrig = value;}
......
......@@ -162,7 +162,7 @@ const int defTriggerSample = 10; // 60-10=50 samples passed to the PSA
//////////////////////// Global Filters /////////////////////////////
/////////////////////////////////////////////////////////////////////
#define NUMGEDETS 20 // max acceptable value of crystal_id is NUMGEDETS-1 (used to limit the size of internal structures and spectra)
#define NUMGEDETS 20 // default number of crystals, used to define size of structures in the global level filters
//////////////////////////
////// GlobalFilter //////
......
......@@ -2604,12 +2604,12 @@ Int_t PSAFilterGridSearch::AlgoSpecificPostProcess(int slot)
int bestDT = int(pSlot->shiftT0*SAMP_STEP + specLenT/2); // the correction, placed in the middle of the spectrum
// CoreT[0] is the original position of the traces, in samples
PsaSpecTZero->incr(36+0, bestDT);
PsaSpecTZero->incr(36+0, pSlot->tCycles); // number of search-tShift cycles (at the beginning of the spectrum)
PsaSpecTZero->incr(36+0, pSlot->tCycles); // number of search-tShift cycles (at the beginning of the spectrum)
PsaSpecTZero->incr(36+1, bestT0);
PsaSpecTZero->incr(36+2, bestT1);
PsaSpecTZero->incr(36+3, bestT2);
if(numNetSegs == 1)
PsaSpecTZero->incr(pOut->netChargeSeg, bestDT);
PsaSpecTZero->incr(pSlot->PsaOut[0].netChargeSeg, bestDT);
}
if(PsaSpecStat && numFit_T0>1)
PsaSpecStat->Incr(0, numNetSegs); // distribution of number of hit segments
......
......@@ -313,7 +313,7 @@ void PostPSAFilter::GetParameters(UInt_t *error_code)
fTrappingFile = data;
KEND;
KPAR("SmearPos", "%f", "smearing parameter for matrix of xyz hits");
KPAR("SmearPos", "%f", "xyz uniform smearing of hits (usually the size of the PSA fine-grid");
ok = 1 == sscanf(data.c_str(), "%f", &fSmearPosHit);
KEND;
......
......@@ -788,7 +788,7 @@ bool PreprocessingFilterPSA::ProcessEvent(int &segMult, float &eSumSG1, float
//// cross-talk correction for the segment energies
// (cross talk of traces is applied to the signal basis in PSAFilterGridSearch)
if(CC.doecalF2) {
if(CC.doEcalXT) {
segMult =
xTalkCorrAll(SegE, sumSegs); // the number of hit segments can change due to the increase of energy produced by the xTalk correction
//xTalkCorrNet(SegE, sumSegs); // BEFORE USING THIS VERSION, CHECK IF IT IS UP-TO-DATE
......
......@@ -130,12 +130,16 @@ TrackingFilter::TrackingFilter() :
fWriteTracked = false;
fNumGeDets = NUMGEDETS;
cryst = new crystdata[fNumGeDets];
fDetEnerMin = 20.f; // should be much wider
fDetEnerMax = 20000.f;
fSumEnerMin = 20.f; // should be much wider
fSumEnerMax = NUMGEDETS*fDetEnerMax;
fSumEnerMax = fNumGeDets*fDetEnerMax;
fDetFoldMin = 0; // i.e. keep events with no gamma detectors
fDetFoldMax = NUMGEDETS;
fDetFoldMax = fNumGeDets;
fEnergyGain = 1.f;
......@@ -186,8 +190,9 @@ TrackingFilter::TrackingFilter() :
TrackMatr_SGSG = NULL;
TrackMatr_Fold = NULL;
for(int ii = 0; ii < 45; ii++)
SpecMap[ii] = ii;
for(int ii = 0; ii < 180; ii++)
SpecMap180[ii] = ii;
#endif //TRF_MULTIHIST
RecoilVc = 0;
......@@ -995,6 +1000,12 @@ void TrackingFilter::GetDataPSA()
if(crystal_status)
return;
// should not happen;
if(crystal_id >= 180) {
crystal_status = 1;
return;
}
ULong64_t timestamp = ((AgataKey *)fFramePSA->GetFrame()->GetKey())->GetTimeStamp();
UInt_t evnumber = ((AgataKey *)fFramePSA->GetFrame()->GetKey())->GetEventNumber();
......@@ -1011,6 +1022,9 @@ void TrackingFilter::GetDataPSA()
if(CoreE[0] < fDetEnerMin || CoreE[0] > fDetEnerMax)
return;
if(number_of_crystals >= fNumGeDets)
return; // no space in cryst
fCoreSumEnergy += CoreE[0];
// loop over number of interaction points in fFramePSA
......@@ -1037,7 +1051,7 @@ void TrackingFilter::GetDataPSA()
CC3D_BefAftRT ->Fill( pt->X, pt->Y, pt->Z);
#endif
Transform(pt, crystal_id);
TransformToGlobalFrame(pt, crystal_id);
#ifdef TRF_FromGRU
CC3D_BefAftRT ->Fill(pt->X, pt->Y, pt->Z);
......@@ -1397,12 +1411,17 @@ void TrackingFilter::GetParameters(UInt_t *error_code)
ok = false;
KEND;
KPAR("NumGeDets", "%d", "Overrides the default number of Ge Crystals [20]");
ok = 1 == sscanf(data.c_str(), "%d", &fNumGeDets);
fNumGeDets = max(1, fNumGeDets);
KEND;
KPAR("SpecMap", "%d %d", "ID SubSpec mapping of crystalID to spectrum SubSpec (given multiple times");
#ifdef TRF_MULTIHIST
int sm1=-1, sm2=-1;
ok = 2 == sscanf(data.c_str(), "%d %d", &sm1, &sm2);
if(sm1 >=0 && sm1 <45 && sm2 >=0 && sm2 < NUMGEDETS)
SpecMap[sm1] = sm2;
if(sm1 >=0 && sm1 <180)
SpecMap180[sm1] = sm2;
else
ok = false;
#else
......@@ -1554,6 +1573,24 @@ void TrackingFilter::InitGenStructures()
{
cout << "\nTrackingFilter::InitGenStructures()" << endl;
if(fNumGeDets != NUMGEDETS) {
if(cryst) {
delete [] cryst;
cryst = new crystdata[fNumGeDets];
}
fSumEnerMax = fNumGeDets*fDetEnerMax;
fDetFoldMax = min(fNumGeDets,fDetFoldMax);
}
#ifdef TRF_MULTIHIST
for(int ii = 0; ii < 180; ii++) {
if(SpecMap180[ii]>=fNumGeDets) {
SpecMap180[ii] = fNumGeDets-1; // all out of range go here ??
}
}
#endif
rotTr = (rotTr3D *)calloc( 180, sizeof(rotTr3D));
pEXYZ = (exyzHit *)calloc(intmax, sizeof(exyzHit));
pGams = (trGamma *)calloc(intmax, sizeof(trGamma));
......@@ -1597,22 +1634,22 @@ void TrackingFilter::InitGenStructures()
#ifdef TRF_MULTIHIST
if(fUseMultiHist) {
TrackSpec_Fold = new MultiHist<unsigned int>(4, NUMGEDETS);
TrackSpec_Fold = new MultiHist<unsigned int>(4, fNumGeDets);
TrackSpec_Fold->setFileName(fOdirPrefix+"Track?Fold.spec");
TrackSpec_Fold->setComment("multiplicity of crystals");
hGroup.add(TrackSpec_Fold);
TrackMatr_Fold = new MultiHist<unsigned int>(NUMGEDETS, NUMGEDETS);
TrackMatr_Fold = new MultiHist<unsigned int>(fNumGeDets, fNumGeDets);
TrackMatr_Fold->setFileName(fOdirPrefix+"Track?Fold.matr");
TrackMatr_Fold->setComment("number of gammas vs number of crystals");
hGroup.add(TrackMatr_Fold);
TrackSpec_EC = new MultiHist<unsigned int>(2, NUMGEDETS, specLenE);
TrackSpec_EC = new MultiHist<unsigned int>(2, fNumGeDets, specLenE);
TrackSpec_EC->setFileName(fOdirPrefix+"Track?EC.spec");
TrackSpec_EC->setComment("core and sumSegs energy sorted by detector");
hGroup.add(TrackSpec_EC);
TrackSpec_ES = new MultiHist<unsigned int>(1, NUMGEDETS, specLenE);
TrackSpec_ES = new MultiHist<unsigned int>(1, fNumGeDets, specLenE);
TrackSpec_ES->setFileName(fOdirPrefix+"Track?ES.spec");
TrackSpec_ES->setComment("sum-energy spectra (gain 1 keV/chan)");
hGroup.add(TrackSpec_ES);
......@@ -1623,22 +1660,22 @@ void TrackingFilter::InitGenStructures()
hGroup.add(TrackSpec_EE);
// no meaning if the integer part of the trigger-point has been moved to timestamp
TrackSpec_TC = new MultiHist<unsigned int>(NUMGEDETS, specLenT1);
TrackSpec_TC = new MultiHist<unsigned int>(fNumGeDets, specLenT1);
TrackSpec_TC->setFileName(fOdirPrefix+"Track?TC.spec");
TrackSpec_TC->setComment("CFD trigger relative to beginning of trace (ns)");
hGroup.add(TrackSpec_TC);
//TrackMatr_ETC = new MultiHist<unsigned short>(1000, NUMGEDETS*1024);
//TrackMatr_ETC = new MultiHist<unsigned short>(1000, fNumGeDets*1024);
//TrackMatr_ETC->setFileName(fOdirPrefix+"Track?ETC.matr");
//TrackMatr_ETC->setComment("gamma-Trigger");
//hGroup.add(TrackMatr_ETC);
TrackSpec_TS = new MultiHist<unsigned int>(NUMGEDETS, NUMGEDETS, specLenT1);
TrackSpec_TS = new MultiHist<unsigned int>(fNumGeDets, fNumGeDets, specLenT1);
TrackSpec_TS->setFileName(fOdirPrefix+"Track?TS.spec");
TrackSpec_TS->setComment("T_gamma-gamma from timestamp (samples)");
hGroup.add(TrackSpec_TS);
TrackSpec_TT = new MultiHist<unsigned int>(NUMGEDETS, NUMGEDETS, specLenT1);
TrackSpec_TT = new MultiHist<unsigned int>(fNumGeDets, fNumGeDets, specLenT1);
TrackSpec_TT->setFileName(fOdirPrefix+"Track?TT.spec");
TrackSpec_TT->setComment("T_gamma-gamma from CFD (ns)");
hGroup.add(TrackSpec_TT);
......@@ -1648,7 +1685,7 @@ void TrackingFilter::InitGenStructures()
TrackSpec_TG->setComment("Tgamma and Tgg of tracked gammas");
hGroup.add(TrackSpec_TG);
//TrackMatr_IEET = new MultiHist<unsigned char>(NUMGEDETS, 100, 100, specLenT2);
//TrackMatr_IEET = new MultiHist<unsigned char>(fNumGeDets, 100, 100, specLenT2);
//TrackMatr_IEET->setFileName(fOdirPrefix+"Track?IEET.matr");
//TrackMatr_IEET->setComment("det_index-energy-energy-time");
//hGroup.add(TrackMatr_IEET);
......@@ -1667,18 +1704,18 @@ void TrackingFilter::InitGenStructures()
hGroup.add(TrackMatr_gg2);
}
TrackMatr_SGSG = new MultiHist<unsigned short>(36*NUMGEDETS, 36*NUMGEDETS);
TrackMatr_SGSG = new MultiHist<unsigned short>(36*fNumGeDets, 36*fNumGeDets);
TrackMatr_SGSG->setFileName(fOdirPrefix+"Track?SGSG.matr");
TrackMatr_SGSG->setComment("segment-segment coincidences for all pairs");
hGroup.add(TrackMatr_SGSG);
if(fAncillary) {
TrackSpec_EA = new MultiHist<unsigned int>(2, NUMGEDETS, specLenE); // 0 alone, 1 Ancillary
TrackSpec_EA = new MultiHist<unsigned int>(2, fNumGeDets, specLenE); // 0 alone, 1 Ancillary
TrackSpec_EA->setFileName(fOdirPrefix+"Track?EA.spec");
TrackSpec_EA->setComment("core energy sorted by detector, in coincidence with ancillary");
hGroup.add(TrackSpec_EA);
TrackSpec_TA = new MultiHist<unsigned int>(2, NUMGEDETS, specLenT1);
TrackSpec_TA = new MultiHist<unsigned int>(2, fNumGeDets, specLenT1);
TrackSpec_TA->setFileName(fOdirPrefix+"Track?TA.spec");
TrackSpec_TA->setComment("Difference of tstamps between crystals and ancillary");
hGroup.add(TrackSpec_TA);
......@@ -1688,12 +1725,12 @@ void TrackingFilter::InitGenStructures()
TrackSpec_Anc->setComment("Projection of all Ancillary Parameters");
hGroup.add(TrackSpec_Anc);
//TrackSpec_TOF = new MultiHist<unsigned int>(NUMGEDETS, 32, specLenT1);
//TrackSpec_TOF = new MultiHist<unsigned int>(fNumGeDets, 32, specLenT1);
//TrackSpec_TOF->setFileName(fOdirPrefix+"Track?TOF.spec");
//TrackSpec_TOF->setComment("Timestamp-corrected time spectra crystals-ancillary");
//hGroup.add(TrackSpec_TOF);
//TrackMatr_TOF = new MultiHist<unsigned short>(NUMGEDETS, 32, 400, 200);
//TrackMatr_TOF = new MultiHist<unsigned short>(fNumGeDets, 32, 400, 200);
//TrackMatr_TOF->setFileName(fOdirPrefix+"Track?TOF.matr");
//TrackMatr_TOF->setComment("Timestamp-corrected Egamma-TOF");
//hGroup.add(TrackMatr_TOF);
......@@ -1803,7 +1840,7 @@ int TrackingFilter::GetRotoTranslations(FILE *fp)
// return false;
// }
//
// for(int ii = 0; ii < NUMGEDETS; ii++) {
// for(int ii = 0; ii < fNumGeDets; ii++) {
// fRecalCoeffs[ii][0] = 0;
// fRecalCoeffs[ii][1] = 1.f;
// }
......@@ -1817,7 +1854,7 @@ int TrackingFilter::GetRotoTranslations(FILE *fp)
// if(fscanf(fp_conf, "%d %f %f", &nd, &xx, &yy) != 3)
// break;
//
// if(nd < 0 || nd >= NUMGEDETS)
// if(nd < 0 || nd >= fNumGeDets)
// return false;
//
// fRecalCoeffs[[nd]0] = xx; // offset
......@@ -1995,6 +2032,8 @@ int TrackingFilter::PostProcessEvent()
trGamma *pg = pGams;
// loop on tracked gammas in the event
float eSumTracked = 0;
float eSumTrackedA = 0;
for(Int_t i = 0; i < number_of_gammas; i++, pg++ ) {
#ifdef TRF_WRITE_TRACKED
if(fWriteTracked) {
......@@ -2007,9 +2046,11 @@ int TrackingFilter::PostProcessEvent()
if(TrackSpec_EE) {
TrackSpec_EE->Incr(0, 2, int(fEnergyGain*pg->E)); // 0-2 Tracked
TrackSpec_EE->Incr(0, 2 + pg->type, int(fEnergyGain*pg->E)); // 0-3 Photo=1 0-4 Compt=2 0-5 pair=3
eSumTracked += pg->E;
if(hasRecoil) {
TrackSpec_EE->Incr(1, 2, int(fEnergyGain*pg->E)); // 1-2 Tracked in coinc with ancillary
TrackSpec_EE->Incr(1, 2 + pg->type, int(fEnergyGain*pg->E)); // 0-3 Photo=1 0-4 Compt=2 0-5 pair=3
eSumTrackedA += pg->E;
}
}
#endif //TRF_MULTIHIST
......@@ -2038,6 +2079,11 @@ int TrackingFilter::PostProcessEvent()
} // loop over the gammas
// sum spectra
TrackSpec_EE->Incr (0, 7, int(fEnergyGain*eSumTracked)); // 0-7 sum energy
if(hasRecoil)
TrackSpec_EE->Incr(1, 7, int(fEnergyGain*eSumTrackedA)); // 1-7 sum energy
#ifdef TRF_MULTIHIST
// time of the individual gammas and time between gammas using only info from the output
if(TrackSpec_TG) {
......@@ -2082,9 +2128,7 @@ bool TrackingFilter::AnalysisOfCrystals()
float esumCC = 0;
for(int nn1 = 0; nn1 < number_of_crystals; nn1++) {
crystdata *pCr1 = cryst + nn1;
int id1 = SM(pCr1->id);
if(id1 < 0 || id1 >= NUMGEDETS) // this should not happen and should be handled as an error ??
continue;
int id1 = SM180(pCr1->id); // off range should return fNumgeDets-1
esumCC += pCr1->CCE0;
if(TrackSpec_EC) {
TrackSpec_EC->Incr(0, id1, int(fEnergyGain*pCr1->CCE0)); // energy of the core sorted by detector
......@@ -2110,9 +2154,7 @@ bool TrackingFilter::AnalysisOfCrystals()
if(nn2 == nn1)
continue;
crystdata *pCr2 = cryst + nn2;
int id2 = SM(pCr2->id);
if(id2 < 0 || id2 >= NUMGEDETS) // this should not happen and should be handled as an error ??
continue;
int id2 = SM180(pCr2->id); // off range should return fNumgeDets-1
dts = (int)(pCr2->tstamp - pCr1->tstamp);
dtf = pCr2->CCT0 - pCr1->CCT0;
dtt = tScale*(dts + dtf) ;
......@@ -2135,7 +2177,7 @@ bool TrackingFilter::AnalysisOfCrystals()
}
}
}
if(TrackSpec_ES && number_of_crystals > 0 && number_of_crystals < NUMGEDETS) {
if(TrackSpec_ES && number_of_crystals > 0 && number_of_crystals < fNumGeDets) {
TrackSpec_ES->Incr(0, 0, int(esumCC)); // total sum energy
TrackSpec_ES->Incr(0, number_of_crystals, int(esumCC)); // sum energy [fold]
}
......@@ -2161,8 +2203,7 @@ bool TrackingFilter::AnalysisGeAncillary()
for(int nn1 = 0; nn1 < number_of_crystals; nn1++) {
crystdata *pCr1 = cryst + nn1;
int ii1 = SM(pCr1->id);
if(ii1 < 0 || ii1 >= NUMGEDETS) continue; // this should not happen
int ii1 = SM180(pCr1->id); // off range should return fNumgeDets-1
int eG = int(fEnergyGain*pCr1->CCE0);
int eS = int(fEnergyGain*pCr1->Esum);
int dT = specLenT1/2; // timestamp difference calculated (and used) with zero at specLenT1/2
......
......@@ -143,10 +143,12 @@ protected:
int number_of_tracked_gammas;
int number_of_crystals;
crystdata cryst[NUMGEDETS];
crystdata *cryst;
Float_t fCoreSumEnergy;
Int_t fNumGeDets; // number of Ge crystals (defaults to NUMGEDETS)
UInt_t *trigType; // statistics of fired triggers
UInt_t trigNums; // number of triggers defined
UInt_t trigUnknown; // triggers fired and not known
......@@ -186,10 +188,10 @@ protected:
MultiHist<unsigned short> *TrackMatr_SGSG; // SG-SG counts for all pairs
MultiHist<unsigned int > *TrackMatr_Fold; // fold distributions
int SpecMap[45]; // SHOULD BE A PARAMETER
int SM(int ii) {return SpecMap[ii];}
int SpecMap180[180]; // in principle there will not be higher ID's than this but it should be done better
int SM180(int ii) {return SpecMap180[ii];}
#else
int SM(int ii) {return ii;}
int SM180(int ii) {return ii;}
#endif //TRF_MULTIHIST
cycleServer cServer; // manager of cyclic operations
......@@ -300,7 +302,7 @@ protected:
void InitGenStructures();
int GetRotoTranslations(FILE *fp_conf);
void Transform(exyzHit *pt, int nb_det) { // transformation to global frame
void TransformToGlobalFrame(exyzHit *pt, int nb_det) {
float Xr = pt->X;
float Yr = pt->Y;
float Zr = pt->Z;
......
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