3 #include "AliRunDigitizer.h"
4 #include "AliTPCDigitizer.h"
11 /// \file AliTPCTestMerge.C
13 /// test of the tpc merging using AliRunDigitizer and
14 /// TPC Hits2Digits, Hits2SDigits and SDigits2Digits macros
16 /// 0. make 2 directorys - ev1 and ev2
17 /// 1.make hits, digits,sdigits and sdigits-digits in both directories
18 /// 1.a galice -b -q grun.C and produce hits
19 /// 1.b. cp galice.root galice.root.hits
20 /// 1.c run AliTPCHits2Digits.C
21 /// 1.d cp galice.root galice.root.digits
22 /// 1.e copy back cp galice.root.hits galice.root
23 /// 1.f run AliTPCSDigits2Digits.C
24 /// 1.g cp galice.root galice.root.sdigits
25 /// 1.h run AliTPCSDigits2Digit.C
26 /// 1.i cp galice.root galice.root.dig2
27 /// 2. cp ev1/galice.root/galice.root.sdigit galice.root
28 /// 3. load this macro and run testmerge()
29 /// 4. run test function bellow to compare merged digits with original one
30 /// 5. to be noticed output ftom testx() function - should be bigger than
31 /// noise because in proces of digitisation we use different seed
32 /// of random numbers for diffusion and gas gain
33 /// -anyway in place where we have the signal should occur the signal in both casses
34 /// - only the amplitude should be different - by factor of sqrt
36 /// \author Marian Ivanov
40 /// merge two example events
42 /// it merge two events -one from current directory -second from directory ev2
44 if(gAlice) delete gAlice;
45 AliRunDigitizer * manager = new AliRunDigitizer(2,1);
46 manager->SetTreeDTPCBaseName("TreeD_75x40_100x60_150x60_");
47 manager->SetInputTreeTPCSBaseName("TreeS_75x40_100x60_150x60_");
48 manager->SetInputStream(0,"galice.root");
49 manager->SetInputStream(1,"ev2/galice.root.sdigits");
50 AliTPCDigitizer *dTPC = new AliTPCDigitizer(manager);
51 manager->SetNrOfEventsToWrite(1);
61 void drawmerged(Int_t sec, Int_t row, Int_t x1=-1, Int_t x2=-1, Int_t y1=-1, Int_t y2=-1)
63 /// if you think that there is memory leak -
64 /// you are tru but othervise graphic doesn't work
67 TFile * f = new TFile("galice.root");
68 TFile * f1= new TFile("ev1/galice.root.digits");
69 TFile * f2= new TFile("ev2/galice.root.digits");
70 TTree * tree = (TTree*)f->Get("TreeD_75x40_100x60_150x60_0");
71 TTree * tree1 = (TTree*)f1->Get("TreeD_75x40_100x60_150x60_0");
72 TTree * tree2 = (TTree*)f2->Get("TreeD_75x40_100x60_150x60_0");
78 tree->GetBranch("Segment")->SetAddress(&dig);
79 tree1->GetBranch("Segment")->SetAddress(&dig1);
80 tree2->GetBranch("Segment")->SetAddress(&dig2);
81 AliTPCParam * param =(AliTPCParam*) f->Get("75x40_100x60");
84 param=new AliTPCParamSR();
86 else param=(AliTPCParam*) f->Get("75x40_100x60_150x60");
87 Int_t index = param->GetIndex(sec,row);
88 tree->GetEvent(index);
89 tree1->GetEvent(index);
90 tree2->GetEvent(index);
92 TCanvas * c = new TCanvas(" Test merged digits", "test",600,900);
96 AliH2F * his = dig->DrawDigits("cont1",x1,x2,y1,y2);
97 his->SetTitle("MergedDigits");
98 his->SetName("Merged Digits");
99 his->GetXaxis()->SetTitle("time");
100 his->GetYaxis()->SetTitle("pad");
101 his->DrawClone("cont1");
104 AliH2F * his1 =dig1->DrawDigits("cont1",x1,x2,y1,y2);
105 his1->SetTitle("background");
106 his1->SetName("background");
107 his1->GetXaxis()->SetTitle("time");
108 his1->GetYaxis()->SetTitle("pad");
109 his1->DrawClone("cont1");
112 AliH2F * his2 =dig2->DrawDigits("cont1",x1,x2,y1,y2);
113 his2->SetTitle("signal");
114 his2->SetName("signal");
115 his2->GetXaxis()->SetTitle("time");
116 his2->GetYaxis()->SetTitle("pad");
117 his2->DrawClone("cont1");
121 void drawd(TFile * f, Int_t amp1, Int_t amp2)
123 TTree * tree = (TTree*)f->Get("TreeD_75x40_100x60_150x60_0");
125 tree->GetBranch("Segment")->SetAddress(&dig);
126 TH1F * his = new TH1F("his","his",amp2-amp1,amp1,amp2);
127 for (Int_t i=0;i<60;i++){
130 Int_t nrows = dig->GetNRows();
131 Int_t ncols = dig->GetNCols();
132 for (Int_t rows=0;rows<nrows; rows++)
133 for (Int_t col=0;col<ncols; col++){
134 Int_t d = dig->GetDigitFast(rows,col);
142 /// test of the merged digits
143 /// compare merged digits with standard digits
145 TFile f("galice.root");
146 TFile f1("ev1/galice.root.digits");
147 TFile f2("ev2/galice.root.digits");
148 TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_150x60_0");
149 TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_150x60_0");
150 TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_150x60_0");
153 AliSimDigits *dig1=0;
154 AliSimDigits *dig2=0;
156 tree->GetBranch("Segment")->SetAddress(&dig);
157 tree1->GetBranch("Segment")->SetAddress(&dig1);
158 tree2->GetBranch("Segment")->SetAddress(&dig2);
160 for (Int_t i=0;i<6000;i++){
165 dig1->ExpandBuffer();
166 dig2->ExpandBuffer();
168 Int_t nrows = dig->GetNRows();
169 Int_t ncols = dig->GetNCols();
170 for (Int_t rows=0;rows<nrows; rows++)
171 for (Int_t col=0;col<ncols; col++){
172 Int_t d = dig->GetDigitFast(rows,col);
173 Int_t d1 = dig1->GetDigitFast(rows,col);
174 Int_t d2 = dig2->GetDigitFast(rows,col);
176 if (abs(d-(d1+d2))>4)
177 printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
183 /// compare merged digits with digits obtained hits2sdig->sdigtodig
185 TFile f("galice.root");
186 TFile f1("ev1/galice.root.dig2");
187 TFile f2("ev2/galice.root.dig2");
188 TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_150x60_0");
189 TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_150x60_0");
190 TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_150x60_0");
193 AliSimDigits *dig1=0;
194 AliSimDigits *dig2=0;
196 tree->GetBranch("Segment")->SetAddress(&dig);
197 tree1->GetBranch("Segment")->SetAddress(&dig1);
198 tree2->GetBranch("Segment")->SetAddress(&dig2);
203 for (Int_t i=0;i<6000;i++){
208 dig1->ExpandBuffer();
209 dig2->ExpandBuffer();
211 Int_t nrows = dig->GetNRows();
212 Int_t ncols = dig->GetNCols();
214 for (Int_t rows=0;rows<nrows; rows++)
215 for (Int_t col=0;col<ncols; col++){
216 Int_t d = dig->GetDigitFast(rows,col);
217 Int_t d1 = dig1->GetDigitFast(rows,col);
218 Int_t d2 = dig2->GetDigitFast(rows,col);
222 printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
228 /// test of the merged digits
230 TFile f("galice.root");
231 TFile f1("ev1/galice.root.sdigits");
232 TFile f2("ev2/galice.root.sdigits");
233 TTree * tree = (TTree*)f.Get("TreeD_75x40_100x60_150x60_0");
234 TTree * tree1 = (TTree*)f1.Get("TreeS_75x40_100x60_150x60_0");
235 TTree * tree2 = (TTree*)f2.Get("TreeS_75x40_100x60_150x60_0");
238 AliSimDigits *dig1=0;
239 AliSimDigits *dig2=0;
241 tree->GetBranch("Segment")->SetAddress(&dig);
242 tree1->GetBranch("Segment")->SetAddress(&dig1);
243 tree2->GetBranch("Segment")->SetAddress(&dig2);
245 for (Int_t i=0;i<6000;i++){
249 if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
250 printf("missed segments\n");
254 dig1->ExpandBuffer();
255 dig2->ExpandBuffer();
257 Int_t nrows = dig->GetNRows();
258 Int_t ncols = dig->GetNCols();
260 for (Int_t rows=0;rows<nrows; rows++)
261 for (Int_t col=0;col<ncols; col++){
262 Int_t d = dig->GetDigitFast(rows,col);
263 Int_t d1 = dig1->GetDigitFast(rows,col)/16;
264 Int_t d2 = dig2->GetDigitFast(rows,col)/16;
267 printf("%d\t%d\t%d\t%d\t%d\t%d\n",i,rows,col,d,d1,d2);
273 void TestSDigitsDig2(){
274 /// test of the digits produced by the Hits2Digits
275 /// and Hits2SDigits - SDigits2Digits chain
277 TFile f1("galice.root.digits");
278 TFile f2("galice.root.dig2");
280 TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_150x60_0");
281 TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_150x60_0");
283 AliSimDigits *dig1=0;
284 AliSimDigits *dig2=0;
286 tree1->GetBranch("Segment")->SetAddress(&dig1);
287 tree2->GetBranch("Segment")->SetAddress(&dig2);
289 for (Int_t i=0;i<6000;i++){
293 //dig->ExpandBuffer();
294 if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
295 printf("miised semgnets\n");
298 dig1->ExpandBuffer();
299 dig2->ExpandBuffer();
300 dig1->ExpandTrackBuffer();
301 dig2->ExpandTrackBuffer();
303 Int_t nrows = dig1->GetNRows();
304 Int_t ncols = dig1->GetNCols();
306 for (Int_t rows=0;rows<nrows; rows++)
307 for (Int_t col=0;col<ncols; col++){
308 Int_t d1 = dig1->GetDigitFast(rows,col);
309 Int_t d2 = dig2->GetDigitFast(rows,col);
310 Int_t t1_1 =dig1->GetTrackIDFast(rows,col,0);
311 Int_t t1_2 =dig1->GetTrackIDFast(rows,col,1);
312 Int_t t1_3 =dig1->GetTrackIDFast(rows,col,2);
314 Int_t t2_1 =dig2->GetTrackIDFast(rows,col,0);
315 Int_t t2_2 =dig2->GetTrackIDFast(rows,col,1);
316 Int_t t2_3 =dig2->GetTrackIDFast(rows,col,2);
318 if ( (d2>0) && (d1>0))
319 if ( ( ( d2>2) || (d1>2)) && (t2_1!=t1_1)) {
320 printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,t1_1,t2_1);
321 printf("\t\t\t\t%d\t%d\n",d1,d2);
327 void TestSDigitsDig1(){
328 /// test of the digits produced by the Hits2Digits
329 /// and Hits2SDigits - SDigits2Digits chain
331 TFile f1("galice.root.digits");
332 TFile f2("galice.root.dig2");
334 TTree * tree1 = (TTree*)f1.Get("TreeD_75x40_100x60_150x60_0");
335 TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_150x60_0");
337 AliSimDigits *dig1=0;
338 AliSimDigits *dig2=0;
340 tree1->GetBranch("Segment")->SetAddress(&dig1);
341 tree2->GetBranch("Segment")->SetAddress(&dig2);
343 for (Int_t i=0;i<6000;i++){
347 //dig->ExpandBuffer();
348 if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
349 printf("miised semgnets\n");
352 dig1->ExpandBuffer();
353 dig2->ExpandBuffer();
354 dig1->ExpandTrackBuffer();
355 dig2->ExpandTrackBuffer();
357 Int_t nrows = dig1->GetNRows();
358 Int_t ncols = dig1->GetNCols();
360 for (Int_t rows=0;rows<nrows; rows++)
361 for (Int_t col=0;col<ncols; col++){
362 // Int_t d = dig->GetDigitFast(rows,col);
363 Int_t d1 = dig1->GetDigitFast(rows,col);
364 Int_t d2 = dig2->GetDigitFast(rows,col);
367 if ( ((d2>4) || (d1>4)) && abs(d1-d2)>4)
368 printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,d1,d2);
376 /// TPC internal test
378 TFile f1("galice.root.sdigits");
379 TFile f2("galice.root.digits");
380 TTree * tree1 = (TTree*)f1.Get("TreeS_75x40_100x60_150x60_0");
381 TTree * tree2 = (TTree*)f2.Get("TreeD_75x40_100x60_150x60_0");
383 AliSimDigits *dig1=0;
384 AliSimDigits *dig2=0;
386 tree1->GetBranch("Segment")->SetAddress(&dig1);
387 tree2->GetBranch("Segment")->SetAddress(&dig2);
389 for (Int_t i=0;i<6000;i++){
393 //dig->ExpandBuffer();
394 if ( (dig1->GetID()!=i) || (dig2->GetID()!=i) ) {
395 printf("miised semgnets\n");
398 dig1->ExpandBuffer();
399 dig2->ExpandBuffer();
401 Int_t nrows = dig1->GetNRows();
402 Int_t ncols = dig1->GetNCols();
404 for (Int_t rows=0;rows<nrows; rows++)
405 for (Int_t col=0;col<ncols; col++){
406 Int_t d1 = dig1->GetDigitFast(rows,col)/16.;
407 Int_t d2 = dig2->GetDigitFast(rows,col);
408 if ((d2>5) &&abs(d1-d2)>2) printf("%d\t%d\t%d\t\t%d\t%d\n",i,rows,col,d1,d2);