Commit 2ff1dea7 authored by dino's avatar dino
Browse files

Added the possibility to write-out the PSA hits in the crystal frame to study...

Added the possibility to write-out the PSA hits in the crystal frame to study the problem of neutron damage using external programs. Consequently the number of spectra/matrices has been reduced.
Added a rescaling of the core energy int GetDataPSA() with values from a file (keyword RescaleHits in getParams()). In this case the energy of the segments is rescaled to the new energy of the core.

git-svn-id: svn://gal-serv.lnl.infn.it/agata/trunk/narval_emulator@900 170316e4-aea8-4b27-aad4-0380ec0519c9
parent a1848222
......@@ -296,9 +296,15 @@ Int_t PSAFilter::SetOutput(int slot)
if(pD->numHits < 0)
return 0;
for(UInt_t nh = 0; nh < pD->numHits; nh++) {
PsaOut_t * pOut = pD->PsaOut;
for(UInt_t nh = 0; nh < pD->numHits; nh++, pOut++) {
ADF::PSAHit *hit = (PSAHit*)fFramePSA->Data()->NewHit();
hit->SetHit(pD->theHits[nh]);
hit->Reset();
hit->SetE(pOut->enerSG);
hit->SetXYZ(pOut->fx, pOut->fy, pOut->fz);
hit->SetID(pOut->netChargeSeg);
hit->SetT(pOut->dTns);
hit->SetDT(pOut->chiSq); // chi2 temporarily stored in the error DT
}
((AgataKey *)fFramePSA->GetKey())->SetEventNumber(pD->evnumber);
......
......@@ -60,16 +60,23 @@ struct CSlot
};
#endif // TCOUNT > 1
struct nDamage_t
// the result of PSA for each hit segment
struct PsaOut_t
{
nDamage_t() : enerCC(0), enerSG(0), fx(0), fy(0), fz(0), fNetSeg(0),
ePath(0), hPath(0), eS_CC(0), hS_CC(0), eS_SG(0), hS_SG(0) {;}
float enerCC, enerSG;
float fx, fy, fz;
float fNetSeg;
float ePath, hPath;
float eS_CC, hS_CC;
float eS_SG, hS_SG;
PsaOut_t() : enerSG(0), fx(0), fy(0), fz(0),
dTns(0), chiSq(0), netChargeSeg(0), bestPoint(0) {;}
void Reset() {
enerSG = 0;
fx = fy = fz = 0;
dTns = chiSq = 0;
netChargeSeg = bestPoint = 0;
}
float enerSG; // segment energy
float fx, fy, fz; // hit position (mm)
float dTns; // time shift (presently unused)
float chiSq; // chi square of the best point
int netChargeSeg; // the index of the segment of this hit
int bestPoint; // used in PostProcess() to extract information stored in the signal basis
};
// Data interface to PSAGridSearch.
......@@ -91,9 +98,8 @@ struct PsaData
Float_t CoreT[ADF::CrystalInterface::kNbCores];
// output data (with this algo there cannot be more hits than segments)
UInt_t numHits;
ADF::PSAHit theHits[ADF::CrystalInterface::kNbSegments];
nDamage_t nDamage[ADF::CrystalInterface::kNbSegments];
UInt_t numHits;
PsaOut_t PsaOut[ADF::CrystalInterface::kNbSegments];
// This is the the resulting "fitted" trace, used to save Psa__xxxx__Traces.samp and to fit T0
// To enable using short or float without including GridSearchParams.h, keep both pointer types
......@@ -112,8 +118,7 @@ struct PsaData
memset(CoreT, 0, sizeof(CoreT));
numHits = 0;
for(UInt_t nh = 0; nh < ADF::CrystalInterface::kNbSegments; nh++)
theHits[nh].Reset();
numHits = 0;
PsaOut[nh].Reset();
sfAmplitude = srAmplitude = NULL;
ffAmplitude = frAmplitude = NULL;
retValue = 0;
......
......@@ -35,9 +35,7 @@ PSAFilterGridSearch::PSAFilterGridSearch()
#ifdef LOCALSPECTRA
specEner = NULL; specTzero = NULL; specSigma = NULL;
matrXYZR = NULL; matrSeg = NULL;
matrRE = NULL; matrZE = NULL; matrZRE = NULL; matrRZE = NULL;
matrDZE = NULL; matrZREcEs = NULL; matrRZEcEs = NULL;
matrXYZR = NULL; matrSeg = NULL; matrRZE = NULL;
#endif //LOCALSPECTRA
bDoSpec = false;
bDoMats = false;
......@@ -166,7 +164,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
}
#ifdef LOCALSPECTRA
specEner = new nDhist<unsigned int>(2, 40, specLenE);
specEner = new nDhist<unsigned int>(40, specLenE);
specEner->setFileName(fOdirPrefix+"Psa?Ener.spec");
specEner->setComment("Energy spectrum of segments and core");
hGroup.add(specEner);
......@@ -181,7 +179,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
//specSigma->setComment ("chisquare");
//hGroup.add(specSigma);
matrXYZR = new nDhist<unsigned short>(2,10, matLen, matLen);
matrXYZR = new nDhist<unsigned short>(20, matLen, matLen);
matrXYZR->setFileName(fOdirPrefix+"Psa?XYZR.matr");
matrXYZR->setComment("XY, XZ, YZ, rZ hits and XY hits as a function of Z");
hGroup.add(matrXYZR);
......@@ -191,43 +189,13 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
matrSeg->setComment("hits in the segments");
hGroup.add(matrSeg);
//matrRE = new nDhist<unsigned short>(3, 50, 2000);
//matrRE->setFileName(fOdirPrefix+"Psa?RE.matr");
//matrRE->setComment("core and segment energy as a function of R {E = (E_keV-1000)*5}");
//hGroup.add(matrRE);
//matrZE = new nDhist<unsigned short>(3, 50, 2000);
//matrZE->setFileName(fOdirPrefix+"Psa?ZE.matr");
//matrZE->setComment("core and segment energy as a function of Z {E = (E_keV-1000)*5}");
//hGroup.add(matrZE);
matrZRE = new nDhist<unsigned short>(3, 50, 50, 2000);
matrZRE->setFileName(fOdirPrefix+"Psa?ZRE.matr");
matrZRE->setComment("core and segment energy as a function of ZR {E = (E_keV-1000)*5}");
hGroup.add(matrZRE);
matrRZE = new nDhist<unsigned short>(3, 50, 50, 2000);
matrRZE = new nDhist<unsigned short>(2, 50, 50, 2000);
matrRZE->setFileName(fOdirPrefix+"Psa?RZE.matr");
matrRZE->setComment("core and segment energy as a function of RZ {E = (E_keV-1000)*5}");
hGroup.add(matrRZE);
matrDZE = new nDhist<unsigned short>(4, 50, 50, 2000);
matrDZE->setFileName(fOdirPrefix+"Psa?DZE.matr");
matrDZE->setComment("core and segment energy as a function of distance to electrodes {E = (E_keV-1000)*5}");
hGroup.add(matrDZE);
//matrZREcEs = new nDhist<unsigned short>(10, 10, 100, 100);
//matrZREcEs->setFileName(fOdirPrefix+"Psa?ZREcEs.matr");
//matrZREcEs->setComment("core-segment energy as a function of ZR {E = (E_keV-1320)*5}");
//hGroup.add(matrZREcEs);
//matrRZEcEs = new nDhist<unsigned short>(10, 10, 100, 100);
//matrRZEcEs->setFileName(fOdirPrefix+"Psa?RZEcEs.matr");
//matrRZEcEs->setComment("core-segment energy as a function of RZ {E = (E_keV-1320)*5}");
//hGroup.add(matrRZEcEs);
bDoSpec = specEner || specTzero || specSigma;
bDoMats = matrXYZR || matrSeg || matrRE || matrZE || matrZRE || matrRZE || matrDZE || matrZREcEs || matrRZEcEs;
bDoMats = matrXYZR || matrSeg || matrRZE;
#endif //LOCALSPECTRA
#ifdef USEADAPTIVE
......@@ -256,7 +224,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
}
}
if(fWritePsaHits) {
if(fWritePsaHits && !fPsaSegCenter) {
std::ostringstream filename;
filename << "Psa__-12-F__Hits.fdat";
fnPsaHits = fOdirPrefix+filename.str();
......@@ -447,15 +415,16 @@ Int_t PSAFilterGridSearch::Process(int slot)
float fy = bestPoint->y;
float fz = bestPoint->z;
// save the results
int numhit = pD->numHits++;
pD->theHits[numhit].Reset();
pD->theHits[numhit].SetE(netChEner);
pD->theHits[numhit].SetXYZ(fx, fy, fz);
pD->theHits[numhit].SetID(netChSeg); // contains the segment number
pD->theHits[numhit].SetT(pS.bestdt*SAMP_STEP); // means that this is saved in ns
pD->theHits[numhit].SetDT(pS.chi2min); // chii2 temporarily stored in the error DT
pD->theHits[numhit].SetStatus((UShort_t)pS.bestPt); // index of the solution temporarily stored in the status
PsaOut_t *pOut = pD->PsaOut + numhit;
// save the results
pOut->enerSG = netChEner;
pOut->fx = fx; pOut->fy = fy; pOut->fz = fz;
pOut->dTns = pS.bestdt*float(SAMP_STEP); // saved in ns
pOut->chiSq = pS.chi2min;
pOut->netChargeSeg = netChSeg;
pOut->bestPoint = pS.bestPt;
// manual selection of events for waves and histograms
if(doSelect) {
......@@ -463,16 +432,6 @@ Int_t PSAFilterGridSearch::Process(int slot)
pD->selectIt = true;
}
if(fWritePsaHits) {
nDamage_t *nDp = pD->nDamage + numhit;
nDp->enerCC = pS.enerGe; nDp->enerSG = netChEner;
nDp->fx = fx; nDp->fy = fy; nDp->fz = fz;
nDp->fNetSeg = (float)netChSeg;
nDp->ePath = bestPoint->ePath; nDp->hPath = bestPoint->hPath;
nDp->eS_CC = bestPoint->eS_CC; nDp->hS_CC = bestPoint->hS_CC;
nDp->eS_SG = bestPoint->eS_SG; nDp->hS_SG = bestPoint->hS_SG;
}
//if(fWriteNumTraces > 0) {
//// The best point is normalized to the experimental amplitude and accumulated in rAmplitude
StorePartialTrace(pS, bestPoint, scaleFact, samp_first, usamp); // always needed to fit T0
......@@ -555,18 +514,17 @@ void PSAFilterGridSearch::PrepareEvent(PsaData *pD, pointExp *pS)
void PSAFilterGridSearch::SetToSegCenter(PsaData *pD, pointExp *pS)
{
int sMult = pS->numNetCharges;
for(int ns = 0; ns < sMult; ns++) {
PsaOut_t *pOut = pD->PsaOut;
for(int ns = 0; ns < sMult; ns++, pOut++) {
int netChSeg = pS->netChargeSegnum[ns];
float netChEner = pS->netChargeEnergy[ns];
float fx = fBasis.averPt[netChSeg].x;
float fy = fBasis.averPt[netChSeg].y;
float fz = fBasis.averPt[netChSeg].z;
pD->theHits[pD->numHits].Reset();
pD->theHits[pD->numHits].SetE(netChEner);
pD->theHits[pD->numHits].SetXYZ(fx, fy, fz);
pD->theHits[pD->numHits].SetT(0); // means that this is saved in ns
pD->theHits[pD->numHits].SetID(netChSeg); // contains the segment number
pD->numHits++;
pOut->enerSG = pS->netChargeEnergy[ns];
pOut->fx = fBasis.averPt[netChSeg].x;
pOut->fy = fBasis.averPt[netChSeg].y;
pOut->fz = fBasis.averPt[netChSeg].z;
pOut->dTns = 0;
pOut->chiSq = 0;
pOut->netChargeSeg = netChSeg;
pOut->bestPoint = -1;
}
}
......@@ -788,14 +746,14 @@ int PSAFilterGridSearch::WriteTraces(PsaData *pD)
pos = 0;
memset(wbuf, 0, sizeof(wbuf));
}
ADF::PSAHit *pH = &pD->theHits[nh];
PsaOut_t *pOut = pD->PsaOut + nh;
pos++;
wbuf[pos++] = pH->GetID();
wbuf[pos++] = (short)(pH->GetE()*sScale);
wbuf[pos++] = (short)(pH->GetX()*dScale);
wbuf[pos++] = (short)(pH->GetY()*dScale);
wbuf[pos++] = (short)(pH->GetZ()*dScale);
wbuf[pos++] = (short)(dScale*sqrt(pH->GetX()*pH->GetX()+pH->GetY()*pH->GetY()));
wbuf[pos++] = pOut->netChargeSeg;;
wbuf[pos++] = (short)(pOut->enerSG*sScale);
wbuf[pos++] = (short)(pOut->fx*dScale);
wbuf[pos++] = (short)(pOut->fy*dScale);
wbuf[pos++] = (short)(pOut->fx*dScale);
wbuf[pos++] = (short)(dScale*sqrt(pOut->fx*pOut->fx+pOut->fy*pOut->fy));
pos += 3;
}
if(pos) {
......@@ -1038,26 +996,19 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
return 0;
float sumSegm = 0;
float corSegm = 0;
if(!pD->selectIt)
goto skipSpectra;
for(int nh = 0; nh < numHits; nh++) {
int netChSeg = pD->theHits[nh].GetID();
float fx = (float)pD->theHits[nh].GetX();
float fy = (float)pD->theHits[nh].GetY();
float fz = (float)pD->theHits[nh].GetZ();
PsaOut_t * pOut = pD->PsaOut;
for(int nh = 0; nh < numHits; nh++, pOut++) {
int netChSeg = pOut->netChargeSeg;
float fx = pOut->fx;
float fy = pOut->fy;
float fz = pOut->fz;
float rr = sqrt(fx*fx + fy*fy);
float fcorr = 1.f;
if(fz > 15.f && rr < 30.f) {
// empirically derived from matrRZE in test4atc/test4atc-100630_HiRate/...Grilled/Data/1B
fcorr += (30.f-rr)*0.00015f;
}
float netChEner = (float)pD->theHits[nh].GetE();
float corChEner = netChEner*fcorr;
float netChEner = pOut->enerSG;
sumSegm += netChEner;
corSegm += corChEner;
#ifdef _FromGRU_
PSAxy_Proje -> Fill(fx,fy);
PSAxz_Proje -> Fill(fx,fz);
......@@ -1067,14 +1018,14 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
if(bDoMats) {
// positions inside the crystal
if(matrXYZR) {
int nz = (int)(fz/10);
int sz = netChSeg%6 + 4;
matrXYZR->Incr(0, 0, int(fx+matOff), int(fy+matOff));
matrXYZR->Incr(0, 1, int(fx+matOff), int(fz ));
matrXYZR->Incr(0, 2, int(fy+matOff), int(fz ));
matrXYZR->Incr(0, 3, int(rr+matOff), int(fz ));
matrXYZR->Incr(0, sz, int(fx+matOff), int(fy+matOff));
matrXYZR->Incr(1, nz, int(fx+matOff), int(fy+matOff));
int nz = (int)(fz/10);
matrXYZR->Incr( 0, int(fx+matOff), int(fy+matOff));
matrXYZR->Incr( 1, int(fx+matOff), int(fz ));
matrXYZR->Incr( 2, int(fy+matOff), int(fz ));
matrXYZR->Incr( 3, int(rr+matOff), int(fz ));
matrXYZR->Incr( sz, int(fx+matOff), int(fy+matOff));
matrXYZR->Incr(10+nz, int(fx+matOff), int(fy+matOff));
}
if(matrSeg) {
matrSeg->Incr(netChSeg, 0, int(fx+matOff), int(fy+matOff));
......@@ -1088,64 +1039,29 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
// 60Co lines
int ienerCC = (int)(5.f*(pD->CoreE[0] - 1000.f));
int ienerSG = (int)(5.f*(netChEner - 1000.f));
int cenerSG = (int)(5.f*(corChEner - 1000.f));
// 2.2 MeV
//int ienerCC = (int)(5.f*(pD->CoreE[0] - 2000.f));
//int ienerSG = (int)(5.f*(netChEner - 2000.f));
//int cenerSG = (int)(5.f*(corChEner - 2000.f));
//if(matrRE) {
// matrRE->Incr(0, irr, ienerCC);
// matrRE->Incr(1, irr, ienerSG);
// matrRE->Incr(2, irr, cenerSG);
//}
//if(matrZE) {
// matrZE->Incr(0, izz, ienerCC);
// matrZE->Incr(1, izz, ienerSG);
// matrZE->Incr(2, izz, cenerSG);
//}
if(matrZRE) {
matrZRE->Incr(0, izz, irr, ienerCC);
matrZRE->Incr(1, izz, irr, ienerSG);
matrZRE->Incr(2, izz, irr, cenerSG);
}
if(matrRZE) {
matrRZE->Incr(0, irr, izz, ienerCC);
matrRZE->Incr(1, irr, izz, ienerSG);
matrRZE->Incr(2, irr, izz, cenerSG);
}
if(matrDZE) {
int ed = (int)pD->nDamage[0].ePath;
int hd = (int)pD->nDamage[0].hPath;
matrDZE->Incr(0, ed, izz, ienerCC);
matrDZE->Incr(1, hd, izz, ienerCC);
matrDZE->Incr(2, ed, izz, ienerSG);
matrDZE->Incr(3, hd, izz, ienerSG);
}
//ienerSG = (int)(5.f*(netChEner - 1320.f));
//ienerCC = (int)(5.f*(pD->CoreE[0] - 1320.f));
//if(matrZREcEs) {
// matrZREcEs->Incr(izz/5, irr/5, ienerCC, ienerSG);
//}
//if(matrRZEcEs) {
// matrRZEcEs->Incr(irr/5, izz/5, ienerCC, ienerSG);
//}
}
}
if(bDoSpec) {
if(specEner) {
// energy of segments
specEner->Incr(0, netChSeg, (int)(netChEner*fEnergyGain));
specEner->Incr(0, 39, netChSeg);
specEner->Incr(1, netChSeg, (int)(corChEner*fEnergyGain));
specEner->Incr(netChSeg, (int)(netChEner*fEnergyGain));
specEner->Incr(39, netChSeg);
}
if(specTzero) {
// t0 after grid-search (possibly not up to date)
int bestdt = (int)(pD->theHits[nh].GetT());
int bestdt = (int)(pOut->dTns);
specTzero->incr(netChSeg, bestdt+specLenT/2);
}
if(specSigma) {
// chisquare of search
int chi2min = (int)(pD->theHits[nh].GetDT());
int chi2min = (int)(pOut->chiSq);
specSigma->Incr(netChSeg, chi2min);
specSigma->Incr(NCHAN-1, chi2min); // OR of all
}
......@@ -1156,10 +1072,9 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
#ifdef LOCALSPECTRA
if(specEner) {
// total energy and some statistics
specEner->Incr(0, NSEGS, (int)(sumSegm*fEnergyGain)); // 36
specEner->Incr(1, NSEGS, (int)(corSegm*fEnergyGain)); // +36
specEner->Incr(0, NCHAN, (int)(pD->CoreE[0]*fEnergyGain)); // 37
specEner->Incr(0, 39, 100 + numHits); // 39
specEner->Incr(NSEGS, (int)(sumSegm*fEnergyGain)); // 36
specEner->Incr(NCHAN, (int)(pD->CoreE[0]*fEnergyGain)); // 37
specEner->Incr( 39, 100 + numHits); // 39
}
// to use this fragment of old code, maxItime, which calculated in
......@@ -1185,8 +1100,28 @@ skipSpectra:
}
if(fpPsaHits && numHits==1) {
size_t nwritten = fwrite(pD->nDamage, sizeof(nDamage_t), 1, fpPsaHits);
if(nwritten != 1)
PsaOut_t *pOut = pD->PsaOut;
int netCharge = pOut->netChargeSeg;
int bestPoint = pOut->bestPoint;
pointPsa *pBpt = &fBasis.Pts[netCharge][bestPoint];
float vv[12];
vv[ 0] = pD->CoreE[0];
vv[ 1] = pOut->enerSG;
vv[ 2] = pOut->fx;
vv[ 3] = pOut->fy;
vv[ 4] = pOut->fz;
vv[ 5] = (float)netCharge;
vv[ 6] = pBpt->ePath;
vv[ 7] = pBpt->hPath;
vv[ 8] = pBpt->eS_CC;
vv[ 9] = pBpt->hS_CC;
vv[10] = pBpt->eS_SG;
vv[11] = pBpt->hS_SG;
size_t nwritten = fwrite(vv, sizeof(float), 12, fpPsaHits);
if(nwritten != 12)
return 1;
}
......
......@@ -53,13 +53,8 @@ private:
nDhist<unsigned int> *specSigma;
nDhist<unsigned short> *matrXYZR;
nDhist<unsigned short> *matrSeg;
nDhist<unsigned short> *matrZE;
nDhist<unsigned short> *matrRE;
nDhist<unsigned short> *matrZRE;
nDhist<unsigned short> *matrRZE;
nDhist<unsigned short> *matrDZE;
nDhist<unsigned short> *matrZREcEs;
nDhist<unsigned short> *matrRZEcEs;
#endif //LOCALSPECTRA
std::string fnPsaTraces;
......
......@@ -91,6 +91,11 @@ fFrameRaw(NULL)
fEnergyGain = 1.f;
fRescaleHits = false;
fRescaleFile.clear();
for(int ii = 0; ii < 20; ii++)
rescaleCoeffs[ii] = 1.f;
rVx = rVy = rVz = rVc = 0;
aPx = aPy = aPz = 0;
......@@ -116,6 +121,7 @@ fFrameRaw(NULL)
OftMatr_AG = NULL;
OftSpec_TTA = NULL;
OftMatr_gg = NULL;
OftSpec_SC = NULL;
#endif //LOCALSPECTRA
RecoilVc = 0;
......@@ -150,7 +156,7 @@ fFrameRaw(NULL)
}
TrackingFilter::~TrackingFilter()
TrackingFilter::~TrackingFilter()
{
cServer.Finish();
......@@ -378,6 +384,23 @@ void TrackingFilter::GetDataPSA()
}
if(fRescaleHits) {
exyzHit *pt1 = pEXYZ + number_of_hits;
CoreE[0] *= rescaleCoeffs[crystal_id];
CoreE[1] *= rescaleCoeffs[crystal_id];
double scaling = CoreE[0]/esumhits;
if(OftSpec_SC) {
int ival = (int)(scaling*1000);
OftSpec_SC->incr( 0, ival);
OftSpec_SC->incr(nbhits, ival);
}
if(scaling>0.9 && scaling<1.2) {
for(UShort_t nh = 0; nh < nbhits; nh++, pt++) {
pt1->E *= scaling; // scale hits to share the whole core eergy
}
}
}
number_of_hits += nbhits;
if(esumhits > 0)
averageT /= esumhits;
......@@ -975,6 +998,15 @@ void TrackingFilter::process_initialise (UInt_t *error_code)
}
fclose(fp_conf);
if(fRescaleFile.size() > 0) {
bool ok = GetRescaleCoeffs();
if(!ok) {
printf("Error reading rescale coefficients %s\n", fRescaleFile.c_str());
*error_code = 100;
return;
}
}
number_of_tracked_gammas = 0;
// version-specific initializations
......@@ -1236,6 +1268,10 @@ void TrackingFilter::GetParameters(UInt_t *error_code)
else if( stringEq(keyw, "EnergyGain") ) {
ok = 1 == sscanf(data.c_str(), "%f", &fEnergyGain);
}
else if( stringEq(keyw, "RescaleHits") ) {
ok = data.size() > 0;
fRescaleFile = data;
}
else if( stringEq(keyw, "TimeWindowGeAnc") || stringEq(keyw, "TimeWindowGeAncillary") ) {
ok = 2 == sscanf(data.c_str(), "%d %d", &fTsGeAnc_min, &fTsGeAnc_max);
}
......@@ -1337,6 +1373,11 @@ void TrackingFilter::InitGenStructures()
OftSpec_TT->setComment("T_gamma-gamma from CFD (ns)");
hGroup.add(OftSpec_TT);
//OftSpec_SC = new nDhist<unsigned int>(10,10000);
//OftSpec_SC->setFileName(fOdirPrefix+"Oft?SC.spec");
//OftSpec_SC->setComment("Renormalization core-segs (ns)");
//hGroup.add(OftSpec_SC);
//OftSpec_TD = new nDhist<unsigned int>(NDETS, specLenT1);
//OftSpec_TD->setFileName(fOdirPrefix+"Oft?TD.spec");
//OftSpec_TD->setComment("time correction from PSA (ns)");
......@@ -1469,6 +1510,40 @@ int TrackingFilter::GetRotoTranslations(FILE *fp)
return nt;
}
bool TrackingFilter::GetRescaleCoeffs()
{
std::string conf_file = fConfPath + fRescaleFile;
printf("Opening file with rescale coefficients for the cores %s\n", conf_file.c_str());
FILE *fp_conf = fopen(conf_file.c_str(), "r");
if(!fp_conf) {
printf("Could not open file %s\n", conf_file.c_str());
return false;
}
for(int ii = 0; ii < 20; ii++)
rescaleCoeffs[ii] = 1.f;
int nt = 0;
while(true) {
int nd;
float xx;
// detector number and rescale coefficient
if(fscanf(fp_conf, "%d %f", &nd, &xx) != 2)
break;
if(nd < 0 || nd >= 20)
return false;
rescaleCoeffs[nd] = xx;
nt++;
}
printf("Found %d coefficients\n", nt);
fRescaleHits = true;
fclose(fp_conf);
return true;
}
// recoil along a fixed direction (beam)
double TrackingFilter::DopplerCorrectFix(trGamma *pg)
{
......
......@@ -87,11 +87,11 @@ protected:
double X;
double Y;
double Z;
double T;
ULong64_t tst;
UInt_t evn;
int Id;
int Sg;
double T; // time info from PSA
ULong64_t tst; // time stamp
UInt_t evn; // event number
int Id; // detector ID
int Sg; // segment number
} exyzHit;
typedef struct {
......@@ -163,6 +163,7 @@ protected:
nDhist<unsigned short> *OftMatr_AG; // ind-anc-gamma for the ancillary detector
nDhist<unsigned short> *OftSpec_TTA; // ind-anc-Tgamma for the ancillary detector
nDhist<unsigned int> *OftMatr_gg; // gamma-gamma
nDhist<unsigned int> *OftSpec_SC; // Scaling factor
#endif //LOCALSPECTRA
// fixed values to correct when there is no Ancillary defining the recoil vector
......@@ -220,6 +221,10 @@ protected:
Float_t fEnergyGain; // gain of energy spectra
Bool_t fRescaleHits;
std::string fRescaleFile;
float rescaleCoeffs[20];
bool fAncillary;
bool fRecoiling;