Subversion Repositories f9daq

// g++  -I`root-config  --incdir` sa02read.C `root-config  --libs` -lz -o sa02read

#ifdef WIN32

#pragma warning( disable : 4996  )

#endif

#include <stdlib.h>

#include <stdio.h>

#include <stdint.h>

#include <iostream>

#include <fstream>

#include <zlib.h>

#include <iostream>

#include <vector>

#include <string>

#include "TString.h"

#include "TFile.h"

#include "TNtuple.h"

#include "TH3F.h"

#include "TH2F.h"

#include "TH1F.h"

#include "TH3D.h"

#include "TH1D.h"

#include "TGraph.h"

#include "TROOT.h"

#include "TObject.h"

#include "SA02_DEF.h"

#include "dataio.h"

#include "H3D.h"

#include "H2D.h"

#include "H1D.h"

#ifdef WIN32

// getopt

//#include "Tchar.h"

//#include "XGetopt.h"

#include <getopt.h>

#else

#include <getopt.h>

#endif

/*

#define ERRREC_ID 0x5

#define SCANREC_ID 0x4

#define POSREC_ID 0x3

#define RUNREC_ID 0x2

#define EVTREC_ID 0x1

*/

int gposx=-1;

int gposy=-1;

char gmapfile[256]="";

int gMapped=0;

int gInitialized=0;

unsigned int gEid2HapdXy[144];

int gBoardNumber=0;

/*

typedef struct {

  unsigned int id,length;

  unsigned int fver,time;

  unsigned int num_events,num_channels,pedestal,xy;

  int nx,x0,dx,ny,y0,dy;

} SCANREC;

*/

MONREC monrec;

RUNREC runrec;

POSREC posrec;

ITERATORREC iteratorrec;

DATREC datrec;

RUNINFO runinfo;

EVTREC evtrec;

H3D h3d;

H2D h2d;

H1D h1d;

int ReadMapFile(char *fname){

    int chip, ch,x,y;

    std::ifstream din;

    din.open (fname);

    for (int i=0; i<144; i++){

      din >> chip >> ch >> x>> y;

      x = 11 -x;

      y = 11 -y;

      gEid2HapdXy[chip*36+ch]=x*12+y;

    }

    din.close();

    return 1;

  }

class readdata {

private:

    TGraph *m_GrTime;

    TGraph *m_GrData;

    TGraph *m_GrEvent;

    TNtuple* m_Ntuple;

    TH2F* m_Scan;

    TH2F* m_Histo;

    TH2F* m_ScanSet;

    TH3F* m_Scan3d;

    TH2F* m_xy[144];

    TH2F* m_xysum;

        TTree * m_ttmon;

        TNtuple *m_ntpos;

        TTree *m_ntrun;

        TTree *m_waveforms;

public:

    uint8_t gCh;

    uint16_t gLen;

        uint16_t gData[10000];

    int m_neve; // number of events to process

    int m_ny;

    int m_dy;

    int m_y0;

    int m_val;

    int m_time;

    int m_cureve;

    int m_y;

    int m_towrite;

public:

    int m_mask;          // maska za status register

    int m_shft;          // shift za status register

    int m_debug;         // debug izpisi

    int m_write;         // ce je ta flag nastavljen, se v root file izpisuje ntuple

    unsigned int m_bitmask ;

    //--------------------------------------------------------------------

