25c5559101bb617b9df50fd62681127cb2f08cf7
[u/mrichter/AliRoot.git] / TPC / AliTPCTestMerge.C
1 #ifndef __CINT__
2 #include "alles.h"
3 #include "AliRunDigitizer.h"
4 #include "AliTPCDigitizer.h"
5 #include "AliH2F.h"
6 #include "TTree.h"
7
8
9 #endif
10
11 // Marian Ivanov
12 // test of the tpc merging using AliRunDigitizer and 
13 // TPC  Hits2Digits, Hits2SDigits and SDigits2Digits macros
14
15 // preparation  
16 // 0. make 2 directorys - ev1 and ev2
17
18 // 1.make hits, digits,sdigits and sdigits-digits in both directories 
19 //  1.a  galice -b -q grun.C and produce hits
20 //  1.b. cp galice.root galice.root.hits
21 //  1.c  run AliTPCHits2Digits.C
22 //  1.d  cp galice.root galice.root.digits
23 //  1.e  copy back cp galice.root.hits galice.root
24 //  1.f  run AliTPCSDigits2Digits.C
25 //  1.g  cp galice.root  galice.root.sdigits
26 //  1.h  run AliTPCSDigits2Digit.C 
27 //  1.i  cp galice.root galice.root.dig2
28
29 // 2. cp ev1/galice.root/galice.root.sdigit galice.root
30 // 3. load this macro and run testmerge()
31
32 // 4. run test function bellow to compare merged digits with original one
33 // 5. to be  noticed output ftom testx() function - should be bigger than
34 //    noise because in proces of digitisation we use different seed
35 //    of random numbers for diffusion and gas gain
36 //    -anyway in place where we have the signal should occur the signal in both casses
37 //    - only the amplitude should be different - by factor of sqrt
38
39       
40 void testmerge()
41 {
42   // merge two example events
43   //
44   //it merge two events -one from current directory -second from directory ev2
45   
46   AliRunDigitizer * manager = new AliRunDigitizer(2,1);
47   manager->SetInputStream(0,"galice.root");
48   manager->SetInputStream(1,"ev2/galice.root.sdigits");
49   AliTPCDigitizer dTPC(manager);
50   manager->SetNrOfEventsToWrite(1); 
51   TStopwatch timer;
52   timer.Start();
53
54   manager->Exec(""); 
55   timer.Stop();
56   timer.Print();
57
58 }
59
60 void drawmerged(Int_t sec, Int_t row, Int_t x1=-1, Int_t x2=-1, Int_t y1=-1, Int_t y2=-1)
61 {
62   //if you think that there is memory leak -
63   //you are tru but othervise graphic doesn't work
64   // sec=0; row =0;
65   TFile * f = new TFile("galice.root");
66   TFile * f1= new TFile("ev1/galice.root.digits");
67   TFile * f2= new TFile("ev2/galice.root.digits");
68   TTree * tree = (TTree*)f->Get("TreeD_75x40_100x60_0");
69   TTree * tree1 = (TTree*)f1->Get("TreeD_75x40_100x60_0");
70   TTree * tree2 = (TTree*)f2->Get("TreeD_75x40_100x60_0");
71   //
72   AliSimDigits *dig=0;
73   AliSimDigits *dig1=0;
74   AliSimDigits *dig2=0;
75   //
76   tree->GetBranch("Segment")->SetAddress(&dig);
77   tree1->GetBranch("Segment")->SetAddress(&dig1);
78   tree2->GetBranch("Segment")->SetAddress(&dig2);
79   AliTPCParam * param =(AliTPCParam*) f->Get("75x40_100x60");
80   Int_t index = param->GetIndex(sec,row);
81   tree->GetEvent(index);
82   tree1->GetEvent(index);
83   tree2->GetEvent(index);
84
85   TCanvas * c = new TCanvas(" Test merged digits", "test",600,900);
86   c->Divide(1,3);
87   //
88   c->cd(1);
89   AliH2F * his = dig->DrawDigits("cont1",x1,x2,y1,y2);
90   his->SetTitle("MergedDigits");
91   his->SetName("Merged Digits");
92   his->GetXaxis()->SetTitle("time");
93   his->GetYaxis()->SetTitle("pad");
94   his->DrawClone("cont1");
95   //
96   c->cd(2);
97   AliH2F * his1 =dig1->DrawDigits("cont1",x1,x2,y1,y2); 
98   his1->SetTitle("background");
99   his1->SetName("background");
100   his1->GetXaxis()->SetTitle("time");
101   his1->GetYaxis()->SetTitle("pad");
102   his1->DrawClone("cont1"); 
103   //
104   c->cd(3);
105   AliH2F * his2 =dig2->DrawDigits("cont1",x1,x2,y1,y2); 
106   his2->SetTitle("signal");
107   his2->SetName("signal");
108   his2->GetXaxis()->SetTitle("time");
109   his2->GetYaxis()->SetTitle("pad");
110   his2->DrawClone("cont1");  
111 }
112
113
114 void drawd(TFile * f, Int_t amp1, Int_t amp2)
115 {
116   TTree * tree = (TTree*)f->Get("TreeD_75x40_100x60_0");
117   AliSimDigits *dig=0;
118   tree->GetBranch("Segment")->SetAddress(&dig); 
119   TH1F * his = new TH1F("his","his",amp2-amp1,amp1,amp2);
120   for (Int_t i=0;i<60;i++){  
121     tree->GetEvent(i);   
122     dig->ExpandBuffer();
123     Int_t nrows = dig->GetNRows();
124     Int_t ncols = dig->GetNCols();
125     for (Int_t rows=0;rows<nrows; rows++)
126       for (Int_t col=0;col<ncols; col++){    
127         Int_t d  = dig->GetDigitFast(rows,col);
128         his->Fill(d);
129       }
130   }
131   his->Draw();
132 }
133
134 void test1(){
135   //test of the merged digits
136   //compare merged digits with standard digits
137   TFile f("galice.root");
138   TFile f1("ev1/galice.root.digits");
139   TFile f2("ev2/galice.root.digits");
140   TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_0");
141   TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_0");
142   TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_0");
143   //
144   AliSimDigits *dig=0;
145   AliSimDigits *dig1=0;
146   AliSimDigits *dig2=0;
147   //
148   tree->GetBranch("Segment")->SetAddress(&dig);
149   tree1->GetBranch("Segment")->SetAddress(&dig1);
150   tree2->GetBranch("Segment")->SetAddress(&dig2);
151   //
152   for (Int_t i=0;i<6000;i++){
153     tree->GetEvent(i);
154     tree1->GetEvent(i);
155     tree2->GetEvent(i);
156     dig->ExpandBuffer();
157     dig1->ExpandBuffer();
158     dig2->ExpandBuffer();
159     //
160     Int_t nrows = dig->GetNRows();
161     Int_t ncols = dig->GetNCols();
162     for (Int_t rows=0;rows<nrows; rows++)
163       for (Int_t col=0;col<ncols; col++){
164         Int_t d  = dig->GetDigitFast(rows,col);
165         Int_t d1 = dig1->GetDigitFast(rows,col);
166         Int_t d2 = dig2->GetDigitFast(rows,col);
167         
168         if (abs(d-(d1+d2))>4)
169           printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
170       }
171   }
172 }
173
174 void test5(){
175   //
176   //compare merged digits with digits obtained hits2sdig->sdigtodig
177   TFile f("galice.root");
178   TFile f1("ev1/galice.root.dig2");
179   TFile f2("ev2/galice.root.dig2");
180   TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_0");
181   TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_0");
182   TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_0");
183
184   AliSimDigits *dig=0;
185   AliSimDigits *dig1=0;
186   AliSimDigits *dig2=0;
187   
188   tree->GetBranch("Segment")->SetAddress(&dig);
189   tree1->GetBranch("Segment")->SetAddress(&dig1);
190   tree2->GetBranch("Segment")->SetAddress(&dig2);
191
192
193
194
195   for (Int_t i=0;i<6000;i++){
196     tree->GetEvent(i);
197     tree1->GetEvent(i);
198     tree2->GetEvent(i);
199     dig->ExpandBuffer();
200     dig1->ExpandBuffer();
201     dig2->ExpandBuffer();
202     
203     Int_t nrows = dig->GetNRows();
204     Int_t ncols = dig->GetNCols();
205
206     for (Int_t rows=0;rows<nrows; rows++)
207       for (Int_t col=0;col<ncols; col++){
208         Int_t d  = dig->GetDigitFast(rows,col);
209         Int_t d1 = dig1->GetDigitFast(rows,col);
210         Int_t d2 = dig2->GetDigitFast(rows,col);
211         
212         if (abs(d-d1-d2)>4)
213         //if (d2>5)
214           printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
215       }
216   }
217 }
218
219 void test3(){
220   //test of the merged digits
221   TFile f("galice.root");
222   TFile f1("ev1/galice.root.sdigits");
223   TFile f2("ev2/galice.root.sdigits");
224   TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_0");
225   TTree * tree1 = (TTree*)f1.Get("TreeS_75x40_100x60_0");
226   TTree * tree2 = (TTree*)f2.Get("TreeS_75x40_100x60_0");
227   //
228   AliSimDigits *dig=0;
229   AliSimDigits *dig1=0;
230   AliSimDigits *dig2=0;
231   //  
232   tree->GetBranch("Segment")->SetAddress(&dig);
233   tree1->GetBranch("Segment")->SetAddress(&dig1);
234   tree2->GetBranch("Segment")->SetAddress(&dig2);
235   //
236   for (Int_t i=0;i<6000;i++){
237     tree->GetEvent(i);
238     tree1->GetEvent(i);
239     tree2->GetEvent(i);
240     if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
241       printf("missed segments\n");
242     }
243     //
244     dig->ExpandBuffer();
245     dig1->ExpandBuffer();
246     dig2->ExpandBuffer();
247     //
248     Int_t nrows = dig->GetNRows();
249     Int_t ncols = dig->GetNCols();
250     //
251     for (Int_t rows=0;rows<nrows; rows++)
252       for (Int_t col=0;col<ncols; col++){
253         Int_t d  = dig->GetDigitFast(rows,col);
254         Int_t d1 = dig1->GetDigitFast(rows,col)/16;
255         Int_t d2 = dig2->GetDigitFast(rows,col)/16;     
256         if (abs(d-d1-d2)>4)
257         //if (d2>5)
258           printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
259       }
260   }
261 }
262
263
264 void TestSDigitsDig2(){
265   //test of the digits produced by the Hits2Digits 
266   //and Hits2SDigits - SDigits2Digits chain
267   TFile f1("galice.root.digits");
268   TFile f2("galice.root.dig2");
269   //
270   TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_0");
271   TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_0");
272   //
273   AliSimDigits *dig1=0;
274   AliSimDigits *dig2=0;
275   //  
276   tree1->GetBranch("Segment")->SetAddress(&dig1);
277   tree2->GetBranch("Segment")->SetAddress(&dig2);
278   //
279   for (Int_t i=0;i<6000;i++){
280     //tree->GetEvent(i);
281     tree1->GetEvent(i);
282     tree2->GetEvent(i);
283     //dig->ExpandBuffer();
284     if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
285       printf("miised semgnets\n");
286     }
287     //        
288     dig1->ExpandBuffer();
289     dig2->ExpandBuffer();
290     dig1->ExpandTrackBuffer();
291     dig2->ExpandTrackBuffer();
292     //
293     Int_t nrows = dig1->GetNRows();
294     Int_t ncols = dig1->GetNCols();
295     //
296     for (Int_t rows=0;rows<nrows; rows++)
297       for (Int_t col=0;col<ncols; col++){
298         Int_t d1 = dig1->GetDigitFast(rows,col);
299         Int_t d2 = dig2->GetDigitFast(rows,col);
300         Int_t t1_1 =dig1->GetTrackIDFast(rows,col,0);
301         Int_t t1_2 =dig1->GetTrackIDFast(rows,col,1);
302         Int_t t1_3 =dig1->GetTrackIDFast(rows,col,2);
303         //
304         Int_t t2_1 =dig2->GetTrackIDFast(rows,col,0);
305         Int_t t2_2 =dig2->GetTrackIDFast(rows,col,1);
306         Int_t t2_3 =dig2->GetTrackIDFast(rows,col,2);
307         //
308         if ( (d2>0) && (d1>0))
309           if ( ( ( d2>2) || (d1>2)) && (t2_1!=t1_1)) {
310             printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,t1_1,t2_1);
311             printf("\t\t\t\t%d\t%d\n",d1,d2);
312           }
313       }
314   }
315 }
316
317 void TestSDigitsDig1(){
318   //test of the digits produced by the Hits2Digits 
319   //and Hits2SDigits - SDigits2Digits chain
320   TFile f1("galice.root.digits");
321   TFile f2("galice.root.dig2");
322   //
323   TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_0");
324   TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_0");
325   //
326   AliSimDigits *dig1=0;
327   AliSimDigits *dig2=0;
328   //  
329   tree1->GetBranch("Segment")->SetAddress(&dig1);
330   tree2->GetBranch("Segment")->SetAddress(&dig2);
331   //
332   for (Int_t i=0;i<6000;i++){
333     //tree->GetEvent(i);
334     tree1->GetEvent(i);
335     tree2->GetEvent(i);
336     //dig->ExpandBuffer();
337     if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
338       printf("miised semgnets\n");
339     }
340     //        
341     dig1->ExpandBuffer();
342     dig2->ExpandBuffer();
343     dig1->ExpandTrackBuffer();
344     dig2->ExpandTrackBuffer();
345     //
346     Int_t nrows = dig1->GetNRows();
347     Int_t ncols = dig1->GetNCols();
348     //
349     for (Int_t rows=0;rows<nrows; rows++)
350       for (Int_t col=0;col<ncols; col++){
351         //      Int_t d  = dig->GetDigitFast(rows,col);
352         Int_t d1 = dig1->GetDigitFast(rows,col);
353         Int_t d2 = dig2->GetDigitFast(rows,col);
354  
355         //      
356         if ( ((d2>4) || (d1>4))  && abs(d1-d2)>4)
357           printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,d1,d2);
358       }
359   }
360 }
361
362
363
364 void test4(){
365   //TPC internal test
366   TFile f1("galice.root.sdigits");
367   TFile f2("galice.root.digits");
368   TTree * tree1 = (TTree*)f1.Get("TreeS_75x40_100x60_0");
369   TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_0");
370   //
371   AliSimDigits *dig1=0;
372   AliSimDigits *dig2=0;
373   //
374   tree1->GetBranch("Segment")->SetAddress(&dig1);
375   tree2->GetBranch("Segment")->SetAddress(&dig2);
376   //
377   for (Int_t i=0;i<6000;i++){
378     //tree->GetEvent(i);
379     tree1->GetEvent(i);
380     tree2->GetEvent(i);
381     //dig->ExpandBuffer();
382     if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
383       printf("miised semgnets\n");
384     }
385     //        
386     dig1->ExpandBuffer();
387     dig2->ExpandBuffer();
388     //
389     Int_t nrows = dig1->GetNRows();
390     Int_t ncols = dig1->GetNCols();
391     //
392     for (Int_t rows=0;rows<nrows; rows++)
393       for (Int_t col=0;col<ncols; col++){
394         Int_t d1 = dig1->GetDigitFast(rows,col)/16.;
395         Int_t d2 = dig2->GetDigitFast(rows,col);                
396         if ((d2>5) &&abs(d1-d2)>2)        printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,d1,d2);
397       }
398   }
399 }
400