CXRadReader.h 8.34 KB
Newer Older
1 2 3
#ifndef CXRadReader_H
#define CXRadReader_H

Jérémie Dudouet's avatar
Jérémie Dudouet committed
4 5 6 7
#include <vector>
#include <signal.h>

#include "Riostream.h"
8 9 10 11 12 13
#include "TString.h"

#define RW_MAXCH 8192   /* max number of channels per dimension */
#define MAXCHS 8192    /*    (must be multiple of 8)    */
#define MAXGAM  1500    /* max number of gammas in scheme */

Jérémie Dudouet's avatar
Jérémie Dudouet committed
14 15 16 17 18
using namespace std;

class TH2;
class TH1;
class TH1D;
19 20 21 22 23 24
class TF1;

class CXRadReader
{
private:

Jérémie Dudouet's avatar
Jérémie Dudouet committed
25 26 27 28 29 30 31 32 33 34 35
    TString fName = "";

    TString fCubeFileName       = "None";
    TString f2DProjFileName     = "None";
    TString fLUTFileName        = "None";
    TString fECalFileName       = "None";
    TString fEffFileName        = "None";
    TString fTotalProjFileName  = "None";
    TString fBackgroundFileName = "None";

    //compression factor
36
    double_t fCompFactor           = 1;
Jérémie Dudouet's avatar
Jérémie Dudouet committed
37

38 39 40 41 42
    // Calibs
    double fGains[6];
    double feff_pars[10];

    // basic histograms
Jérémie Dudouet's avatar
Jérémie Dudouet committed
43 44 45 46
    TH1 *hLUTHist      = nullptr;
    TH1D *hTotalProj   = nullptr;
    TH1*hBackground    = nullptr;
    TH2 *h2DProj       = nullptr;
47 48 49 50 51

    // calib functions
    TF1 *fECalibFunc    = nullptr;
    TF1 *fEffFunc       = nullptr;

Jérémie Dudouet's avatar
Jérémie Dudouet committed
52 53 54
    // Gated spectra
    TH1 *fGatedSpectrum   = nullptr;
    TH1 *fGatedSpectrumNoBG   = nullptr;
55 56

    // Binning of histograms
Jérémie Dudouet's avatar
Jérémie Dudouet committed
57 58
    double_t edges[RW_MAXCH+1];

59 60 61
    double_t fBGLenghtFactor = 1.;
    double_t fDefaultWidth = 6.;
    double_t fFWHMPars[3] = {9.0,0.004,0.0};
62 63 64 65 66

public:
    CXRadReader();
    ~CXRadReader();

Jérémie Dudouet's avatar
Jérémie Dudouet committed
67 68 69 70 71 72 73 74 75
    TString GetName(){return fName;}

    TString GetCubeFileName(){return fCubeFileName;}
    TString Get2DProjFileName(){return f2DProjFileName;}
    TString GetLUTFileName(){return fLUTFileName;}
    TString GetECalFileName(){return fECalFileName;}
    TString GetEffFileName(){return fEffFileName;}
    TString GetTotalProjFileName(){return fTotalProjFileName;}
    TString GetBackgroundFileName(){return fBackgroundFileName;}
76
    double_t   GetCompFactor(){return fCompFactor;}
Jérémie Dudouet's avatar
Jérémie Dudouet committed
77 78


79
    // Read Cube conf
80 81 82 83 84 85
    Int_t ReadCube(const TString &FileName);
    Int_t ReadTabFile(const TString &FileName);
    Int_t ReadTotProj(const TString &FileName);
    Int_t ReadBackground(const TString &FileName);
    Int_t Read2DProj(const TString &FileName);
    Int_t ReadCalibs(const TString &EnergyFileName, const TString &EfficiencyFileName, double_t CompFactor);
Jérémie Dudouet's avatar
Jérémie Dudouet committed
86 87 88 89

    Int_t ReadGG(TH2 *h);

    Int_t AutoBuildProj(TString FileName, Int_t Mode);
90 91 92