  int FillHistograms(unsigned int *data, double yset, double y){

    for (int ch=0;ch<144;ch++){

      int id=(17-ch/8)+gBoardNumber*18;

      int shft=(ch%8)*4;

      for (int i=0;i<4;i++){

        if (data[id] & (1<<(i+shft)) )  {

         if ( m_Scan3d ) m_Scan3d->Fill(ch,i,y);

         if ( m_Histo  ) m_Histo->Fill(ch,i);

        }

      }

      //const unsigned int mask = 0xFF;

      if (data[id] & (m_bitmask<<shft) ) {

        if (m_Scan)    m_Scan->Fill(ch,y);

        if (m_ScanSet) m_ScanSet->Fill(ch,yset);

        unsigned int ii= gEid2HapdXy[ch];

        if (m_xy[ii]) {

          double xs;

          double ys;

          if (runrec.fver){

            xs =  posrec.ix*runrec.dx+runrec.x0;

            ys =  posrec.iy*runrec.dy+runrec.y0;

          } else {

            xs =  posrec.xset;

            ys =  posrec.yset;

          }

          if (m_xysum) m_xysum->Fill( xs, ys);

          m_xy[ii]->Fill(xs, ys);

          //printf("%3d c=%d ch=%2d x=%6d y=%6d\n",ch, ii/36,ii%36,posrec.set_pos_x, posrec.set_pos_y);

        } else {

          printf("error m_xy[%d] does not exist\n",ii);

        }

      }

    }

    return 1;

  }

    //--------------------------------------------------------------------

short

#define MAXCH 144

    int Initialize(int nx, int ny, double ny0, double ny1) {

        printf("Initialize %d %d %2.2f %2.2f\n", nx,ny,ny0 ,ny1);

        if (m_write) {

            m_Ntuple = new TNtuple("nt","SA02 RAW data","time:eve:x:y");

            m_GrTime = new TGraph(ny);

            m_GrTime->SetTitle("grtime");

            m_GrEvent= new TGraph(ny);

            m_GrEvent->SetTitle("grevent");

            m_GrData = new TGraph(ny);

            m_GrData->SetTitle("grdata");

        }

        m_Scan  = new TH2F("Scan","Scan;SA02 ch;Threshold value",nx,-0.5,nx-0.5, ny,ny0,ny1 );

       // m_Scan  = new TH2F("Scan","Scan;SA02 ch;Threshold value",nx,-0.5,nx-0.5, (ny+1.2452)/0.0024,ny0,ny1 );

        m_Scan3d  = new TH3F("Scan3d","Scan3d;SA02 ch;bit;Threshold value",nx,-0.5,nx-0.5,  8,-0.5,7.5 ,ny,ny0,ny1 );

       // m_Scan3d  = new TH3F("Scan3d","Scan3d;SA02 ch;bit;Threshold value",nx,-0.5,nx-0.5,  8,-0.5,7.5 ,(ny+1.2452)/0.0024,ny0,ny1 );

        m_Histo = new TH2F("Histo","Histo;SA02 ch;bit",nx,-0.5,nx-0.5, 8,-0.5,7.5 );

        m_ScanSet  = new TH2F("ScanSet","Scan;SA02 ch;set value",nx,-0.5,nx-0.5, ny,ny0,ny1 );

        return 1;

    }

    //--------------------------------------------------------------------

