quadrant finding moved to singlePhotnDetector class

git-svn-id: file:///afs/psi.ch/project/sls_det_software/svn/slsDetectorCalibration@30 113b152e-814d-439b-b186-022a431db7b5
This commit is contained in:
bergamaschi 2014-02-03 13:58:07 +00:00
parent fea8670994
commit f3e58f95ea
2 changed files with 86 additions and 110 deletions

View File

@ -83,7 +83,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
THStack *hs=new THStack("hs",fformat);
int iev=0;
TH2F *h1=new TH2F("h1",tit,nbins,hmin-0.5,hmax-0.5,NC*NR,-0.5,NC*NR-0.5);
hs->Add(h1);
@ -113,7 +112,6 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
tall=filter->initEventTree(tit, &iFrame);
double tot, tl, tr, bl, br, v;
#ifdef MY_DEBUG
@ -152,16 +150,16 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
if (hitfinder) {
filter->newFrame();
filter->newFrame();
//calculate pedestals and common modes
if (cmsub) {
if (cmsub) {
// cout << "cm" << endl;
for (ix=xmin-1; ix<xmax+1; ix++)
for (iy=ymin-1; iy<ymax+1; iy++) {
thisEvent=filter->getEventType(buff, ix, iy,0);
}
}
for (ix=xmin-1; ix<xmax+1; ix++)
for (iy=ymin-1; iy<ymax+1; iy++) {
thisEvent=filter->getEventType(buff, ix, iy,0);
}
}
}
// cout << "new frame " << endl;
@ -179,87 +177,14 @@ THStack *moenchReadData(char *fformat, char *tit, int runmin, int runmax, int nb
#endif
if (nf>1000) {
h1->Fill(filter->getClusterElement(0,0), iy+NR*ix);
h1->Fill(filter->getClusterTotal(1), iy+NR*ix);
if (hitfinder) {
tot=0;
tl=0;
tr=0;
bl=0;
br=0;
// if (nf%1000==0 && ix==20 && iy==20) cout << " val="<< decoder->getClusterElement(0,0)<< endl;
if (thisEvent==PHOTON_MAX ) {
nph++;
for (ir=-1; ir<2; ir++) {
for (ic=-1; ic<2; ic++) {
v=filter->getClusterElement(ic,ir);
// data[ic+1][ir+1]=v;
tot+=v;
if (ir<1) {
if (ic<1) {
bl+=v;
}
if (ic>-1) {
br+=v;
}
}
if (ir>-1) {
if (ic<1) {
tl+=v;
}
if (ic>-1) {
tr+=v;
}
}
}
}
if (bl>br && bl>tl && bl>tr) {
h2->Fill(bl, iy+NR*ix);
if (bl>0) {
hetaX->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(0,-1))/bl,iy+NR*ix);
hetaY->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(-1,0))/bl,iy+NR*ix);
hetaX->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(0,-1))/bl,iy+NR*(ix-1));
hetaY->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(-1,0))/bl,(iy-1)+NR*ix);
}
} else if (br>bl && br>tl && br>tr) {
h2->Fill(br, iy+NR*ix);
if (br>0) {
hetaX->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(0,-1))/br,iy+NR*ix);
hetaY->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(1,0))/br,iy+NR*ix);
hetaX->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(0,-1))/br,iy+NR*(ix+1));
hetaY->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(1,0))/br,iy-1+NR*ix);
}
} else if (tl>br && tl>bl && tl>tr) {
h2->Fill(tl, iy+NR*ix);
if (tl>0) {
hetaX->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(0,1))/tl,iy+NR*ix);
hetaY->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(-1,0))/tl,iy+NR*ix);
hetaX->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(0,1))/tl,iy+NR*(ix-1));
hetaY->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(-1,0))/tl,iy+1+NR*ix);
}
} else if (tr>bl && tr>tl && tr>br) {
h2->Fill(tr, iy+NR*ix);
if (tr>0) {
hetaX->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(0,1))/tr,iy+NR*ix);
hetaY->Fill((filter->getClusterElement(0,0)+filter->getClusterElement(1,0))/tr,iy+NR*ix);
hetaX->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(0,1))/tr,iy+NR*(ix+1));
hetaY->Fill(1.-(filter->getClusterElement(0,0)+filter->getClusterElement(1,0))/tr,iy+1+NR*ix);
}
}
h3->Fill(tot, iy+NR*ix);
h2->Fill(filter->getClusterTotal(2), iy+NR*ix);
h3->Fill(filter->getClusterTotal(3), iy+NR*ix);
iFrame=decoder->getFrameNumber(buff);
tall->Fill();

