diff --git a/slsDetectorCalibration/Makefile.moench_eta b/slsDetectorCalibration/Makefile.moench_eta new file mode 100644 index 000000000..5f9b8ed56 --- /dev/null +++ b/slsDetectorCalibration/Makefile.moench_eta @@ -0,0 +1,22 @@ + +CBFLIBDIR= /home/l_msdetect/CBFlib-0.9.5 +ZMQLIB=../slsReceiverSoftware/include +LIBRARYCBF=$(CBFLIBDIR)/lib/*.o +INCDIR=-IslsDetectorCalibration -I../slsReceiverSoftware/include -I$(CBFLIBDIR)/include/ -I. -IetaVEL +LIBHDF5= +LDFLAG= -L/usr/lib64/ -lpthread -lm -lstdc++ -L. -lzmq -pthread -lrt -L$(CBFLIBDIR)/lib/ -lhdf5 -ltiff -L$(ZMQLIB) +#-L../../bin +MAIN=moench03ZmqOnTheFlyEta.cpp + +#DESTDIR?=../bin + +all: moench03ZmqOnTheFlyEta + + + +moench03ZmqOnTheFlyEta: $(MAIN) $(INCS) clean + g++ -o moench03ZmqOnTheFlyEta $(MAIN) $(LDFLAG) $(INCDIR) $(LIBHDF5) $(LIBRARYCBF) tiffIO.cpp -DSAVE_ALL + +clean: + rm -f moench03ZmqOnTheFlyEta + diff --git a/slsDetectorCalibration/Makefile.moench_zmq b/slsDetectorCalibration/Makefile.moench_zmq new file mode 100644 index 000000000..f2ed06a32 --- /dev/null +++ b/slsDetectorCalibration/Makefile.moench_zmq @@ -0,0 +1,22 @@ + +CBFLIBDIR=/afs/psi.ch/project/sls_det_software/CBFlib-0.9.5 +ZMQLIB=../slsReceiverSoftware/include +LIBRARYCBF=$(CBFLIBDIR)/lib/*.o +INCDIR=-IslsDetectorCalibration -I../slsReceiverSoftware/include -I$(CBFLIBDIR)/include/ -I. -IetaVEL +LIBHDF5= +LDFLAG= -L/usr/lib64/ -lpthread -lm -lstdc++ -L. -lzmq -pthread -lrt -L$(CBFLIBDIR)/lib/ -lhdf5 -ltiff -L$(ZMQLIB) +#-L../../bin +MAIN=moench03ZmqOnTheFly.cpp + +#DESTDIR?=../bin + +all: moench03ZmqOnTheFly + + + +moench03ZmqOnTheFly: $(MAIN) $(INCS) clean + g++ -o moench03ZmqOnTheFly $(MAIN) $(LDFLAG) $(INCDIR) $(LIBHDF5) $(LIBRARYCBF) tiffIO.cpp -DSAVE_ALL + +clean: + rm -f moench03ZmqOnTheFlyEta + diff --git a/slsDetectorCalibration/analogDetector.h b/slsDetectorCalibration/analogDetector.h index 3fb9e7f5f..9d2fad96b 100644 --- a/slsDetectorCalibration/analogDetector.h +++ b/slsDetectorCalibration/analogDetector.h @@ -6,8 +6,12 @@ #include "pedestalSubtraction.h" #include "commonModeSubtraction.h" #include "tiffIO.h" + #ifndef FRAMEMODE_DEF #define FRAMEMODE_DEF +/** +enum to define the flags of the data set, which are needed to seect the type of processing it should undergo +*/ enum frameMode { eFrame, ePedestal, eFlat }; #endif @@ -23,20 +27,20 @@ template class analogDetector { /** Constructor (no error checking if datasize and offsets are compatible!) - \param d detector data structure to be used - \param csize cluster size (should be an odd number). Defaults to 3 - \param nsigma number of rms to discriminate from the noise. Defaults to 5 - \param sign 1 if photons are positive, -1 if negative - \param cm common mode subtraction algorithm, if any. Defaults to NULL i.e. none + \param d detector data structure to be used - if null it is assumed that the data are in ordered ip=iy*nx+ix + \param sign is the sign of the data \param nped number of samples for pedestal averaging - \param nd number of dark frames to average as pedestals without photon discrimination at the beginning of the measurement + \param cm common mode subtraction algorithm, if any. Defaults to NULL i.e. none + \param nnx detector size in x - must be specified if no data structure is defined, otherwise defaults to the size of the data structure. + \param nny detector size in y - must be specified if no data structure is defined, otherwise defaults to the size of the data structure. + \param gm pointer to tha gain map matrix */ analogDetector(slsDetectorData *d, int sign=1, - commonModeSubtraction *cm=NULL, int nped=1000, int nnx=-1, int nny=-1, double *gm=NULL) : det(d), nx(nnx), ny(nny), stat(NULL), cmSub(cm), iframe(-1), dataSign(sign), gmap(gm) { + commonModeSubtraction *cm=NULL, int nped=1000, int nnx=-1, int nny=-1, double *gm=NULL) : det(d), nx(nnx), ny(nny), stat(NULL), cmSub(cm), iframe(-1), dataSign(sign), gmap(gm), id(0) { if (det) det->getDetectorSize(nx,ny); @@ -48,16 +52,17 @@ template class analogDetector { stat[i][ix].SetNPedestals(nped); } } - + image=new int[nx*ny]; }; /** destructor. Deletes the cluster structure and the pdestalSubtraction array */ - virtual ~analogDetector() {for (int i=0; i class analogDetector { iframe=orig->iframe; gmap=orig->gmap; cmSub=orig->cmSub; - + id=orig->id; + stat=new pedestalSubtraction*[ny]; for (int i=0; i class analogDetector { setPedestal(ix,iy,orig->getPedestal(ix,iy),orig->getPedestalRMS(ix,iy),orig->GetNPedestals(ix,iy)); } } + image=new int[nx*ny]; + } + /** + clone. Must be virtual! + */ virtual analogDetector *Clone() { return new analogDetector(this); } + /** + Gives an id to the structure. For debugging purposes in case of multithreading. + \param i is to be set + \returns current id + */ + int setId(int i){id=i; return id;}; + + /** + Returns id of the structure. For debugging purposes in case of multithreading. + \returns current id + */ + int getId() {return id; }; + /** + Returns data size of the detector data structure + \returns data size of the detector data structurein bytes + */ int getDataSize(){return det->getDataSize();}; + /** + Returns data size of the detector image matrix + \param nnx reference to image size in x + \param nny reference to image size in y + \param nns reference to number of subpixels for interpolating detector, will always be 1 in this case + \returns number of pixels of the detector image + */ + virtual int getImageSize(int &nnx, int &nny, int &nns){nnx=nx; nny=ny; nns=1; return nx*ny;}; + /** + Returns data size of the detector image matrix + \param nnx reference to pixel size in x + \param nny reference to pixel size in y + \returns number of pixels of the detector image + */ + virtual int getDetectorSize(int &nnx, int &nny){nnx=nx; nny=ny; return nx*ny;}; /** set gain map + \param gm pointer to gain map matrix to be set - NULL unsets + \returns pointer to current gain map */ double *setGainMap(double *gm) {gmap=gm; return gmap;}; /** return gain map + \returns pointer to current gain map */ double *getGainMap() {return gmap;}; - + /** + reads a 32 bit tiff file of the size of the detector and sets its values as gain map for the detector. If file does not exist returns NULL, but does not change gainmap compared to previous settings. + \param imgname complete name of the file containing the gain map data + \returns pointer to current gain map is file reading succeeded, NULL is file reading didn't work. + */ double *readGainMap(const char * imgname) { uint32 nnx, nny; float *gm=ReadFromTiff( imgname, nny, nnx); @@ -115,7 +164,11 @@ template class analogDetector { } return NULL; } - + /** + wrties a 32 bit tiff file of the size of the detector and contaning the gain map value, if any. If file doesn'e exist or gainmap is undefined, does not do anything. + \param imgname complete name of the file to be written + \returns NULL + */ void *writeGainMap(const char * imgname) { float *gm=NULL; if (gmap) { @@ -125,7 +178,8 @@ template class analogDetector { gm[iy*nx+ix]=gmap[iy*nx+ix]; } } - return WriteToTiff(gm, imgname, ny, nx); + WriteToTiff(gm, imgname, ny, nx); + delete [] gm; } return NULL; } @@ -175,6 +229,7 @@ template class analogDetector { \param ix pixel x coordinate \param iy pixel y coordinate \param cm 0 (default) without common mode subtraction, 1 with common mode subtraction (if defined) + \returns pedestal value */ virtual double getPedestal(int ix, int iy, int cm=0){if (ix>=0 && ix=0 && iy0) return stat[iy][ix].getPedestal()-cmSub->getCommonMode(); else return stat[iy][ix].getPedestal(); else return -1;}; @@ -182,8 +237,57 @@ template class analogDetector { gets pedestal rms (i.e. noise) \param ix pixel x coordinate \param iy pixel y coordinate + \returns pedestal rms */ - double getPedestalRMS(int ix, int iy){if (ix>=0 && ix=0 && iy=0 && ix=0 && iy class analogDetector { \param ix pixel x coordinate \param iy pixel y coordinate \param val value to set + \param rms rms to be set if any, defaults to 0 + \param m number of pedestal samples to be set or the moving stat structure is any, defaults to 0 */ virtual void setPedestal(int ix, int iy, double val, double rms=0, int m=-1){if (ix>=0 && ix=0 && iy=0 && ix=0 && iy class analogDetector { return NULL; } + int readImage(const char * imgname) { + uint32 nnx, nny; + float *gm=ReadFromTiff( imgname, nny, nnx); + if (nnx>nx) nnx=nx; + if (nny>ny) nny=ny; + + if (gm) { + for (int ix=0; ix class analogDetector { virtual int getNPhotons(char *data, int ix, int iy=0, int thr=-1) { double g=1.; + int v; if (ix>=0 && ix=0 && iy class analogDetector { if (thr<=0) thr=-1*thr*getPedestalRMS(ix,iy)/g; - if (det) - return (dataSign*det->getValue(data, ix, iy)-getPedestal(ix,iy))/g/thr; - else - return (((double*)data)[(iy)*nx+ix]-getPedestal(ix,iy))/g/thr; - + /* if (det) */ + /* return (dataSign*det->getValue(data, ix, iy)-getPedestal(ix,iy))/g/thr; */ + /* else { */ + //if (ix==30 && iy==50 ) + v=subtractPedestal(data,ix,iy)/g/thr; + if (v>0) + return v; + //cout << v << endl; + return 0; + // } } return 0; @@ -376,7 +560,7 @@ template class analogDetector { double val; if (nph==NULL) - nph=new int[nx*ny]; + nph=image;//image=new int[nx*ny]; double tthr=thr; newFrame(); for (int ix=0; ix class analogDetector { } - + virtual void clearImage(){ for (int ix=0; ix class analogDetector { } }; + virtual char *getInterpolation(){return NULL;}; protected: @@ -443,6 +629,9 @@ template class analogDetector { int dataSign; /**< sign of the data i.e. 1 if photon is positive, -1 if negative */ int iframe; /**< frame number (not from file but incremented within the dataset every time newFrame is called */ double *gmap; + int *image; + int id; + //int xmin, xmax, ymin, ymax; }; #endif diff --git a/slsDetectorCalibration/etaVEL/etaInterpolationBase.h b/slsDetectorCalibration/etaVEL/etaInterpolationBase.h index aac4b409c..4a3517151 100644 --- a/slsDetectorCalibration/etaVEL/etaInterpolationBase.h +++ b/slsDetectorCalibration/etaVEL/etaInterpolationBase.h @@ -18,8 +18,8 @@ class etaInterpolationBase : public slsInterpolation { if (nb<=0) nbeta=nSubPixels*10; if (etamin>=etamax) { - etamin=-0.1; - etamax=1.1; + etamin=-1; + etamax=2; } etastep=(etamax-etamin)/nbeta; #ifdef MYROOT1 @@ -79,6 +79,11 @@ class etaInterpolationBase : public slsInterpolation { } return heta; }; + TH2D *setFlatField(TH2D *h, int nb=-1, double emin=1, double emax=0) + { + return setEta(h, nb, emin, emax); + }; + TH2D *getFlatField(){return setEta(NULL);}; #endif @@ -91,17 +96,27 @@ class etaInterpolationBase : public slsInterpolation { etamin=emin; etamax=emax; if (etamin>=etamax) { - etamin=-0.1; - etamax=1.1; + etamin=-1; + etamax=2; } etastep=(etamax-etamin)/nbeta; } return heta; }; + int *setFlatField(int *h, int nb=-1, double emin=1, double emax=0) + { + return setEta(h, nb, emin, emax); + }; int *getFlatField(){return setEta(NULL);}; - + int *getFlatField(int &nb, double &emin, double &emax){ + nb=nbeta; + //cout << "igff* ff has " << nb << " bins " << endl; + emin=etamin; + emax=etamax; + return getFlatField(); + }; void *writeFlatField(const char * imgname) { @@ -122,8 +137,8 @@ class etaInterpolationBase : public slsInterpolation { if (emin<=0) etamin=emin; if (etamin>=etamax) { - etamin=-0.1; - etamax=1.1; + etamin=-1; + etamax=2; } etastep=(etamax-etamin)/nbeta; @@ -207,44 +222,92 @@ int *gethhx() int corner; corner=calcQuad(data, tot, totquad, sDum); - - calcEta(totquad, sDum, etax, etay); - getInterpolatedPosition(x,y,etax,etay,corner,int_x,int_y); + if (nSubPixels>2) + calcEta(totquad, sDum, etax, etay); + getInterpolatedPosition(x,y,etax,etay,corner,int_x,int_y); + return; }; + virtual void getInterpolatedPosition(int x, int y, double totquad,int quad,double *cl,double &int_x, double &int_y) { + + double cc[2][2]; + double *cluster[3]; + int xoff, yoff; + cluster[0]=cl; + cluster[1]=cl+3; + cluster[2]=cl+6; + + switch (quad) { + case BOTTOM_LEFT: + xoff=0; + yoff=0; + break; + case BOTTOM_RIGHT: + xoff=1; + yoff=0; + break; + case TOP_LEFT: + xoff=0; + yoff=1; + break; + case TOP_RIGHT: + xoff=1; + yoff=1; + break; + default: + ; + } + double etax, etay; + if (nSubPixels>2) { + cc[0][0]=cluster[yoff][xoff]; + cc[1][0]=cluster[yoff+1][xoff]; + cc[0][1]=cluster[yoff][xoff+1]; + cc[1][1]=cluster[yoff+1][xoff+1]; + calcEta(totquad,cc,etax,etay); + } + return getInterpolatedPosition(x,y,etax, etay,quad,int_x,int_y); + + } + + + + + virtual void getInterpolatedPosition(int x, int y, double etax, double etay, int corner, double &int_x, double &int_y) { - double xpos_eta,ypos_eta; + double xpos_eta=0,ypos_eta=0; double dX,dY; int ex,ey; switch (corner) { case TOP_LEFT: - dX=-1.; - dY=+1.; + dX=-.99; + dY=+.99; break; case TOP_RIGHT: - dX=+1.; - dY=+1.; + dX=+.99; + dY=+.99; break; case BOTTOM_LEFT: - dX=-1.; - dY=-1.; + dX=-.99; + dY=-.99; break; case BOTTOM_RIGHT: - dX=+1.; - dY=-1.; + dX=+.99; + dY=-.99; break; default: + cout << "bad quadrant" << endl; dX=0.; dY=0.; } + if (nSubPixels>2) { #ifdef MYROOT1 xpos_eta=(hhx->GetBinContent(hhx->GetXaxis()->FindBin(etax),hhy->GetYaxis()->FindBin(etay)))/((double)nSubPixels); @@ -265,9 +328,14 @@ int *gethhx() //return 0; #endif - + } else { + xpos_eta=-dX*0.25; + ypos_eta=-0.25*dY; + } + int_x=((double)x) + 0.5*dX + xpos_eta; int_y=((double)y) + 0.5*dY + ypos_eta; + // cout << "***"<< x <<" " << y << " " << int_x << " " << int_y << endl; // cout << etax << " " << ex << " " << etay << " " << ey << " " << xpos_eta << " " << int_x << " " << ypos_eta << " " << int_y << endl; //return 1; @@ -315,6 +383,56 @@ int *gethhx() return; }; + virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) { + double cc[2][2]; + double *cluster[3]; + int xoff, yoff; + cluster[0]=cl; + cluster[1]=cl+3; + cluster[2]=cl+6; + + switch (quad) { + case BOTTOM_LEFT: + xoff=0; + yoff=0; + break; + case BOTTOM_RIGHT: + xoff=1; + yoff=0; + break; + case TOP_LEFT: + xoff=0; + yoff=1; + break; + case TOP_RIGHT: + xoff=1; + yoff=1; + break; + default: + ; + } + cc[0][0]=cluster[yoff][xoff]; + cc[1][0]=cluster[yoff+1][xoff]; + cc[0][1]=cluster[yoff][xoff+1]; + cc[1][1]=cluster[yoff+1][xoff+1]; + + /* cout << cl[0] << " " << cl[1] << " " << cl[2] << endl; */ + /* cout << cl[3] << " " << cl[4] << " " << cl[5] << endl; */ + /* cout << cl[6] << " " << cl[7] << " " << cl[8] << endl; */ + /* cout <<"******"<=0 && ey>=0) - heta[ey*nbeta+ex]++; + heta[ey*nbeta+ex]++; + // cout << "*"<< etax << " " << etay << endl; + /* cout << etax << " " << etay << " " << ex << " " << ey << " " << ey*nbeta+ex << endl; */ + /* cout <<"********"<< endl << endl ; */ #endif return 0; }; diff --git a/slsDetectorCalibration/etaVEL/linearInterpolation.h b/slsDetectorCalibration/etaVEL/linearInterpolation.h index 21c323686..ae7f52e3a 100644 --- a/slsDetectorCalibration/etaVEL/linearInterpolation.h +++ b/slsDetectorCalibration/etaVEL/linearInterpolation.h @@ -32,13 +32,62 @@ class linearInterpolation : public slsInterpolation{ int corner; corner=calcQuad(data, tot, totquad, sDum); - - calcEta(totquad, sDum, etax, etay); + if (ns>2) { + calcEta(totquad, sDum, etax, etay); + } getInterpolatedPosition(x, y, etax,etay, corner, int_x, int_y); return; }; - + + + virtual int getInterpolatedPosition(int x, int y, double totquad,int quad,double *cl,double &etax, double &etay) { + + if (ns>2) { + double cc[2][2]; + double *cluster[3]; + cluster[0]=cl; + cluster[1]=cl+3; + cluster[2]=cl+6; + + switch (quad) { + case BOTTOM_LEFT: + xoff=0; + yoff=0; + break; + case BOTTOM_RIGHT: + xoff=1; + yoff=0; + break; + case TOP_LEFT: + xoff=0; + yoff=1; + break; + case TOP_RIGHT: + xoff=1; + yoff=1; + break; + default: + ; + } + cc[0][0]=cluster[yoff][xoff]; + cc[1][0]=cluster[yoff+1][xoff]; + cc[0][1]=cluster[yoff][xoff+1]; + cc[1][1]=cluster[yoff+1][xoff+1]; + double eta_x, eta_y; + calcEta(quadTot,cc,eta_x,eta_y); + } + return getInterpolatedPosition(x,y,etax, etay,quad,int_x,int_y); + + + + + + + + } + + virtual void getInterpolatedPosition(int x, int y, double etax, double etay, int corner, double &int_x, double &int_y) { @@ -69,9 +118,13 @@ class linearInterpolation : public slsInterpolation{ } - xpos_eta=(etax); - ypos_eta=(etay); - + if (ns>2) { + xpos_eta=(etax); + ypos_eta=(etay); + } else { + xpos_eta=0; + xpos_eta=0; + } int_x=((double)x) + 0.5*dX + xpos_eta; int_y=((double)y) + 0.5*dY + ypos_eta; @@ -108,7 +161,7 @@ class linearInterpolation : public slsInterpolation{ virtual int addToFlatField(double *cluster, double &etax, double &etay){}; virtual int addToFlatField(double etax, double etay){}; - + virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) {}; protected: ; diff --git a/slsDetectorCalibration/etaVEL/noInterpolation.h b/slsDetectorCalibration/etaVEL/noInterpolation.h index 936360a67..3f47d6c71 100644 --- a/slsDetectorCalibration/etaVEL/noInterpolation.h +++ b/slsDetectorCalibration/etaVEL/noInterpolation.h @@ -42,6 +42,9 @@ class noInterpolation : public slsInterpolation{ { return getInterpolatedPosition(x, y, NULL, int_x, int_y); }; + virtual int getInterpolatedPosition(int x, int y, double totquad,int quad,double *cl,double &etax, double &etay){ + return getInterpolatedPosition(x, y, NULL, int_x, int_y); + } ////////////////////////////////////////////////////////////////////////////////////// virtual void getPositionETA3(int x, int y, double *data, double &int_x, double &int_y) @@ -56,6 +59,9 @@ class noInterpolation : public slsInterpolation{ virtual int addToFlatField(double *cluster, double &etax, double &etay){return 0;}; virtual int addToFlatField(double etax, double etay){return 0;}; + virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay){return 0;}; + + protected: ; // TRandom *eventGenerator; diff --git a/slsDetectorCalibration/etaVEL/slsInterpolation.h b/slsDetectorCalibration/etaVEL/slsInterpolation.h index 7bdd1e488..506bc1ae5 100644 --- a/slsDetectorCalibration/etaVEL/slsInterpolation.h +++ b/slsDetectorCalibration/etaVEL/slsInterpolation.h @@ -102,19 +102,30 @@ hint=new TH2F("hint","hint",ns*nx, 0, nx, ns*ny, 0, ny); - - - - - - - - - //return position inside the pixel for the given photon virtual void getInterpolatedPosition(int x, int y, double *data, double &int_x, double &int_y)=0; //return position inside the pixel for the given photon virtual void getInterpolatedPosition(int x, int y, double etax, double etay, int quad, double &int_x, double &int_y)=0; + virtual void getInterpolatedPosition(int x, int y, double totquad,int quad,double *cluster,double &etax, double &etay)=0; + + + //return position inside the pixel for the given photon + virtual void clearInterpolatedImage() { + + +#ifdef MYROOT1 + hint->Reset(); +#endif +#ifndef MYROOT1 + for (int ix=0; ix=0 && ix<(nPixelsX*nSubPixels) && iy<(nSubPixels*nPixelsY) && iy>=0 )(*(hint+ix+iy*nPixelsX*nSubPixels))+=1; return hint; }; @@ -133,22 +145,26 @@ hint=new TH2F("hint","hint",ns*nx, 0, nx, ns*ny, 0, ny); virtual int addToFlatField(double *cluster, double &etax, double &etay)=0; virtual int addToFlatField(double etax, double etay)=0; - + + virtual int addToFlatField(double totquad,int quad,double *cluster,double &etax, double &etay)=0; + #ifdef MYROOT1 virtual TH2D *getFlatField(){return NULL;}; - virtual TH2D *setFlatField(TH2D *h){return NULL;}; + virtual TH2D *setFlatField(TH2D *h, int nb=-1, double emin=-1, double emax=-1){return NULL;}; + virtual TH2D *getFlatField(int &nb, double &emin, double &emax){nb=0; emin=0; emax=0; return getFlatField();}; #endif - + #ifndef MYROOT1 virtual int *getFlatField(){return NULL;}; - virtual int *setFlatField(int *h){return NULL;}; + virtual int *setFlatField(int *h, int nb=-1, double emin=-1, double emax=-1){return NULL;}; void *writeFlatField(const char * imgname){return NULL;}; void *readFlatField(const char * imgname, int nb=-1, double emin=1, double emax=0){return NULL;}; + virtual int *getFlatField(int &nb, double &emin, double &emax){nb=0; emin=0; emax=0; return getFlatField();}; #endif //virtual void Streamer(TBuffer &b); - + static int calcQuad(double *cl, double &sum, double &totquad, double sDum[2][2]){ int corner = UNDEFINED_QUADRANT; @@ -303,6 +319,35 @@ hint=new TH2F("hint","hint",ns*nx, 0, nx, ns*ny, 0, ny); } + static int calcMyEta(double totquad, int quad, double *cl, double &etax, double &etay) { + double l,r,t,b, sum; + int yoff; + switch (quad) { + case BOTTOM_LEFT: + case BOTTOM_RIGHT: + yoff=0; + break; + case TOP_LEFT: + case TOP_RIGHT: + yoff=1; + break; + default: + ; + } + l=cl[0+yoff*3]+cl[0+yoff*3+3]; + r=cl[2+yoff*3]+cl[2+yoff*3+3]; + b=cl[0+yoff*3]+cl[1+yoff*3]*cl[2+yoff*3]; + t=cl[0+yoff*3+3]+cl[1+yoff*3+3]*cl[0+yoff*3+3]; + sum=t+b; + if (sum>0) { + etax=(-l+r)/sum; + etay=(+t)/sum; + } + + return -1; + } + + static int calcEta3X(double *cl, double &etax, double &etay, double &sum) { double l,r,t,b; diff --git a/slsDetectorCalibration/interpolatingDetector.h b/slsDetectorCalibration/interpolatingDetector.h index 0f92ea618..4135bb35b 100644 --- a/slsDetectorCalibration/interpolatingDetector.h +++ b/slsDetectorCalibration/interpolatingDetector.h @@ -6,6 +6,7 @@ #include "slsInterpolation.h" +#define M015 #ifdef MYROOT1 #include @@ -16,7 +17,10 @@ #include using namespace std; - +#define XMIN 350/2 +#define XMAX 600/2 +#define YMIN 0 +#define YMAX 400 class interpolatingDetector : public singlePhotonDetector { @@ -46,13 +50,21 @@ class interpolatingDetector : public singlePhotonDetector { commonModeSubtraction *cm=NULL, int nped=1000, int nd=100, int nnx=-1, int nny=-1) : - singlePhotonDetector(d, 3,nsigma,sign, cm, nped, nd, nnx, nny) , interp(inte), id(0) {}; + singlePhotonDetector(d, 3,nsigma,sign, cm, nped, nd, nnx, nny) , interp(inte), id(0), xmin(XMIN), xmax(XMAX), ymin(YMIN), ymax(YMAX) { + cout << "**"<< xmin << " " << xmax << " " << ymin << " " << ymax << endl; + +}; interpolatingDetector(interpolatingDetector *orig) : singlePhotonDetector(orig) { interp=(orig->interp)->Clone(); id=orig->id; + xmin=orig->xmin; + xmax=orig->xmax; + ymin=orig->ymin; + ymax=orig->ymax; + } @@ -75,19 +87,27 @@ class interpolatingDetector : public singlePhotonDetector { writeGainMap(tit); } #endif - interp->prepareInterpolation(ok); + if (interp) + interp->prepareInterpolation(ok); } - int getImageSize(int &nnx, int &nny, int &ns) {return interp->getImageSize(nnx, nny, ns);}; + + void clearImage() {if (interp) interp->clearInterpolatedImage(); else singlePhotonDetector::clearImage();}; + + int getImageSize(int &nnx, int &nny, int &ns) {if (interp) return interp->getImageSize(nnx, nny, ns); else return analogDetector::getImageSize(nnx, nny, ns);}; #ifdef MYROOT1 - virtual TH2F *getInterpolatedImage() + virtual TH2F *getImage() #endif #ifndef MYROOT1 - virtual int *getInterpolatedImage() + virtual int *getImage() #endif { + // cout << "image " << endl; if (interp) return interp->getInterpolatedImage(); + else + return analogDetector::getImage(); + //cout << "null " << endl; } #ifdef MYROOT1 @@ -101,7 +121,8 @@ class interpolatingDetector : public singlePhotonDetector { if (interp) return interp->getInterpolatedImage(); else - return NULL; + singlePhotonDetector::getImage(); + //return NULL; }; @@ -119,15 +140,19 @@ class interpolatingDetector : public singlePhotonDetector { return NULL; }; - void *writeInterpolatedImage(const char * imgname) { - cout << id << "=" << imgname<< endl; - interp->writeInterpolatedImage(imgname); + void *writeImage(const char * imgname) { + // cout << id << "=" << imgname<< endl; + if (interp) + interp->writeInterpolatedImage(imgname); + else + analogDetector::writeImage(imgname); return NULL; } int addFrame(char *data, single_photon_hit *clusters=NULL, int ff=0) { + double g=1; single_photon_hit *cl; single_photon_hit clust; @@ -136,16 +161,22 @@ class interpolatingDetector : public singlePhotonDetector { else cl=&clust; + + int ccs=clusterSize; + int ccy=clusterSizeY; + double int_x,int_y, eta_x, eta_y; int nph=0; - double val[ny][nx]; + double rest[ny][nx]; int cy=(clusterSizeY+1)/2; int cs=(clusterSize+1)/2; int ir, ic; double cc[2][2]; - double max=0, tl=0, tr=0, bl=0,br=0, *v, vv; + double max=0, tl=0, tr=0, bl=0,br=0, v, vv; int xoff,yoff; int skip=0; + // cout <<"fr"<< endl; + double tthr; if (iframerms=getPedestalRMS(ix,iy); - //(clusters+nph)->rms=getPedestalRMS(ix,iy); - // cout << iframe << " " << nph << " " << ix << " " << iy << endl; - - 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)=0 && (ix+ic)1 && ir<=0)) { - ; - } else if (ix>1 && ic<=0) { - ; - } else { - val[iy+ir][ix+ic]=subtractPedestal(data,ix+ic,iy+ir); - eventMask[iy+ir][ix+ic]=PEDESTAL; - } - - /* if ((iy+ir)>=iy && (ix+ic)>=ix) { */ - /* val[iy+ir][ix+ic]=subtractPedestal(data,ix+ic,iy+ir); */ - /* eventMask[iy+ir][ix+ic]=PEDESTAL; */ - /* } */ - - // cout << ir << " " << ic << " " << val[iy+ir][ix+ic] << endl; - v=&(val[iy+ir][ix+ic]); - - if (ir==0 && ic==0) { - if (*v<-nSigma*cl->rms) { - eventMask[iy][ix]=NEGATIVE_PEDESTAL; - // cout << "neg ped" << endl; - } - } - - // if (skip==0) { - 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 (*v>max) { - max=*v; - } - // } - // } else skip=1; - } - } - - if (bl>=br && bl>=tl && bl>=tr) { - cl->quad=BOTTOM_LEFT; - cl->quadTot=bl; - xoff=0; - yoff=0; - } else if (br>=bl && br>=tl && br>=tr) { - cl->quad=BOTTOM_RIGHT; - cl->quadTot=br; - xoff=1; - yoff=0; - } else if (tl>=br && tl>=bl && tl>=tr) { - cl->quad=TOP_LEFT; - cl->quadTot=tl; - xoff=0; - yoff=1; - } else if (tr>=bl && tr>=tl && tr>=br) { - cl->quad=TOP_RIGHT; - cl->quadTot=tr; - xoff=1; - yoff=1; - } - - if (max>nSigma*cl->rms || tot>sqrt(clusterSizeY*clusterSize)*nSigma*cl->rms || (cl->quadTot)>sqrt(cy*cs)*nSigma*cl->rms) { - if (val[iy][ix]>=max) { - // cout << "max" << endl; - eventMask[iy][ix]=PHOTON_MAX; - cl->tot=tot; - cl->x=ix; - cl->y=iy; - cl->ped=getPedestal(ix,iy, 0); - // cout << iframe << " " << ix << " " << iy << " "<< (clusters+nph)->tot << " " << (clusters+nph)->quadTot << " " << (clusters+nph)->ped<< " " << (clusters+nph)->rms << endl; + tthr=nSigma*getPedestalRMS(ix,iy)/g; + if (ix==xmin || iy==ymin) + rest[iy][ix]=subtractPedestal(data,ix,iy); + + + max=0; + tl=0; + tr=0; + bl=0; + br=0; + tot=0; + quadTot=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)=0 && (ix+ic)set_data(val[iy+ir][ix+ic],ic,ir); + if ((iy+ir)>=ymin && (iy+ir)=xmin && (ix+ic)set_data(rest[iy+ir][ix+ic], ic, ir); - - // cout << val[iy+ir][ix+ic] << " " ; + if (ir>=0 && ic>=0 ) + rest[iy+ir][ix+ic]=subtractPedestal(data,ix+ic,iy+ir); + + v=rest[iy+ir][ix+ic];//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 (v>max) { + max=v; + } + // if (ir==0 && ic==0) { + if (v>tthr) { + eventMask[iy][ix]=NEIGHBOUR; + } + //} + } + } + } + if (rest[iy][ix]<=-tthr) { + eventMask[iy][ix]=NEGATIVE_PEDESTAL; + //if (cluster->get_data(0,0)>=max) { + } else if (max>tthr || tot>sqrt(ccy*ccs)*tthr || quadTot>sqrt(cy*cs)*tthr) { + if (rest[iy][ix]>=max) { + if (bl>=br && bl>=tl && bl>=tr) { + cl->quad=BOTTOM_LEFT; + cl->quadTot=bl; + } else if (br>=bl && br>=tl && br>=tr) { + cl->quad=BOTTOM_RIGHT; + cl->quadTot=br; + } else if (tl>=br && tl>=bl && tl>=tr) { + cl->quad=TOP_LEFT; + cl->quadTot=tl; + } else if (tr>=bl && tr>=tl && tr>=br) { + cl->quad=TOP_RIGHT; + cl->quadTot=tr; + } + + + + /* cout << ix << " " << iy << " " << rest[iy][ix] <<" " << tot << " " << quadTot << endl; */ + /* for (int ir=-(clusterSizeY/2); ir<(clusterSizeY/2)+1; ir++) { */ + /* for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) */ + /* cout << rest[iy+ir][ix+ic] << " " ; */ + /* cout << endl; */ + /* } */ + eventMask[iy][ix]=PHOTON_MAX; + cl->tot=tot; + cl->x=ix; + cl->y=iy; + cl->ped=getPedestal(ix,iy, 0); + cl->rms=getPedestalRMS(ix,iy); + + for (int ir=-(clusterSizeY/2); ir<(clusterSizeY/2)+1; ir++) { + for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) { + if ((iy+ir)>=ymin && (iy+ir)=xmin && (ix+ic)set_data(rest[iy+ir][ix+ic],ic,ir); } } - // cout << endl << " " ; } - // cout << endl << " " ; - cc[0][0]=cl->get_data(-1+xoff,-1+yoff); - cc[1][0]=cl->get_data(-1+xoff,0+yoff); - cc[0][1]=cl->get_data(0+xoff,-1+yoff); - cc[1][1]=cl->get_data(0+xoff,0+yoff); - - // cout << cc[0][0] << " " << cc[0][1] << endl; - // cout << cc[1][0] << " " << cc[1][1] << endl; + if (interp) { + if (ff) { +#ifdef M015 + if (iy>100) +#endif + interp->addToFlatField(cl->quadTot,cl->quad,cl->get_cluster(),eta_x, eta_y); + // if ((eta_x<0.1 || eta_x>0.9)&&(eta_y<0.1 || eta_y>0.9)) + // cout << ix << " " << iy << " " << eta_x <<" " << eta_y << endl; + + } else { + interp->getInterpolatedPosition(ix, iy, cl->quadTot,cl->quad,cl->get_cluster(),int_x, int_y); + interp->addToImage(int_x, int_y); + } + + } else + image[ix+nx*iy]++; + + nph++; + + if (clusters) cl=(clusters+nph); - // cout << endl << " " ; + // rest[iy][ix]-=tthr; + } else + eventMask[iy][ix]=PHOTON; + //else if (thr<=0 ) { + //addToPedestal(data,ix,iy); + // } + } + if (eventMask[iy][ix]==PEDESTAL) { + addToPedestal(data,ix,iy); + } + } + } + return nph; + } + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + /* for (int ix=0; ixrms=getPedestalRMS(ix,iy); */ + /* //(clusters+nph)->rms=getPedestalRMS(ix,iy); */ + + /* // cout << iframe << " " << nph << " " << ix << " " << iy << endl; */ + /* if (ix==0 || iy==0) */ + /* val[iy][ix]=subtractPedestal(data,ix,iy); */ + + /* if (val[iy][ix]<-nSigma*cl->rms) { */ + /* eventMask[iy][ix]=NEGATIVE_PEDESTAL; */ + /* // cout << "neg ped" << endl; */ + /* } else { */ + /* for (int ir=-(clusterSizeY/2); ir<(clusterSizeY/2)+1; ir++) { */ + /* for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) { */ - if (interp) { - //interp->calcEta((clusters+nph)->quadTot,cc,eta_x, eta_y); - interp->calcEtaL(cl->quadTot,cl->quad,cc,eta_x, eta_y); - // cout << eta_x << " " << eta_y << endl; - // cout << "eta" << endl; - if (ff) { - interp->addToFlatField(eta_x,eta_y); - // cout << "**************************************************************************"<< endl; - } else { - // cout << "interp" << endl; - interp->getInterpolatedPosition(ix,iy,eta_x,eta_y,cl->quad,int_x,int_y); - // cout << "add" << endl; - interp->addToImage(int_x, int_y); - } - // cout << "done" << endl; - } - if (clusters) cl=(clusters+nph); + + /* if ((iy+ir)>=0 && (iy+ir)=0 && (ix+ic)=0 && ic>=0) { */ + /* val[iy+ir][ix+ic]=subtractPedestal(data,ix+ic,iy+ir); */ + /* eventMask[iy+ir][ix+ic]=PEDESTAL; */ + /* } */ + + /* // cout << ir << " " << ic << " " << val[iy+ir][ix+ic] << endl; */ + /* v=&(val[iy+ir][ix+ic]); */ + /* // if (skip==0) { */ + /* 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 (*v>max) { */ + /* max=*v; */ + /* } */ + + /* } */ + /* } */ + /* } */ - nph++; - } else { - // cout << "ph" << endl; - eventMask[iy][ix]=PHOTON; - } - } else if (eventMask[iy][ix]==PEDESTAL) { - // cout << "ped" << endl; - addToPedestal(data,ix,iy); - } + /* if (bl>=br && bl>=tl && bl>=tr) { */ + /* cl->quad=BOTTOM_LEFT; */ + /* cl->quadTot=bl; */ + /* } else if (br>=bl && br>=tl && br>=tr) { */ + /* cl->quad=BOTTOM_RIGHT; */ + /* cl->quadTot=br; */ + /* } else if (tl>=br && tl>=bl && tl>=tr) { */ + /* cl->quad=TOP_LEFT; */ + /* cl->quadTot=tl; */ + /* } else if (tr>=bl && tr>=tl && tr>=br) { */ + /* cl->quad=TOP_RIGHT; */ + /* cl->quadTot=tr; */ + /* } */ - } - } + /* if (max>nSigma*cl->rms || tot>sqrt(clusterSizeY*clusterSize)*nSigma*cl->rms || (cl->quadTot)>sqrt(cy*cs)*nSigma*cl->rms) { */ + /* if (val[iy][ix]>=max) { */ + /* eventMask[iy][ix]=PHOTON_MAX; */ + /* cl->tot=tot; */ + /* cl->x=ix; */ + /* cl->y=iy; */ + /* cl->ped=getPedestal(ix,iy, 0); */ + + /* if (interp) { */ + /* if (ff) { */ + /* interp->addToFlatField(cl->quadTot,cl->quad,cl->get_cluster(),eta_x, eta_y); */ + /* // if ((eta_x<0.1 || eta_x>0.9)&&(eta_y<0.1 || eta_y>0.9)) */ + /* // cout << ix << " " << iy << " " << eta_x <<" " << eta_y << endl; */ + /* } else { */ + /* interp->getInterpolatedPosition(ix, iy, cl->quadTot,cl->quad,cl->get_cluster(),int_x, int_y); */ + /* interp->addToImage(int_x, int_y); */ + /* } */ + /* } else */ + /* image[ix+nx*iy]++; */ + + /* if (clusters) cl=(clusters+nph); */ + + /* nph++; */ + + /* } else { */ + /* eventMask[iy][ix]=PHOTON; */ + /* } */ + /* } else if (eventMask[iy][ix]==PEDESTAL) { */ + /* addToPedestal(data,ix,iy); */ + /* } */ + + /* } */ + /* } */ + /* } */ + /* return nph; */ - return nph; - - }; + /* }; */ + + + + + + + virtual void processData(char *data, frameMode i=eFrame, int *val=NULL) { + if (interp){ switch(i) { case ePedestal: addToPedestal(data); break; case eFlat: - addFrame(data,NULL,1); + addFrame(data,NULL,1); break; default: - addFrame(data); + addFrame(data,NULL,0); } + } else + singlePhotonDetector::processData(data,i,val); + + }; + virtual char *getInterpolation(){return (char*)interp;}; protected: slsInterpolation *interp; int id; + //should put it to analogDetector + int xmin, xmax, ymin, ymax; }; diff --git a/slsDetectorCalibration/moench03T1ZmqData.h b/slsDetectorCalibration/moench03T1ZmqData.h new file mode 100644 index 000000000..5c01eb3bc --- /dev/null +++ b/slsDetectorCalibration/moench03T1ZmqData.h @@ -0,0 +1,268 @@ +#ifndef MOENCH03T1ZMQDATA_H +#define MOENCH03T1ZMQDATA_H +#include "slsDetectorData.h" + + + +class moench03T1ZmqData : public slsDetectorData { + + private: + + int iframe; + int nadc; + int sc_width; + int sc_height; + const int nPackets; /**(400, 400, ps*npackets), packetSize(ps), nPackets(npackets) { + + int nadc=32; + int sc_width=25; + int sc_height=200; + + int adc_nr[32]={300,325,350,375,300,325,350,375, \ + 200,225,250,275,200,225,250,275,\ + 100,125,150,175,100,125,150,175,\ + 0,25,50,75,0,25,50,75}; + + int row, col; + + int isample; + int iadc; + int ix, iy; + + // int npackets=40; + int i; + int adc4(0); + + for (int ip=0; ip=8192*40) + cout << "Error: pointer " << dataMap[row][col] << " out of range "<< endl; + } + } + } + } + + int ipacket; + int ibyte; + int ii=0; + for (int ipacket=0; ipacket0) { */ +/* iframe++; */ +/* // cout << ib << "-" << endl; */ +/* return (char*)afifo_cont; */ +/* } else { */ +/* delete [] afifo_cont; */ +/* return NULL; */ +/* } */ +/* } */ +/* return NULL; */ +/* }; */ + + + virtual char *readNextFrame(ifstream &filebin) { + int ff=-1, np=-1; + return readNextFrame(filebin, ff, np); + }; + + virtual char *readNextFrame(ifstream &filebin, int &ff) { + int np=-1; + return readNextFrame(filebin, ff, np); + }; + + virtual char *readNextFrame(ifstream &filebin, int& ff, int &np) { + char *data=new char[packetSize*nPackets]; + char *d=readNextFrame(filebin, ff, np, data); + if (d==NULL) {delete [] data; data=NULL;} + return data; + } + + + + + virtual char *readNextFrame(ifstream &filebin, int& ff, int &np, char *data) { + char *retval=0; + int nd; + int fnum = -1; + np=0; + int pn; + + + if (ff>=0) + fnum=ff; + + if (filebin.is_open()) { + if (filebin.read(data, packetSize*nPackets) ){ + iframe++; + ff=iframe; + return data; + } + } + return NULL; + + + + }; + + + + /** + + Loops over a memory slot until a complete frame is found (i.e. all packets 0 to nPackets, same frame number). purely virtual func + \param data pointer to the memory to be analyzed + \param ndata reference to the amount of data found for the frame, in case the frame is incomplete at the end of the memory slot + \param dsize size of the memory slot to be analyzed + \returns pointer to the beginning of the last good frame (might be incomplete if ndata smaller than dataSize), or NULL if no frame is found + + */ + virtual char *findNextFrame(char *data, int &ndata, int dsize){ + if (dsize -#include "interpolatingDetector.h" +#include "analogDetector.h" #include "circularFifo.h" - +include "etaVEL/slsInterpolation.h" @@ -21,7 +21,7 @@ class threadedDetector { public: - threadedDetector(interpolatingDetector *d, int fs=10000) { + threadedDetector(analogDetector *d, int fs=10000) { char *mem, *mm; det=d; fifoFree=new CircularFifo(fs); @@ -31,47 +31,54 @@ public: mm=mem+i*det->getDataSize(); fifoFree->push(mm); } + busy=0; stop=1; fMode=eFrame; } - int setFrameMode(int fm) {if (fMode>=0) fMode=fm; return fMode;} + virtual int setFrameMode(int fm) {if (fMode>=0) fMode=fm; return fMode;} - void prepareInterpolation(int &ok) { - cout << "-" << endl; - det->prepareInterpolation(ok); - }; + /* void prepareInterpolation(int &ok) { */ + /* cout << "-" << endl; */ + /* det->prepareInterpolation(ok); */ + /* }; */ - int *getInterpolatedImage() { + virtual int *getImage() { return det->getInterpolatedImage(); } - int getImageSize(int &nnx, int &nny, int &ns) {return det->getImageSize(nnx, nny, ns);}; - virtual ~threadedDetector() {StopThread(); free(mem); delete fifoFree; delete fifoData;} + virtual int getImageSize(int &nnx, int &nny, int &ns) {return det->getImageSize(nnx, nny, ns);}; + + ~threadedDetector() {StopThread(); free(mem); delete fifoFree; delete fifoData;} /** Returns true if the thread was successfully started, false if there was an error starting the thread */ - bool StartThread() - { stop=0; + virtual bool StartThread() + { stop=0; return (pthread_create(&_thread, NULL, processData, this) == 0); } - void StopThread() + virtual void StopThread() { stop=1; (void) pthread_join(_thread, NULL); } - - bool pushData(char* &ptr) { + + virtual bool pushData(char* &ptr) { fifoData->push(ptr); } - bool popFree(char* &ptr) { + virtual bool popFree(char* &ptr) { fifoFree->pop(ptr); } + virtual int isBusy() {return busy;} + //protected: /** Implement this method in your subclass with the code you want your thread to run. */ //virtual void InternalThreadEntry() = 0; - void *writeInterpolatedImage(const char * imgname) {cout << "a" <writeInterpolatedImage(imgname);}; + virtual void *writeImage(const char * imgname) {cout << "a" <writeImage(imgname);}; + + virtual void clearImage(){det->clearImage();}; + private: interpolatingDetector *det; int fMode; @@ -81,18 +88,21 @@ private: CircularFifo *fifoFree; CircularFifo *fifoData; int stop; + int busy; char *data; static void * processData(void * ptr) { threadedDetector *This=((threadedDetector *)ptr); return This->processData(); } - void * processData() { + busy=1; while (!stop) { if (fifoData->isEmpty()) { + busy=0; usleep(100); } else { + busy=1; fifoData->pop(data); //blocking! det->processData(data,(frameMode)fMode); fifoFree->push(data); @@ -103,9 +113,6 @@ private: - - - }; @@ -113,7 +120,7 @@ private: class multiThreadedDetector { public: - multiThreadedDetector(interpolatingDetector *d, int n, int fs=1000) : nThreads(n), ithread(0) { + multiThreadedDetector(analogDetector *d, int n, int fs=1000) : nThreads(n), ithread(0) { dd[0]=d; if (nThreads==1) dd[0]->setId(100); @@ -146,24 +153,15 @@ public: int setFrameMode(int fm) { int ret; for (int i=0; isetFrameMode(fm); return ret;}; - void prepareInterpolation(int &ok) { - int oo; - ok=1; - for (int i=0; iprepareInterpolation(oo); - //if (oo<1) ok=0; - } - }; - - int *getInterpolatedImage() { + + int *getImage() { int *img; int nnx, nny, ns; int nn=dets[0]->getImageSize(nnx, nny, ns); //for (i=0; igetInterpolatedImage(); + img=dets[ii]->getImage(); for (int i=0; iclearImage(); + } + + } + + + + + + + + void *writeImage(const char * imgname) { #ifdef SAVE_ALL for (int ii=0; iiwriteInterpolatedImage(tit); + dets[ii]->writeImage(tit); } #endif - getInterpolatedImage(); + getImage(); int nnx, nny, ns; int nn=dets[0]->getImageSize(nnx, nny, ns); float *gm=new float[ nn]; @@ -214,6 +227,17 @@ public: } + + int isBusy() { + int ret=0, ret1; + for (int i=0; iisBusy(); + ret|=ret1; + if (ret1) cout << "thread " << i <<" still busy " << endl; + } + return ret; + } + bool pushData(char* &ptr) { dets[ithread]->pushData(ptr); @@ -229,6 +253,14 @@ public: return ithread; } + virtual void prepareInterpolation(int &ok){ + slsInterpolation * + +}; + + + + private: bool stop; const int nThreads; diff --git a/slsDetectorCalibration/singlePhotonDetector.h b/slsDetectorCalibration/singlePhotonDetector.h index c053c0ad2..4e9f81123 100644 --- a/slsDetectorCalibration/singlePhotonDetector.h +++ b/slsDetectorCalibration/singlePhotonDetector.h @@ -94,7 +94,6 @@ public analogDetector { eventMask=new eventType*[ny]; for (int i=0; i { clusterSize=orig->clusterSize; clusterSizeY=orig->clusterSizeY; cluster=new single_photon_hit(clusterSize,clusterSizeY); + quad=UNDEFINED_QUADRANT; tot=0; quadTot=0; + gmap=orig->gmap; } @@ -139,119 +140,166 @@ public analogDetector { virtual int *getNPhotons(char *data, double thr=-1, int *nph=NULL) { - + double val; if (nph==NULL) - nph=new int[nx*ny]; + nph=image; + //nph=new int[nx*ny]; + double rest[ny][nx]; int cy=(clusterSizeY+1)/2; int cs=(clusterSize+1)/2; - double g=-1.; - int ccs=clusterSize; int ccy=clusterSizeY; + + double g=1.; - + double tthr=thr; int nn=0; double max=0, tl=0, tr=0, bl=0,br=0, v; - - + + if (thr>=0) { cy=1; cs=1; ccs=1; ccy=1; } - - - for (int ix=0; ix::getNPhotons(data,ix,iy,tthr); - nph[ix+nx*iy]+=nn; - - rest[iy][ix]=(val-nn*tthr); - - - } - } - for (int ix=clusterSize/2; ix=0 && (iy+ir)=0 && (ix+ic)set_data(rest[iy+ir][ix+ic], ic, ir); - v=rest[iy+ir][ix+ic];//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 (v>max) { - max=v; - } - if (ir==0 && ic==0) { - if (v>tthr) { - eventMask[iy][ix]=PHOTON; + if (iframe0) { + cout << "threshold" << endl; + for (int ix=0; ix::getNPhotons(data,ix,iy,tthr); + nph[ix+nx*iy]+=nn; + rest[iy][ix]=(val-nn*tthr); + + } + } + } + for (int ix=0; ix=0 && (iy+ir)=0 && (ix+ic)set_data(rest[iy+ir][ix+ic], ic, ir); + + if (thr<=0 && ir>=0 && ic>=0 ) + rest[iy+ir][ix+ic]=subtractPedestal(data,ix+ic,iy+ir); + + v=rest[iy+ir][ix+ic];//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 (v>max) { + max=v; } + // if (ir==0 && ic==0) { + if (v>tthr) { + eventMask[iy][ix]=NEIGHBOUR; + } + //} } } } - } - //if (cluster->get_data(0,0)>=max) { - if (rest[iy][ix]>=max) { - 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; + if (rest[iy][ix]<=-tthr) { + eventMask[iy][ix]=NEGATIVE_PEDESTAL; + //if (cluster->get_data(0,0)>=max) { + } else if (rest[iy][ix]>=max) { + 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; + } + if (max>tthr || tot>sqrt(ccy*ccs)*tthr || quadTot>sqrt(cy*cs)*tthr) { + /* cout << ix << " " << iy << " " << rest[iy][ix] <<" " << tot << " " << quadTot << endl; */ + /* for (int ir=-(clusterSizeY/2); ir<(clusterSizeY/2)+1; ir++) { */ + /* for (int ic=-(clusterSize/2); ic<(clusterSize/2)+1; ic++) */ + /* cout << rest[iy+ir][ix+ic] << " " ; */ + /* cout << endl; */ + /* } */ + eventMask[iy][ix]=PHOTON; + nph[ix+nx*iy]++; + // rest[iy][ix]-=tthr; + } //else if (thr<=0 ) { + //addToPedestal(data,ix,iy); + // } } - if (rest[iy][ix]>tthr || tot>sqrt(ccy*ccs)*tthr || quadTot>sqrt(cy*cs)*tthr) { - nph[ix+nx*iy]++; - rest[iy][ix]-=tthr; + if (thr<0 && eventMask[iy][ix]==PEDESTAL) { + addToPedestal(data,ix,iy); } } } - } + + // } + // } + + + return nph; } @@ -551,6 +599,7 @@ int getClusters(char *data, single_photon_hit *clusters) { #endif virtual void processData(char *data, frameMode i=eFrame, int *val=NULL) { + // cout << "sp" << endl; switch(i) { case ePedestal: addToPedestal(data); @@ -558,6 +607,8 @@ int getClusters(char *data, single_photon_hit *clusters) { default: getNPhotons(data,-1,val); } + iframe++; + // cout << "done" << endl; }; diff --git a/slsDetectorCalibration/single_photon_hit.h b/slsDetectorCalibration/single_photon_hit.h index 54bd4060f..ba06c363a 100644 --- a/slsDetectorCalibration/single_photon_hit.h +++ b/slsDetectorCalibration/single_photon_hit.h @@ -57,7 +57,7 @@ class single_photon_hit { \returns value of the cluster element */ double get_data(int ix, int iy=0){return data[(iy+dy/2)*dx+ix+dx/2];}; - + double *get_cluster() {return data;}; int x; /**< x-coordinate of the center of hit */ int y; /**< x-coordinate of the center of hit */ double rms; /**< noise of central pixel l -- at some point it can be removed*/ diff --git a/slsDetectorCalibration/slsDetectorCalibration.doxy b/slsDetectorCalibration/slsDetectorCalibration.doxy new file mode 100644 index 000000000..a1d50fe5b --- /dev/null +++ b/slsDetectorCalibration/slsDetectorCalibration.doxy @@ -0,0 +1,87 @@ +# If the EXTRACT_ALL tag is set to YES doxygen will assume all entities in +# documentation are documented, even if no documentation was available. +# Private class members and static file members will be hidden unless +# the EXTRACT_PRIVATE and EXTRACT_STATIC tags are set to YES + +EXTRACT_ALL = YES + +# If the EXTRACT_PRIVATE tag is set to YES all private members of a class +# will be included in the documentation. + +EXTRACT_PRIVATE = NO + + + +# If the EXTRACT_STATIC tag is set to YES all static members of a file +# will be included in the documentation. + +EXTRACT_STATIC = YES + +# If the EXTRACT_LOCAL_CLASSES tag is set to YES classes (and structs) +# defined locally in source files will be included in the documentation. +# If set to NO only classes defined in header files are included. + +EXTRACT_LOCAL_CLASSES = YES + +# This flag is only useful for Objective-C code. When set to YES local +# methods, which are defined in the implementation section but not in +# the interface are included in the documentation. +# If set to NO (the default) only methods in the interface are included. + +EXTRACT_LOCAL_METHODS = YES + +# If this flag is set to YES, the members of anonymous namespaces will be +# extracted and appear in the documentation as a namespace called +# 'anonymous_namespace{file}', where file will be replaced with the base +# name of the file that contains the anonymous namespace. By default +# anonymous namespace are hidden. + +EXTRACT_ANON_NSPACES = NO + +# If the HIDE_UNDOC_MEMBERS tag is set to YES, Doxygen will hide all +# undocumented members of documented classes, files or namespaces. +# If set to NO (the default) these members will be included in the +# various overviews, but no documentation section is generated. +# This option has no effect if EXTRACT_ALL is enabled. + +HIDE_UNDOC_MEMBERS = NO + +# If the HIDE_UNDOC_CLASSES tag is set to YES, Doxygen will hide all +# undocumented classes that are normally visible in the class hierarchy. +# If set to NO (the default) these classes will be included in the various +# overviews. This option has no effect if EXTRACT_ALL is enabled. + +HIDE_UNDOC_CLASSES = NO + +# If the HIDE_FRIEND_COMPOUNDS tag is set to YES, Doxygen will hide all +# friend (class|struct|union) declarations. +# If set to NO (the default) these declarations will be included in the +# documentation. + +HIDE_FRIEND_COMPOUNDS = NO + +INTERNAL_DOCS = NO + +SHOW_INCLUDE_FILES = NO + +SHOW_FILES = NO + +SHOW_NAMESPACES = NO + +COMPACT_LATEX = YES + +PAPER_TYPE = a4 + +PDF_HYPERLINKS = YES + +USE_PDFLATEX = YES + +LATEX_HIDE_INDICES = YES + +PREDEFINED = __cplusplus + +INPUT = analogDetector.h pedestalSubtraction.h MovingStat.h singlePhotonDetector.h interpolatingDetector.h tiffIO.h single_photon_hit.h slsDetectorData.h moench03Ctb10GbT1Data.h moench03TCtbData.h moench03T1CtbData.h moench03CtbData.h moench03Ctb10GbData.h moench03TCtb10GbData.h Mythen3_01_jctbData.h adcSar2_jctbData.h eigerHalfModuleData.h energyCalibration.h slsReceiverData.h gotthardModuleData.h gotthardShortModuleData.h jungfrau02Data.h jungfrau10ModuleData.h moench02Ctb10GbData.h moench02CtbData.h moench02ModuleData.h moench03CommonMode.h moenchCommonMode.h chiptestBoardData.h commonModeSubtraction.h RunningStat.h etaVEL/etaInterpolationBase.h etaVEL/slsInterpolation.h etaVEL/etaInterpolationPosXY.h etaVEL/linearInterpolation.h etaVEL/noInterpolation.h etaVEL/etaInterpolationGlobal.h etaVEL/interpolation_EtaVEL.h etaVEL/EtaVEL.h etaVEL/iterativeEtaInterpolation.h multiThreadedDetector.h moench03T1ZmqData.h + + +OUTPUT_DIRECTORY = slsDetectorCalibrationDocs + diff --git a/slsDetectorCalibration/slsDetectorData.h b/slsDetectorCalibration/slsDetectorData.h index cd6da7cab..bc52e76a6 100644 --- a/slsDetectorCalibration/slsDetectorData.h +++ b/slsDetectorCalibration/slsDetectorData.h @@ -333,9 +333,6 @@ class slsDetectorData { - - - }; diff --git a/slsDetectorCalibration/tiffIO.h b/slsDetectorCalibration/tiffIO.h index 3076d641d..b27fdd938 100644 --- a/slsDetectorCalibration/tiffIO.h +++ b/slsDetectorCalibration/tiffIO.h @@ -1,5 +1,6 @@ -#ifndef TIFF_IO_H -#define TIFF_IO_H +#ifndef MY_TIFF_IO_H +#define MY_TIFF_IO_H + #include #include @@ -28,65 +29,8 @@ } \ } -void *WriteToTiff(float * imgData, const char * imgname, int nrow, int ncol){ - int sampleperpixel=1; - // unsigned char * buff=NULL; - tsize_t linebytes; - cout << "--" <