    int ProcessData(unsigned int *rdata, int count,int neve) {

      int nc=0;

        for (int j=0; j<count; j++) {

          unsigned int recid   = rdata[j++];

            unsigned int len     = rdata[j++];

            if (m_debug ) printf(" recid=0x%0x len=%d pos=%d(maxpos %d) val=%d\n",recid, len, j, count, rdata[j]);

            if (m_debug>2 && len>2) printf("\n");

            if (recid!=3) {

                    m_val=(rdata[j] >> m_shft) & m_mask ;

                if (m_write) {

                    m_Ntuple->Fill(m_time,m_cureve,m_y,m_val);

                    m_GrTime->SetPoint(m_GrTime->GetN(),m_time,m_val);

                    m_GrEvent->SetPoint(m_GrEvent->GetN(),m_cureve,m_val);

                    m_GrData->SetPoint(m_GrData->GetN(),m_y,m_val);

                        }

            } else {

                if ( len !=36 ) printf(" recid=0x%0x len=%u pos=%d(maxpos %d) val=%d\n",recid, len, j, count, rdata[j]);

                nc+=FillHistograms(&rdata[j], m_y, m_val);

                if ( nc == neve ) return nc;

            }

            if(len>20) return nc;

            for (unsigned int i=0; i<len; i++) {

                if (m_debug>2) if (j< count) printf("0x%0x\t", rdata[j]);

                if (m_debug>2 &&i%4==3) printf("\n");

                j++;

            }

            if (m_debug>2)      printf("\n");

            if (len) j--;

            if(j+4>=count) return nc;

        }

        return nc;

    }

double H3DGetBinContentFrom( H3D *h ,unsigned int atx, unsigned int aty, unsigned int atz){

        unsigned int idx;

        if (!h) return 0;

        if (h->nx <= atx)  return 0;

        if (h->ny <= aty)  return 0;

        if (h->nz <= atz)  return 0;

        idx = atx+(aty+atz*h->ny)*h->nx;

        if (idx*sizeof(double)  < h->size ) {

                //if (h->data[idx]>0)  printf("xy %d %d %g\n",atx,aty,h->data[idx]);

                return h->data[idx];

        }

        return 0;

}

double H2DGetBinContentFrom( H2D *h ,unsigned int atx, unsigned int aty){

        unsigned int idx;

        if (!h) return 0;

        if (h->nx <= atx)  return 0;

        if (h->ny <= aty)  return 0;

        idx = atx+aty*h->nx;

        if (idx*sizeof(double)  < h->size ) {

                //if (h->data[idx]>0)  printf("xy %d %d %g\n",atx,aty,h->data[idx]);

                return h->data[idx];

        }

        return 0;

}

double H1DGetBinContentFrom( H1D *h ,unsigned int atx){

        if (!h) return 0;

        if (h->nx <= atx)  return 0;

        if (atx*sizeof(double)  < h->size ) {

                return h->data[atx];

        }

        return 0;

}

TH2D *hwf[4];

int V729_initialized = 0;

int V729_init(){

const double maxy=512;

//const double maxy=512;

const int ny=100;

//const int h0=10;

int len=100;

        TObject *hobj = gROOT->FindObject("gwf0"); if (hobj !=NULL) hobj->Write();

        hobj = gROOT->FindObject("gwf1"); if (hobj !=NULL) hobj->Write();

        hobj = gROOT->FindObject("gwf2"); if (hobj !=NULL) hobj->Write();

        hobj = gROOT->FindObject("gwf3"); if (hobj !=NULL) hobj->Write();

        hwf[0]= new TH2D("gwf0","CAENV_729 ch0",len,-0.5,len-0.5,ny,-maxy,maxy);

        hwf[1]= new TH2D("gwf1","CAENV_729 ch1",len,-0.5,len-0.5,ny,-maxy,maxy);

        hwf[2]= new TH2D("gwf2","CAENV_729 ch2",len,-0.5,len-0.5,ny,-maxy,maxy);

        hwf[3]= new TH2D("gwf3","CAENV_729 ch3",len,-0.5,len-0.5,ny,-maxy,maxy);

        // h1d=new TH1F("h1d","Charge;signal charge(a.u.);N",500,-maxy,5*maxy);

       return 0;

}

int V729_plot(uint16_t *buf ){

if (!V729_initialized) {

        V729_init();

        V729_initialized=1;

}

int start=310;

int end=335;

//const double maxy=2048;

//const double maxy=512;

//const double maxy=512;

//const int ny=100;

//const int h0=10;

int i;

uint32_t *hdr = (uint32_t *)buf;

int nb= hdr[1]-2*sizeof(uint32_t);

int nitems = nb /sizeof(uint16_t);

unsigned int recid=hdr[0];

uint16_t *data = &buf[4];

double x,y;

//double baseline=0;

double charge=0;

while (nitems>0){

int id = data[0]-1;

int len= data[1]/sizeof(uint16_t)-2;

if (m_debug)  printf("\t->> recid=%d nb=%d nitems=%d id=%d len=%d\n", recid, hdr[1],nitems,id,len);

if (id>3) return 0;

data+=2;

//for (i=0; i<200;i++) baseline += data[i];

//      baseline /=200.;

gCh = id;

gLen=len;

for (i=0; i<len;i++) {

        gData[i]=data[i];

        if (m_debug >2) printf("ADC %d %d\n", i, data[i]);

        x=i;

        y=data[i];

        //y-=baseline;

        if (i>start && i<end) charge += y;

        if (hwf[id]) hwf[id]->Fill(x,y-2000);

                else {

                  printf("gwf%d not initialized\n", id);

                  break;

                }

}

if (m_write) m_waveforms->Fill();

data+=len;

nitems=nitems - len-2;

//if (id==2) h1d->Fill(charge);

}

return 0;

}