View File

@ -28,7 +28,15 @@ using namespace std;
PHOTON=2,
PHOTON_MAX=3,
NEGATIVE_PEDESTAL=4,
UNDEFINED=-1
UNDEFINED_EVENT=-1
};
enum quadrant {
TOP_LEFT=0,
TOP_RIGHT=1,
BOTTOM_LEFT=2,
BOTTOM_RIGHT=3,
UNDEFINED_QUADRANT=-1
};
@ -62,7 +70,7 @@ class singlePhotonDetector {
int sign=1,
commonModeSubtraction *cm=NULL,
int nped=1000,
int nd=100) : det(d), nx(0), ny(0), stat(NULL), cmSub(cm), nDark(nd), eventMask(NULL),nSigma (nsigma), clusterSize(csize), clusterSizeY(csize), cluster(NULL), iframe(-1), dataSign(sign) {
int nd=100) : det(d), nx(0), ny(0), stat(NULL), cmSub(cm), nDark(nd), eventMask(NULL),nSigma (nsigma), clusterSize(csize), clusterSizeY(csize), cluster(NULL), iframe(-1), dataSign(sign), quad(UNDEFINED_QUADRANT), tot(0), quadTot(0) {
det->getDetectorSize(nx,ny);
@ -94,7 +102,7 @@ class singlePhotonDetector {
void newDataSet(){iframe=-1; for (int iy=0; iy<ny; iy++) for (int ix=0; ix<nx; ix++) stat[iy][ix].Clear(); if (cmSub) cmSub->Clear(); };
/** resets the eventMask to undefined and the commonModeSubtraction */
void newFrame(){iframe++; for (int iy=0; iy<ny; iy++) for (int ix=0; ix<nx; ix++) eventMask[iy][ix]=UNDEFINED; if (cmSub) cmSub->newFrame();};
void newFrame(){iframe++; for (int iy=0; iy<ny; iy++) for (int ix=0; ix<nx; ix++) eventMask[iy][ix]=UNDEFINED_EVENT; if (cmSub) cmSub->newFrame();};
/** sets the commonModeSubtraction algorithm to be used
@ -184,8 +192,13 @@ class singlePhotonDetector {
eventType getEventType(char *data, int ix, int iy, int cm=0) {
// eventType ret=PEDESTAL;
double tot=0, max=0;
double max=0, tl=0, tr=0, bl=0,br=0, v;
// cout << iframe << endl;
tot=0;
quadTot=0;
quad=UNDEFINED_QUADRANT;
if (iframe<nDark) {
if (cm==0) {
@ -194,7 +207,7 @@ class singlePhotonDetector {
cluster->set_data(dataSign*(det->getValue(data, ix, iy)-getPedestal(ix,iy,cm)), 0,0 );
// cout << "=" << endl;
}
return UNDEFINED;
return UNDEFINED_EVENT;
}
@ -214,12 +227,22 @@ class singlePhotonDetector {
for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) {
if ((iy+ir)>=0 && (iy+ir)<ny && (ix+ic)>=0 && (ix+ic)<nx) {
cluster->set_data(dataSign*(det->getValue(data, ix+ic, iy+ir)-getPedestal(ix+ic,iy+ir,cm)), ic, ir );
tot+=cluster->get_data(ic,ir);
v=cluster->get_data(ic,ir);
tot+=v;
if (ir<=0 && ic<=0)
bl+=v;
if (ir<=0 && ic>=0)
br+=v;
if (ir>=0 && ic<=0)
tl+=v;
if (ir>=0 && ic>=0)
tr+=v;
if (cluster->get_data(ic,ir)>max) {
max=cluster->get_data(ic,ir);
max=v;
}
if (ir==0 && ic==0) {
if (cluster->get_data(ic,ir)>nSigma*cluster->rms) {
if (v>nSigma*cluster->rms) {
eventMask[iy][ix]=PHOTON;
} else if (cluster->get_data(ic,ir)<-nSigma*cluster->rms)
eventMask[iy][ix]=NEGATIVE_PEDESTAL;
@ -233,31 +256,54 @@ class singlePhotonDetector {
} else if (eventMask[iy][ix]==PHOTON) {
if (cluster->get_data(0,0)>=max) {
eventMask[iy][ix]=PHOTON_MAX;
/* if (iframe%1000==0) { */
/* for (int ir=-(clusterSizeY/2); ir<(clusterSizeY/2)+1; ir++) { */
/* for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) { */
/* if ((iy+ir)>=0 && (iy+ir)<ny && (ix+ic)>=0 && (ix+ic)<nx) { */
/* /\* if (eventMask[iy+ir][ix+ic]==UNDEFINED) *\/ */
/* /\* eventMask[iy+ir][ix+ic]=NEIGHBOUR; *\/ */
/* cout << cluster->get_data(ic,ir) << " "; */
/* } */
/* } */
/* } */
/* cout << endl;; */
/* } */
}
} else if (eventMask[iy][ix]==PEDESTAL) {
if (cm==0)
addToPedestal(det->getValue(data, ix, iy),ix,iy);
}
// }
if (bl>=br && bl>=tl && bl>=tr) {
quad=BOTTOM_LEFT;
quadTot=bl;
} else if (br>=bl && br>=tl && br>=tr) {
quad=BOTTOM_RIGHT;
quadTot=br;
} else if (tl>=br && tl>=bl && tl>=tr) {
quad=TOP_LEFT;
quadTot=tl;
} else if (tr>=bl && tr>=tl && tr>=br) {
quad=TOP_RIGHT;
quadTot=tr;
}
return eventMask[iy][ix];
};
/**<
retrurns the total signal in a cluster
\param size cluser size should be 1,2 or 3
\returns cluster center if size=1, sum of the maximum quadrant if size=2, total of the cluster if size=3 or anything else
*/
double getClusterTotal(int size) {
switch (size) {
case 1:
return getClusterElement(0,0);
case 2:
return quadTot;
default:
return tot;
};
};
/**<
retrurns the quadrant with maximum signal
\returns quadrant where the cluster is located */
quadrant getQuadrant() {return quad;};
/** sets/gets number of samples for moving average pedestal calculation
\param i number of samples to be set (0 or negative gets)
@ -279,6 +325,7 @@ class singlePhotonDetector {
*/
eventType getEventMask(int ic, int ir=0){return eventMask[ir][ic];};
#ifdef MYROOT1
/** generates a tree and maps the branches
\param tname name for the tree
@ -298,8 +345,8 @@ class singlePhotonDetector {
char tit[100];
sprintf(tit,"data[%d]/D",clusterSize*clusterSizeY);
tall->Branch("data",cluster->data,tit);
tall->Branch("pedestal",&(cluster->ped),"pedestal/D");
tall->Branch("rms",&(cluster->rms),"rms/D");
// tall->Branch("pedestal",&(cluster->ped),"pedestal/D");
// tall->Branch("rms",&(cluster->rms),"rms/D");
return tall;
};
#endif
@ -322,7 +369,11 @@ class singlePhotonDetector {
single_photon_hit *cluster; /**< single photon hit data structure */
int iframe; /**< frame number (not from file but incremented within the dataset every time newFrame is called */
int dataSign; /**< sign of the data i.e. 1 if photon is positive, -1 if negative */
quadrant quad; /**< quadrant where the photon is located */
double tot; /**< sum of the 3x3 cluster */
double quadTot; /**< sum of the maximum 2x2cluster */
};