Subversion Repositories f9daq

Rev

Blame | Last modification | View Log | RSS feed

  1. // g++  -I`root-config  --incdir` sa02read.C `root-config  --libs` -lz -o sa02read
  2. #ifdef WIN32
  3. #pragma warning( disable : 4996  )
  4. #endif
  5. #include <stdlib.h>
  6. #include <stdio.h>
  7. #include <stdint.h>
  8. #include <iostream>
  9. #include <fstream>
  10. #include <zlib.h>
  11. #include <iostream>
  12. #include <vector>
  13. #include <string>
  14. #include "TString.h"
  15. #include "TFile.h"
  16. #include "TNtuple.h"
  17. #include "TH3F.h"
  18. #include "TH2F.h"
  19. #include "TH1F.h"
  20. #include "TH3D.h"
  21. #include "TH1D.h"
  22. #include "TGraph.h"
  23. #include "TROOT.h"
  24. #include "TObject.h"
  25.  
  26. #include "SA02_DEF.h"
  27. #include "dataio.h"
  28. #include "H3D.h"
  29. #include "H2D.h"
  30. #include "H1D.h"
  31. #ifdef WIN32
  32. // getopt
  33. //#include "Tchar.h"
  34. //#include "XGetopt.h"
  35. #include <getopt.h>
  36. #else
  37. #include <getopt.h>
  38. #endif
  39.  
  40. /*
  41. #define ERRREC_ID 0x5
  42. #define SCANREC_ID 0x4
  43. #define POSREC_ID 0x3
  44.  
  45. #define RUNREC_ID 0x2
  46. #define EVTREC_ID 0x1
  47. */
  48.  
  49. int gposx=-1;
  50. int gposy=-1;
  51. char gmapfile[256]="";
  52. int gMapped=0;
  53. int gInitialized=0;
  54. unsigned int gEid2HapdXy[144];
  55. int gBoardNumber=0;
  56.  
  57. /*
  58. typedef struct {
  59.   unsigned int id,length;
  60.   unsigned int fver,time;
  61.   unsigned int num_events,num_channels,pedestal,xy;
  62.   int nx,x0,dx,ny,y0,dy;
  63. } SCANREC;
  64.  
  65. */
  66.  
  67. MONREC monrec;
  68. RUNREC runrec;
  69. POSREC posrec;
  70. ITERATORREC iteratorrec;
  71. DATREC datrec;
  72. RUNINFO runinfo;
  73. EVTREC evtrec;
  74. H3D h3d;
  75. H2D h2d;
  76. H1D h1d;
  77.  
  78. int ReadMapFile(char *fname){
  79.     int chip, ch,x,y;
  80.     std::ifstream din;
  81.     din.open (fname);
  82.     for (int i=0; i<144; i++){
  83.       din >> chip >> ch >> x>> y;
  84.       x = 11 -x;
  85.       y = 11 -y;
  86.       gEid2HapdXy[chip*36+ch]=x*12+y;
  87.     }
  88.     din.close();
  89.     return 1;
  90.   }
  91.  
  92. class readdata {
  93. private:
  94.     TGraph *m_GrTime;
  95.     TGraph *m_GrData;
  96.     TGraph *m_GrEvent;
  97.     TNtuple* m_Ntuple;
  98.     TH2F* m_Scan;
  99.     TH2F* m_Histo;
  100.     TH2F* m_ScanSet;
  101.     TH3F* m_Scan3d;
  102.     TH2F* m_xy[144];
  103.     TH2F* m_xysum;
  104.         TTree * m_ttmon;
  105.         TNtuple *m_ntpos;
  106.         TTree *m_ntrun;
  107.         TTree *m_waveforms;
  108.  
  109. public:
  110.     uint8_t gCh;
  111.     uint16_t gLen;
  112.         uint16_t gData[10000];
  113.     int m_neve; // number of events to process
  114.     int m_ny;
  115.     int m_dy;
  116.     int m_y0;
  117.     int m_val;
  118.     int m_time;
  119.     int m_cureve;
  120.     int m_y;
  121.     int m_towrite;
  122.  
  123. public:
  124.     int m_mask;          // maska za status register
  125.     int m_shft;          // shift za status register
  126.     int m_debug;         // debug izpisi
  127.     int m_write;         // ce je ta flag nastavljen, se v root file izpisuje ntuple
  128.     unsigned int m_bitmask ;
  129.  
  130.     //--------------------------------------------------------------------
  131.  
  132.  
  133.  
  134.   int FillHistograms(unsigned int *data, double yset, double y){
  135.  
  136.     for (int ch=0;ch<144;ch++){
  137.       int id=(17-ch/8)+gBoardNumber*18;
  138.       int shft=(ch%8)*4;
  139.       for (int i=0;i<4;i++){
  140.         if (data[id] & (1<<(i+shft)) )  {
  141.          if ( m_Scan3d ) m_Scan3d->Fill(ch,i,y);
  142.          if ( m_Histo  ) m_Histo->Fill(ch,i);
  143.         }
  144.       }
  145.  
  146.       //const unsigned int mask = 0xFF;
  147.       if (data[id] & (m_bitmask<<shft) ) {
  148.         if (m_Scan)    m_Scan->Fill(ch,y);
  149.         if (m_ScanSet) m_ScanSet->Fill(ch,yset);
  150.  
  151.         unsigned int ii= gEid2HapdXy[ch];
  152.         if (m_xy[ii]) {
  153.  
  154.           double xs;
  155.           double ys;
  156.           if (runrec.fver){
  157.             xs =  posrec.ix*runrec.dx+runrec.x0;
  158.             ys =  posrec.iy*runrec.dy+runrec.y0;
  159.           } else {
  160.             xs =  posrec.xset;
  161.             ys =  posrec.yset;
  162.           }
  163.           if (m_xysum) m_xysum->Fill( xs, ys);
  164.           m_xy[ii]->Fill(xs, ys);
  165.           //printf("%3d c=%d ch=%2d x=%6d y=%6d\n",ch, ii/36,ii%36,posrec.set_pos_x, posrec.set_pos_y);
  166.         } else {
  167.           printf("error m_xy[%d] does not exist\n",ii);
  168.         }
  169.  
  170.       }
  171.  
  172.     }
  173.     return 1;
  174.   }
  175.  
  176.  
  177.  
  178.     //--------------------------------------------------------------------
  179. short
  180. #define MAXCH 144
  181.     int Initialize(int nx, int ny, double ny0, double ny1) {
  182.         printf("Initialize %d %d %2.2f %2.2f\n", nx,ny,ny0 ,ny1);
  183.         if (m_write) {
  184.             m_Ntuple = new TNtuple("nt","SA02 RAW data","time:eve:x:y");
  185.             m_GrTime = new TGraph(ny);
  186.             m_GrTime->SetTitle("grtime");
  187.             m_GrEvent= new TGraph(ny);
  188.             m_GrEvent->SetTitle("grevent");
  189.             m_GrData = new TGraph(ny);
  190.             m_GrData->SetTitle("grdata");
  191.  
  192.         }
  193.         m_Scan  = new TH2F("Scan","Scan;SA02 ch;Threshold value",nx,-0.5,nx-0.5, ny,ny0,ny1 );
  194.        // m_Scan  = new TH2F("Scan","Scan;SA02 ch;Threshold value",nx,-0.5,nx-0.5, (ny+1.2452)/0.0024,ny0,ny1 );
  195.         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 );
  196.        // 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 );
  197.         m_Histo = new TH2F("Histo","Histo;SA02 ch;bit",nx,-0.5,nx-0.5, 8,-0.5,7.5 );
  198.         m_ScanSet  = new TH2F("ScanSet","Scan;SA02 ch;set value",nx,-0.5,nx-0.5, ny,ny0,ny1 );
  199.  
  200.  
  201.         return 1;
  202.     }
  203.     //--------------------------------------------------------------------
  204.  
  205.     int ProcessData(unsigned int *rdata, int count,int neve) {
  206.       int nc=0;
  207.         for (int j=0; j<count; j++) {
  208.           unsigned int recid   = rdata[j++];
  209.             unsigned int len     = rdata[j++];
  210.             if (m_debug ) printf(" recid=0x%0x len=%d pos=%d(maxpos %d) val=%d\n",recid, len, j, count, rdata[j]);
  211.             if (m_debug>2 && len>2) printf("\n");
  212.  
  213.             if (recid!=3) {
  214.                     m_val=(rdata[j] >> m_shft) & m_mask ;
  215.                 if (m_write) {
  216.                     m_Ntuple->Fill(m_time,m_cureve,m_y,m_val);
  217.                     m_GrTime->SetPoint(m_GrTime->GetN(),m_time,m_val);
  218.                     m_GrEvent->SetPoint(m_GrEvent->GetN(),m_cureve,m_val);
  219.                     m_GrData->SetPoint(m_GrData->GetN(),m_y,m_val);
  220.                         }
  221.             } else {
  222.                 if ( len !=36 ) printf(" recid=0x%0x len=%u pos=%d(maxpos %d) val=%d\n",recid, len, j, count, rdata[j]);
  223.                 nc+=FillHistograms(&rdata[j], m_y, m_val);
  224.                 if ( nc == neve ) return nc;
  225.             }
  226.             if(len>20) return nc;
  227.             for (unsigned int i=0; i<len; i++) {
  228.                 if (m_debug>2) if (j< count) printf("0x%0x\t", rdata[j]);
  229.                 if (m_debug>2 &&i%4==3) printf("\n");
  230.                 j++;
  231.             }
  232.             if (m_debug>2)      printf("\n");
  233.             if (len) j--;
  234.             if(j+4>=count) return nc;
  235.         }
  236.  
  237.         return nc;
  238.     }
  239.  
  240.  
  241. double H3DGetBinContentFrom( H3D *h ,unsigned int atx, unsigned int aty, unsigned int atz){
  242.  
  243.         unsigned int idx;
  244.         if (!h) return 0;
  245.         if (h->nx <= atx)  return 0;
  246.         if (h->ny <= aty)  return 0;
  247.         if (h->nz <= atz)  return 0;
  248.         idx = atx+(aty+atz*h->ny)*h->nx;
  249.         if (idx*sizeof(double)  < h->size ) {
  250.                 //if (h->data[idx]>0)  printf("xy %d %d %g\n",atx,aty,h->data[idx]);
  251.                 return h->data[idx];
  252.         }
  253.  
  254.         return 0;
  255. }
  256.  
  257. double H2DGetBinContentFrom( H2D *h ,unsigned int atx, unsigned int aty){
  258.  
  259.         unsigned int idx;
  260.         if (!h) return 0;
  261.         if (h->nx <= atx)  return 0;
  262.         if (h->ny <= aty)  return 0;
  263.         idx = atx+aty*h->nx;
  264.         if (idx*sizeof(double)  < h->size ) {
  265.                 //if (h->data[idx]>0)  printf("xy %d %d %g\n",atx,aty,h->data[idx]);
  266.                 return h->data[idx];
  267.         }
  268.  
  269.  
  270.         return 0;
  271. }
  272.  
  273.  
  274. double H1DGetBinContentFrom( H1D *h ,unsigned int atx){
  275.  
  276.         if (!h) return 0;
  277.         if (h->nx <= atx)  return 0;
  278.  
  279.  
  280.         if (atx*sizeof(double)  < h->size ) {
  281.                 return h->data[atx];
  282.         }
  283.  
  284.  
  285.         return 0;
  286. }
  287.  
  288.  
  289. TH2D *hwf[4];
  290.  
  291. int V729_initialized = 0;
  292. int V729_init(){
  293. const double maxy=512;
  294. //const double maxy=512;
  295. const int ny=100;
  296. //const int h0=10;
  297. int len=100;
  298.  
  299.         TObject *hobj = gROOT->FindObject("gwf0"); if (hobj !=NULL) hobj->Write();
  300.         hobj = gROOT->FindObject("gwf1"); if (hobj !=NULL) hobj->Write();
  301.         hobj = gROOT->FindObject("gwf2"); if (hobj !=NULL) hobj->Write();
  302.         hobj = gROOT->FindObject("gwf3"); if (hobj !=NULL) hobj->Write();
  303.  
  304.         hwf[0]= new TH2D("gwf0","CAENV_729 ch0",len,-0.5,len-0.5,ny,-maxy,maxy);
  305.         hwf[1]= new TH2D("gwf1","CAENV_729 ch1",len,-0.5,len-0.5,ny,-maxy,maxy);
  306.         hwf[2]= new TH2D("gwf2","CAENV_729 ch2",len,-0.5,len-0.5,ny,-maxy,maxy);
  307.         hwf[3]= new TH2D("gwf3","CAENV_729 ch3",len,-0.5,len-0.5,ny,-maxy,maxy);
  308.  
  309.         // h1d=new TH1F("h1d","Charge;signal charge(a.u.);N",500,-maxy,5*maxy);
  310.        return 0;
  311. }
  312.  
  313. int V729_plot(uint16_t *buf ){
  314. if (!V729_initialized) {
  315.         V729_init();
  316.         V729_initialized=1;
  317. }
  318. int start=310;
  319. int end=335;
  320.  
  321. //const double maxy=2048;
  322. //const double maxy=512;
  323. //const double maxy=512;
  324. //const int ny=100;
  325. //const int h0=10;
  326.  
  327. int i;
  328. uint32_t *hdr = (uint32_t *)buf;
  329. int nb= hdr[1]-2*sizeof(uint32_t);
  330. int nitems = nb /sizeof(uint16_t);
  331. unsigned int recid=hdr[0];
  332. uint16_t *data = &buf[4];
  333. double x,y;
  334. //double baseline=0;
  335. double charge=0;
  336.  
  337. while (nitems>0){
  338. int id = data[0]-1;
  339. int len= data[1]/sizeof(uint16_t)-2;
  340.  
  341. if (m_debug)  printf("\t->> recid=%d nb=%d nitems=%d id=%d len=%d\n", recid, hdr[1],nitems,id,len);
  342. if (id>3) return 0;
  343.  
  344. data+=2;
  345.  
  346.  
  347.  
  348. //for (i=0; i<200;i++) baseline += data[i];
  349. //      baseline /=200.;
  350.  
  351. gCh = id;
  352. gLen=len;
  353. for (i=0; i<len;i++) {
  354.         gData[i]=data[i];
  355.         if (m_debug >2) printf("ADC %d %d\n", i, data[i]);
  356.         x=i;
  357.         y=data[i];
  358.         //y-=baseline;
  359.         if (i>start && i<end) charge += y;
  360.  
  361.         if (hwf[id]) hwf[id]->Fill(x,y-2000);
  362.                 else {
  363.                   printf("gwf%d not initialized\n", id);
  364.                   break;
  365.                 }
  366. }
  367. if (m_write) m_waveforms->Fill();
  368. data+=len;
  369. nitems=nitems - len-2;
  370. //if (id==2) h1d->Fill(charge);
  371. }
  372.  
  373. return 0;
  374. }
  375.  
  376.     //--------------------------------------------------------------------
  377.     int Process( const char *fnamein, int nevetoprocess, int histoonly) {
  378.         gzFile fp=gzopen(fnamein,"r");
  379.         printf("Process....\n");
  380.         if (!fp) {
  381.             printf("File %s cannot be opened!\n", fnamein);
  382.             return 0;
  383.         } else {
  384.             printf("Processing file %s\n", fnamein);
  385.         }
  386.         unsigned int bsize=1000000; // initial buffer size
  387.  
  388.         unsigned int *buf = new unsigned int[bsize];
  389.         unsigned int *data= buf+4; // data shift
  390.         int nr=0;
  391.  
  392.         HEADERREC *hdr = (HEADERREC *) buf;
  393.  
  394. //        int nread=0;
  395.         while (!gzeof(fp) ) {
  396.                         int nitems  = gzread(fp,buf,   sizeof(HEADERREC) );
  397.                 if (nitems<=0) break;
  398.             if (nitems !=sizeof(HEADERREC)) break;
  399.             int recid   = hdr->id;
  400.             int nb      = hdr->len - sizeof(HEADERREC);
  401.  
  402.             if ( hdr->len > bsize*sizeof(int)) {
  403.                 printf ("[INFO] Requested data size=%lu bigger than buffer size %u\n",hdr->len/sizeof(uint32_t), bsize);
  404.                 bsize=2*hdr->len/sizeof(int);
  405.                 printf("[INFO] Increasing data buffer to %d elements\n", bsize);
  406.                 delete buf;
  407.                 buf = new unsigned int[bsize];
  408.                 buf[0] = recid;
  409.                 buf[1] = nb + sizeof(HEADERREC);
  410.                 hdr = (HEADERREC *) buf;
  411.                 data = buf +4;
  412.             }
  413.  
  414.                     switch (recid){
  415.                                 case H3D_ID :
  416.                                 case H2D_ID :
  417.         case H1D_ID :
  418.                                 //  gzseek(fp, nb , SEEK_CUR); /* add one zero byte */
  419.                                 case POSREC_ID :
  420.                                 case MONREC_ID :
  421.                                 case RUNREC_ID :
  422.                                 case RUNINFO_ID:
  423.                                 case EVTREC_ID :
  424.                                 case ITERATORREC_ID :
  425.                                 case DATREC_ID :
  426.                                 case CAENV729REC_ID :
  427.                                   nitems=gzread(fp, ((unsigned char *) buf)+sizeof(HEADERREC),  nb);
  428.                                   break;
  429.                                 default:
  430.                                   do { // scan the data until the next MONREC_ID
  431.  
  432.                                         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));
  433.                                         gzseek(fp, -sizeof(HEADERREC)+1 , SEEK_CUR); /* skip 1  byte */
  434.                                         nitems  = gzread(fp,buf,   sizeof(HEADERREC) );
  435.                                         if (nitems<=0) break;
  436.                                         if (nitems !=sizeof(HEADERREC)) break;
  437.                                         recid   = hdr->id;
  438.                                         nb      = hdr->len - sizeof(HEADERREC);
  439.                                   } while ( recid != MONREC_ID && nb!=492);
  440.                   nitems=gzread(fp, ((unsigned char *) buf)+sizeof(HEADERREC),  nb);nr=0;
  441.                                   break;
  442.                         }
  443.  
  444.  
  445.  
  446.                         if (nitems<=0) break;
  447.             nr++;
  448.  
  449.             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));
  450.             if (nitems!=nb) {
  451.                 printf("Read Error nb=%d nitems=%d. Exiting ...\n",nb,nitems);
  452.                 break;
  453.                 //exit(0);
  454.             }
  455.             switch (recid) {
  456.                     //case 0x3a98 : {
  457.  
  458.                         //}
  459.             case POSREC_ID : {
  460.                                 V729_initialized = 0;
  461.                 memcpy(&posrec.id, buf, sizeof(posrec));
  462.                 printf("POSRECID = %u pos= %d %d\n",posrec.id, posrec.ix,posrec.iy);
  463.  
  464.                 if (m_debug || nr<20) {
  465.                                         printf("Length = %u\n",posrec.len);
  466.                                         printf("Time = %u\n",posrec.time);
  467.                                         printf("Iteration X = %d\n",posrec.ix);
  468.                                         printf("MIKRO Position X = %d\n",posrec.x);
  469.                                         printf("Set position X = %d\n",posrec.xset);
  470.                                         printf("Iteration Y = %d\n",posrec.iy);
  471.                                         printf("MIKRO Position Y = %d\n",posrec.y);
  472.                                         printf("Set position Y = %d\n",posrec.yset);
  473.                 }
  474.                                 m_ntpos->Fill(posrec.time-runrec.time, posrec.ix,posrec.x,posrec.xset, posrec.iy,posrec.y,posrec.yset);
  475.                 break;
  476.             }
  477.  
  478.                         case MONREC_ID : {
  479.                                 memcpy(&monrec.id, buf, sizeof(monrec));
  480.  
  481.                                 //m_ttmon->Branch("monrec",&monrec,"id/i:len:time:imon[6]/I:vmon[6]:status[6]");
  482.  
  483.         if (!m_ttmon){
  484.           m_ttmon = new TTree("Monitor","Voltage, Current and Status Monitor");
  485.           m_ttmon->Branch("id",&monrec.id,"id/i");
  486.           m_ttmon->Branch("len",&monrec.len,"len/i");
  487.           m_ttmon->Branch("time",&monrec.time,"time/i");
  488.           m_ttmon->Branch("imon",&monrec.imon[0],"imon[24]/I");
  489.           m_ttmon->Branch("iset",&monrec.iset[0],"iset[24]/I");
  490.           m_ttmon->Branch("vmon",&monrec.vmon[0],"vmon[24]/I");
  491.           m_ttmon->Branch("vset",&monrec.vset[0],"vset[24]/I");
  492.           m_ttmon->Branch("status",&monrec.status[0],"status[24]/I");
  493.         }
  494.  
  495.                                 m_ttmon->Fill();
  496.  
  497.                                 break;
  498.                         }
  499.  
  500.  
  501.             case RUNREC_ID :{
  502.  
  503.                 memcpy(&runrec.id, buf, sizeof(runrec));
  504.  
  505.                                 m_ntrun->Branch("time",&runrec.time,"time/i");
  506.                                 m_ntrun->Branch("fver",&runrec.fver,"fver/i");
  507.                                 m_ntrun->Branch("nev",&runrec.nev,"nev/i");
  508.                                 m_ntrun->Branch("ped",&runrec.ped,"ped/i");
  509.                                 m_ntrun->Branch("nx",&runrec.nx,"nx/I");
  510.                                 m_ntrun->Branch("x0",&runrec.x0,"x0/I");
  511.                                 m_ntrun->Branch("dx",&runrec.dx,"dx/I");
  512.                                 m_ntrun->Branch("ny",&runrec.ny,"ny/I");
  513.                                 m_ntrun->Branch("y0",&runrec.y0,"y0/I");
  514.                                 m_ntrun->Branch("dy",&runrec.dy,"dy/I");
  515.                                 m_ntrun->Fill();
  516.  
  517.                         printf("***************RECID = %u\n",runrec.id);
  518.                                 printf("Length = %u\n",runrec.len);
  519.                                 printf("File version = %u\n",runrec.fver);
  520.                                 printf("Time = %u\n",runrec.time);
  521.                                 printf("Number of events per step = %u\n",runrec.nev);
  522.                                 printf("Pedestal = %u\n",runrec.ped);
  523.  
  524.                 //              printf("Scan type = %u :: 0 -> single data scan :: 1 -> XY position scan\n",scanrec.xy);
  525.                                 printf("Number of steps in X = %d\n",runrec.nx);
  526.                                 printf("Start position X = %d\n",runrec.x0);
  527.                                 printf("Step size direction X = %d\n",runrec.dx);
  528.                                 printf("Number of steps in Y = %d\n",runrec.ny);
  529.                                 printf("Start position Y = %d\n",runrec.y0);
  530.                                 printf("Step size direction Y = %d\n",runrec.dy);
  531. /*
  532.                 scanrec.nx = (scanrec.nx-1) /scanrec.dx;
  533.                 scanrec.ny = (scanrec.ny-1) / scanrec.dy;
  534.                 scanrec.nx ++;
  535.                 scanrec.ny ++;
  536.                 printf("Corrected Steps direction X = %d\n",scanrec.nx);
  537.                 printf("Corrected Steps direction Y = %d\n",scanrec.ny);
  538. */
  539.                 for (int i=0;i<144;i++){
  540.                    char hname[256];
  541.                    if (gMapped) sprintf(hname,"ch_%02d_%02d",i%12,i/12);
  542.                    else sprintf(hname,"xy_%02d_%02d",i/36,i%36);
  543.                    m_xy[i]= new TH2F(hname,hname, runrec.nx, runrec.x0-runrec.dx*0.5,  runrec.x0+runrec.dx*(runrec.nx-0.5),
  544.                                                   runrec.ny, runrec.y0-runrec.dy*0.5,  runrec.y0+runrec.dy*(runrec.ny-0.5) );
  545.  
  546.                 }
  547.                 m_xysum= (TH2F *) (m_xy[0]->Clone("xy"));
  548.                 m_xysum->SetTitle("xy all channels");
  549.                 break;
  550.             }
  551.  
  552.             case RUNINFO_ID: {
  553.                 memcpy(&runinfo.id, buf, sizeof(runinfo));
  554.                 if (!gInitialized){
  555.                   m_y0  = runinfo.x0;
  556.                   m_dy  = runinfo.dx;
  557.                   m_ny  = runinfo.nx;
  558.  
  559.                   m_towrite = runinfo.writemode;
  560.                   if (m_towrite ==3) m_y0=runinfo.chip*36+runinfo.ch;
  561.                   printf("RUNINFO x0=%d nx=%d dx=%d\n", runinfo.x0,runinfo.dx,runinfo.nx);
  562.                   gInitialized=Initialize(144,m_ny,m_y0-0.5,m_y0-0.5+m_dy*m_ny);
  563.  
  564.                 }
  565.                 break;
  566.             }
  567.  
  568.             case EVTREC_ID:
  569.                 memcpy(&evtrec.id, buf, sizeof(evtrec));
  570.                 m_time   = evtrec.time;
  571.                 m_cureve = evtrec.nev;
  572.                 m_y = m_y0 +m_cureve * m_dy;
  573.                 break;
  574.             case ITERATORREC_ID:
  575.                 memcpy(&iteratorrec.id, buf, sizeof(iteratorrec));
  576.                 printf("ITERATORREC level=%d n=%d value=%f step=%f\n", iteratorrec.level,iteratorrec.n,iteratorrec.value, iteratorrec.step);
  577.                 break;
  578.             case DATREC_ID:
  579.                 memcpy(&datrec.id, buf, sizeof(datrec));
  580.                 printf("DATREC chip=%d channel=%d cmd=%d data=%d retval=%d\n", datrec.chip,datrec.channel,datrec.cmd, datrec.data, datrec.retval);
  581.                 break;
  582.             case CAENV729REC_ID:{
  583.                  if (!histoonly) {
  584.                                          V729_plot( (uint16_t*)buf );
  585.                                  }
  586.  
  587.             }
  588.                 break;
  589.              case H3D_ID: {
  590.                 int hlen= sizeof(h3d)-sizeof(double *);
  591.                 memcpy(&h3d.id, buf, hlen);
  592.                 h3d.data = new double[h3d.size/sizeof(double)];
  593.                 int align=0;
  594.                 if (sizeof(double *)==4) align=0;
  595.                 memcpy(h3d.data, ((char *) buf)+hlen-align, h3d.size);
  596. //                memcpy(h2d.data, ((char *) buf)+hlen, h2d.size);
  597.  
  598.  
  599.                 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 );
  600.  
  601.                 //printf("H3D sizeof(H2D)=%lu len-datasize=%d len=%lu datasize=%d\t",hlen,h3d.len-h3d.size,h3d.len,h3d.size);
  602.                 //printf("H3D sz=%d\n",sizeof(double *));
  603.                 H3D  *h = &h3d;
  604.  
  605.                 char hname[100];
  606.  
  607.                 sprintf(hname,"%s_%d",h->name,datrec.data);
  608.  
  609.  
  610.                 TObject *hobj = gROOT->FindObject(hname);
  611.  
  612.                 if (hobj !=NULL) {
  613.                   printf("TH3D duplicate name %s!\n",hname);
  614.                   hobj->Write();
  615.                   delete hobj;
  616.                   //break;
  617.                 }
  618.  
  619.                 if (h->stepx==0) h->stepx=1;
  620.                 if (h->stepy==0) h->stepy=1;
  621.                 if (h->stepz==0) h->stepz=1;
  622.                 TH3D *h3 = new TH3D(hname, h->title,
  623.                      h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,
  624.                      h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy,
  625.                      h->nz, h->minz-0.5*h->stepz,h->minz+(h->nz-0.5)*h->stepz);
  626.  
  627.                 h3->SetTitle(h->title);
  628.                 h3->GetXaxis()->SetTitle(h->titlex);
  629.                 h3->GetYaxis()->SetTitle(h->titley);
  630.                 h3->GetZaxis()->SetTitle(h->titlez);
  631.                 if (m_debug){
  632.                   printf("TH3D name=%s %s %s %s %s",hname, h->title,h->titlex, h->titley, h->titlez);
  633.                   printf("\tnx=%d %f %f \tny=%d %f %f\tnz=%d %f %f\n",
  634.                      h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,
  635.                      h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy,
  636.                      h->nz, h->minz-0.5*h->stepz,h->minz+(h->nz-0.5)*h->stepz);
  637.                 }
  638.  
  639.                 for (unsigned int ix=0;ix<h->nx;ix++){
  640.                   for (unsigned int iy=0;iy<h->ny;iy++){
  641.                      for (unsigned int iz=0;iz<h->nz;iz++){
  642.                         double val = H3DGetBinContentFrom(h,ix,iy,iz);
  643.                         if (m_debug>2) printf("%d %d %d %g\n", ix,iy,iz,val);
  644.                         h3->SetBinContent(ix+1,iy+1,iz+1,val);
  645.                      }
  646.                   }
  647.                 }
  648.                 if (m_debug){
  649.                   h3->ls();
  650.                   h3->Print();
  651.                   printf("TH3D Sum %f\n", h3->GetMaximum());
  652.                 }
  653.  
  654.                 }
  655.                 break;
  656.             case H2D_ID: {
  657.                 int hlen= sizeof(h2d)-sizeof(double *);
  658.                 memcpy(&h2d.id, buf, hlen);
  659.                 h2d.data = new double[h2d.size/sizeof(double)];
  660.                 int align=0;
  661.                 if (sizeof(double *)==4) align=0;
  662.                 memcpy(h2d.data, ((char *) buf)+hlen-align, h2d.size);
  663. //                memcpy(h2d.data, ((char *) buf)+hlen, h2d.size);
  664.  
  665.  
  666.                 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 );
  667.  
  668.                 //printf("H2D sizeof(H2D)=%lu len-datasize=%d len=%lu datasize=%d\t",hlen,h2d.len-h2d.size,h2d.len,h2d.size);
  669.                 //printf("H2D sz=%d\n",sizeof(double *));
  670.                 H2D  *h = &h2d;
  671.  
  672.                 char hname[100];
  673.  
  674.                 sprintf(hname,"%s_%d",h->name,datrec.data);
  675.                 TObject *hobj = gROOT->FindObject(hname);
  676.                 if (hobj !=NULL) {
  677.                   printf("TH2D duplicate name %s!\n",hname);
  678.                   hobj->Write();
  679.                   delete hobj;
  680.                   //break;
  681.                 }
  682.                 if (h->stepx==0) h->stepx=1;
  683.                 if (h->stepy==0) h->stepy=1;
  684.                 TH2D *h2 = new TH2D(hname, h->title,
  685.                      h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,
  686.                      h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy);
  687.                 h2->SetTitle(h->title);
  688.                 h2->GetXaxis()->SetTitle(h->titlex);
  689.                 h2->GetYaxis()->SetTitle(h->titley);
  690.                 if (m_debug) printf("TH2D name=%s %s %s %s\tnx=%d %g %g \tny=%d %g %g\n",hname, h->title,h->titlex, h->titley,
  691.                      h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx,
  692.                      h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy);
  693.  
  694.  
  695.  
  696.                 for (unsigned int ix=0;ix<h->nx;ix++){
  697.                   for (unsigned int iy=0;iy<h->ny;iy++){
  698.                     double val = H2DGetBinContentFrom(h,ix,iy);
  699.                     if (m_debug>2) printf("%d %d %g\n", ix,iy,val);
  700.                     h2->SetBinContent(ix+1,iy+1,val);
  701.                   }
  702.                 }
  703.                 if (m_debug){
  704.                   h2->ls();
  705.                   h2->Print();
  706.                   printf("TH2D Sum %f\n", h2->GetMaximum());
  707.                 }
  708.  
  709.                 }
  710.                 break;
  711. case H1D_ID: {
  712.                 int hlen= sizeof(h1d)-sizeof(double *);
  713.                 memcpy(&h1d.id, buf, hlen);
  714.                 h1d.data = new double[h1d.size/sizeof(double)];
  715.                 int align=0;
  716.                 if (sizeof(double *)==4) align=0;
  717.                 memcpy(h1d.data, ((char *) buf)+hlen-align, h1d.size);
  718.  
  719.                 if (m_debug) printf("H1D nx=%d minx=%f stepx=%f min=%f max=%f size=%d title='%s' x='%s' y='%s' \n", h1d.nx, h1d.minx, h1d.stepx, h1d.min, h1d.max, h1d.size, h1d.title, h1d.titlex, h1d.titley );
  720.  
  721.  
  722.                 H1D  *h = &h1d;
  723.  
  724.                 char hname[100];
  725.  
  726.                 sprintf(hname,"%s_%d",h->name,datrec.data);
  727.                 TObject *hobj = gROOT->FindObject(hname);
  728.                 if (hobj !=NULL) {
  729.                   printf("TH1D duplicate name %s!\n",hname);
  730.                   hobj->Write();
  731.                   delete hobj;
  732.                   //break;
  733.                 }
  734.  
  735.                 if (h->stepx==0) h->stepx=1;
  736.                 TH1D *h1 = new TH1D(hname, h->title, h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx);
  737.                 h1->SetTitle(h->title);
  738.                 h1->GetXaxis()->SetTitle(h->titlex);
  739.                 h1->GetYaxis()->SetTitle(h->titley);
  740.  
  741.  
  742.                 h1->GetXaxis()->SetTitle("timebin");
  743.                 h1->GetYaxis()->SetTitle(h->titlex);
  744.  
  745.                 if (m_debug) printf("TH1D name=%s %s %s %s\tnx=%d %g %g\n",hname, h->title,h->titlex, h->titley,
  746.                      h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx);
  747.  
  748.  
  749.  
  750.                 for (unsigned int ix=0;ix<h->nx;ix++){
  751.  
  752.                     double val = H1DGetBinContentFrom(h,ix);
  753.                     if (m_debug>2) printf("%d %g\n", ix,val);
  754.                     h1->SetBinContent(ix+1,val);
  755.  
  756.                 }
  757. if (m_debug){
  758. h1->ls();
  759. h1->Print();
  760. printf("TH1D Sum %f\n", h1->GetMaximum());
  761. }
  762.  
  763.                 } break;
  764.  
  765.             default:
  766.                 printf("Unknown record 0x%x \n", recid);
  767.                                 gzclose(fp);
  768.                                 delete buf;
  769.                                 return nr;
  770.  
  771.             }
  772.  
  773.             if (recid!=EVTREC_ID)  continue;
  774.             if (histoonly) continue;
  775.             if (m_debug) printf("EVTREC %d histo=%d\n",recid,histoonly);
  776.             //if (nread++==nevetoprocess) break;
  777.  
  778.             if (gposx>=0 && gposx!=posrec.ix ) continue;
  779.             if (gposy>=0 && gposy!=posrec.iy ) continue;
  780.             int nc=ProcessData(data,nb/sizeof(unsigned int)-2 , nevetoprocess);
  781.             if (m_debug) printf("EVTREC Events %d nb=%ld\n",nc,nb/sizeof(unsigned int)-2);
  782.         } //  while (!gzeof(fp) )
  783.         gzclose(fp);
  784.         delete buf;
  785.         return nr;
  786.     }
  787.  
  788.  
  789.  
  790.     readdata( std::vector<TString> &fnamein, const char *fnameout, int write2root=0, int nevetoread=-1, unsigned int bmask=0xff,  int mask=0xFFFFFFFF , int shift=0,int neve=-1,int debug=0 , int histoonly=0) {
  791.         m_debug = debug;
  792.         m_write = write2root;
  793.         m_neve  = nevetoread;
  794.         m_mask  = mask;
  795.         m_shft =  shift;
  796.         m_towrite = 0;
  797.         m_bitmask = bmask;
  798.  
  799.         char rootname[0xFF];
  800.         sprintf(rootname,"%s",fnameout);
  801.         TFile *rootFile = new TFile(rootname,"RECREATE");
  802.                 m_ntpos = new TNtuple("pos","Position record", "time:ix:x:xset:iy:y:yset");
  803.                 m_ntrun = new TTree("run","Run record");
  804.  
  805.  
  806.                 if (m_write){
  807.                   m_GrTime=NULL;
  808.       m_GrData=NULL;
  809.       m_GrEvent=NULL;
  810.                   m_waveforms = new TTree("v729","Waveform collection") ;
  811.             m_waveforms->Branch("x",&posrec.ix,"x/I") ;
  812.             m_waveforms->Branch("y",&posrec.iy,"y/I") ;
  813.             m_waveforms->Branch("ch",&gCh,"ch/b") ;
  814.             m_waveforms->Branch("time",&posrec.time,"time/i") ;
  815.                   m_waveforms->Branch("len",&gLen,"len/s") ;
  816.             m_waveforms->Branch("adc",gData,"adc[len]/s") ;
  817.     }
  818.  
  819.                 //Initialize();
  820.         printf("Konec inicializacije ...\n");
  821.  
  822.          if (m_debug) printf("readdata histo=%d\n",histoonly);
  823.  
  824.         int nr = 0;
  825.         for (unsigned int i=0; i<  fnamein.size() ; i++) {
  826.             int neve = m_neve -nr;
  827.             if (m_neve == -1)  nr += Process(fnamein[i].Data(), m_neve, histoonly);
  828.             else if (neve>0)   nr += Process(fnamein[i].Data(), neve, histoonly);
  829.         }
  830.         fprintf(stderr, "End...\nreaddata::Number of events read=%d required=%d\n",nr,m_neve);
  831.  
  832.         if (m_write) {
  833.             if (m_GrTime) m_GrTime->Write();
  834.             if (m_GrData) m_GrData->Write();
  835.             if (m_GrEvent) m_GrEvent->Write();
  836.                         if (m_waveforms) m_waveforms->Write();
  837.         }
  838.         rootFile->Write();
  839.         if (m_debug) rootFile->ls();
  840.  
  841.         rootFile->Close();
  842.     }
  843.  
  844.     ~readdata() { };
  845.  
  846. };
  847.  
  848.  
  849. #ifdef MAIN
  850.  
  851.  
  852. int main(int argc, char **argv) {
  853.  
  854.     std::cout << "Usage:" << argv[0] <<  "  -i input.dat -o output -n nevents -w writentuple -d debuglevel -s shift -m mask -l boardnumber -h histoonly -f electronicmap.map "<< std::endl;
  855.     //-----------------------------------------------------------------
  856.     char outputfile[256]="/tmp/sa02asic.root";
  857.  
  858.     int nevetoread      =-1;
  859.     int writentuple     =0;
  860.     int verbose         =0;
  861.     int shift=0;
  862.     int histoonly=0;
  863.     Int_t mask =0xFFFFFFFF;
  864.     unsigned int bitmask=0xf;
  865.     char c;
  866.  
  867.         //char getopt_string[256];
  868.         //sprintf(getopt_string, "l:i:o:n:w:d:m:s:x:y:f:b:h:");
  869.         //char *getopt_cp = getopt_string;
  870.  
  871.     for (int i=0;i<144;i++) gEid2HapdXy[i]=i;
  872.  
  873.     std::vector<TString> inputfilelist;
  874.  
  875.         int option_index = 0;
  876.     static struct option long_options[] = {
  877. //                  {"xyz",     no_argument, 0,  'z'},
  878.             {"board",     required_argument, 0,  'l'},
  879.             {"input",     required_argument, 0,  'i' },
  880.             {"output",    required_argument, 0,  'o' },
  881.             {"neve", required_argument,       0,  'n' },
  882.             {"ntuple",  required_argument, 0, 'w'},
  883.             {"verbose",    required_argument, 0,  'd' },
  884.                         {"shift",    required_argument, 0,  's' },
  885.                         {"mask",    required_argument, 0,  'm' },
  886.                         {"bitmask",    required_argument, 0,  'b' },
  887.                         {"xposition",    required_argument, 0,  'x' },
  888.                         {"yposition",    required_argument, 0,  'y' },
  889.                         {"histogramsonly",    required_argument, 0,  'h' },
  890.                         {"mapfile",    required_argument, 0,  'm' },
  891.             {0,         0,                 0,  0 }
  892.     };
  893.  
  894.  
  895.         while (1){
  896.                 //c = getopt(argc, argv, "l:i:o:n:w:d:m:s:x:y:f:b:h:");
  897.                 c = getopt_long(argc, argv, "l:i:o:n:w:d:m:s:x:y:f:b:h:", long_options, &option_index);
  898.                 if (c == -1) break;
  899.         switch (c) {
  900.         case 0:
  901.                    printf("option %s", long_options[option_index].name);
  902.             if (optarg)
  903.                 printf(" with arg %s", optarg);
  904.             printf("\n");
  905.             break;
  906.         case 'i':
  907.             inputfilelist.push_back(TString(optarg));
  908.             break;
  909.         case 'o':
  910.             sprintf(outputfile,"%s",optarg);
  911.             break;
  912.         case 'n':
  913.             nevetoread = atoi(optarg)    ;
  914.             break;
  915.         case 'w':
  916.             writentuple = atoi(optarg);
  917.             break;
  918.         case 'd':
  919.             verbose = atoi(optarg);
  920.             break;
  921.         case 'l':
  922.             gBoardNumber = atoi(optarg);
  923.             break;
  924.         case 's':
  925.             shift = strtoul(optarg,NULL,0);
  926.             break;
  927.         case 'm':
  928.             mask  = strtoul(optarg,NULL,0);
  929.             break;
  930.         case 'b':
  931.             bitmask  = strtoul(optarg,NULL,0);
  932.             break;
  933.         case 'x':
  934.             gposx = atoi(optarg)    ;
  935.             break;
  936.         case 'y':
  937.             gposy = atoi(optarg)    ;
  938.             break;
  939.         case 'h':
  940.             histoonly = atoi(optarg)    ;
  941.             break;
  942.         case 'f':
  943.             sprintf(gmapfile,"%s", optarg)    ;
  944.             gMapped=ReadMapFile(gmapfile);
  945.             break;
  946.  
  947.  
  948.         default:
  949.             abort();
  950.         }
  951.     }
  952.     std::cout << "Used: " << argv[0]  << " -o " << outputfile
  953.               << " -n " << nevetoread << " -w " << writentuple
  954.               << " -d  " << verbose  << " -m  0x" << std::hex << mask
  955.               << " -s  0x" << std::hex << shift  << " -l  " << gBoardNumber << " -h  " << histoonly ;
  956.     for (unsigned int i=0; i<  inputfilelist.size() ; i++) std::cout << " -i " << inputfilelist[i];
  957.     std::cout << std::endl;
  958.     if (!inputfilelist.size()) {
  959.         std::cout << "No input file specified" << std::endl;
  960.         return 0;
  961.     }
  962.     //-----------------------------------------------------------------
  963.  
  964.     new readdata(inputfilelist,outputfile, writentuple, nevetoread,bitmask, mask, shift,-1,verbose,histoonly);
  965.  
  966.     std::cout << "Output data files " << outputfile <<  std::endl;
  967.     return 0;
  968. }
  969.  
  970. #endif
  971.