Subversion Repositories f9daq

Rev

Rev 216 | Go to most recent revision | Details | Last modification | View Log | RSS feed

Rev Author Line No. Line
215 f9daq 1
#include <stdlib.h>
2
#include <stdio.h>
3
#include <TFile.h>
4
#include <TH1F.h>
5
#include <TH1D.h>
6
#include <TH2F.h>
7
#include <TStyle.h>
8
#include <TCanvas.h>
9
#include <TLegend.h>
10
// Program za analizo podatkov zajetih z AITSiPMDAQ 
11
 
12
int aitana( const char *fname="a1.dat", int range=10, float cut=500, int *offset = NULL)
13
{
14
        int debug=0;
15
 
16
        FILE *fin;
17
 
18
        unsigned short hdr[2];
19
 
20
        int off[4]={507,521,497,459};
21
                if (offset) for (int i=0;i<4;i++){off[i]=offset[i];};
22
    //    gStyle->SetOptStat(1111111);
23
                gStyle->SetOptStat(0);
24
//  gStyle->SetPalette(52);
25
 
26
        char rootname[0xFF];
27
                sprintf(rootname, "%s.root", fname);
28
        TFile *froot = new  TFile(rootname,"RECREATE");
29
                char histname[128];
30
            TLegend * leg = new TLegend(0.8, 0.8, 1, 1);
31
                TH1F *hadc[4];
32
                TH1F *hadcsum;
33
                TH1F *hadcsumx;
34
                TH1F *hadcsumy;
35
                TH2F *hxy;
36
                TH2F *h2d=NULL;
37
                TH2F *hxcorr=NULL;
38
                TH2F *hycorr=NULL;
39
                unsigned short buf[10000];
40
 
41
            for (int i=0; i<4; i++) {
42
         sprintf(histname,"adc%d",i);
43
         hadc[i]=new TH1F(histname,histname,200,0.5,200*range+0.5);
44
                }
45
 
46
                sprintf(histname,"adcsum");
47
                hadcsum=new TH1F(histname,histname,200,0.5,8*200*range+0.5);
48
 
49
                sprintf(histname,"adcsumx");
50
                hadcsumx=new TH1F(histname,histname,200,0.5,200*range+0.5);
51
 
52
                sprintf(histname,"adcsumy");
53
                hadcsumy=new TH1F(histname,histname,200,0.5,200*range+0.5);
54
 
55
        sprintf(histname,"hxy");
56
        hxy=new TH2F(histname,"Center of gravity;ycog(a.u.);xcog(a.u.)",200,0,1,200,0,1);
57
 
58
 
59
 
60
        fin=fopen(fname,"rb");
61
        if (fin==NULL) {
62
                printf("Error opening file %s\n",fname);
63
                return 0;
64
        }
65
 
66
        int pos[7]={0,0,0,0,0,0,0};
67
                int poshdr[9]={0,0,0,0,0,0,0,0,0};
68
        float b[4], sum=0;
69
                unsigned short a[4];
70
 
71
        while(!feof(fin)) {
72
 
73
                int stat=fread(hdr,1,4,fin);
74
 
75
                //if (hdr[0]!= 1 ) printf("#%d %d\n",hdr[0], hdr[1] );
76
 
77
                                stat=fread(buf,1,hdr[1],fin);
78
                                if (stat!=hdr[1]){
79
 
80
                                        printf("READ ERROR!\n");
81
                                        break;
82
 
83
                                }
84
                switch (hdr[0]) {
85
 
86
                                        case 0x1:
87
 
88
                                                   for (int i=0;i<hdr[1]/4/sizeof(unsigned short);i++){
89
                                                           unsigned short *adc=&buf[i*4];
90
                                                           sum=0;                                                          
91
                                                           for (int j=0;j<4;j++){
92
                                                                 a[j]=adc[j]&0xFFF;
93
                                                                 sum+=a[j];
94
                                                                 hadc[j]->Fill(a[j]);
95
                                 //const int off[4]={87,125,95,80};
96
 
97
 
98
                                 b[j]=a[j]-off[j]+10;
99
                                 sum+=b[j];                                                              
100
                                                           }
101
                                                           hadcsum->Fill(sum);
102
                                                           hadcsumx->Fill(a[0]+a[1]);
103
                                                           hadcsumy->Fill(a[2]+a[3]);
104
                               float y = (b[0]+b[1])? b[0]/(b[0]+b[1]) : 0;                                                    
105
                               float x = (b[2]+b[3])? b[2]/(b[2]+b[3]) : 0;
106
                                                           if (h2d) h2d->Fill(pos[0],pos[1],sum);
107
                                                           if (sum>cut) {
108
                                                                   hxy->Fill(y,x);
109
                                                                   if (hxcorr) hxcorr->Fill(pos[0],y);
110
                                                                   if (hycorr) hycorr->Fill(pos[1],x);
111
                                                           }       
112
 
113
                                                   }
114
                                                   break;
115
                                                case 0x2:{
116
                           int *p= (int *) &buf[0];    
117
                                                   for (int k=0;k<7;k++) pos[k]=p[k];
118
                           printf("Position: %d %d\n",pos[3], pos[4] );
119
                        }                                                  
120
                                                   break;
121
                                                case 0x3:{
122
                                                int *p= (int *) &buf[0];
123
                        for (int k=0;k<9;k++) poshdr[k]=p[k];                                          
124
                                                sprintf(histname,"h2d");    
125
                                h2d=new TH2F(histname,"Position dependency of sum;x(step);y(step)",poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5),
126
                                                                                                                   poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5)   );
127
 
128
                                sprintf(histname,"hxcorr");    
129
                                hxcorr=new TH2F(histname,"Correlation stage position vs cog;x(step);xcog(a.u.)",
130
                                                                                                                   poshdr[6],poshdr[0]-poshdr[3]*0.5 ,poshdr[0]+poshdr[6]*(poshdr[3]-0.5), 100 , 0, 1);
131
 
132
                                sprintf(histname,"hycorr");    
133
                                hycorr=new TH2F(histname,"Correlation stage position vs cog;y(step);ycog(a.u.)",
134
                                                                                                                   poshdr[7],poshdr[1]-poshdr[4]*0.5 ,poshdr[1]+poshdr[7]*(poshdr[4]-0.5), 100 , 0, 1);
135
                                                }
136
                                                break;
137
                                                default: break;
138
                                }
