CXRadReader.h 5.92 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 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 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195
#ifndef CXRadReader_H
#define CXRadReader_H

#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 */

class TH2F;
class TH1F;
class TF1;

class CXRadReader
{

private:

    // Calibs
    double fGains[6];
    double feff_pars[10];

    // compression factor
    Int_t fCompFactor   = 1;

    // basic histograms
    TH1F *hLUTHist      = nullptr;
    TH1F *hTotalProj    = nullptr;
    TH1F *hBackground   = nullptr;
    TH2F *h2DProj       = nullptr;

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


    // Binning of histograms
    Float_t edges[RW_MAXCH+1];

public:
    CXRadReader();
    ~CXRadReader();

    // Read Cube conf
    Int_t ReadCube(TString FileName);
    Int_t ReadTabFile(TString FileName);
    Int_t ReadTotProj(TString FileName);
    Int_t ReadBackground(TString FileName);
    Int_t Read2DProj(TString FileName);
    Int_t ReadCalibs(TString EnergyFileName, TString EfficiencyFileName, Int_t CompFactor);

    // Handle histos
    void BuildHistos();
    TH1F *GetLUTSpectrum(){return hLUTHist;}
    TH1F *GetTotalProjection(){return hTotalProj;}
    TH1F *GetBackground(){return hBackground;}
    TH2F *Get2DProj(){return h2DProj;}

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

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





    //******************************* 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];

    int gate_label_h_offset = 0;

    /* 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;
        float spec[6][MAXCHS], old_spec[6][MAXCHS], bspec[6][MAXCHS], background[MAXCHS];
        float v_depth[MAXCHS], v_width;
        char  name_gat[80], old_name_gat[80];  /* changed from 20 or 30 */

        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];

        char  cubnam[80], levfile[80], fdname[80], fdgfile[80];
        char  progname[8];

        float eff_sp[MAXCHS], energy_sp[MAXCHS];
        float elo_sp[MAXCHS], ehi_sp[MAXCHS], ewid_sp[MAXCHS];

        int   luch[MAXCHS+1], matchs ,le2pro2d;
        float pro2d[MAXCHS*(MAXCHS+1)/2];
        float e2pro2d[MAXCHS*(MAXCHS+1)/2], e2e2spec[MAXCHS], e2e2e2sum;

        float bf1, bf2, bf4, bf5;

        float pk_shape[MAXGAM][15];
        float pk_deriv[MAXGAM][15];
        float w_deriv [MAXGAM][15];

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

        /* stuff added for linear combinations of cubes: */
        float dpro2d[MAXCHS*(MAXCHS+1)/2];
        char  cubenam1[5][80];
        float cubefact[5];
        int   many_cubes;
    } xxgd;

private:

    int read_tab_file(char *infilnam, int *outnclook, int *outlookmin, int *outlookmax, short *outlooktab, int inmaxchs);
    int setext(char *filnam, const char *cext, int filnam_len);
    FILE *open_readonly(char *filnam);
    int get_file_rec(FILE *fd, void *data, int maxbytes, int swap_bytes);
    int file_error(const char *inmessage, char *infilnam);
    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);
    int read_spe_file(char *filnam, float *sp, char *namesp, int *numch, int idimsp);
    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);

};

#endif