    //--------------------------------------------------------------------

    int Process( const char *fnamein, int nevetoprocess, int histoonly) {

        gzFile fp=gzopen(fnamein,"r");

        printf("Process....\n");

        if (!fp) {

            printf("File %s cannot be opened!\n", fnamein);

            return 0;

        } else {

            printf("Processing file %s\n", fnamein);

        }

        unsigned int bsize=1000000; // initial buffer size

        unsigned int *buf = new unsigned int[bsize];

        unsigned int *data= buf+4; // data shift

        int nr=0;

        HEADERREC *hdr = (HEADERREC *) buf;

//        int nread=0;

        while (!gzeof(fp) ) {

                        int nitems  = gzread(fp,buf,   sizeof(HEADERREC) );

                if (nitems<=0) break;

            if (nitems !=sizeof(HEADERREC)) break;

            int recid   = hdr->id;

            int nb      = hdr->len - sizeof(HEADERREC);

            if ( hdr->len > bsize*sizeof(int)) {

                printf ("[INFO] Requested data size=%lu bigger than buffer size %u\n",hdr->len/sizeof(uint32_t), bsize);

                bsize=2*hdr->len/sizeof(int);

                printf("[INFO] Increasing data buffer to %d elements\n", bsize);

                delete buf;

                buf = new unsigned int[bsize];

                buf[0] = recid;

                buf[1] = nb + sizeof(HEADERREC);

                hdr = (HEADERREC *) buf;

                data = buf +4;

            }

                    switch (recid){

                                case H3D_ID :

                                case H2D_ID :

        case H1D_ID :

                                //  gzseek(fp, nb , SEEK_CUR); /* add one zero byte */

                                case POSREC_ID :

                                case MONREC_ID :

                                case RUNREC_ID :

                                case RUNINFO_ID:

                                case EVTREC_ID :

                                case ITERATORREC_ID :

                                case DATREC_ID :

                                case CAENV729REC_ID :

                                  nitems=gzread(fp, ((unsigned char *) buf)+sizeof(HEADERREC),  nb);

                                  break;

                                default:

                                  do { // scan the data until the next MONREC_ID

                                        printf("*0x%06x [unknown RECID=%u 0x%x] Length       = %u items=%lu\n",gztell (fp), hdr->id,hdr->id,hdr->len,hdr->len/sizeof(uint32_t));

                                        gzseek(fp, -sizeof(HEADERREC)+1 , SEEK_CUR); /* skip 1  byte */

                                        nitems  = gzread(fp,buf,   sizeof(HEADERREC) );

                                        if (nitems<=0) break;

                                        if (nitems !=sizeof(HEADERREC)) break;

                                        recid   = hdr->id;

                                        nb      = hdr->len - sizeof(HEADERREC);

                                  } while ( recid != MONREC_ID && nb!=492);

                  nitems=gzread(fp, ((unsigned char *) buf)+sizeof(HEADERREC),  nb);nr=0;

                                  break;

                        }

                        if (nitems<=0) break;

            nr++;

            if (m_debug || nr<30 || recid==0x3a98) printf("0x%06x [RECID=%u 0x%x] Length       = %u 0x%x items=%lu\n",gztell (fp),hdr->id,hdr->id,hdr->len,hdr->len,hdr->len/sizeof(uint32_t));

            if (nitems!=nb) {

                printf("Read Error nb=%d nitems=%d. Exiting ...\n",nb,nitems);

                break;

                //exit(0);

            }

            switch (recid) {

                    //case 0x3a98 : {

                        //}

            case POSREC_ID : {

                                V729_initialized = 0;

                memcpy(&posrec.id, buf, sizeof(posrec));

                printf("POSRECID = %u pos= %d %d\n",posrec.id, posrec.ix,posrec.iy);

                if (m_debug || nr<20) {

                                        printf("Length = %u\n",posrec.len);

                                        printf("Time = %u\n",posrec.time);

                                        printf("Iteration X = %d\n",posrec.ix);

                                        printf("MIKRO Position X = %d\n",posrec.x);

                                        printf("Set position X = %d\n",posrec.xset);

                                        printf("Iteration Y = %d\n",posrec.iy);

                                        printf("MIKRO Position Y = %d\n",posrec.y);

                                        printf("Set position Y = %d\n",posrec.yset);

                }

                                m_ntpos->Fill(posrec.time-runrec.time, posrec.ix,posrec.x,posrec.xset, posrec.iy,posrec.y,posrec.yset);

                break;

            }

                        case MONREC_ID : {

                                memcpy(&monrec.id, buf, sizeof(monrec));

                                //m_ttmon->Branch("monrec",&monrec,"id/i:len:time:imon[6]/I:vmon[6]:status[6]");

        if (!m_ttmon){

          m_ttmon = new TTree("Monitor","Voltage, Current and Status Monitor");

          m_ttmon->Branch("id",&monrec.id,"id/i");

          m_ttmon->Branch("len",&monrec.len,"len/i");

          m_ttmon->Branch("time",&monrec.time,"time/i");

          m_ttmon->Branch("imon",&monrec.imon[0],"imon[24]/I");

          m_ttmon->Branch("iset",&monrec.iset[0],"iset[24]/I");

          m_ttmon->Branch("vmon",&monrec.vmon[0],"vmon[24]/I");

          m_ttmon->Branch("vset",&monrec.vset[0],"vset[24]/I");

          m_ttmon->Branch("status",&monrec.status[0],"status[24]/I");

        }

                                m_ttmon->Fill();

                                break;

                        }

            case RUNREC_ID :{

                memcpy(&runrec.id, buf, sizeof(runrec));

                                m_ntrun->Branch("time",&runrec.time,"time/i");

                                m_ntrun->Branch("fver",&runrec.fver,"fver/i");

                                m_ntrun->Branch("nev",&runrec.nev,"nev/i");

                                m_ntrun->Branch("ped",&runrec.ped,"ped/i");

                                m_ntrun->Branch("nx",&runrec.nx,"nx/I");

                                m_ntrun->Branch("x0",&runrec.x0,"x0/I");

                                m_ntrun->Branch("dx",&runrec.dx,"dx/I");

                                m_ntrun->Branch("ny",&runrec.ny,"ny/I");

                                m_ntrun->Branch("y0",&runrec.y0,"y0/I");

                                m_ntrun->Branch("dy",&runrec.dy,"dy/I");

                                m_ntrun->Fill();

                        printf("***************RECID = %u\n",runrec.id);

                                printf("Length = %u\n",runrec.len);

                                printf("File version = %u\n",runrec.fver);

                                printf("Time = %u\n",runrec.time);

                                printf("Number of events per step = %u\n",runrec.nev);

                                printf("Pedestal = %u\n",runrec.ped);

                //              printf("Scan type = %u :: 0 -> single data scan :: 1 -> XY position scan\n",scanrec.xy);

                                printf("Number of steps in X = %d\n",runrec.nx);

                                printf("Start position X = %d\n",runrec.x0);

                                printf("Step size direction X = %d\n",runrec.dx);

                                printf("Number of steps in Y = %d\n",runrec.ny);

                                printf("Start position Y = %d\n",runrec.y0);

                                printf("Step size direction Y = %d\n",runrec.dy);

/*

                scanrec.nx = (scanrec.nx-1) /scanrec.dx;

                scanrec.ny = (scanrec.ny-1) / scanrec.dy;

                scanrec.nx ++;

                scanrec.ny ++;

                printf("Corrected Steps direction X = %d\n",scanrec.nx);

                printf("Corrected Steps direction Y = %d\n",scanrec.ny);

*/

                for (int i=0;i<144;i++){

                   char hname[256];

                   if (gMapped) sprintf(hname,"ch_%02d_%02d",i%12,i/12);

                   else sprintf(hname,"xy_%02d_%02d",i/36,i%36);

                   m_xy[i]= new TH2F(hname,hname, runrec.nx, runrec.x0-runrec.dx*0.5,  runrec.x0+runrec.dx*(runrec.nx-0.5),

                                                  runrec.ny, runrec.y0-runrec.dy*0.5,  runrec.y0+runrec.dy*(runrec.ny-0.5) );

                }

                m_xysum= (TH2F *) (m_xy[0]->Clone("xy"));

                m_xysum->SetTitle("xy all channels");

                break;

            }

            case RUNINFO_ID: {

                memcpy(&runinfo.id, buf, sizeof(runinfo));

                if (!gInitialized){

                  m_y0  = runinfo.x0;

                  m_dy  = runinfo.dx;

                  m_ny  = runinfo.nx;

                  m_towrite = runinfo.writemode;

                  if (m_towrite ==3) m_y0=runinfo.chip*36+runinfo.ch;

                  printf("RUNINFO x0=%d nx=%d dx=%d\n", runinfo.x0,runinfo.dx,runinfo.nx);

                  gInitialized=Initialize(144,m_ny,m_y0-0.5,m_y0-0.5+m_dy*m_ny);

                }

                break;

            }

            case EVTREC_ID:

                memcpy(&evtrec.id, buf, sizeof(evtrec));

                m_time   = evtrec.time;

                m_cureve = evtrec.nev;

                m_y = m_y0 +m_cureve * m_dy;

                break;

            case ITERATORREC_ID:

                memcpy(&iteratorrec.id, buf, sizeof(iteratorrec));

                printf("ITERATORREC level=%d n=%d value=%f step=%f\n", iteratorrec.level,iteratorrec.n,iteratorrec.value, iteratorrec.step);

                break;

            case DATREC_ID:

                memcpy(&datrec.id, buf, sizeof(datrec));

                printf("DATREC chip=%d channel=%d cmd=%d data=%d retval=%d\n", datrec.chip,datrec.channel,datrec.cmd, datrec.data, datrec.retval);

                break;

            case CAENV729REC_ID:{

                 if (!histoonly) {

                                         V729_plot( (uint16_t*)buf );

                                 }

            }

                break;

             case H3D_ID: {

                int hlen= sizeof(h3d)-sizeof(double *);

                memcpy(&h3d.id, buf, hlen);

                h3d.data = new double[h3d.size/sizeof(double)];

                int align=0;

                if (sizeof(double *)==4) align=0;

                memcpy(h3d.data, ((char *) buf)+hlen-align, h3d.size);

//                memcpy(h2d.data, ((char *) buf)+hlen, h2d.size);

                if (m_debug) printf("H3D %s nx=%d minx=%f stepx=%f ny=%d miny=%f stepy=%f  nz=%d minz=%f stepz=%f  min=%f max=%f size=%d title='%s' x='%s' y='%s' z='%s'  \n", h3d.name, h3d.nx, h3d.minx, h3d.stepx, h3d.ny, h3d.miny, h3d.stepy,h3d.nz, h3d.minz, h3d.stepz, h3d.min, h3d.max, h3d.size, h3d.title, h3d.titlex, h3d.titley, h3d.titlez );

                //printf("H3D sizeof(H2D)=%lu len-datasize=%d len=%lu datasize=%d\t",hlen,h3d.len-h3d.size,h3d.len,h3d.size);

                //printf("H3D sz=%d\n",sizeof(double *));

                H3D  *h = &h3d;

                char hname[100];

                sprintf(hname,"%s_%d",h->name,datrec.data);

                TObject *hobj = gROOT->FindObject(hname);

                if (hobj !=NULL) {

                  printf("TH3D duplicate name %s!\n",hname);

                  hobj->Write();

                  delete hobj;

                  //break;

                }

                if (h->stepx==0) h->stepx=1;

                if (h->stepy==0) h->stepy=1;

                if (h->stepz==0) h->stepz=1;

                TH3D *h3 = new TH3D(hname, h->title,

                     h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,

                     h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy,

                     h->nz, h->minz-0.5*h->stepz,h->minz+(h->nz-0.5)*h->stepz);

                h3->SetTitle(h->title);

                h3->GetXaxis()->SetTitle(h->titlex);

                h3->GetYaxis()->SetTitle(h->titley);

                h3->GetZaxis()->SetTitle(h->titlez);

                if (m_debug){

                  printf("TH3D name=%s %s %s %s %s",hname, h->title,h->titlex, h->titley, h->titlez);

                  printf("\tnx=%d %f %f \tny=%d %f %f\tnz=%d %f %f\n",

                     h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,

                     h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy,

                     h->nz, h->minz-0.5*h->stepz,h->minz+(h->nz-0.5)*h->stepz);

                }

                for (unsigned int ix=0;ix<h->nx;ix++){

                  for (unsigned int iy=0;iy<h->ny;iy++){

                     for (unsigned int iz=0;iz<h->nz;iz++){

                        double val = H3DGetBinContentFrom(h,ix,iy,iz);

                        if (m_debug>2) printf("%d %d %d %g\n", ix,iy,iz,val);

                        h3->SetBinContent(ix+1,iy+1,iz+1,val);

                     }

                  }

                }

                if (m_debug){

                  h3->ls();

                  h3->Print();

                  printf("TH3D Sum %f\n", h3->GetMaximum());

                }

                }

                break;

            case H2D_ID: {

                int hlen= sizeof(h2d)-sizeof(double *);

                memcpy(&h2d.id, buf, hlen);

                h2d.data = new double[h2d.size/sizeof(double)];

                int align=0;

                if (sizeof(double *)==4) align=0;

                memcpy(h2d.data, ((char *) buf)+hlen-align, h2d.size);

//                memcpy(h2d.data, ((char *) buf)+hlen, h2d.size);

                if (m_debug) printf("H2D nx=%d minx=%f stepx=%f ny=%d miny=%f stepy=%f min=%f max=%f size=%d title='%s' x='%s' y='%s' \n", h2d.nx, h2d.minx, h2d.stepx, h2d.ny, h2d.miny, h2d.stepy, h2d.min, h2d.max, h2d.size, h2d.title, h2d.titlex, h2d.titley );

                //printf("H2D sizeof(H2D)=%lu len-datasize=%d len=%lu datasize=%d\t",hlen,h2d.len-h2d.size,h2d.len,h2d.size);

                //printf("H2D sz=%d\n",sizeof(double *));

                H2D  *h = &h2d;

                char hname[100];

                sprintf(hname,"%s_%d",h->name,datrec.data);

                TObject *hobj = gROOT->FindObject(hname);

                if (hobj !=NULL) {

                  printf("TH2D duplicate name %s!\n",hname);

                  hobj->Write();

                  delete hobj;

                  //break;

                }

                if (h->stepx==0) h->stepx=1;

                if (h->stepy==0) h->stepy=1;

                TH2D *h2 = new TH2D(hname, h->title,

                     h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,