Commit cdcba22e authored by dino's avatar dino
Browse files

Removed simplification of grid-search for very small energies

git-svn-id: svn://gal-serv.lnl.infn.it/agata/trunk/narval_emulator@1054 170316e4-aea8-4b27-aad4-0380ec0519c9
parent 288c9174
......@@ -25,7 +25,7 @@ const int numXYZ = 10; // sub matrXYZ for up to numXYZ-1 hits
const size_t nPsaHitsValues = 16; // number of parameters written in Psa__-nn-F__Hits.fdat
const bool doFit_T0 = true; // enables the GridSearch-FitT0AfterPSA looping
const bool doFit_T0 = false /*true*/ ; // enables the GridSearch-FitT0AfterPSA looping
const float tauDefault = 35.f; // default decay constant of the preamp response function for segments and core (in nanoseconds)
......@@ -39,7 +39,7 @@ const float lowSegEnergy = 50.f; // segment-energies smaller than this are deco
//#define TRACESXYZ // WriteTraces() sorts on position files xx_yy_zz_xyz.dat
//#define MATRCHI2 // 3D matrix with distribution of chi2; meaningful only for segmult=1; done only in SearchFullGrid
// If enabled, these macros are used to write the waves at various stages of the processing of one event.
// If enabled, the following macros are used to write the waves at various stages of the processing of one event.
// The output file debug__0-37-56-F__waves.fsamp is written in pwd and may or may not be rewound at the beginning of every event.
// They produce a large amoutn of output and are normally used in Debug mode as a graphical inspection tool.
// They are not thread-safe and, due to the fixed-name ouput file, incompatible with multiple crystal analyses.
......@@ -50,7 +50,6 @@ const float lowSegEnergy = 50.f; // segment-energies smaller than this are deco
# define WRITEWORKINGWAVE(p1,p2)
# define WRITEPARTIALWAVE(p1,p2,p3)
#endif
bool firstEventPSA = true;
PSAFilterGridSearch::PSAFilterGridSearch()
{
......@@ -68,6 +67,8 @@ PSAFilterGridSearch::PSAFilterGridSearch()
gPsaSegCenter = false;
gFullSearch = false;
firstEventPSA = true;
memset(hmask, ' ', sizeof(hmask));
for(int nn = 0; nn < NSEGS; nn++) { // segments and core
hmask[nn][INDCC+1] = 0; // +1 to close each line as a string
......@@ -78,7 +79,7 @@ PSAFilterGridSearch::PSAFilterGridSearch()
#ifdef PSA_LOCALSPECTRA
baseZXY = NULL; specEner = NULL; specStat = NULL; specTzero = NULL; specSigma = NULL; specMetrics = NULL;
matrXYZR = NULL; matrSeg = NULL; matrRZE = NULL;
matrXYZR = NULL; matrXYZ = NULL; matrSeg = NULL; matrRZE = NULL;
#endif //PSA_LOCALSPECTRA
fHitSegThreshold = minSegEnergy; // keV
......@@ -225,9 +226,11 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
}
tauSGCC[36] = fTauSlice[6];
}
fBasis.ExponentialResponse(tauSGCC); // convolution with an exponential function
//fBasis.ExperimentalResponse(); // convolution with the experimental response derived from the pulser
// convolution with an exponential function
fBasis.ExponentialResponse(tauSGCC);
// convolution with the experimental response derived from the pulser
//fBasis.ExperimentalResponse();
// apply the cross-talk correction to the basis signals
fBasis.XtalkCorrection(fDeadSegment);
......@@ -248,14 +251,15 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
#if 0 && defined(PSA_LOCALSPECTRA)
{
// The whole basis (2 mm packed). Normalization is 30000
// Write the complete basis as a <short> spectrum. 2 mm packed. Normalization is 30000.
// write only NSEGS+CC as they are going to be used (i.e. with differentiated netChargeSeg and CC)
// To get the non-differentiated version uncomment //#define USE_NETCHARGES in FridSearchParameters.h
// To get the whole internal basis replace NCHAN with TCHAN in the next row and in the copy to pd loop
// To get the whole internal basis replace NCHAN with TCHAN in the definition of numTraces
const int fineStep = fstep;
const int detRadius = 40;
const int detLength = 90;
nDhist<short> *baseXYZ = new nDhist<short>((2*detRadius)/fineStep, (2*detRadius)/fineStep, detLength/fineStep, NCHAN, NTIME);
const int numTraces = NCHAN;
nDhist<short> *baseXYZ = new nDhist<short>((2*detRadius)/fineStep, (2*detRadius)/fineStep, detLength/fineStep, numTraces, NTIME);
baseXYZ->setFileName(fOdirPrefix+"Psa?BaseXYZ.spec");
baseXYZ->setComment("the whole PSA base");
int errcount = 0;
......@@ -279,7 +283,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
int iz = int(bpt->z )/fineStep;
short *pd = baseXYZ->getData(ix, iy, iz, 0);
if(pd) {
for (int is = 0; is < NCHAN; is++) {
for (int is = 0; is < numTraces; is++) {
for(int it = 0; it < NTIME; it++) {
*pd++ = short(bb[is][it]*fact);
}
......@@ -386,7 +390,7 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
largest++; // this line is just to be able to set a breakpoint
#endif // #else // #if 0
}
#if 1 && defined(PSA_LOCALSPECTRA)
#if 1 && defined(PSA_LOCALSPECTRA) && !defined(USE_SSE_VERSION)
{
// The precalculated metrics
nDhist<float> *vMetrics = new nDhist<float>(NMETRIC);
......@@ -418,6 +422,11 @@ Int_t PSAFilterGridSearch::AlgoSpecificInitialise()
matrXYZR->setComment("XY, XZ, YZ, rZ hits and XY hits as a function of Z");
hGroup.add(matrXYZR);
matrXYZ = new nDhist<unsigned short>(3, matLen/matPack, matLen/matPack, matLen/matPack);
matrXYZ->setFileName(fOdirPrefix+"Psa?XYZ.matr");
matrXYZ->setComment("zXY, yXZ, xYZ distributions of hits");
hGroup.add(matrXYZ);
specStat = new nDhist<unsigned int>(40, 100);
specStat->setFileName(fOdirPrefix+"Psa?Stat.spec");
specStat->setComment ("Distribution of hits");
......@@ -753,7 +762,7 @@ Int_t PSAFilterGridSearch::Process(int slot, int nslots)
// stop looping if shift too small or total-shift too large
if(abs(dT0)<dT0min || abs(pD->shiftT0)>shiftT0max)
break;
// at the momentim looping through the one-hit search is not compatible with the two-hits solutions
// at the moment, looping through the one-hit search is not compatible with the two-hits solutions
if(gTryTwoHits)
break;
......@@ -819,10 +828,10 @@ int PSAFilterGridSearch::ProcessTheEvent(pointFull &pS)
// this must become a setup parameter from the configuration file
const bool onlyFirst = false; // just to see the effect of placing the results at the center of the segment
pS.samp_used = SearchAdaptive1(&pS, netChSeg, onlyCoarse, onlyFirst); // coarse-fine or coarse-only with one hit
if(simplified || gCoarseOnly) {
int tmp = pS.bestPt1; pS.bestPt1 = pS.randPt1; pS.randPt1 = tmp; // just take a random point
// REMEMBER TO DO IT ALSO IN ProcessTwoHits
}
//if(simplified || gCoarseOnly) {
// int tmp = pS.bestPt1; pS.bestPt1 = pS.randPt1; pS.randPt1 = tmp; // just take a random point
// // REMEMBER TO DO IT ALSO IN ProcessTwoHits
//}
pS.bestFac = 1.f; // second point not present
pS.bestPt2 = -1;
}
......@@ -2052,11 +2061,6 @@ void PSAFilterGridSearch::PreSearchCoarse(pointFull &pS)
Int_t PSAFilterGridSearch::PostProcess(int slot)
{
//#define xAddr(fx) ( int((fx+matOff)/matPack) )
//#define yAddr(fy) ( int((fy+matOff)/matPack) )
//#define zAddr(fz) ( int((fz )/matPack) )
//#define rAddr(rr) ( int((rr+matOff)/matPack) )
PsaData *pD = DD + slot;
int numNetSegs = pD->numNetSegs;
......@@ -2154,7 +2158,8 @@ Int_t PSAFilterGridSearch::PostProcess(int slot)
PSAyz_Proje -> Fill(fy,fz);
#endif
// positions inside the crystal
if(matrXYZR || matrSeg) {
if(matrXYZR || matrSeg ||matrXYZ) {
incrPosMats(-1, 0, fx, fy, fz, rr); // matrXYZ
if(pOut->enerSG >= lowSegEnergy) { // not too-low energy
incrPosMats(netChSeg, 0, fx, fy, fz, rr); // 0 --> all hits
if(nh < numXYZ-3)
......@@ -2297,7 +2302,15 @@ void PSAFilterGridSearch::incrPosMats(int netChSeg, int nhi, float fx, float fy,
#define zAddr(fz) ( int((fz )/matPack) )
#define rAddr(rr) ( int((rr+matOff)/matPack) )
if(matrXYZR) {
if(netChSeg < 0) {
if(matrXYZR) {
matrXYZ->Incr(0, zAddr(fz), xAddr(fx), yAddr(fy) );
matrXYZ->Incr(1, yAddr(fy), xAddr(fx), zAddr(fz) );
matrXYZ->Incr(2, xAddr(fx), yAddr(fy), zAddr(fz) );
}
return;
}
if(matrXYZ) {
int sz = netChSeg%6;
int nz = (int)(fz/10);
int nw = (int)(fz/2 );
......@@ -2552,12 +2565,12 @@ Float_t PSAFilterGridSearch::FitT0AfterPSA(pointFull &pS)
rfit = tRfit + SAMP0/2;
const int nsteps = 3;
const int nsamps = SAMP0/2+20; // calculate chi2 using SAMP0/2 of pretrigger and 20 of risetime
float x2[2*nsteps+1]; memset(x2, 0, sizeof(x2));
for(int ll = 0; ll <= 2*nsteps; ll++) {
float *data = tData + SAMP0/2 - nsteps + ll;
float x2t = 0;
// calculate chi2 using SAMP0/2 of pretrigger and 20 of risetime
for(int ii = nsteps; ii < SAMP0/2+20; ii++) {
for(int ii = nsteps; ii < nsamps; ii++) {
float diff = float(data[ii] - rfit[ii]);
x2t += diff*diff;
}
......@@ -2583,7 +2596,7 @@ Float_t PSAFilterGridSearch::FitT0AfterPSA(pointFull &pS)
float ys = x2[ssmin ];
float yp = x2[ssmin+1];
float y0 = -0.5f*(yp-ym)/(yp+ym-2.f*ys);
float y0 = -0.5f*(yp-ym)/(yp+ym-2.f*ys); // can produce a divide by 0 !!
return ssmin+y0-nsteps;
}
......@@ -93,6 +93,8 @@ private:
bool gPsaSegCenter;
bool gFullSearch;
bool firstEventPSA;
float tauSGCC[37]; // the 36 segments and the CC
#ifdef PSA_LOCALSPECTRA
......@@ -102,6 +104,7 @@ private:
nDhist<unsigned int > *specTzero;
nDhist<unsigned int > *specSigma;
nDhist<unsigned int > *specMetrics;
nDhist<unsigned short> *matrXYZ;
nDhist<unsigned short> *matrXYZR;
nDhist<unsigned short> *matrSeg;
nDhist<unsigned short> *matrRZE;
......
......@@ -22,7 +22,7 @@ const float MAXNORM = 1000.f; // normalization of the signal basis (1 MeV,
#define USE_NETCHARGES // to use the net charge signals in the grid search
#define USE_CORETRACE // to use also the core signal in the grid search
#define SIMPLIFY_SMALL_ENERGIES // simplified search for segments with low energy if multiplicity > 2
//#define SIMPLIFY_SMALL_ENERGIES // simplified search for low energy hits if segment multiplicity > 2
#ifdef USE_NETCHARGES
const int DIFFLAG = 0; // use the net-charge signals as they are
......
......@@ -54,7 +54,6 @@ const int speTTlen = int(1000*tScale);
# define WRITEWORKINGWAVE1(p1,p2)
# define WRITEWORKINGWAVE2(p1,p2)
#endif
bool firstEventPP = true;
PreprocessingFilterPSA::PreprocessingFilterPSA() :
#ifdef PPF_LOCALSPECTRA
......@@ -77,6 +76,8 @@ PreprocessingFilterPSA::PreprocessingFilterPSA() :
fTC = fTS = fTT = fWW = NULL;
fTFA = fCFD = NULL;
firstEventPP = true;
#ifdef READ_TSHIFT_CORRECTION
tstmpFileI = NULL;
fTimestampOld = 0;
......@@ -564,8 +565,7 @@ bool PreprocessingFilterPSA::initPresort()
}
// The new version with triggering on the sum of core + net-charge segments
// - the timing is done in FindTimeCFD by a CFD
// - the timing is done in FindTriggerTime by fitting a stright line to the first samples of the trace
// Done in FindTriggerTime by fitting a stright line to the first samples of the trace
bool PreprocessingFilterPSA::ProcessEvent(int &smult, float &eSumSG1, float &eSumSG2)
{
......
......@@ -111,6 +111,8 @@ private:
int segF1; // gate on segment multiplicity
int segF2; // gate on segment multiplicity
bool firstEventPP;
float deadSegFactor; // dawn-scaling of core trace to consider xTalk correction
std::string fnPrepTraces;
......
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