    // Handle histos
    void BuildHistos();
Jérémie Dudouet's avatar
Jérémie Dudouet committed
93 94 95 96
    TH1 *GetLUTSpectrum(){return hLUTHist;}
    TH1D *GetTotalProjection(){return hTotalProj;}
    TH1 *GetBackground(){return hBackground;}
    TH2  *Get2DProj(){return h2DProj;}
97 98 99 100

    TF1 *GetECalibFunc(){return fECalibFunc;}
    TF1 *GetEEffFunc(){return fEffFunc;}

101 102
    TH1 *GetHistFromSpeFile(const TString &filename);
    void SaveBackground(const TString &filename);
Jérémie Dudouet's avatar
Jérémie Dudouet committed
103

104 105
    double EffFuncFormula(double *x, double *p);

106
    void ReEvalBackground(double_t Xmin, double_t XMax,double_t LenghtFactor, double_t Width, double_t FWHM_0, double_t FWHM_1, double_t FWHM_2);
107

108
    TH1 *Project(const vector<pair<double, double> > &gatesX, const vector<pair<double, double> > &gatesY, Bool_t BGSubtract = true);
Jérémie Dudouet's avatar
Jérémie Dudouet committed
109 110
    int Project3D(vector< pair<double, double> > gatesX, vector< pair<double, double> > gatesY);
    int Project2D(vector< pair<double, double> > gates);
111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180



    //******************************* RW util.c Stuffs ***********************************//

private:

    /************************************************
     *  COMPRESSED 3D CUBE FILE FORMAT
     *  1024byte file header
     *  4096byte data records
     *
     *      each record contains:
     *      compressed 8x8x4 mini-cubes, originally 4bytes/channel
     */

    /* 3d cube file header */
    typedef struct {
        char id[16];            /* "Incub8r3/Pro4d  " */
        int  numch;             /* number of channels on axis */
        int  bpc;               /* bytes per channel, = 4 */
        int  cps;               /* 1/cps symmetry compression, = 6 */
        int  numrecs;           /* number of 4kB data records in the file */
        char resv[992];         /* FUTURE flags */
    } FHead3D;

    /* 3Drecord header */
    typedef struct {
        int minmc;               /* start minicube number, starts at 0 */
        unsigned short nummc;    /* number of minicubes stored in here */
        unsigned short offset;   /* offset in bytes to first full minicube */
    } RHead3D;

    typedef struct {
        RHead3D h;
        unsigned char d[4088];  /* the bit compressed data */
    } Record3D;               /* see the compression alg for details */

    /* Some global data */

    int lumx3d[RW_MAXCH+16]; /* look-up table, maps 3d ch to linear minicube */
    int lumy3d[RW_MAXCH+16];
    int lumz3d[RW_MAXCH+16];

    /* look up minicube addr */
#define LUMC3D(x,y,z) lumx3d[x] + lumy3d[y] + lumz3d[z]

#define MCLEN(mcptr) (mcptr[0] + (mcptr[0] >= 7*32 ? mcptr[1]*32+2 : 1))

    typedef struct {
        FILE *CubeFile3d;
        int length;     /* length of axis on cube */
        int nummc;      /* number of minicubes */
        int *mcnum;     /* SuperRecord minicube number table */
        int numrecs;    /* number of 4kB data records in the file */
        int cubetype;   /* 0 for compressed 1/6 cube,
                 1 for compressed 1/2 cube,
                 2 for uncompressed 1/6 cube, or
                 3 for uncompressed 1/2 cube  */
        int numsr;
        int swapbytes;
    } CubeData3D;

    /* stuff added for linear combinations of cubes: */
    CubeData3D gd3dn[6], *gd3d = &gd3dn[0];
    int recloaded;

    /* Common Block Declarations */
    struct {
        int   numchs;
181
        double spec[6][MAXCHS], old_spec[6][MAXCHS], bspec[6][MAXCHS], background[MAXCHS], tmpspec[MAXCHS];
Jérémie Dudouet's avatar
Jérémie Dudouet committed
182 183
        double v_depth[MAXCHS], v_width;
        char  name_gat[800], old_name_gat[800];  /* changed from 20 or 30 */
184 185 186 187 188 189

