Commit 217c1dd9 authored by dino's avatar dino
Browse files

Preliminary implementation neutron damage treatment.

The latest Bart's signal bases contain the path length and sensitivity to charge trapping for electrons and holes.
To develop the correction methods the PSA writes out a binary file Psa__-12-F__Hits.fdat with the needed information (--> should become a Root Tree!).
Be aware that at the moment this is done in PSAFilterGridSearch:Process() (and not in PostProcess) )and therefore is not thread-safe.

git-svn-id: svn://gal-serv.lnl.infn.it/agata/trunk/narval_emulator@892 170316e4-aea8-4b27-aad4-0380ec0519c9
parent b06ff59e
......@@ -15,7 +15,7 @@ using namespace std;
//
const float dScale = 10.f; // mm
const float sScale = 1.f; // chan/keV of saved traces traces
const float sScale = 1.f; // chan/keV of saved traces
const int specLenE = 16*1024;
const int specLenT = 2*1000;
const int matLen = 100;
......@@ -26,6 +26,7 @@ const bool doSelect = false; // the selection has to be programmed manually i
PSAFilterGridSearch::PSAFilterGridSearch()
{
fpPsaTraces = NULL;
fpPsaHits = NULL;
for(int nn = 0; nn < NCHAN; nn++) {
memset(hmask[nn], ' ', NCHAN);
......@@ -35,7 +36,8 @@ PSAFilterGridSearch::PSAFilterGridSearch()
#ifdef LOCALSPECTRA
specEner = NULL; specTzero = NULL; specSigma = NULL;
matrXYZR = NULL; matrSeg = NULL;
matrRE = NULL; matrZE = NULL; matrZRE = NULL; matrRZE = NULL; matrZREcEs = NULL; matrRZEcEs = NULL;
matrRE = NULL; matrZE = NULL; matrZRE = NULL; matrRZE = NULL;
matrDZE = NULL; matrZREcEs = NULL; matrRZEcEs = NULL;
#endif //LOCALSPECTRA
bDoSpec = false;
bDoMats = false;
......@@ -144,7 +146,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
fBasis.TshiftCorrection();
// the netcharge segment is duplicated into netcharge to allow for the use of its delayed-differentiated version
// pay attention to the fact that the time shifts should be done before this
// pay attention to the fact that cross-talk corrections, signal filtering and time shifts should be done before this
fBasis.StoreNetCharges();
if(DIFFLAG > 0) {
......@@ -189,15 +191,15 @@ 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);
//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);
//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");
......@@ -209,6 +211,11 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
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}");
......@@ -220,7 +227,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
//hGroup.add(matrRZEcEs);
bDoSpec = specEner || specTzero || specSigma;
bDoMats = matrXYZR || matrSeg || matrRE || matrZE || matrZRE || matrRZE || matrZREcEs || matrRZEcEs;
bDoMats = matrXYZR || matrSeg || matrRE || matrZE || matrZRE || matrRZE || matrDZE || matrZREcEs || matrRZEcEs;
#endif //LOCALSPECTRA
#ifdef USEADAPTIVE
......@@ -249,6 +256,21 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
}
}
#if 1
// be aware that the present implementation does not work with threads
{
std::ostringstream filename;
filename << "Psa__-12-F__Hits.fdat";
fnPsaHits = fOdirPrefix+filename.str();
cout << "Opening PSA Hits file " << fnPsaHits << endl;
fpPsaHits = fopen(fnPsaHits.c_str(), "wb");
if(!fpPsaHits) {
cout << "Error: could not open " << fnPsaHits << endl;
return 1;
}
}
#endif
return 0;
}
......@@ -437,6 +459,8 @@ Int_t PSAFilterGridSearch::Process(int slot)
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
pD->theHits[numhit].SetDX(bestPoint->ePath); // electron path-length temporarily stored in the error Dx
pD->theHits[numhit].SetDY(bestPoint->hPath); // hole path-length temporarily stored in the error Dy
// manual selection of events for waves and histograms
if(doSelect) {
......@@ -444,6 +468,20 @@ Int_t PSAFilterGridSearch::Process(int slot)
pD->selectIt = true;
}
if(fpPsaHits && sMult==1) {
// be aware that the present implementation does not work with threads
float wtmp[12];
wtmp[ 0] = pS.enerGe; wtmp[1] = netChEner;
wtmp[ 2] = fx; wtmp[3] = fy; wtmp[4] = fz;
wtmp[ 5] = (float)netChSeg;
wtmp[ 6] = bestPoint->ePath; wtmp[ 7] = bestPoint->hPath;
wtmp[ 8] = bestPoint->eS_CC; wtmp[ 9] = bestPoint->hS_CC;
wtmp[10] = bestPoint->eS_SG; wtmp[11] = bestPoint->hS_SG;
size_t nwritten = fwrite(wtmp, sizeof(float), 12, fpPsaHits);
if(nwritten != 12)
return 1;
}
//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
......@@ -1056,19 +1094,24 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
if(numHits==1) {
int irr = (int)rr;
int izz = (int)(0.5*fz);
// 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));
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);
}
// 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);
......@@ -1079,14 +1122,22 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
matrRZE->Incr(1, irr, izz, ienerSG);
matrRZE->Incr(2, irr, izz, cenerSG);
}
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(matrDZE) {
int ed = (int)pD->theHits[nh].GetDX();
int hd = (int)pD->theHits[nh].GetDY();
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) {
......
#ifndef ADF_PSAFILTERGRID_SEARCH_H
#define ADF_PSAFILTERGRID_SEARCH_H
#include "NarvalInterface.h"
#include "AgataFrameFactory.h"
#include "AgataKeyFactory.h"
#include "Trigger.h"
#include "commonDefs.h"
//#include "NarvalInterface.h"
//#include "AgataFrameFactory.h"
//#include "AgataKeyFactory.h"
//#include "Trigger.h"
//
//#include "commonDefs.h"
#include "PSAFilter.h"
#include "SignalBasis.h"
......@@ -57,6 +57,7 @@ private:
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
......@@ -64,6 +65,9 @@ private:
std::string fnPsaTraces;
FILE *fpPsaTraces;
std::string fnPsaHits;
FILE *fpPsaHits;
#ifdef _FromGRU_
GNetServerRoot *PSANetworkRoot;
GSpectra *PSASpectraDB;
......
......@@ -24,14 +24,22 @@ public:
char masksignals[NCHAN+1]; // usage of segments (0 added at the end to close the string)
bool fullChargeCollected; // reached full amplitude
float ePath; // total distance travelled by the electrons
float hPath; // total distance travelled by the holes
float eS_CC; // Bart's sensitivity for electrons to the anode
float hS_CC; // Bart's sensitivity for holes to the Anode
float eS_SG; // Bart's sensitivity for electrons to the cathode
float hS_SG; // Bart's sensitivity for holes to the cathode
pointPsa() {
Init();
}
void Init() {
x = 0;
y = 0;
z = 0;
x = y = z = 0;
ePath = hPath= 0;
eS_CC = hS_CC = 0;
eS_SG = hS_SG = 0;
netChargeSegment = NCHAN-1; // set to CC ??
addr_first = TSHIFT;
addr_last = NTIME-1;
......
......@@ -147,15 +147,22 @@ int SignalBasis::ReadBasisFormatVentu(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
if(ntseg != ii)
nerr++;
}
if(!keep) {
delete [] Pts[ii];
Pts[ii] = NULL;
}
//if(!keep) {
// delete [] Pts[ii];
// Pts[ii] = NULL;
//}
}
fclose(ofp1);
CalculateSegmentCenters();
if(!keep) {
for(int ns = 0; ns < NSEGS; ns++) {
delete [] Pts[ns];
Pts[ns] = NULL;
}
}
bCenters = true;
bFullGrid = keep;
......@@ -288,9 +295,11 @@ int SignalBasis::ReadBasisFormatBartB(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
int jj = numPts[ii];
Pts[ii][jj].fullChargeCollected = true;
Pts[ii][jj].netChargeSegment = ii;
Pts[ii][jj].x = bbhit.Pos[0];
Pts[ii][jj].y = bbhit.Pos[1];
Pts[ii][jj].z = bbhit.Pos[2] - 1.f; // remove the 1 mm z-offset due to SimIon
Pts[ii][jj].ClearAmplitude(); // zero all samples
// move the SAMP_STEP-spaced samples to pts.amplitude
......@@ -306,6 +315,12 @@ int SignalBasis::ReadBasisFormatBartB(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
}
Pts[ii][jj].amplitude[ns][ann] = (gs_type)last;
}
if(ns==ii) { // the net charge
Pts[ii][jj].ePath = ptf[1]; // distance travelled buy the electrons
Pts[ii][jj].eS_SG = ptf[2]; // sensitivity for electrons
Pts[ii][jj].hPath = ptf[5]; // distance travelled buy the holes
Pts[ii][jj].hS_SG = ptf[6]; // sensitivity for holes
}
}
// in Bart's files, CC is in the first position
{
......@@ -319,6 +334,10 @@ int SignalBasis::ReadBasisFormatBartB(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
}
Pts[ii][jj].amplitude[NSEGS][ann] = (gs_type)last;
}
Pts[ii][jj].ePath = ptf[1]; // distance travelled buy the electrons (to check if it is the same)
Pts[ii][jj].eS_CC = ptf[2]; // sensitivity for electrons
Pts[ii][jj].hPath = ptf[5]; // distance travelled buy the holes
Pts[ii][jj].hS_CC = ptf[6]; // sensitivity for holes
}
// insert mask provided by the caller
......@@ -330,6 +349,8 @@ int SignalBasis::ReadBasisFormatBartB(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
}
CalculateSegmentCenters();
if(!keep) {
for(int ns = 0; ns < NSEGS; ns++) {
delete [] Pts[ns];
......@@ -337,8 +358,6 @@ int SignalBasis::ReadBasisFormatBartB(FILE * ofp1, char hmask[NCHAN][NCHAN+1], b
}
}
CalculateSegmentCenters();
bCenters = true;
bFullGrid = true;
......
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