139
        }                              
140
 
141
        TCanvas* c= new TCanvas("c","ADC",750,50,700,700);
142
                sprintf(rootname, "%s.pdf", fname);
143
                c->Print(TString(rootname) + "[","pdf");
144
 
145
        c->Divide(2,3);
146
        for (int i=0; i<4; i++) {
147
                  c->cd(i+1)->SetLogy();
148
                  hadc[i]->DrawCopy();
149
            }
150
                c->cd(5)->SetLogy();hadcsumx->DrawCopy();
151
                c->cd(6)->SetLogy();hadcsumy->DrawCopy();
152
 
153
                /*
154
                c->cd(7);hadcsumx->DrawCopy();
155
                hadc[0]->DrawCopy("same");
156
                hadc[1]->DrawCopy("same");
157
                c->cd(8);hadcsumy->DrawCopy();
158
                hadc[2]->DrawCopy("same");
159
                hadc[3]->DrawCopy("same");
160
                */
161
                c->Print(rootname,"pdf");
162
 
163
 
164
                c= new TCanvas("c1","ait reconstruction",750,50,700,700);
165
                c->Divide(2,3);
166
                if (!h2d){
167
              c->Clear();
168
                  c->Divide(1,3);
169
                }      
170
        c->cd(3)->SetLogy(); hadcsum->DrawCopy();
171
 
172
                hxy->SetMinimum(-1);
173
                hxy->SetMaximum(1000);
174
                c->cd(1)->SetLogz(); hxy->DrawCopy("colz");
175
 
176
                if (h2d){
177
                  hxcorr->SetMinimum(-1);
178
                  hycorr->SetMinimum(-1);
179
                  c->cd(4); hxcorr->DrawCopy("colz");
180
                  c->cd(5); hycorr->DrawCopy("colz");
181
                  c->cd(6); h2d->DrawCopy("colz");
182
                }  
183
                TH1D *px = hxy->ProjectionX("_px",0,-1);
184
                c->cd(2); px->DrawCopy();
185
                c->Modified();
186
                c->Update();
187
                c->Print(rootname,"pdf");
188
                c->Print(TString(rootname) + "]","pdf");
189
                froot->Write();
190
                froot->Close();
191
                if (fin) fclose(fin);
192
 
193
        return 0;
194
}