        int   lo_ch[MAXGAM], hi_ch[MAXGAM];
        int   npks_ch[MAXCHS];                     /* no. of gammas for each ch. */
        short pks_ch[MAXCHS][60];                  /* gammas w/ counts in ch. */
        int   fitchx[MAXCHS], fitchy[MAXCHS];

Jérémie Dudouet's avatar
Jérémie Dudouet committed
190
        char  cubnam[800], levfile[800], fdname[800], fdgfile[800];
191 192
        char  progname[8];

Jérémie Dudouet's avatar
Jérémie Dudouet committed
193 194
        double eff_sp[MAXCHS], energy_sp[MAXCHS];
        double elo_sp[MAXCHS], ehi_sp[MAXCHS], ewid_sp[MAXCHS];
195 196

        int   luch[MAXCHS+1], matchs ,le2pro2d;
Jérémie Dudouet's avatar
Jérémie Dudouet committed
197 198
        double pro2d[MAXCHS*(MAXCHS+1)/2];
        double e2pro2d[MAXCHS*(MAXCHS+1)/2], e2e2spec[MAXCHS], e2e2e2sum;
199

200
        double bf1, bf2, bf4, bf5;
201

Jérémie Dudouet's avatar
Jérémie Dudouet committed
202 203 204
        double pk_shape[MAXGAM][15];
        double pk_deriv[MAXGAM][15];
        double w_deriv [MAXGAM][15];
205 206 207 208 209

        short looktab[16384];
        int   nclook, lookmin, lookmax;

        /* stuff added for linear combinations of cubes: */
Jérémie Dudouet's avatar
Jérémie Dudouet committed
210 211 212
        double dpro2d[MAXCHS*(MAXCHS+1)/2];
        char  cubenam1[5][800];
        double cubefact[5];
213 214 215
        int   many_cubes;
    } xxgd;

Jérémie Dudouet's avatar
Jérémie Dudouet committed
216 217
public:

218
    int autobkgnd(double_t Xmin=-1, double_t XMax=-1);
Jérémie Dudouet's avatar
Jérémie Dudouet committed
219

220 221
private:

222
    int read_tab_file(char *infilnam, int *outnclook, int *outlookmin, int *outlookmax, short *outlooktab, int dimlook);
223
    int setext(char *filnam, const char *cext, int filnam_len);
224
    FILE *open_readonly(char *filename);
225
    int get_file_rec(FILE *fd, void *data, int maxbytes, int swap_bytes);
226
    int file_error(const char *error_type, char *filename);
227 228 229 230 231 232 233
    void swapb8(char *buf);
    void swapb4(char *buf);
    void swapb2(char *buf);
    int open3dfile(char *name, int *numchs_ret); /* returns 1 on error */
    int close3dfile(void);
    void swap4(int *in);
    void swap2(unsigned short *in);
234
    int read_spe_file(char *filnam, double *sp, char *namesp, int *numch, int idimsp);
235 236 237 238 239 240
    int inq_file(char *filename, int *reclen);

    int read_cal_file(char *filnam, char *title, int *iorder, double *gain);
    int read_eff_file(char *filnam, char *title, double *pars);
    void setcubenum(int cubenum);

Jérémie Dudouet's avatar
Jérémie Dudouet committed
241 242 243 244 245 246
    int add_gate(double elo, double ehi);
    int add_trip_gate(double elo1, double ehi1, double elo2, double ehi2);
    void getmc3d(int num, unsigned int mc[256]);
    int read3dcube(int lo1, int hi1, int lo2, int hi2, int *spec_ret);
    int read_cube(int ilo1, int ihi1, int ilo2, int ihi2);
    void decompress3d (unsigned char in[1024], unsigned int out[256]);
247 248
    int find_bg2(int loch, int hich, double *x1, double *y1, double *x2, double *y2);
    int wspec(char *filnam, double *spec, int idim);
Jérémie Dudouet's avatar
Jérémie Dudouet committed
249 250
    FILE *open_new_file(char *filename, int force_open);
    int put_file_rec(FILE *fd, void *data, int numbytes);
251 252 253
};

#endif