Blame |
Last modification |
View Log
| RSS feed
// 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,
h->ny, h->miny-0.5*h->stepy,h->miny+(h->ny-0.5)*h->stepy);
h2->SetTitle(h->title);
h2->GetXaxis()->SetTitle(h->titlex);
h2->GetYaxis()->SetTitle(h->titley);
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,
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);
for (unsigned int ix=0;ix<h->nx;ix++){
for (unsigned int iy=0;iy<h->ny;iy++){
double val = H2DGetBinContentFrom(h,ix,iy);
if (m_debug>2) printf("%d %d %g\n", ix,iy,val);
h2->SetBinContent(ix+1,iy+1,val);
}
}
if (m_debug){
h2->ls();
h2->Print();
printf("TH2D Sum %f\n", h2->GetMaximum());
}
}
break;
case H1D_ID: {
int hlen= sizeof(h1d)-sizeof(double *);
memcpy(&h1d.id, buf, hlen);
h1d.data = new double[h1d.size/sizeof(double)];
int align=0;
if (sizeof(double *)==4) align=0;
memcpy(h1d.data, ((char *) buf)+hlen-align, h1d.size);
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 );
H1D *h = &h1d;
char hname[100];
sprintf(hname,"%s_%d",h->name,datrec.data);
TObject *hobj = gROOT->FindObject(hname);
if (hobj !=NULL) {
printf("TH1D duplicate name %s!\n",hname);
hobj->Write();
delete hobj;
//break;
}
if (h->stepx==0) h->stepx=1;
TH1D *h1 = new TH1D(hname, h->title, h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx);
h1->SetTitle(h->title);
h1->GetXaxis()->SetTitle(h->titlex);
h1->GetYaxis()->SetTitle(h->titley);
h1->GetXaxis()->SetTitle("timebin");
h1->GetYaxis()->SetTitle(h->titlex);
if (m_debug) printf("TH1D name=%s %s %s %s\tnx=%d %g %g\n",hname, h->title,h->titlex, h->titley,
h->nx, h->minx-0.5*h->stepx,h->minx+(h->nx-0.5)*h->stepx);
for (unsigned int ix=0;ix<h->nx;ix++){
double val = H1DGetBinContentFrom(h,ix);
if (m_debug>2) printf("%d %g\n", ix,val);
h1->SetBinContent(ix+1,val);
}
if (m_debug){
h1->ls();
h1->Print();
printf("TH1D Sum %f\n", h1->GetMaximum());
}
} break;
default:
printf("Unknown record 0x%x \n", recid);
gzclose(fp);
delete buf;
return nr;
}
if (recid!=EVTREC_ID) continue;
if (histoonly) continue;
if (m_debug) printf("EVTREC %d histo=%d\n",recid,histoonly);
//if (nread++==nevetoprocess) break;
if (gposx>=0 && gposx!=posrec.ix ) continue;
if (gposy>=0 && gposy!=posrec.iy ) continue;
int nc=ProcessData(data,nb/sizeof(unsigned int)-2 , nevetoprocess);
if (m_debug) printf("EVTREC Events %d nb=%ld\n",nc,nb/sizeof(unsigned int)-2);
} // while (!gzeof(fp) )
gzclose(fp);
delete buf;
return nr;
}
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) {
m_debug = debug;
m_write = write2root;
m_neve = nevetoread;
m_mask = mask;
m_shft = shift;
m_towrite = 0;
m_bitmask = bmask;
char rootname[0xFF];
sprintf(rootname,"%s",fnameout);
TFile *rootFile = new TFile(rootname,"RECREATE");
m_ntpos = new TNtuple("pos","Position record", "time:ix:x:xset:iy:y:yset");
m_ntrun = new TTree("run","Run record");
if (m_write){
m_GrTime=NULL;
m_GrData=NULL;
m_GrEvent=NULL;
m_waveforms = new TTree("v729","Waveform collection") ;
m_waveforms->Branch("x",&posrec.ix,"x/I") ;
m_waveforms->Branch("y",&posrec.iy,"y/I") ;
m_waveforms->Branch("ch",&gCh,"ch/b") ;
m_waveforms->Branch("time",&posrec.time,"time/i") ;
m_waveforms->Branch("len",&gLen,"len/s") ;
m_waveforms->Branch("adc",gData,"adc[len]/s") ;
}
//Initialize();
printf("Konec inicializacije ...\n");
if (m_debug) printf("readdata histo=%d\n",histoonly);
int nr = 0;
for (unsigned int i=0; i< fnamein.size() ; i++) {
int neve = m_neve -nr;
if (m_neve == -1) nr += Process(fnamein[i].Data(), m_neve, histoonly);
else if (neve>0) nr += Process(fnamein[i].Data(), neve, histoonly);
}
fprintf(stderr, "End...\nreaddata::Number of events read=%d required=%d\n",nr,m_neve);
if (m_write) {
if (m_GrTime) m_GrTime->Write();
if (m_GrData) m_GrData->Write();
if (m_GrEvent) m_GrEvent->Write();
if (m_waveforms) m_waveforms->Write();
}
rootFile->Write();
if (m_debug) rootFile->ls();
rootFile->Close();
}
~readdata() { };
};
#ifdef MAIN
int main(int argc, char **argv) {
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;
//-----------------------------------------------------------------
char outputfile[256]="/tmp/sa02asic.root";
int nevetoread =-1;
int writentuple =0;
int verbose =0;
int shift=0;
int histoonly=0;
Int_t mask =0xFFFFFFFF;
unsigned int bitmask=0xf;
char c;
//char getopt_string[256];
//sprintf(getopt_string, "l:i:o:n:w:d:m:s:x:y:f:b:h:");
//char *getopt_cp = getopt_string;
for (int i=0;i<144;i++) gEid2HapdXy[i]=i;
std::vector<TString> inputfilelist;
int option_index = 0;
static struct option long_options[] = {
// {"xyz", no_argument, 0, 'z'},
{"board", required_argument, 0, 'l'},
{"input", required_argument, 0, 'i' },
{"output", required_argument, 0, 'o' },
{"neve", required_argument, 0, 'n' },
{"ntuple", required_argument, 0, 'w'},
{"verbose", required_argument, 0, 'd' },
{"shift", required_argument, 0, 's' },
{"mask", required_argument, 0, 'm' },
{"bitmask", required_argument, 0, 'b' },
{"xposition", required_argument, 0, 'x' },
{"yposition", required_argument, 0, 'y' },
{"histogramsonly", required_argument, 0, 'h' },
{"mapfile", required_argument, 0, 'm' },
{0, 0, 0, 0 }
};
while (1){
//c = getopt(argc, argv, "l:i:o:n:w:d:m:s:x:y:f:b:h:");
c = getopt_long(argc, argv, "l:i:o:n:w:d:m:s:x:y:f:b:h:", long_options, &option_index);
if (c == -1) break;
switch (c) {
case 0:
printf("option %s", long_options[option_index].name);
if (optarg)
printf(" with arg %s", optarg);
printf("\n");
break;
case 'i':
inputfilelist.push_back(TString(optarg));
break;
case 'o':
sprintf(outputfile,"%s",optarg);
break;
case 'n':
nevetoread = atoi(optarg) ;
break;
case 'w':
writentuple = atoi(optarg);
break;
case 'd':
verbose = atoi(optarg);
break;
case 'l':
gBoardNumber = atoi(optarg);
break;
case 's':
shift = strtoul(optarg,NULL,0);
break;
case 'm':
mask = strtoul(optarg,NULL,0);
break;
case 'b':
bitmask = strtoul(optarg,NULL,0);
break;
case 'x':
gposx = atoi(optarg) ;
break;
case 'y':
gposy = atoi(optarg) ;
break;
case 'h':
histoonly = atoi(optarg) ;
break;
case 'f':
sprintf(gmapfile,"%s", optarg) ;
gMapped=ReadMapFile(gmapfile);
break;
default:
abort();
}
}
std::cout << "Used: " << argv[0] << " -o " << outputfile
<< " -n " << nevetoread << " -w " << writentuple
<< " -d " << verbose << " -m 0x" << std::hex << mask
<< " -s 0x" << std::hex << shift << " -l " << gBoardNumber << " -h " << histoonly ;
for (unsigned int i=0; i< inputfilelist.size() ; i++) std::cout << " -i " << inputfilelist[i];
std::cout << std::endl;
if (!inputfilelist.size()) {
std::cout << "No input file specified" << std::endl;
return 0;
}
//-----------------------------------------------------------------
new readdata(inputfilelist,outputfile, writentuple, nevetoread,bitmask, mask, shift,-1,verbose,histoonly);
std::cout << "Output data files " << outputfile << std::endl;
return 0;
}
#endif