mirror of
https://github.com/slsdetectorgroup/slsDetectorPackage.git
synced 2025-04-23 06:50:02 +02:00
included files for data compression, saving only hits
git-svn-id: file:///afs/psi.ch/project/sls_det_software/svn/slsDetectorSoftware@671 951219d9-93cf-4727-9268-0efd64621fa3
This commit is contained in:
l_maliakal_d
parent
a86cf645ed
commit
4811001680
@ -11,7 +11,7 @@ INCLUDES?= -IcommonFiles -IslsDetector -IMySocketTCP -IusersFunctions -ImultiSls
|
||||
CC=g++
|
||||
|
||||
|
||||
SRC_CLNT= slsDetectorAnalysis/fileIO.cpp MySocketTCP/MySocketTCP.cpp usersFunctions/usersFunctions.cpp slsDetector/slsDetectorUtils.cpp slsDetector/slsDetectorCommand.cpp slsDetectorAnalysis/angularConversion.cpp slsDetectorAnalysis/angularConversionStatic.cpp slsDetectorAnalysis/energyConversion.cpp slsDetector/slsDetectorActions.cpp slsDetectorAnalysis/postProcessing.cpp slsDetector/slsDetector.cpp multiSlsDetector/multiSlsDetector.cpp slsDetector/slsDetectorUsers.cpp slsDetectorAnalysis/postProcessingFuncs.cpp slsReceiverInterface/receiverInterface.cpp slsReceiver/slsReceiverFunctionList.cpp slsReceiver/slsReceiver_funcs.cpp slsReceiver/slsReceiverUsers.cpp
|
||||
SRC_CLNT= slsDetectorAnalysis/fileIO.cpp MySocketTCP/MySocketTCP.cpp usersFunctions/usersFunctions.cpp slsDetector/slsDetectorUtils.cpp slsDetector/slsDetectorCommand.cpp slsDetectorAnalysis/angularConversion.cpp slsDetectorAnalysis/angularConversionStatic.cpp slsDetectorAnalysis/energyConversion.cpp slsDetector/slsDetectorActions.cpp slsDetectorAnalysis/postProcessing.cpp slsDetector/slsDetector.cpp multiSlsDetector/multiSlsDetector.cpp slsDetector/slsDetectorUsers.cpp slsDetectorAnalysis/postProcessingFuncs.cpp slsReceiverInterface/receiverInterface.cpp slsReceiver/slsReceiverFunctionList.cpp slsReceiver/slsReceiver_funcs.cpp slsReceiver/slsReceiverUsers.cpp slsDetectorAnalysis/singlePhotonFilter.cpp
|
||||
|
||||
OBJS = $(SRC_CLNT:.cpp=.o)
|
||||
|
||||
|
||||
140
slsDetectorSoftware/slsDetectorAnalysis/movingStat.h
Normal file
140
slsDetectorSoftware/slsDetectorAnalysis/movingStat.h
Normal file
@ -0,0 +1,140 @@
|
||||
/********************************************//**
|
||||
* @file movingStat.h
|
||||
* @short handles pedestal data and moves accordingly
|
||||
***********************************************/
|
||||
#ifndef MOVINGSTAT_H
|
||||
#define MOVINGSTAT_H
|
||||
|
||||
#include "sls_detector_defs.h"
|
||||
|
||||
#include <iostream>
|
||||
#include <fstream>
|
||||
#include <iomanip>
|
||||
#include <cstring>
|
||||
#include <string>
|
||||
#include <sstream>
|
||||
#include <queue>
|
||||
#include <math.h>
|
||||
|
||||
using namespace std;
|
||||
/**
|
||||
@short class handling pedestal data and moves according to data
|
||||
*/
|
||||
|
||||
class movingStat{
|
||||
|
||||
public:
|
||||
/**
|
||||
* Constructor
|
||||
*/
|
||||
movingStat() : n(0), m_n(0), m_oldM(0), m_newM(0), m_oldM2(0), m_newM2(0){}
|
||||
|
||||
/**
|
||||
* Clear number of data values
|
||||
*/
|
||||
void Clear(){
|
||||
m_n = 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Set Pedestal
|
||||
*/
|
||||
void SetN(int i) {n=i;};
|
||||
|
||||
/**
|
||||
* Get Pedestal
|
||||
*/
|
||||
int GetN() {return n;};
|
||||
|
||||
/**
|
||||
* Calculate Pedestal
|
||||
*/
|
||||
void Calc(double x) {
|
||||
if (m_n<n) Add(x);
|
||||
else Push(x);
|
||||
}
|
||||
|
||||
/**
|
||||
* Adding Pedestal
|
||||
*/
|
||||
void Add(double x){
|
||||
m_n++;
|
||||
|
||||
// See Knuth TAOCP vol 2, 3rd edition, page 232
|
||||
if (m_n == 1){
|
||||
m_oldM = m_newM = x;
|
||||
m_oldM2 = x*x;
|
||||
}else{
|
||||
m_newM = m_oldM + x;
|
||||
m_newM2 = m_oldM2 + x*x;
|
||||
//m_newM2 = m_oldM2 + (x*x - m_oldM*m_oldM)/m_n;
|
||||
// set up for next iteration
|
||||
m_oldM = m_newM;
|
||||
m_oldM2 = m_newM2;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Push Pedestal
|
||||
*/
|
||||
void Push(double x){
|
||||
if (m_n == 1){
|
||||
m_oldM = m_newM = x;
|
||||
m_oldM2 = x*x;
|
||||
}else{
|
||||
m_newM = m_oldM + (x - m_oldM/m_n);
|
||||
m_newM2 = m_oldM2 + (x*x - m_oldM2/m_n);
|
||||
//m_newM2 = m_oldM2 + (x*x/m_n - m_oldM*m_oldM/(m_n*m_n));
|
||||
// set up for next iteration
|
||||
m_oldM = m_newM;
|
||||
m_oldM2 = m_newM2;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Get number of data values
|
||||
*/
|
||||
int NumDataValues() const{
|
||||
return m_n;
|
||||
}
|
||||
|
||||
/**
|
||||
* Get mean
|
||||
*/
|
||||
double Mean() const{
|
||||
return (m_n > 0) ? m_newM/m_n : 0.0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Get mean 2
|
||||
*/
|
||||
double M2() const{
|
||||
return ( (m_n > 1) ? m_newM2/m_n : 0.0 );
|
||||
}
|
||||
|
||||
/**
|
||||
* Get variance
|
||||
*/
|
||||
double Variance() const{
|
||||
return ( (m_n > 1) ? (M2()-Mean()*Mean()) : 0.0 );
|
||||
}
|
||||
|
||||
/**
|
||||
* Get standard deviation
|
||||
*/
|
||||
double StandardDeviation() const{
|
||||
return ( (Variance() > 0) ? sqrt( Variance() ) : -1 );
|
||||
}
|
||||
|
||||
private:
|
||||
/** pedestal */
|
||||
int n;
|
||||
/** number of data values */
|
||||
int m_n;
|
||||
|
||||
/** old and new mean */
|
||||
double m_oldM, m_newM, m_oldM2, m_newM2;
|
||||
};
|
||||
|
||||
|
||||
#endif
|
||||
89
slsDetectorSoftware/slsDetectorAnalysis/runningStat.h
Normal file
89
slsDetectorSoftware/slsDetectorAnalysis/runningStat.h
Normal file
@ -0,0 +1,89 @@
|
||||
/********************************************//**
|
||||
* @file runningStat.h
|
||||
* @short handles pedestal data and doesnt move
|
||||
***********************************************/
|
||||
#ifndef RUNNINGSTAT_H
|
||||
#define RUNNINGSTAT_H
|
||||
|
||||
#include <math.h>
|
||||
|
||||
#include "sls_detector_defs.h"
|
||||
|
||||
|
||||
/**
|
||||
@short class handling pedestal data that is static
|
||||
*/
|
||||
|
||||
class runningStat{
|
||||
|
||||
public:
|
||||
/**
|
||||
* Constructor
|
||||
*/
|
||||
runningStat() : m_n(0),m_oldM(0),m_newM(0),m_oldS(0),m_newS(0) {}
|
||||
|
||||
/**
|
||||
* Clear number of data values
|
||||
*/
|
||||
void Clear(){
|
||||
m_n = 0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Push Pedestal
|
||||
*/
|
||||
void Push(double x){
|
||||
m_n++;
|
||||
|
||||
// See Knuth TAOCP vol 2, 3rd edition, page 232
|
||||
if (m_n == 1){
|
||||
m_oldM = m_newM = x;
|
||||
m_oldS = 0.0;
|
||||
}else{
|
||||
m_newM = m_oldM + (x - m_oldM)/m_n;
|
||||
m_newS = m_oldS + (x - m_oldM)*(x - m_newM);
|
||||
|
||||
// set up for next iteration
|
||||
m_oldM = m_newM;
|
||||
m_oldS = m_newS;
|
||||
}
|
||||
}
|
||||
|
||||
/**
|
||||
* Get number of data values
|
||||
*/
|
||||
int NumDataValues() const{
|
||||
return m_n;
|
||||
}
|
||||
|
||||
/**
|
||||
* Get mean
|
||||
*/
|
||||
double Mean() const{
|
||||
return (m_n > 0) ? m_newM : 0.0;
|
||||
}
|
||||
|
||||
/**
|
||||
* Get variance
|
||||
*/
|
||||
double Variance() const{
|
||||
return ( (m_n > 1) ? m_newS/(m_n - 1) : 0.0 );
|
||||
}
|
||||
|
||||
/**
|
||||
* Get standard deviation
|
||||
*/
|
||||
double StandardDeviation() const{
|
||||
return sqrt( Variance() );
|
||||
}
|
||||
|
||||
private:
|
||||
/** number of data values */
|
||||
int m_n;
|
||||
|
||||
/** old and new mean */
|
||||
double m_oldM, m_newM, m_oldS, m_newS;
|
||||
|
||||
};
|
||||
|
||||
#endif
|
||||
544
slsDetectorSoftware/slsDetectorAnalysis/singlePhotonFilter.cpp
Normal file
544
slsDetectorSoftware/slsDetectorAnalysis/singlePhotonFilter.cpp
Normal file
@ -0,0 +1,544 @@
|
||||
/********************************************//**
|
||||
* @file singlePhotonFilter.cpp
|
||||
* @short single photon filter using trees
|
||||
***********************************************/
|
||||
#include "singlePhotonFilter.h"
|
||||
|
||||
|
||||
/*#define MOENCH_FRAME_INDEX_MASK 0xFFFFFF00
|
||||
#define MOENCH_PACKET_INDEX_MASK 0xFF
|
||||
#define MOENCH_FRAME_INDEX_OFFSET 8
|
||||
*/
|
||||
|
||||
singlePhotonFilter::singlePhotonFilter(int x, int y, vector <vector<int16_t> >m, vector <vector<int16_t> >s, int d):
|
||||
#ifdef MYROOT
|
||||
myTree(NULL),
|
||||
myFile(NULL),
|
||||
#else
|
||||
myFile(NULL),
|
||||
#endif
|
||||
nChannelsX(x),
|
||||
nChannelsY(y),
|
||||
nClusterX(CLUSTER_SIZE),
|
||||
nClusterY(CLUSTER_SIZE),
|
||||
map(m),
|
||||
dataSize(d),
|
||||
mask(s),
|
||||
nped(DEFAULT_NUM_PEDESTAL),
|
||||
nsigma(DEFAULT_SIGMA),
|
||||
nbackground(DEFAULT_BACKGROUND),
|
||||
outcorr(DEFAULT_CORRECTION),
|
||||
fnum(0),
|
||||
pnum(0),
|
||||
ptot(0),
|
||||
f0(0),
|
||||
frame_index_mask(0x00000000),
|
||||
packet_index_mask(0x00000000),
|
||||
frame_index_offset(0),
|
||||
packet_index_offset(0),
|
||||
packets_per_frame(0){
|
||||
|
||||
|
||||
if (x == 1)
|
||||
nClusterX = 1;
|
||||
|
||||
stat.resize(x);
|
||||
for(int i=0; i<x; i++)
|
||||
stat[i].resize(y);
|
||||
|
||||
//struct
|
||||
myPhotonHit.data.resize(nClusterX);
|
||||
for(int i=0; i<nClusterX; i++)
|
||||
myPhotonHit.data[i].resize(nClusterY);
|
||||
myPhotonHit.x = 0;
|
||||
myPhotonHit.y = 0;
|
||||
myPhotonHit.rms = 0;
|
||||
myPhotonHit.ped = 0;
|
||||
myPhotonHit.iframe = -1;
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
void singlePhotonFilter::initTree(char *outdir, char *fname){
|
||||
//path
|
||||
char outfname[MAX_STR_LENGTH];
|
||||
sprintf(outfname, "%s/%s.root", outdir, fname);
|
||||
|
||||
#ifdef MYROOT
|
||||
char c1[10],c2[10],cdata[100];
|
||||
sprintf(c1,"%d",nClusterX);
|
||||
sprintf(c2,"%d",nClusterY);
|
||||
sprintf(cdata,"data[%s][%s]/D",c1,c2);
|
||||
|
||||
//file
|
||||
myFile = new TFile(outfname, "RECREATE"); /** later return error if it exists */
|
||||
//tree
|
||||
myTree = new TTree(fname, fname);
|
||||
myTree->Branch("x",&myPhotonHit.x,"x/I");
|
||||
myTree->Branch("y",&myPhotonHit.y,"y/I");
|
||||
myTree->Branch("data",myPhotonHit.data,cdata);
|
||||
myTree->Branch("pedestal",&myPhotonHit.ped,"pedestal/D");
|
||||
myTree->Branch("rms",&myPhotonHit.rms,"rms/D");
|
||||
#else
|
||||
myFile = fopen(outfname, "w");
|
||||
#endif
|
||||
//initialize
|
||||
for (int ir=0; ir<nChannelsX; ir++){
|
||||
for (int ic=0; ic<nChannelsY; ic++){
|
||||
stat[ir][ic].Clear();
|
||||
stat[ir][ic].SetN(nbackground);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
||||
int singlePhotonFilter::writeToFile(){
|
||||
#ifdef MYROOT
|
||||
if((myTree) && (myFile)){
|
||||
myTree->Write();
|
||||
|
||||
myFile = myTree->GetCurrentFile();
|
||||
myFile->Close();
|
||||
delete myFile;
|
||||
}
|
||||
#else
|
||||
if(myFile){ /*&& (number of structs?)*/
|
||||
fwrite((void*)(&myPhotonHit), 1, sizeof(myPhotonHit), myFile);
|
||||
fclose(myFile);
|
||||
delete myFile;
|
||||
return OK;
|
||||
}
|
||||
#endif
|
||||
return FAIL;
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
int singlePhotonFilter::findHits(int16_t *myData, int myDataSize){
|
||||
|
||||
int nHits = 0;
|
||||
int hits[nChannelsX][nChannelsY];
|
||||
|
||||
|
||||
int ir,ic; // for indexing row, column
|
||||
int dum;
|
||||
double tot; // total value of pixel
|
||||
|
||||
//initialize to 0
|
||||
for (ir=0; ir<nChannelsX; ir++)
|
||||
for (ic=0; ic<nChannelsY;ic++)
|
||||
hits[ir][ic]=0;
|
||||
|
||||
|
||||
//for each pixel
|
||||
for (ir=0; ir<160; ir++){
|
||||
for (ic=0; ic<160;ic++){
|
||||
|
||||
//validate mapping
|
||||
if ((map[ir][ic] < 0) || (map[ir][ic] >= (myDataSize))){
|
||||
cout << "Bad Channel Mapping index: " << map[ir][ic] << endl;
|
||||
continue;
|
||||
}
|
||||
|
||||
|
||||
//if frame within pedestal number
|
||||
if (myPhotonHit.iframe < nped)
|
||||
stat[ir][ic].Calc((double)(mask[ir][ic] ^ myData[map[ir][ic]]));
|
||||
|
||||
// frame outside pedestal number
|
||||
else{
|
||||
|
||||
//if hit wasnt registered
|
||||
if (hits[ir][ic] == 0){
|
||||
|
||||
myPhotonHit.rms = stat[ir][ic].StandardDeviation();
|
||||
myPhotonHit.ped = stat[ir][ic].Mean();
|
||||
|
||||
dum = 1;
|
||||
tot = 0;
|
||||
|
||||
//for 1d and 2d
|
||||
int c = 1;
|
||||
int d = 1;
|
||||
if (nChannelsX == 1)
|
||||
c = 0;
|
||||
|
||||
//center pixel
|
||||
myPhotonHit.data[c][d] = ((double)(mask[ir][ic] ^ myData[map[ir][ic]])) - myPhotonHit.ped;
|
||||
|
||||
//check if neighbours are involved
|
||||
for (int ih = -c; ih <= c ; ih++){
|
||||
for (int iv = -d; iv <= d; iv++){
|
||||
//validate neighbouring pixels (not really needed)
|
||||
if (((ir+ih) >= 0) && ((ir+ih) < nChannelsX) && ((ic+iv) >= 0) && ((ic+iv) < nChannelsY)) {
|
||||
//validate mapping
|
||||
if ((map[ir+ih][ic+iv] < 0) || (map[ir+ih][ic+iv] >= myDataSize)){
|
||||
cout << "Bad Channel Mapping index: " << map[ir][ic] << endl;
|
||||
continue;
|
||||
}
|
||||
|
||||
myPhotonHit.data[iv+d][ih+c] = (double)(mask[ir][ic] ^ myData[map[ir+ih][ic+iv]])-stat[ir+ih][ic+iv].Mean();
|
||||
tot += myPhotonHit.data[iv+c][ih+d];
|
||||
if (myPhotonHit.data[iv+c][ih+d] > myPhotonHit.data[c][d])
|
||||
dum = 2;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
if (tot < CLUSTER_SIZE * nsigma * myPhotonHit.rms)
|
||||
dum = 0;
|
||||
|
||||
if (myPhotonHit.data[c][d] < nsigma * myPhotonHit.rms && dum != 0)
|
||||
dum = 3;
|
||||
|
||||
if (myPhotonHit.data[c][d] > -nsigma * myPhotonHit.rms &&
|
||||
myPhotonHit.data[c][d] < nsigma * myPhotonHit.rms &&
|
||||
dum == 0){
|
||||
//Appriximated running average
|
||||
stat[ir][ic].Calc((double)(mask[ir][ic]^myData[map[ir][ic]]));
|
||||
}
|
||||
else if (dum == 1){
|
||||
myPhotonHit.x = ic;
|
||||
myPhotonHit.y = ir;
|
||||
#ifdef MYROOT
|
||||
myTree->Fill();
|
||||
#endif
|
||||
hits[ir][ic] = 1;
|
||||
nHits++;
|
||||
}
|
||||
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (myPhotonHit.iframe%1000 == 0)
|
||||
cout << "Frame: " << myPhotonHit.iframe << " Hits: " << nHits << endl;
|
||||
|
||||
return nHits;
|
||||
}
|
||||
|
||||
|
||||
/*
|
||||
//inside the receiver
|
||||
//pop fifo
|
||||
|
||||
int16_t myData[160*160];
|
||||
if (singlePhotonFilter::verifyFrame(wbuf,rc, myData,firstTimeHere))
|
||||
if(firstTimeHere){
|
||||
firsttimeHere = 0;
|
||||
//f0 = fnum;
|
||||
}
|
||||
//singlePhotonFilter:SetFrameNumber(fnum - f0); singlePhotonFilter:SetFrameNumber(int i){myPhotonHit.iframe = i;};
|
||||
singlePhotonFilter.findhits(myData,1286*40);
|
||||
*/
|
||||
|
||||
|
||||
void singlePhotonFilter::initialize(int fmask, int pmask, int foffset, int poffset, int pperf){
|
||||
fnum = 0; pnum = 0; ptot = 0; f0 = 0;
|
||||
frame_index_mask = fmask;
|
||||
packet_index_mask = pmask;
|
||||
frame_index_offset = foffset;
|
||||
packet_index_offset = poffset;
|
||||
packets_per_frame = pperf;
|
||||
}
|
||||
|
||||
|
||||
|
||||
int singlePhotonFilter::verifyFrame(char *inData, int inDataSize, int16_t* myData, int firstTime){
|
||||
int offset = 0;
|
||||
int fIndex = 0;
|
||||
int pIndex = 0;
|
||||
|
||||
|
||||
while(offset < inDataSize){
|
||||
fIndex = ((((int)(*((int*)inData))) & (frame_index_mask)) >> frame_index_offset);
|
||||
pIndex = ((((int)(*((int*)inData))) & (packet_index_mask)) >> packet_index_offset);
|
||||
|
||||
offset += 4;
|
||||
|
||||
//check validity of packet index
|
||||
if ((pIndex < 0) && (pIndex >= packets_per_frame)){
|
||||
cout << "cannot decode packet index:" << pIndex << endl;
|
||||
//its already dealt with cuz this frame will be discarded in the end
|
||||
}
|
||||
|
||||
//put first packet last
|
||||
if (pIndex == 0)
|
||||
pIndex = packets_per_frame;
|
||||
|
||||
|
||||
//copy packet to correct spot in outData
|
||||
memcpy(((char*)(myData+(pIndex-1)*640)), (inData + offset), 1280);
|
||||
offset += 1280 + 2;
|
||||
|
||||
//firsttime
|
||||
if (firstTime){
|
||||
firstTime = 0;
|
||||
f0 = fIndex;
|
||||
fnum = fIndex;
|
||||
pnum = pIndex - 1; //should b 0 at first
|
||||
ptot = 0;
|
||||
}
|
||||
|
||||
//if it is not matching withthe frame number
|
||||
if (fIndex != fnum){
|
||||
cout << "Missing Packet! " << fnum << " "
|
||||
"Expected f" << fnum << " p " << pnum + 1 << " received f " << fnum << " p " << pIndex << endl;
|
||||
fnum = fIndex;
|
||||
pnum = pIndex - 1;
|
||||
ptot = 0;
|
||||
}
|
||||
|
||||
//if missing a packet, discard
|
||||
if (pIndex != pnum + 1){
|
||||
cout << "Missing Packet! " << fnum << " "
|
||||
"Expected f" << fnum << " p " << pnum + 1 << " received f " << fnum << " p " << pIndex << endl;
|
||||
pnum = pIndex;
|
||||
ptot = 0;/** should it be 0 or just not incrememnted */
|
||||
}
|
||||
//no missing packet
|
||||
else{
|
||||
pnum++;
|
||||
ptot++;
|
||||
}
|
||||
|
||||
|
||||
//if its the last index
|
||||
if (pIndex == packets_per_frame){
|
||||
//got all packets
|
||||
if (ptot == packets_per_frame){
|
||||
myPhotonHit.iframe = fnum - f0;/** ?? */
|
||||
findHits(myData,nChannelsX*nChannelsY*2);
|
||||
return 1;
|
||||
}else{
|
||||
cout << "* Some packets have been missed!" << fnum << endl;
|
||||
ptot = 0;
|
||||
pnum = 0;
|
||||
fnum = fIndex + 1;
|
||||
}
|
||||
}
|
||||
|
||||
//index not 40, but total is 40.. strange
|
||||
else if (ptot == packets_per_frame){
|
||||
cout << "** Some packets have been missed! " << fnum << endl;
|
||||
ptot = 0;
|
||||
pnum = pIndex - 1;
|
||||
fnum = fIndex;
|
||||
}
|
||||
}
|
||||
|
||||
delete myData;
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
/*
|
||||
void singlePhotonFilter::makeRunTree(char *fformat, int runmin, int runmax, char *tname, int off) {
|
||||
|
||||
char buff[1286];
|
||||
int *header=(int*)buff;
|
||||
char *footer=buff+1284;
|
||||
Short_t mydata[160*160];//=(Short_t*)(buff+4);
|
||||
int pnum;
|
||||
int fnum;
|
||||
int ifr;
|
||||
int ipa;
|
||||
Double_t dd[3][3];
|
||||
int ncol=colwidth*4;//=160;
|
||||
int nch=nrow*ncol;
|
||||
char fname[1000];
|
||||
Int_t chan[160][160];
|
||||
Int_t nHits=0;
|
||||
char vtype[100];
|
||||
int afifo_length;
|
||||
ifstream filebin;
|
||||
int ptot=0;
|
||||
|
||||
|
||||
for (int irun=runmin; irun<=runmax; irun++) {
|
||||
sprintf(fname,fformat,irun);
|
||||
cout << "file name " << fname << " " << tall->GetEntries() << endl;
|
||||
filebin.open((const char *)(fname), ios::in | ios::binary);
|
||||
if (filebin.is_open()) {
|
||||
while (readNextFrame(filebin, mydata, fnum)) {
|
||||
//to get rid of the first frame number
|
||||
if (iii==0) {
|
||||
f0=fnum;
|
||||
iii = 1;
|
||||
}
|
||||
iframe=fnum-f0;
|
||||
findHits(mydata, tall);
|
||||
}
|
||||
if (filebin.is_open())
|
||||
filebin.close();
|
||||
} else
|
||||
cout << "could not open file " << fname << endl;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
||||
|
||||
int readNextFrame(ifstream &filebin, int16_t *mydata, int &ifr){
|
||||
|
||||
char buff[1286];
|
||||
int *header=(int*)buff;
|
||||
char *footer=buff+1284;
|
||||
// Short_t mydata[160*160];//=(Short_t*)(buff+4);
|
||||
int pnum;
|
||||
int fnum;
|
||||
// int ifr;
|
||||
int ipa;
|
||||
int iii=0;
|
||||
int nrow=160;
|
||||
int colwidth=40;
|
||||
// int nsigma=5;
|
||||
int ncol=colwidth*4;//=160;
|
||||
int nch=nrow*ncol;
|
||||
int f0;
|
||||
char fname[1000];
|
||||
int ptot=0;
|
||||
if (filebin.is_open()) {
|
||||
while (filebin.read(buff,sizeof(4))) {
|
||||
// while (filebin.read(buff,sizeof(buff))) {
|
||||
pnum=(*header)&0xff;//packet num
|
||||
fnum=((*header)&(0xffffff00))>>8;//frame num
|
||||
if (pnum==0)
|
||||
pnum=40;
|
||||
filebin.read((char*)(mydata+(pnum-1)*640),1280);
|
||||
filebin.read(footer,2);
|
||||
if (iii==0) {
|
||||
// cout << "fnum " << fnum << " pnum " << pnum << coff << endl;
|
||||
ifr=fnum;
|
||||
ipa=pnum-1;
|
||||
f0=fnum;//not used
|
||||
iii=1;
|
||||
ptot=0;
|
||||
}
|
||||
if (fnum!=ifr) {
|
||||
cout << "Missing packet! "<< ifr << " ";
|
||||
cout << "Expected f " <<ifr << " p " << ipa+1 << " received f " << fnum << " p " << pnum << endl;
|
||||
ifr=fnum;
|
||||
ipa=pnum-1;
|
||||
ptot=0;
|
||||
}
|
||||
if (pnum!=ipa+1) {
|
||||
cout << "Missing packet! "<< ifr << " ";
|
||||
cout << "Expected f " <<ifr << " p " << ipa+1 << " received f " << fnum << " p " << pnum << endl;
|
||||
ipa=pnum;
|
||||
ptot=0;
|
||||
} else {
|
||||
ipa++;
|
||||
ptot++;
|
||||
}
|
||||
|
||||
if (pnum==40) {
|
||||
if (ptot==40) {
|
||||
return 1;
|
||||
} else {
|
||||
cout << "* Some packets have been missed! " << ifr << endl;
|
||||
ptot=0;
|
||||
ipa=0;
|
||||
ifr=fnum+1;
|
||||
// return 0;
|
||||
}
|
||||
} else if (ptot==40) {
|
||||
cout << "** Some packets have been missed! " << ifr << endl;
|
||||
ptot=0;
|
||||
ipa=pnum-1;
|
||||
ifr=fnum;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
|
||||
|
||||
*/
|
||||
/*
|
||||
int nch=nrow*ncol;
|
||||
|
||||
int f0;
|
||||
|
||||
char fname[1000];
|
||||
|
||||
int ptot=0;
|
||||
|
||||
if (filebin.is_open()) {
|
||||
|
||||
while (filebin.read(buff,sizeof(4))) {
|
||||
// while (filebin.read(buff,sizeof(buff))) {
|
||||
|
||||
pnum=(*header)&0xff;
|
||||
fnum=((*header)&(0xffffff00))>>8;
|
||||
|
||||
if (pnum==0)
|
||||
pnum=40;
|
||||
|
||||
|
||||
filebin.read((char*)(mydata+(pnum-1)*640),1280);
|
||||
filebin.read(footer,2);
|
||||
|
||||
if (iii==0) {
|
||||
// cout << "fnum " << fnum << " pnum " << pnum << coff << endl;
|
||||
ifr=fnum;
|
||||
ipa=pnum-1;
|
||||
f0=fnum;
|
||||
iii=1;
|
||||
ptot=0;
|
||||
|
||||
}
|
||||
|
||||
if (fnum!=ifr) {
|
||||
cout << "Missing packet! "<< ifr << " ";
|
||||
cout << "Expected f " <<ifr << " p " << ipa+1 << " received f " << fnum << " p " << pnum << endl;
|
||||
ifr=fnum;
|
||||
ipa=pnum-1;
|
||||
ptot=0;
|
||||
}
|
||||
|
||||
|
||||
if (pnum!=ipa+1) {
|
||||
cout << "Missing packet! "<< ifr << " ";
|
||||
cout << "Expected f " <<ifr << " p " << ipa+1 << " received f " << fnum << " p " << pnum << endl;
|
||||
ipa=pnum;
|
||||
ptot=0;
|
||||
} else {
|
||||
ipa++;
|
||||
ptot++;
|
||||
}
|
||||
|
||||
if (pnum==40) {
|
||||
if (ptot==40) {
|
||||
return 1;
|
||||
} else {
|
||||
cout << "* Some packets have been missed! " << ifr << endl;
|
||||
ptot=0;
|
||||
ipa=0;
|
||||
ifr=fnum+1;
|
||||
// return 0;
|
||||
}
|
||||
} else if (ptot==40) {
|
||||
cout << "** Some packets have been missed! " << ifr << endl;
|
||||
ptot=0;
|
||||
ipa=pnum-1;
|
||||
ifr=fnum;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
return 0;
|
||||
}
|
||||
|
||||
*/
|
||||
|
||||
260
slsDetectorSoftware/slsDetectorAnalysis/singlePhotonFilter.h
Normal file
260
slsDetectorSoftware/slsDetectorAnalysis/singlePhotonFilter.h
Normal file
@ -0,0 +1,260 @@
|
||||
/********************************************//**
|
||||
* @file singlePhotonFilter.h
|
||||
* @short single photon filter using trees
|
||||
***********************************************/
|
||||
#ifndef SINGLEPHOTONFILTER_H
|
||||
#define SINGLEPHOTONFILTER_H
|
||||
|
||||
//#include "sls_detector_defs.h"
|
||||
|
||||
|
||||
#ifdef __CINT
|
||||
#define MYROOT
|
||||
#endif
|
||||
|
||||
|
||||
#ifdef MYROOT
|
||||
#include <TFile.h>
|
||||
#include <TMath.h>
|
||||
#include <TTree.h>
|
||||
#include <TChain.h>
|
||||
#endif
|
||||
|
||||
#include <stdio.h>
|
||||
#include <iostream>
|
||||
#include <deque>
|
||||
#include <list>
|
||||
#include <queue>
|
||||
#include <fstream>
|
||||
|
||||
|
||||
#include "runningStat.h"
|
||||
#include "movingStat.h"
|
||||
using namespace std;
|
||||
|
||||
|
||||
typedef double double32_t;
|
||||
typedef float float32_t;
|
||||
typedef int int32_t;
|
||||
|
||||
|
||||
/**
|
||||
return values
|
||||
*/
|
||||
enum {
|
||||
OK, /**< function succeeded */
|
||||
FAIL /**< function failed */
|
||||
};
|
||||
/**
|
||||
@short structure for a single photon hit
|
||||
*/
|
||||
typedef struct{
|
||||
vector < vector<double> > data; /**< data size */
|
||||
int x; /**< x-coordinate of the center of hit */
|
||||
int y; /**< x-coordinate of the center of hit */
|
||||
double rms; /**< noise of central pixel */
|
||||
double ped; /**< pedestal of the central pixel */
|
||||
int iframe; /**< frame number */
|
||||
}single_photon_hit;
|
||||
|
||||
|
||||
|
||||
|
||||
/**
|
||||
@short class handling trees and its data file
|
||||
*/
|
||||
|
||||
class singlePhotonFilter{
|
||||
public:
|
||||
/**
|
||||
* Constructor
|
||||
* @param nChannelsX Number of Channels in X direction
|
||||
* @param nChannelsY Number of Channels in Y direction
|
||||
* @param m Map to data without headers
|
||||
* @param s mask as to which adcs are inverted
|
||||
* @param d Size of data with the headers
|
||||
*/
|
||||
|
||||
/** why is the datasize -1, you need to know the datasize with headers so that you dont go over the limits */
|
||||
singlePhotonFilter(int x, int y, vector <vector<int16_t> >m, vector <vector<int16_t> >s, int d = -1);
|
||||
/*map[56][63]=656; data[map[56][63]] ^ 0x7fff*/
|
||||
|
||||
/** virtual destructor */
|
||||
virtual ~singlePhotonFilter(){};
|
||||
|
||||
/**
|
||||
* Construct a tree, populate struct for the single photon hit
|
||||
* @param outdir Output file directory
|
||||
* @param fname Output file name
|
||||
*/
|
||||
void initTree(char *outdir, char *fname);
|
||||
|
||||
/**
|
||||
* Reset Indices before starting acquisition and provide all the masks and offsets
|
||||
* @param fmask frame index mask
|
||||
* @param pmask packet index mask
|
||||
* @param foffset frame index offset
|
||||
* @param poffset packet index offset
|
||||
* @param pperf packets per frame
|
||||
*/
|
||||
void initialize(int fmask, int pmask, int foffset, int poffset, int pperf);
|
||||
|
||||
/** Verify if all packets exist for the frame
|
||||
* @param inData the data from socket to be verified
|
||||
* @param inDataSize datasize of packet
|
||||
* @param myData frame with all the packets
|
||||
* @param firstTime the first frame received from socket
|
||||
* */
|
||||
int verifyFrame(char *inData, int inDataSize, int16_t* myData, int firstTime);
|
||||
|
||||
/**
|
||||
* Writes tree/struct to file
|
||||
* returns OK if successful, else FAIL
|
||||
*/
|
||||
int writeToFile();
|
||||
|
||||
/**
|
||||
* Find Hits frame by frame
|
||||
* @param myData data for one frame
|
||||
* @param myDataSize data size for one frame with headers
|
||||
* returns number of hits
|
||||
*/
|
||||
int findHits(int16_t *myData, int myDataSize);
|
||||
|
||||
|
||||
|
||||
#ifdef MYROOT
|
||||
/**
|
||||
* returns tree
|
||||
*/
|
||||
TTree *getTree(){ return myTree; };
|
||||
#endif
|
||||
|
||||
/**
|
||||
* returns struct
|
||||
*/
|
||||
single_photon_hit getStructure(){ return myPhotonHit; };
|
||||
|
||||
/** Set number of frames to calculate pedestal at beginning */
|
||||
void setNPed(int n){ nped = n; };
|
||||
|
||||
/** Get number of frames to calculate pedestal at beginning */
|
||||
int getNPed(){return nped;};
|
||||
|
||||
/** Set Distance from pedestal to detect a hit */
|
||||
void setNSigma(int n){ nsigma = n; };
|
||||
|
||||
/** Get Distance from pedestal to detect a hit */
|
||||
int getNSigma(){return nsigma;};
|
||||
|
||||
/** Set background */
|
||||
void setNBackground(int n){ nbackground = n; };
|
||||
|
||||
/** Get background */
|
||||
int getNBackground(){return nbackground;};
|
||||
|
||||
/** Set correction */
|
||||
void setOutCorr(double d){ outcorr = d; };
|
||||
|
||||
/** Get correction */
|
||||
double getOutCorr(){return outcorr;};
|
||||
|
||||
|
||||
|
||||
|
||||
private:
|
||||
|
||||
#ifdef MYROOT
|
||||
/** Tree where the hits are stored */
|
||||
TTree *myTree;
|
||||
|
||||
/** File where the tree is saved */
|
||||
TFile *myFile;
|
||||
#else
|
||||
FILE *myFile;
|
||||
#endif
|
||||
|
||||
/** Number of Channels in X direction */
|
||||
int nChannelsX;
|
||||
|
||||
/** Number of Channels in Y direction */
|
||||
int nChannelsY;
|
||||
|
||||
/** Cluster size in X direction */
|
||||
int nClusterX;
|
||||
|
||||
/** Cluster size in Y direction */
|
||||
int nClusterY;
|
||||
|
||||
/** map to the data without headers */
|
||||
vector <vector<int16_t> >map;
|
||||
|
||||
/** Size of data with headers */
|
||||
int dataSize;
|
||||
|
||||
/** mask as to which adcs are inverted */
|
||||
vector <vector<int16_t> >mask;
|
||||
|
||||
/** movingStat object */
|
||||
vector <vector<movingStat> > stat;
|
||||
|
||||
/** single Photon Hit structure */
|
||||
single_photon_hit myPhotonHit;
|
||||
|
||||
/** Cluster size */
|
||||
const static int CLUSTER_SIZE = 3;
|
||||
|
||||
/** Default Number of frames at the beginning to calculate pedestal */
|
||||
const static int DEFAULT_NUM_PEDESTAL = 500;
|
||||
|
||||
/** Default Distance from pedestal to detect a hit */
|
||||
const static int DEFAULT_SIGMA = 5;
|
||||
|
||||
/** Default Background */
|
||||
const static int DEFAULT_BACKGROUND = 1000;
|
||||
|
||||
/** Default Correction Percentage */
|
||||
const static double DEFAULT_CORRECTION = 1.0;
|
||||
|
||||
/** Number of frames at the beginning to calculate pedestal */
|
||||
int nped;
|
||||
|
||||
/** Distance from pedestal to detect a hit */
|
||||
int nsigma;
|
||||
|
||||
/** background */
|
||||
int nbackground;
|
||||
|
||||
/** correction */
|
||||
double outcorr;
|
||||
|
||||
/** previous frame index */
|
||||
unsigned int fnum;
|
||||
|
||||
/** previous packet index */
|
||||
unsigned int pnum;
|
||||
|
||||
/** total packets received */
|
||||
unsigned int ptot;
|
||||
|
||||
/** first frame number */
|
||||
unsigned int f0;
|
||||
|
||||
/** frame index mask */
|
||||
int frame_index_mask;
|
||||
|
||||
/** packet index mask */
|
||||
int packet_index_mask;
|
||||
|
||||
/** frame index offset */
|
||||
int frame_index_offset;
|
||||
|
||||
/** packet index offset */
|
||||
int packet_index_offset;
|
||||
|
||||
/** number of packets per frame */
|
||||
int packets_per_frame;
|
||||
|
||||
};
|
||||
|
||||
#endif
|
||||
Loading…
x
Reference in New Issue
Block a user