diff --git a/slsDetectorCalibration/interpolations/etaInterpolationRosenblatt.h b/slsDetectorCalibration/interpolations/etaInterpolationRosenblatt.h new file mode 100644 index 000000000..24156064e --- /dev/null +++ b/slsDetectorCalibration/interpolations/etaInterpolationRosenblatt.h @@ -0,0 +1,230 @@ +// SPDX-License-Identifier: LGPL-3.0-or-other +// Copyright (C) 2021 Contributors to the SLS Detector Package +#ifndef ETA_INTERPOLATION_POSXY_H +#define ETA_INTERPOLATION_POSXY_H + +//#include "sls/tiffIO.h" +#include "eta2InterpolationBase.h" +#include "eta3InterpolationBase.h" +#include "etaInterpolationBase.h" + +class etaInterpolationRosenblatt : public virtual etaInterpolationBase { + public: + etaInterpolationRosenblatt(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){ + // std::cout << "epxy " << nb << " " << emin << " " << emax << std::endl; + // std::cout << nbeta << " " << etamin << " " << etamax << std::endl; + }; + + etaInterpolationRosenblatt(etaInterpolationRosenblatt *orig) + : etaInterpolationBase(orig){}; + + virtual etaInterpolationRosenblatt *Clone() = 0; /* { */ + + /* return new etaInterpolationPosXY(this); */ + + /* }; */ + + virtual void prepareInterpolation(int &ok) { + ok = 1; + + ///*Eta Distribution Rebinning*/// + // double bsize=1./nSubPixels; //precision + // std::cout<<"nPixelsX = "<= 0 && etay <= 1) + hy[iby] = heta[ib + iby * nbetaX]; + else + hy[iby] = 0; + // tot_eta_y+=hy[iby]; + } + + hiy[0] = hy[0]; + + for (int iby = 1; iby < nbetaY; iby++) { + hiy[iby] = hiy[iby - 1] + hy[iby]; + } + + tot_eta_y = hiy[nbetaY - 1] + 1; + + for (int iby = 0; iby < nbetaY; iby++) { + if (tot_eta_y <= 0) { + hhy[ib + iby * nbetaX] = -1; + // ii=(ibx*nSubPixels)/nbeta; + } else { + // if (hiy[ibx]>tot_eta_y*(ii+1)/nSubPixels) ii++; + hhy[ib + iby * nbetaX] = hiy[iby] / tot_eta_y; + } + } + } + + for (int ib = 0; ib < nbetaY; ib++) { + + for (int ibx = 0; ibx < nbetaX; ibx++) { + etax = etamin + ibx * etastepX; + // std::cout << etax << std::endl; + if (etax >= 0 && etax <= 1) + hx[ibx] = heta[ibx + ib * nbetaX]; + else { + hx[ibx] = 0; + } + } + hix[0] = hx[0]; + + for (int ibx = 1; ibx < nbetaX; ibx++) { + hix[ibx] = hix[ibx - 1] + hx[ibx]; + } + + // tot_eta_x = hix[nbetaX - 1] + 1; + + for (int ibx = 0; ibx < nbetaX; ibx++) { + //if (tot_eta_x <= 0) { + // hhx[ibx + ib * nbetaX] = -1; + // } else { + hhx[ibx + ib * nbetaX] = hix[ibx];// / tot_eta_x; + //} + } + + + } + + + for (int ib = 0; ib < nbetaY; ib++) { + int val=0; + for (int ibx = 0; ibx < nbetaX; ibx++) { + val+=hhx[ibx + ib * nbetaX]; + } + for (int ibx = 0; ibx < nbetaX; ibx++) { + hhx[ibx + ib * nbetaX]=val; + } + } + tot_eta_x = hhx[nbetaX - 1] + 1; + for (int ib = 0; ib < nbetaY; ib++) { + + //tot_eta_x = hix[nbetaX - 1] + 1; + + for (int ibx = 0; ibx < nbetaX; ibx++) { + if (tot_eta_x <= 0) { + hhx[ibx + ib * nbetaX] = -1; + } else { + hhx[ibx + ib * nbetaX] = hhx[ibx + ib * nbetaX] / tot_eta_x; + } + } + } + + int ibx, iby, ib; + + iby = 0; + while (hhx[iby * nbetaY + nbetaY / 2] < 0) + iby++; + for (ib = 0; ib < iby; ib++) { + for (ibx = 0; ibx < nbetaX; ibx++) + hhx[ibx + nbetaX * ib] = hhx[ibx + nbetaX * iby]; + } + iby = nbetaY - 1; + + while (hhx[iby * nbetaY + nbetaY / 2] < 0) + iby--; + for (ib = iby + 1; ib < nbetaY; ib++) { + for (ibx = 0; ibx < nbetaX; ibx++) + hhx[ibx + nbetaX * ib] = hhx[ibx + nbetaX * iby]; + } + + iby = 0; + while (hhy[nbetaX / 2 * nbetaX + iby] < 0) + iby++; + for (ib = 0; ib < iby; ib++) { + for (ibx = 0; ibx < nbetaY; ibx++) + hhy[ib + nbetaX * ibx] = hhy[iby + nbetaX * ibx]; + } + iby = nbetaX - 1; + + while (hhy[nbetaX / 2 * nbetaX + iby] < 0) + iby--; + for (ib = iby + 1; ib < nbetaX; ib++) { + for (ibx = 0; ibx < nbetaY; ibx++) + hhy[ib + nbetaX * ibx] = hhy[iby + nbetaX * ibx]; + } + +#ifdef SAVE_ALL + debugSaveAll(); +#endif + delete[] hx; + delete[] hy; + delete[] hix; + delete[] hiy; + + return; + } +}; + +class eta2InterpolationRosenblatt : public virtual eta2InterpolationBase, + public virtual etaInterpolationRosenblatt { + public: + eta2InterpolationRosenblatt(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(nx, ny, ns, nsy, nb, nby, emin, emax), + etaInterpolationRosenblatt(nx, ny, ns, nsy, nb, nby, emin, emax){ + // std::cout << "e2pxy " << nb << " " << emin << " " << emax << std::endl; + }; + + eta2InterpolationRosenblatt(eta2InterpolationRosenblatt *orig) + : etaInterpolationBase(orig), etaInterpolationRosenblatt(orig){}; + + virtual eta2InterpolationRosenblatt *Clone() { + return new eta2InterpolationRosenblatt(this); + }; +}; + +class eta3InterpolationRosenblatt : public virtual eta3InterpolationBase, + public virtual etaInterpolationRosenblatt { + public: + eta3InterpolationRosenblatt(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), + eta3InterpolationBase(nx, ny, ns, nsy, nb, nby, emin, emax), + etaInterpolationRosenblatt(nx, ny, ns, nsy, nb, nby, emin, emax){ + // std::cout << "e3pxy " << nbeta << " " << etamin << " " << etamax + // << " " << nSubPixels<< std::endl; + }; + + eta3InterpolationRosenblatt(eta3InterpolationRosenblatt *orig) + : etaInterpolationBase(orig), etaInterpolationRosenblatt(orig){}; + + virtual eta3InterpolationRosenblatt *Clone() { + return new eta3InterpolationRosenblatt(this); + }; +}; + +#endif