diff --git a/examples/moench03_T1_lab.config b/examples/moench03_T1_lab.config index 41bcfd872..a9788f2ef 100644 --- a/examples/moench03_T1_lab.config +++ b/examples/moench03_T1_lab.config @@ -1,4 +1,4 @@ -hostname bchip064+ +hostname bchip020+ patword 0000 0000000000000000 patword 0001 0000000000000000 @@ -429,6 +429,7 @@ rx_zmqport 50003 0:rx_hostname mpc2011 + ####pcmoench01 #0:rx_tcpport 1977 #0:rx_udpip 10.1.1.100 diff --git a/slsDetectorCalibration/interpolations/eta2InterpolationBase.h b/slsDetectorCalibration/interpolations/eta2InterpolationBase.h deleted file mode 100644 index 26431ed92..000000000 --- a/slsDetectorCalibration/interpolations/eta2InterpolationBase.h +++ /dev/null @@ -1,406 +0,0 @@ -#ifndef ETA2_INTERPOLATION_BASE_H -#define ETA2_INTERPOLATION_BASE_H - -#ifdef MYROOT1 -#include -#include -#include -#include -#endif - -#include "etaInterpolationBase.h" - -class eta2InterpolationBase : public virtual etaInterpolationBase { - - public: - eta2InterpolationBase(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : etaInterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax) { - - /* if (etamin>=etamax) { */ - /* etamin=-1; */ - /* etamax=2; */ - /* // cout << ":" <2 || nSubPixelsY>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 *data, double &int_x, double &int_y) - { - double sDum[2][2]; - double tot, totquad; - double etax=0,etay=0; - - int corner; - corner=calcQuad(data, tot, totquad, sDum); - if (nSubPixelsX>2 || nSubPixelsY>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]; - int xoff=0, yoff=0; - 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=0, etay=0; - if (nSubPixelsX>2 || nSubPixelsY>2) { - cc[0][0]=cl[xoff+3*yoff]; - cc[1][0]=cl[xoff+3*(yoff+1)]; - cc[0][1]=cl[xoff+1+3*yoff]; - cc[1][1]=cl[xoff+1+3*(yoff+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 totquad,int quad,int *cl,double &int_x, double &int_y) { - - double cc[2][2]; - int xoff=0, yoff=0; - - 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=0, etay=0; - if (nSubPixelsX>2 || nSubPixelsY>2) { - cc[0][0]=cl[xoff+3*yoff]; - cc[1][0]=cl[xoff+3*(yoff+1)]; - cc[0][1]=cl[xoff+1+3*yoff]; - cc[1][1]=cl[xoff+1+3*(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=0,ypos_eta=0; - double dX,dY; - int ex,ey; - switch (corner) - { - case TOP_LEFT: - dX=-1.; - dY=0; - break; - case TOP_RIGHT: - ; - dX=0; - dY=0; - break; - case BOTTOM_LEFT: - dX=-1.; - dY=-1.; - break; - case BOTTOM_RIGHT: - dX=0; - dY=-1.; - break; - default: - cout << "bad quadrant" << endl; - dX=0.; - dY=0.; - } - - - if (nSubPixelsX>2 || nSubPixelsY>2 ) { - - ex=(etax-etamin)/etastepX; - ey=(etay-etamin)/etastepY; - if (ex<0) { - cout << "x*"<< ex << endl; - ex=0; - } - if (ex>=nbetaX) { - cout << "x?"<< ex << endl; - ex=nbetaX-1; - } - if (ey<0) { - cout << "y*"<< ey << " " << nbetaY << endl; - ey=0; - } - if (ey>=nbetaY) { - cout << "y?"<< ey << " " << nbetaY << endl; - ey=nbetaY-1; - } - - - - xpos_eta=(((double)hhx[(ey*nbetaX+ex)]))+dX ;///((double)nSubPixels); - ypos_eta=(((double)hhy[(ey*nbetaX+ex)]))+dY ;///((double)nSubPixels); - - } else { - xpos_eta=0.5*dX+0.25; - ypos_eta=0.5*dY+0.25; - } - - int_x=((double)x) + xpos_eta+0.5; - int_y=((double)y) + ypos_eta+0.5; - - - } - - - - virtual int addToFlatField(double totquad,int quad,int *cl,double &etax, double &etay) { - double cc[2][2]; - int xoff=0, yoff=0; - - 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]=cl[xoff+3*yoff]; - cc[1][0]=cl[xoff+3*(yoff+1)]; - cc[0][1]=cl[xoff+1+3*yoff]; - cc[1][1]=cl[xoff+1+3*(yoff+1)]; - - //calcMyEta(totquad,quad,cl,etax, etay); - calcEta(totquad, cc,etax, etay); - - // cout <<"******"<< etax << " " << etay << endl; - - - return addToFlatField(etax,etay); - } - - virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) { - double cc[2][2]; - int xoff=0, yoff=0; - - 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]=cl[xoff+3*yoff]; - cc[1][0]=cl[(yoff+1)*3+xoff]; - cc[0][1]=cl[yoff*3+xoff+1]; - cc[1][1]=cl[(yoff+1)*3+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 <<"******"<Fill(etax,etay); -#endif -#ifndef MYROOT1 - int ex,ey; - ex=(etax-etamin)/etastepX; - ey=(etay-etamin)/etastepY; - if (ey=0 && ey>=0) - heta[ey*nbetaX+ex]++; -#endif - return 0; - }; - - virtual int *getInterpolatedImage(){ - int ipx, ipy; - // cout << "ff" << endl; - calcDiff(1, hhx, hhy); //get flat - double avg=0; - for (ipx=0; ipx0) - hintcorr[ibx+iby*nSubPixelsX*nPixelsX]=hint[ibx+iby*nSubPixelsX*nPixelsX]*(avg/flat[ipx+ipy*nSubPixelsX]); - else - hintcorr[ibx+iby*nSubPixelsX*nPixelsX]=hint[ibx+iby*nSubPixelsX*nPixelsX]; - - - } - } - - - - return hintcorr; - }; - -/* protected: */ - -/* #ifdef MYROOT1 */ -/* TH2D *heta; */ -/* TH2D *hhx; */ -/* TH2D *hhy; */ -/* #endif */ -/* #ifndef MYROOT1 */ -/* int *heta; */ -/* float *hhx; */ -/* float *hhy; */ -/* #endif */ -/* int nbeta; */ -/* double etamin, etamax, etastep; */ - -}; - -#endif diff --git a/slsDetectorCalibration/interpolations/eta3InterpolationBase.h b/slsDetectorCalibration/interpolations/eta3InterpolationBase.h deleted file mode 100644 index a00969bb2..000000000 --- a/slsDetectorCalibration/interpolations/eta3InterpolationBase.h +++ /dev/null @@ -1,271 +0,0 @@ -#ifndef ETA3_INTERPOLATION_BASE_H -#define ETA3_INTERPOLATION_BASE_H - -#ifdef MYROOT1 -#include -#include -#include -#include -#endif - -#include "etaInterpolationBase.h" - -class eta3InterpolationBase : public virtual etaInterpolationBase { - - public: - eta3InterpolationBase(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : etaInterpolationBase(nx, ny, ns, nsy, nb, nby, emin, emax) { - // cout << "e3ib " << nb << " " << emin << " " << emax << endl; - /* if (nbeta<=0) { */ - /* nbeta=nSubPixels*10; */ - /* } */ - - // cout << nbeta << " " << etamin << " " << etamax << endl; -}; - - eta3InterpolationBase(eta3InterpolationBase *orig): etaInterpolationBase(orig){ }; - - /* virtual eta3InterpolationBase* Clone()=0; */ - - - - - // virtual void prepareInterpolation(int &ok){}; - - - ////////////////////////////////////////////////////////////////////////////// - //////////// /*It return position hit for the event in input */ ////////////// - virtual void getInterpolatedPosition(int x, int y, int *data, double &int_x, double &int_y) - { - double tot; - double etax,etay; - - int corner=calcEta3(data,etax,etay, tot); - - getInterpolatedPosition(x,y,etax,etay,corner,int_x,int_y); - - return; - }; - - - virtual void getInterpolatedPosition(int x, int y, double *data, double &int_x, double &int_y) - { - //double sDum[2][2]; - double tot; - double etax,etay; - - int corner=calcEta3(data,etax,etay, tot); - - 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 etax, etay; - if (nSubPixelsX>2 || nSubPixelsY>2 ) { - calcEta3(cl,etax,etay, totquad); - } - return getInterpolatedPosition(x,y,etax, etay,quad,int_x,int_y); - - } - - - - virtual void getInterpolatedPosition(int x, int y, double totquad,int quad,int *cl,double &int_x, double &int_y) { - - - double etax, etay; - if (nSubPixelsX>2 || nSubPixelsY>2 ) { - calcEta3(cl,etax,etay, totquad); - } - 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=0,ypos_eta=0; - int ex,ey; - - if (nSubPixelsX>2 || nSubPixelsY>2 ) { - -#ifdef MYROOT1 - xpos_eta=(hhx->GetBinContent(hhx->GetXaxis()->FindBin(etax),hhy->GetYaxis()->FindBin(etay)))/((double)nSubPixelsX); - ypos_eta=(hhy->GetBinContent(hhx->GetXaxis()->FindBin(etax),hhy->GetYaxis()->FindBin(etay)))/((double)nSubPixelsY); -#endif -#ifndef MYROOT1 - ex=(etax-etamin)/etastepX; - ey=(etay-etamin)/etastepY; - if (ex<0) { - /* cout << etax << " " << etamin << " "; */ - /* cout << "3x*"<< ex << endl; */ - ex=0; - } - if (ex>=nbetaX) { - /* cout << etax << " " << etamin << " "; */ - /* cout << "3x?"<< ex << endl; */ - ex=nbetaX-1; - } - if (ey<0) { - /* cout << etay << " " << etamin << " "; */ - /* cout << "3y*"<< ey << endl; */ - ey=0; - } - if (ey>=nbetaY) { - /* cout << etay << " " << etamin << " "; */ - /* cout << "3y?"<< ey << endl; */ - ey=nbetaY-1; - - } - xpos_eta=(((double)hhx[(ey*nbetaX+ex)]));///((double)nSubPixels); - ypos_eta=(((double)hhy[(ey*nbetaX+ex)]));///((double)nSubPixels); - -#endif - - } else { - switch (corner) { - case BOTTOM_LEFT: - xpos_eta=-0.25; - ypos_eta=-0.25; - break; - case BOTTOM_RIGHT: - xpos_eta=0.25; - ypos_eta=-0.25; - break; - case TOP_LEFT: - xpos_eta=-0.25; - ypos_eta=0.25; - break; - case TOP_RIGHT: - xpos_eta=0.25; - ypos_eta=0.25; - break; - default: - xpos_eta=0; - ypos_eta=0; - } - - } - int_x=((double)x) + xpos_eta; - int_y=((double)y) + ypos_eta; - // int_x=5. + xpos_eta; - // int_y=5. + ypos_eta; - - - } - - - - - - -/* ///////////////////////////////////////////////////////////////////////////////////////////////// */ -/* virtual void getPositionETA3(int x, int y, double *data, double &int_x, double &int_y) */ -/* { */ -/* double sDum[2][2]; */ -/* double tot, totquad; */ -/* double eta3x,eta3y; */ -/* double ex,ey; */ - -/* calcQuad(data, tot, totquad, sDum); */ -/* calcEta3(data,eta3x, eta3y,tot); */ - -/* double xpos_eta,ypos_eta; */ - -/* #ifdef MYROOT1 */ -/* xpos_eta=((hhx->GetBinContent(hhx->GetXaxis()->FindBin(eta3x),hhy->GetYaxis()->FindBin(eta3y))))/((double)nSubPixels); */ -/* ypos_eta=((hhy->GetBinContent(hhx->GetXaxis()->FindBin(eta3x),hhy->GetYaxis()->FindBin(eta3y))))/((double)nSubPixels); */ - -/* #endif */ -/* #ifndef MYROOT1 */ -/* ex=(eta3x-etamin)/etastep; */ -/* ey=(eta3y-etamin)/etastep; */ - -/* if (ex<0) ex=0; */ -/* if (ex>=nbeta) ex=nbeta-1; */ -/* if (ey<0) ey=0; */ -/* if (ey>=nbeta) ey=nbeta-1; */ - -/* xpos_eta=(((double)hhx[(int)(ey*nbeta+ex)]))/((double)nSubPixels); */ -/* ypos_eta=(((double)hhy[(int)(ey*nbeta+ex)]))/((double)nSubPixels); */ -/* #endif */ - -/* int_x=((double)x) + xpos_eta; */ -/* int_y=((double)y) + ypos_eta; */ - -/* return; */ -/* }; */ - - virtual int addToFlatField(double totquad,int quad,int *cl,double &etax, double &etay) { - - calcEta3(cl, etax, etay, totquad); - return addToFlatField(etax,etay); - } - - virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) { - - - calcEta3(cl, etax, etay, totquad); - return addToFlatField(etax,etay); - } - - - - ////////////////////////////////////////////////////////////////////////////////////// - virtual int addToFlatField(double *cluster, double &etax, double &etay){ - double totquad; - calcEta3(cluster, etax, etay, totquad); - return addToFlatField(etax,etay); - - }; - - virtual int addToFlatField(int *cluster, double &etax, double &etay){ - - double totquad; - - calcEta3(cluster, etax, etay, totquad); - return addToFlatField(etax,etay); - - }; - - - virtual int addToFlatField(double etax, double etay){ -#ifdef MYROOT1 - heta->Fill(etax,etay); -#endif -#ifndef MYROOT1 - int ex,ey; - ex=(etax-etamin)/etastepX; - ey=(etay-etamin)/etastepY; - if (ey=0 && ey>=0) - heta[ey*nbetaX+ex]++; -#endif - return 0; - }; - - -/* protected: */ - -/* #ifdef MYROOT1 */ -/* TH2D *heta; */ -/* TH2D *hhx; */ -/* TH2D *hhy; */ -/* #endif */ -/* #ifndef MYROOT1 */ -/* int *heta; */ -/* float *hhx; */ -/* float *hhy; */ -/* #endif */ -/* int nbeta; */ -/* double etamin, etamax, etastep; */ - -}; - -#endif diff --git a/slsDetectorCalibration/interpolations/etaInterpolationBase.h b/slsDetectorCalibration/interpolations/etaInterpolationBase.h index 52273aab0..df3679aba 100644 --- a/slsDetectorCalibration/interpolations/etaInterpolationBase.h +++ b/slsDetectorCalibration/interpolations/etaInterpolationBase.h @@ -12,6 +12,19 @@ #include "tiffIO.h" class etaInterpolationBase : public slsInterpolation { + + protected: + float *hhx; + float *hhy; + int *heta; + int nbetaX, nbetaY; + double etamin, etamax, etastepX, etastepY; + double rangeMin, rangeMax; + + + + double *flat; + int *hintcorr; public: @@ -190,19 +203,10 @@ float *gethhx() // hhy->Scale((double)nSubPixels); return hhy; }; - virtual int addToFlatField(double etax, double etay){ - int ex,ey; - ex=(etax-etamin)/etastepX; - ey=(etay-etamin)/etastepY; - if (ey=0 && ey>=0) - heta[ey*nbetaX+ex]++; - return 0; - }; - - - // virtual void prepareInterpolation(int &ok)=0; + + void debugSaveAll(int ind=0) { int ibx, iby; char tit[10000]; @@ -223,15 +227,15 @@ float *gethhx() for (int ii=0; ii0) */ + /* hintcorr[ibx+iby*nSubPixelsX*nPixelsX]=hint[ibx+iby*nSubPixelsX*nPixelsX]*(avg/flat[ipx+ipy*nSubPixelsX]); */ + /* else */ + /* hintcorr[ibx+iby*nSubPixelsX*nPixelsX]=hint[ibx+iby*nSubPixelsX*nPixelsX]; */ + /* } */ + /* } */ + + + return hint; + /* return hintcorr; */ + }; + + protected: @@ -363,19 +398,575 @@ float *gethhx() return sqrt(diff); } - float *hhx; - float *hhy; - int *heta; - int nbetaX, nbetaY; - double etamin, etamax, etastepX, etastepY; - double rangeMin, rangeMax; - - - - double *flat; - int *hintcorr; - }; + + + +class eta2InterpolationBase : public virtual etaInterpolationBase { + + public: + eta2InterpolationBase(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : etaInterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax) { + + + + }; + + eta2InterpolationBase(eta2InterpolationBase *orig): etaInterpolationBase(orig){ }; + + + ////////////////////////////////////////////////////////////////////////////// + //////////// /*It return position hit for the event in input */ ////////////// + virtual void getInterpolatedPosition(int x, int y, int *data, double &int_x, double &int_y) + { + double sDum[2][2]; + double tot, totquad; + double etax=0,etay=0; + + int corner; + corner=calcQuad(data, tot, totquad, sDum); + if (nSubPixelsX>2 || nSubPixelsY>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 *data, double &int_x, double &int_y) + { + double sDum[2][2]; + double tot, totquad; + double etax=0,etay=0; + + int corner; + corner=calcQuad(data, tot, totquad, sDum); + if (nSubPixelsX>2 || nSubPixelsY>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]; + int xoff=0, yoff=0; + 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=0, etay=0; + if (nSubPixelsX>2 || nSubPixelsY>2) { + cc[0][0]=cl[xoff+3*yoff]; + cc[1][0]=cl[xoff+3*(yoff+1)]; + cc[0][1]=cl[xoff+1+3*yoff]; + cc[1][1]=cl[xoff+1+3*(yoff+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 totquad,int quad,int *cl,double &int_x, double &int_y) { + + double cc[2][2]; + int xoff=0, yoff=0; + + 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=0, etay=0; + if (nSubPixelsX>2 || nSubPixelsY>2) { + cc[0][0]=cl[xoff+3*yoff]; + cc[1][0]=cl[xoff+3*(yoff+1)]; + cc[0][1]=cl[xoff+1+3*yoff]; + cc[1][1]=cl[xoff+1+3*(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=0,ypos_eta=0; + double dX,dY; + int ex,ey; + switch (corner) + { + case TOP_LEFT: + dX=-1.; + dY=0; + break; + case TOP_RIGHT: + ; + dX=0; + dY=0; + break; + case BOTTOM_LEFT: + dX=-1.; + dY=-1.; + break; + case BOTTOM_RIGHT: + dX=0; + dY=-1.; + break; + default: + cout << "bad quadrant" << endl; + dX=0.; + dY=0.; + } + + + if (nSubPixelsX>2 || nSubPixelsY>2 ) { + + ex=(etax-etamin)/etastepX; + ey=(etay-etamin)/etastepY; + if (ex<0) { + cout << "x*"<< ex << endl; + ex=0; + } + if (ex>=nbetaX) { + cout << "x?"<< ex << endl; + ex=nbetaX-1; + } + if (ey<0) { + cout << "y*"<< ey << " " << nbetaY << endl; + ey=0; + } + if (ey>=nbetaY) { + cout << "y?"<< ey << " " << nbetaY << endl; + ey=nbetaY-1; + } + + + + xpos_eta=(((double)hhx[(ey*nbetaX+ex)])); + ypos_eta=(((double)hhy[(ey*nbetaX+ex)])); + + + if (xpos_eta<0 || xpos_eta>1) + xpos_eta=-100; + else + xpos_eta+=dX ;///((double)nSubPixels); + + if (ypos_eta<0 || ypos_eta>1) + ypos_eta=-100; + else + ypos_eta+=dY ;///((double)nSubPixels); + + } else { + xpos_eta=0.5*dX+0.25; + ypos_eta=0.5*dY+0.25; + } + + + if (xpos_eta<-1) + int_x=-100; + else + int_x=((double)x) + xpos_eta+0.5; + + if (ypos_eta<-1) + int_y=-100; + else + int_y=((double)y) + ypos_eta+0.5; + + + + + } + + + + virtual int addToFlatField(double totquad,int quad,int *cl,double &etax, double &etay) { + double cc[2][2]; + int xoff=0, yoff=0; + + 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]=cl[xoff+3*yoff]; + cc[1][0]=cl[xoff+3*(yoff+1)]; + cc[0][1]=cl[xoff+1+3*yoff]; + cc[1][1]=cl[xoff+1+3*(yoff+1)]; + + //calcMyEta(totquad,quad,cl,etax, etay); + calcEta(totquad, cc,etax, etay); + + // cout <<"******"<< etax << " " << etay << endl; + + + return addToFlatField(etax,etay); + } + + virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) { + double cc[2][2]; + int xoff=0, yoff=0; + + 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]=cl[xoff+3*yoff]; + cc[1][0]=cl[(yoff+1)*3+xoff]; + cc[0][1]=cl[yoff*3+xoff+1]; + cc[1][1]=cl[(yoff+1)*3+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 <<"******"<Fill(etax,etay); +#endif +#ifndef MYROOT1 + int ex,ey; + ex=(etax-etamin)/etastepX; + ey=(etay-etamin)/etastepY; + if (ey=0 && ey>=0) + heta[ey*nbetaX+ex]++; +#endif + return 0; + }; + +}; + + + + + + + + + + + + + + + + + +class eta3InterpolationBase : public virtual etaInterpolationBase { + + public: + eta3InterpolationBase(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : etaInterpolationBase(nx, ny, ns, nsy, nb, nby, emin, emax) { + cout << "e3ib " << nb << " " << emin << " " << emax << endl; + +}; + + eta3InterpolationBase(eta3InterpolationBase *orig): etaInterpolationBase(orig){ }; + + + ////////////////////////////////////////////////////////////////////////////// + //////////// /*It return position hit for the event in input */ ////////////// + virtual void getInterpolatedPosition(int x, int y, int *data, double &int_x, double &int_y) + { + double tot; + double etax,etay; + + int corner=calcEta3(data,etax,etay, tot); + + getInterpolatedPosition(x,y,etax,etay,corner,int_x,int_y); + + return; + }; + + + virtual void getInterpolatedPosition(int x, int y, double *data, double &int_x, double &int_y) + { + //double sDum[2][2]; + double tot; + double etax,etay; + + int corner=calcEta3(data,etax,etay, tot); + + 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 etax, etay; + if (nSubPixelsX>2 || nSubPixelsY>2 ) { + calcEta3(cl,etax,etay, totquad); + } + return getInterpolatedPosition(x,y,etax, etay,quad,int_x,int_y); + + } + + + + virtual void getInterpolatedPosition(int x, int y, double totquad,int quad,int *cl,double &int_x, double &int_y) { + + + double etax, etay; + if (nSubPixelsX>2 || nSubPixelsY>2 ) { + calcEta3(cl,etax,etay, totquad); + } + 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) + { + + + //cout << "**"; + double xpos_eta=0,ypos_eta=0; + int ex,ey; + + if (nSubPixelsX>2 || nSubPixelsY>2 ) { + +#ifdef MYROOT1 + xpos_eta=(hhx->GetBinContent(hhx->GetXaxis()->FindBin(etax),hhy->GetYaxis()->FindBin(etay)))/((double)nSubPixelsX); + ypos_eta=(hhy->GetBinContent(hhx->GetXaxis()->FindBin(etax),hhy->GetYaxis()->FindBin(etay)))/((double)nSubPixelsY); +#endif +#ifndef MYROOT1 + ex=(etax-etamin)/etastepX; + ey=(etay-etamin)/etastepY; + if (ex<0) { + /* cout << etax << " " << etamin << " "; */ + /* cout << "3x*"<< ex << endl; */ + ex=0; + } + if (ex>=nbetaX) { + /* cout << etax << " " << etamin << " "; */ + /* cout << "3x?"<< ex << endl; */ + ex=nbetaX-1; + } + if (ey<0) { + /* cout << etay << " " << etamin << " "; */ + /* cout << "3y*"<< ey << endl; */ + ey=0; + } + if (ey>=nbetaY) { + /* cout << etay << " " << etamin << " "; */ + /* cout << "3y?"<< ey << endl; */ + ey=nbetaY-1; + + } + xpos_eta=(((double)hhx[(ey*nbetaX+ex)]));///((double)nSubPixels); + ypos_eta=(((double)hhy[(ey*nbetaX+ex)]));///((double)nSubPixels); + +#endif + + } else { + switch (corner) { + case BOTTOM_LEFT: + xpos_eta=-0.25; + ypos_eta=-0.25; + break; + case BOTTOM_RIGHT: + xpos_eta=0.25; + ypos_eta=-0.25; + break; + case TOP_LEFT: + xpos_eta=-0.25; + ypos_eta=0.25; + break; + case TOP_RIGHT: + xpos_eta=0.25; + ypos_eta=0.25; + break; + default: + xpos_eta=0; + ypos_eta=0; + } + + } + /* if (xpos_eta<0 || xpos_eta>1) */ + /* int_x=-1; */ + /* else */ + int_x=((double)x) + xpos_eta; + /* if (ypos_eta<0 || ypos_eta>1) */ + /* int_y=-1; */ + /* else */ + int_y=((double)y) + ypos_eta; + // int_x=5. + xpos_eta; + // int_y=5. + ypos_eta; + + + } + + + virtual int addToFlatField(double totquad,int quad,int *cl,double &etax, double &etay) { + + calcEta3(cl, etax, etay, totquad); + return addToFlatField(etax,etay); + } + + virtual int addToFlatField(double totquad,int quad,double *cl,double &etax, double &etay) { + + + calcEta3(cl, etax, etay, totquad); + return addToFlatField(etax,etay); + } + + + + ////////////////////////////////////////////////////////////////////////////////////// + virtual int addToFlatField(double *cluster, double &etax, double &etay){ + double totquad; + calcEta3(cluster, etax, etay, totquad); + return addToFlatField(etax,etay); + + }; + + virtual int addToFlatField(int *cluster, double &etax, double &etay){ + + double totquad; + + calcEta3(cluster, etax, etay, totquad); + return addToFlatField(etax,etay); + + }; + + + + virtual int addToFlatField(double etax, double etay){ +#ifdef MYROOT1 + heta->Fill(etax,etay); +#endif +#ifndef MYROOT1 + int ex,ey; + ex=(etax-etamin)/etastepX; + ey=(etay-etamin)/etastepY; + if (ey=0 && ey>=0) + heta[ey*nbetaX+ex]++; +#endif + return 0; + }; + + +}; + #endif diff --git a/slsDetectorCalibration/interpolations/etaInterpolationGlobal.h b/slsDetectorCalibration/interpolations/etaInterpolationGlobal.h index 615b7c04b..5eec2eac3 100644 --- a/slsDetectorCalibration/interpolations/etaInterpolationGlobal.h +++ b/slsDetectorCalibration/interpolations/etaInterpolationGlobal.h @@ -4,11 +4,11 @@ #include "etaInterpolationBase.h" -class etaInterpolationGlobal : public etaInterpolationBase{ +class etaInterpolationGlobal : public virtual etaInterpolationBase { public: - globalEtaInterpolation(int nx=400, int ny=400, int ns=25, int nb=-1, double emin=1, double emax=0) : etaInterpolationBase(nx,ny,ns, nb, emin,emax){}; + etaInterpolationGlobal(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : etaInterpolationBase(nx,ny,ns,nsy, nb, nby,emin,emax){}; + - virtual void prepareInterpolation(int &ok) { @@ -24,62 +24,106 @@ class etaInterpolationGlobal : public etaInterpolationBase{ ///*Eta Distribution Rebinning*/// - double bsize=1./nSubPixels; //precision + // double bsizeX=1./nSubPixelsX; //precision // cout<<"nPixelsX = "<(ib+1)*tot_eta*bsize) ib++; - for (int iby=0; ibySetBinContent(ibx+1,iby+1,ib); -#endif -#ifndef MYROOT1 - hhx[ibx+iby*nbeta]=ib; -#endif - } - } - ib=0; - for (int iby=0; iby(ib+1)*tot_eta*bsize) ib++; - for (int ibx=0; ibxSetBinContent(ibx+1,iby+1,ib); -#endif -#ifndef MYROOT1 - hhy[ibx+iby*nbeta]=ib; -#endif + + for (int ibx=0; ibx(ib+1)*tot_eta*bsize) ib++; + for (int iby=0; iby0) + hhx[ibx+iby*nbetaX]=hix[ibx]/tot_x; + else + hhx[ibx+iby*nbetaX]=-1; + + if (tot_y>0) + hhy[ibx+iby*nbetaX]=hiy[iby]/tot_y; + else + hhy[ibx+iby*nbetaX]=-1; } } +/* ib=0; */ +/* for (int iby=0; iby(ib+1)*tot_eta*bsize) ib++; */ +/* for (int ibx=0; ibxSetBinContent(ibx+1,iby+1,ib); */ +/* #endif */ +/* #ifndef MYROOT1 */ +/* hhy[ibx+iby*nbeta]=ib; */ +/* #endif */ +/* } */ +/* } */ return ; }; + + etaInterpolationGlobal(etaInterpolationGlobal *orig): etaInterpolationBase(orig) {}; }; + +//etaInterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax), + +class eta2InterpolationGlobal : public virtual eta2InterpolationBase, public virtual etaInterpolationGlobal { + + public: + eta2InterpolationGlobal(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : eta2InterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax), etaInterpolationGlobal(nx,ny, ns, nsy, nb, nby, emin, emax){ + cout << "e2pxy " << nbetaX << " " << nbetaY << etamin << " " << etamax << " " << nSubPixelsX<< " " << nSubPixelsY << endl; + }; + + eta2InterpolationGlobal(eta2InterpolationGlobal *orig): etaInterpolationBase(orig), etaInterpolationGlobal(orig) {}; + + virtual eta2InterpolationGlobal* Clone() { return new eta2InterpolationGlobal(this);}; + +}; + +//etaInterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax), + +class eta3InterpolationGlobal : public virtual eta3InterpolationBase, public virtual etaInterpolationGlobal { + public: + eta3InterpolationGlobal(int nx=400, int ny=400, int ns=25, int nsy=25, int nb=-1, int nby=-1, double emin=1, double emax=0) : eta3InterpolationBase(nx,ny, ns, nsy, nb, nby, emin, emax), etaInterpolationGlobal(nx,ny, ns, nsy, nb, nby, emin, emax){ + cout << "e3pxy " << nbetaX << " " << nbetaY << etamin << " " << etamax << " " << nSubPixelsX<< " " << nSubPixelsY << endl; + }; + + eta3InterpolationGlobal(eta3InterpolationGlobal *orig): etaInterpolationBase(orig), etaInterpolationGlobal(orig) {}; + + virtual eta3InterpolationGlobal* Clone() {return new eta3InterpolationGlobal(this);}; +}; + #endif diff --git a/slsDetectorCalibration/interpolations/etaInterpolationPosXY.h b/slsDetectorCalibration/interpolations/etaInterpolationPosXY.h index 0d29c1658..d6b46d619 100644 --- a/slsDetectorCalibration/interpolations/etaInterpolationPosXY.h +++ b/slsDetectorCalibration/interpolations/etaInterpolationPosXY.h @@ -4,8 +4,6 @@ //#include "tiffIO.h" #include "etaInterpolationBase.h" -#include "eta2InterpolationBase.h" -#include "eta3InterpolationBase.h" class etaInterpolationPosXY : public virtual etaInterpolationBase{ public: @@ -44,9 +42,8 @@ class etaInterpolationPosXY : public virtual etaInterpolationBase{ double hiy[nbetaY]; //integral of projection y // int ii=0; double etax, etay; - for (int ib=0; ib=0 && etay<=1) + // if (etay>=0 && etay<=1) hy[iby]=heta[ib+iby*nbetaX]; - else - hy[iby]=0; + // else + // hy[iby]=0; // tot_eta_y+=hy[iby]; } @@ -71,7 +68,7 @@ class etaInterpolationPosXY : public virtual etaInterpolationBase{ tot_eta_y=hiy[nbetaY-1]+1; for (int iby=0; iby=0 && etax<=1) + // if (etax>=0 && etax<=1) hx[ibx]=heta[ibx+ib*nbetaX]; - else { - hx[ibx]=0; - } + // else { + // hx[ibx]=0; + // } } hix[0]=hx[0]; @@ -104,7 +101,7 @@ class etaInterpolationPosXY : public virtual etaInterpolationBase{ tot_eta_x=hix[nbetaX-1]+1; for (int ibx=0; ibx0) iby--; */ + + /* for (ib=iby+1; ib) iby--; */ + /* for (ib=iby+1; ib0) { diff --git a/slsDetectorCalibration/moenchExecutables/Makefile.cluster_finder_rh7 b/slsDetectorCalibration/moenchExecutables/Makefile.cluster_finder_rh7 index 64cd5d73e..afe472590 100644 --- a/slsDetectorCalibration/moenchExecutables/Makefile.cluster_finder_rh7 +++ b/slsDetectorCalibration/moenchExecutables/Makefile.cluster_finder_rh7 @@ -17,22 +17,22 @@ all: moenchClusterFinder moenchMakeEta moenchInterpolation moenchNoInterpolati moenchClusterFinder : moench03ClusterFinder.cpp $(INCS) clean - g++ -o moenchClusterFinder moench03ClusterFinder.cpp $(LDFLAG) $(INCDIR) -DSAVE_ALL -DNEWRECEIVER + g++ -o moenchClusterFinder moench03ClusterFinder.cpp $(LDFLAG) $(INCDIR) -DSAVE_ALL -DNEWRECEIVER -DWRITE_QUAD moenchMakeEta : moench03Interpolation.cpp $(INCS) clean - g++ -o moenchMakeEta moench03Interpolation.cpp $(LDFLAG) $(INCDIR) $(LIBHDF5) $(LIBRARYCBF) -DFF + g++ -o moenchMakeEta moench03Interpolation.cpp $(LDFLAG) $(INCDIR) $(LIBHDF5) $(LIBRARYCBF) -DFF -DWRITE_QUAD moenchInterpolation : moench03Interpolation.cpp $(INCS) clean - g++ -o moenchInterpolation moench03Interpolation.cpp $(LDFLAG) $(INCDIR) + g++ -o moenchInterpolation moench03Interpolation.cpp $(LDFLAG) $(INCDIR) -DWRITE_QUAD moenchNoInterpolation : moench03NoInterpolation.cpp $(INCS) clean - g++ -o moenchNoInterpolation moench03NoInterpolation.cpp $(LDFLAG) $(INCDIR) + g++ -o moenchNoInterpolation moench03NoInterpolation.cpp $(LDFLAG) $(INCDIR) -DWRITE_QUAD moenchPhotonCounter : moenchPhotonCounter.cpp $(INCS) clean - g++ -o moenchPhotonCounter moenchPhotonCounter.cpp $(LDFLAG) $(INCDIR) $ -DNEWRECEIVER + g++ -o moenchPhotonCounter moenchPhotonCounter.cpp $(LDFLAG) $(INCDIR) $ -DNEWRECEIVER -DWRITE_QUAD moenchAnalog : moenchPhotonCounter.cpp $(INCS) clean - g++ -o moenchAnalog moenchPhotonCounter.cpp $(LDFLAG) $(INCDIR) -DNEWRECEIVER -DANALOG + g++ -o moenchAnalog moenchPhotonCounter.cpp $(LDFLAG) $(INCDIR) -DNEWRECEIVER -DANALOG -DWRITE_QUAD clean: rm -f moenchClusterFinder moenchMakeEta moenchInterpolation moenchNoInterpolation moenchPhotonCounter moenchAnalog diff --git a/slsDetectorCalibration/moenchExecutables/moench03Interpolation.cpp b/slsDetectorCalibration/moenchExecutables/moench03Interpolation.cpp index ca6bf5d4e..5c82988c3 100644 --- a/slsDetectorCalibration/moenchExecutables/moench03Interpolation.cpp +++ b/slsDetectorCalibration/moenchExecutables/moench03Interpolation.cpp @@ -14,7 +14,12 @@ #include "single_photon_hit.h" #endif +#ifdef RECT +#include "etaInterpolationGlobal.h" +#endif +#ifndef RECT #include "etaInterpolationPosXY.h" +#endif #include "noInterpolation.h" //#include "etaInterpolationCleverAdaptiveBins.h" //#include "etaInterpolationRandomBins.h" @@ -36,6 +41,7 @@ using namespace std; int main(int argc, char *argv[]) { #ifndef FF + cout << "INTERPOLATING" << endl; if (argc<9) { cout << "Wrong usage! Should be: "<< argv[0] << " infile etafile outfile runmin runmax ns cmin cmax" << endl; return 1; @@ -43,6 +49,7 @@ int main(int argc, char *argv[]) { #endif #ifdef FF + cout << "CALC ETA" << endl; if (argc<7) { cout << "Wrong usage! Should be: "<< argv[0] << " infile etafile runmin runmax cmin cmax" << endl; return 1; @@ -74,11 +81,11 @@ int main(int argc, char *argv[]) { cout << "Energy min: " << cmin << endl; cout << "Energy max: " << cmax << endl; //int etabins=500; - int etabins=1000;//nsubpix*2*100; + int etabins=999;//nsubpix*2*100; int etabinsY=etabins;//nsubpix*2*100; - double etamin=-1, etamax=2; + double etamin=0, etamax=1; //double etamin=-0.1, etamax=1.1; - double eta3min=-2, eta3max=2; + double eta3min=-1, eta3max=1; int quad; double sum, totquad; double sDum[2][2]; @@ -92,6 +99,9 @@ int main(int argc, char *argv[]) { int nph=0, badph=0, totph=0; FILE *f=NULL; + + char fff[10000]; + #ifdef DOUBLE_SPH single_photon_hit_double cl(3,3); #endif @@ -102,12 +112,23 @@ int main(int argc, char *argv[]) { int nSubPixels=nsubpix; int nSubPixelsY=nsubpix; #ifdef RECT - nSubPixels=1; //might make more sense using 2? - etabins=2; + nSubPixels=5; //might make more sense using 2? + // etabins=5; #endif #ifndef NOINTERPOLATION - eta2InterpolationPosXY *interp=new eta2InterpolationPosXY(NC, NR, nSubPixels, nSubPixelsY, etabins, etabinsY, etamin, etamax); + etaInterpolationBase *interp; +#ifndef RECT + interp=new eta2InterpolationPosXY(NC, NR, nSubPixels, nSubPixelsY, etabins, etabinsY, etamin, etamax); //eta2InterpolationCleverAdaptiveBins *interp=new eta2InterpolationCleverAdaptiveBins(NC, NR, nsubpix, etabins, etamin, etamax); +#endif +#ifdef RECT + interp=new eta2InterpolationGlobal(NC, NR, nSubPixels, nSubPixelsY, etabins, etabinsY, eta3min, eta3max); + + // eta3InterpolationPosXY *interpOdd=new eta3InterpolationPosXY(NC, NR, nSubPixels, nSubPixelsY, etabins, etabinsY, eta3min, eta3max); + //eta2InterpolationCleverAdaptiveBins *interp=new eta2InterpolationCleverAdaptiveBins(NC, NR, nsubpix, etabins, etamin, etamax); + +#endif + #endif #ifdef NOINTERPOLATION noInterpolation *interp=new noInterpolation(NC, NR, nsubpix, nSubPixelsY); @@ -119,8 +140,21 @@ int main(int argc, char *argv[]) { #ifndef NOINTERPOLATION cout << "read ff " << argv[2] << endl; sprintf(fname,"%s",argv[2]); + //#ifndef RECT interp->readFlatField(fname, etamin, etamax); interp->prepareInterpolation(ok);//, MAX_ITERATIONS); + interp->debugSaveAll(0); + //#endif +// #ifdef RECT +// sprintf(fff,"%s.even",fname); +// interpEven->readFlatField(fff, etamin, etamax); +// interpEven->prepareInterpolation(ok);, MAX_ITERATIONS); +// interpEven->debugSaveAll(1); +// sprintf(fff,"%s.odd",fname); +// interpOdd->readFlatField(fff, etamin, etamax); +// interpOdd->prepareInterpolation(ok);, MAX_ITERATIONS); +// interpOdd->debugSaveAll(2); +// #endif #endif // return 0; #endif @@ -143,6 +177,9 @@ int main(int argc, char *argv[]) { #ifdef FF sprintf(outfname,argv[2]); #endif + + + int tl=0, tr=0, bl=0, br=0; int irun; for (irun=runmin; iruncalcEta3(cl.get_cluster(), etax, etay, sum); quad=interp->calcEta(cl.get_cluster(), etax, etay, sum, totquad, sDum); - if (sum>cmin && totquad/sum>0.8 && totquad/sum<1.2 && sum0) + // interp=interpEven; + // else + // interp=interpOdd; + // else + // if (etay>0) + // interp=interpOdd; + // else + // interp=interpEven; + +#endif + +#ifndef RECT + quad=interp->calcEta(cl.get_cluster(), etax, etay, sum, totquad, sDum); + + +#endif + switch(quad) { + case TOP_LEFT: + tl++; + break; + case TOP_RIGHT: + tr++; + break; + case BOTTOM_LEFT: + bl++; + break; + case BOTTOM_RIGHT: + br++; + break; + } + // if (sum>cmin && totquad/sum>0.8 && totquad/sum<1.2 && sumcmin && sum200 && sum<580) { - // interp->getInterpolatedPosition(cl.x,cl.y, totquad,quad,cl.get_cluster(),int_x, int_y); -// #ifdef SOLEIL -// if (cl.x>210 && cl.x<240 && cl.y>210 && cl.y<240) { -// #endif -#ifndef FF - // interp->getInterpolatedPosition(cl.x,cl.y, cl.get_cluster(),int_x, int_y); +#ifndef FF interp->getInterpolatedPosition(cl.x,cl.y, etax, etay, quad,int_x, int_y); - // cout <<"**************"<< endl; - // cout << cl.x << " " << cl.y << " " << sum << endl; - // cl.print(); - // cout << int_x << " " << int_y << endl; - // cout <<"**************"<< endl; - // if (etax!=0 && etay!=0 && etax!=1 && etay!=1) interp->addToImage(int_x, int_y); - if (int_x<0 || int_y<0 || int_x>NC || int_y>NR) { - cout <<"**************"<< endl; - cout << cl.x << " " << cl.y << " " << sum << endl; - cl.print(); - cout << int_x << " " << int_y << endl; - cout <<"**************"<< endl; - } #endif #ifdef FF - // interp->addToFlatField(cl.get_cluster(), etax, etay); -// #ifdef UCL -// if (cl.x>50) -// #endif -// if (etax!=0 && etay!=0 && etax!=1 && etay!=1) interp->addToFlatField(etax, etay); - // if (etax==0 || etay==0) cout << cl.x << " " << cl.y << endl; - #endif -// #ifdef SOLEIL -// } -// #endif if (nph%1000000==0) cout << nph << endl; if (nph%10000000==0) { #ifndef FF - interp->writeInterpolatedImage(outfname); + + cout << "writing interpolated iamge" << endl; + //#ifndef RECT + interp->writeInterpolatedImage(outfname); + //#endif +// #ifdef RECT +// sprintf(fff,"%s.even",outfname); +// interpEven->writeInterpolatedImage(fff); +// sprintf(fff,"%s.odd",outfname); +// interpOdd->writeInterpolatedImage(fff); +// #endif + cout << "done" << endl; #endif #ifdef FF - interp->writeFlatField(outfname); + cout << "writing eta" << endl; + //#ifndef RECT + interp->writeFlatField(outfname); + //#endif +// #ifdef RECT +// sprintf(fff,"%s.even",outfname); +// interpEven->writeFlatField(fff); +// sprintf(fff,"%s.odd",outfname); +// interpOdd->writeFlatField(fff); +// #endif + cout << "done" << endl; #endif } @@ -221,24 +285,64 @@ int main(int argc, char *argv[]) { fclose(f); #ifdef FF + //#ifndef RECT interp->writeFlatField(outfname); +// #endif +// #ifdef RECT +// sprintf(fff,"%s.even",outfname); +// interpEven->writeFlatField(outfname); +// sprintf(fff,"%s.odd",outfname); +// interpOdd->writeFlatField(outfname); +// #endif #endif #ifndef FF + //#ifndef RECT interp->writeInterpolatedImage(outfname); +// #endif +// #ifdef RECT +// sprintf(fff,"%s.even",outfname); +// interpEven->writeInterpolatedImage(fff); +// sprintf(fff,"%s.odd",outfname); +// interpOdd->writeInterpolatedImage(fff); +// #endif +//#ifdef RECT + //interp=interpEven; + //#endif img=interp->getInterpolatedImage(); for (ix=0; ixclearInterpolatedImage(); +// #endif + + #endif } else @@ -250,7 +354,15 @@ int main(int argc, char *argv[]) { #endif #ifdef FF + //#ifndef RECT interp->writeFlatField(outfname); + //#endif +// #ifdef RECT +// sprintf(fff,"%s.even",outfname); +// interpEven->writeFlatField(fff); +// sprintf(fff,"%s.odd",outfname); +// interpOdd->writeFlatField(fff); +// #endif #endif cout << "Filled " << nph << " (/"<< totph <<") " << endl;