Rev 218 |
Blame |
Compare with Previous |
Last modification |
View Log
| RSS feed
#include <stdlib.h>
#include <stdio.h>
#include <TFile.h>
#include <TH1F.h>
#include <TH1D.h>
#include <TH2F.h>
#include <TStyle.h>
#include <TCanvas.h>
#include <TLegend.h>
// Program za analizo podatkov zajetih z AITSiPMDAQ
int aitana( const char *fname="a1.dat", int range=5, float cut=300, int *offset = NULL)
{
int debug=0;
//int off[4]={507,521,497,459};
//int off[4]={198,203,190,150}; /// run 1
int off[4]={120,120,120,120}; /// run 6
if (offset) for (int i=0;i<4;i++){off[i]=offset[i];};
gStyle->SetOptStat(0);
char rootname[0xFF];
sprintf(rootname, "%s.root", fname);
TFile *froot = new TFile(rootname,"RECREATE");
char histname[128];
TLegend * leg = new TLegend(0.8, 0.8, 1, 1);
TH1F *hadc[4];
TH1F *hadcsum;
TH1F *hadcsumx;
TH1F *hadcsumy;
TH2F *hxy,*hxy1,*hxy2;
TH2F *h2d=NULL;
TH2F *h2dx=NULL;
TH2F *h2dy=NULL;
TH2F *hxcorr=NULL;
TH2F *hycorr=NULL;
for (int i=0; i<4; i++) {
sprintf(histname,"adc%d",i);
hadc[i]=new TH1F(histname,histname,200,0.5,200*range+0.5);
}
hadcsum=new TH1F("adcsum","adcsum",200,0.5,8*200*range+0.5);
hadcsumx=new TH1F("adcsumx","adcsumx",200,0.5,200*range+0.5);
hadcsumy=new TH1F("adcsumy","adcsumy",200,0.5,200*range+0.5);
hxy=new TH2F("hxy","Center of gravity;ycog(a.u.);xcog(a.u.)",200,-0.1,1.1,200,-0.1,1.1);
hxy1=new TH2F("hxy1","Center of gravity for position 1;ycog(a.u.);xcog(a.u.)",200,-0.1,1.1,200,-0.1,1.1);
hxy2=new TH2F("hxy2","Center of gravity for position 2;ycog(a.u.);xcog(a.u.)",200,-0.1,1.1,200,-0.1,1.1);
FILE *fin=fopen(fname,"rb");
if (fin==NULL) {
printf("Error opening file %s\n",fname);
return 0;
}
unsigned short hdr[2];
int pos[7]={0,0,0,0,0,0,0};
int poshdr[9]={0,0,0,0,0,0,0,0,0};
float b[4], sum=0;
unsigned short a[4];
unsigned short buf[10000];
float posx=0, posy=0;
while(!feof(fin)) {
int stat=fread(hdr,1,4,fin);
//if (hdr[0]!= 1 ) printf("#%d %d\n",hdr[0], hdr[1] );
stat=fread(buf,1,hdr[1],fin);
if (stat!=hdr[1]){
printf("READ ERROR!\n");
break;
}
switch (hdr[0]) {
case 0x1:
for (int i=0;i<hdr[1]/4/sizeof(unsigned short);i++){
unsigned short *adc=&buf[i*4];
sum=0;
for (int j=0;j<4;j++){
a[j]=adc[j]&0xFFF;
//sum+=a[j];
hadc[j]->Fill(a[j]);
b[j]=a[j]-off[j];
if (b[j]< 0) b[j]=0;
sum+=b[j];
}
hadcsum->Fill(sum);
hadcsumx->Fill(a[0]+a[1]);
hadcsumy->Fill(a[2]+a[3]);
float y = (b[0]+b[1])? b[0]/(b[0]+b[1]) : 0;
float x = (b[2]+b[3])? b[2]/(b[2]+b[3]) : 0;
if (h2d) h2d->Fill(posx,posy,sum);
if (h2dx) h2dx->Fill(posx,posy,b[0]+b[1]);
if (h2dy) h2dy->Fill(posx,posy,b[2]+b[3]);
if (sum>cut) {
hxy->Fill(y,x);
if (pos[3]==poshdr[6]/4 && pos[4]==poshdr[7]/4) hxy1->Fill(y,x);
if (pos[3]==poshdr[6]/4*3 && pos[4]==poshdr[7]/4*3) hxy2->Fill(y,x);
if (hxcorr) hxcorr->Fill(pos[0],y);
if (hycorr) hycorr->Fill(pos[1],x);
}
}
break;
case 0x2:{
int *p= (int *) &buf[0];
for (int k=0;k<7;k++) pos[k]=p[k];
printf("Position: %d %d\n",pos[3], pos[4] );
posx = poshdr[0]+pos[3]* poshdr[3];
posy = poshdr[1]+pos[4]* poshdr[4];
}
break;
case 0x3:{
int *p= (int *) &buf[0];
for (int k=0;k<9;k++) poshdr[k]=p[k];
h2d=new TH2F("h2d","Position dependency of sum;x(step);y(step)",
poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5),
poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5) );
h2dx=new TH2F("h2dx","Position dependency of xsum;x(step);y(step)",
poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5),
poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5) );
h2dy=new TH2F("h2dy","Position dependency of ysum;x(step);y(step)",
poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5),
poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5) );
hxcorr=new TH2F("hxcorr","Correlation stage position vs cog;x(step);xcog(a.u.)",
poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5), 100 , 0, 1);
hycorr=new TH2F("hycorr","Correlation stage position vs cog;y(step);ycog(a.u.)",
poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5), 100 , 0, 1);
}
break;
default: break;
}
}
TCanvas* c= new TCanvas("c","ADC",750,50,700,700);
sprintf(rootname, "%s.pdf", fname);
c->Print(TString(rootname) + "[","pdf");
c->Divide(2,3);
for (int i=0; i<4; i++) {
c->cd(i+1)->SetLogy();
hadc[i]->DrawCopy();
}
c->cd(5)->SetLogy();hadcsumx->DrawCopy();
c->cd(6)->SetLogy();hadcsumy->DrawCopy();
/*
c->cd(7);hadcsumx->DrawCopy();
hadc[0]->DrawCopy("same");
hadc[1]->DrawCopy("same");
c->cd(8);hadcsumy->DrawCopy();
hadc[2]->DrawCopy("same");
hadc[3]->DrawCopy("same");
*/
c->Print(rootname,"pdf");
c= new TCanvas("c1","ait reconstruction",750,50,700,700);
c->Divide(2,4);
if (!h2d){
c->Clear();
c->Divide(2,3);
}
hxy->SetMinimum(-1);
hxy->SetMaximum(1000);
c->cd(1)->SetLogz(); hxy->DrawCopy("colz");
c->cd(3)->SetLogz(); hxy1->DrawCopy("colz");
c->cd(4)->SetLogz(); hxy2->DrawCopy("colz");
c->cd(2)->SetLogy(); hadcsum->DrawCopy();
if (h2d){
hxcorr->SetMinimum(-1);
hycorr->SetMinimum(-1);
c->cd(7); hxcorr->DrawCopy("colz");
c->cd(8); hycorr->DrawCopy("colz");
}
TH1D *px = hxy->ProjectionX("_px",0,-1);
c->cd(5); px->DrawCopy();
TH1D *py = hxy->ProjectionY("_py",0,-1);
c->cd(6); py->DrawCopy();
c->Modified();
c->Update();
c->Print(rootname,"pdf");
if (h2d){
c= new TCanvas("c2","position dependence sum, xsum, ysum",750,50,700,700);
c->Divide(1,3);
c->cd(1); h2d->DrawCopy("colz");
c->cd(2); h2dx->DrawCopy("colz");
c->cd(3); h2dy->DrawCopy("colz");
c->Modified();
c->Update();
c->Print(rootname,"pdf");
}
c->Print(TString(rootname) + "]","pdf");
froot->Write();
froot->Close();
if (fin) fclose(fin);
return 0;
}