Added virtual function Init for setting the slice, patch and n_eta_seg information...
[u/mrichter/AliRoot.git] / HLT / programs / MakePileup.cxx
1 #ifndef __CINT__
2 #include "AliL3FileHandler.h"
3 #include "AliL3DigitData.h"
4 #include "AliL3Transform.h"
5 #include "MakePileup.h"
6 #include <stdio.h>
7 #include <iostream.h>
8
9 void QSort(AliL3DigitData **a,Int_t first,Int_t last);
10 Int_t CompareDigits(AliL3DigitData *a,AliL3DigitData *b);
11 #endif
12
13 void MakePileup(Int_t triggerevent)
14 {
15   Int_t srow[2] = {0,175};
16   UInt_t size;
17   const Int_t NEvents = 25;
18   //Int_t Offset[25] = {0,10,20,30,40,50,60,70,80,90,100,110,120,130,140,150,160,170,180,190,200,210,220,230,240};
19   //Int_t Offset[25] = {0,20,40,60,80,100,120,140,160,180,200,220,240,-20,-40,-60,-80,-100,-120,-140,-160,-180,-200,-220,-240};
20   Int_t Offset[25]   = {30,60,90,120,150,180,210,240,270,300,330,360,-30,-60,-90,-120,-150,-180,-210,-240,-270,-300,-330,-360,-390};
21   Offset[triggerevent]=0;
22   AliL3FileHandler *hand[NEvents];
23   AliL3DigitRowData *data[NEvents];
24   AliL3DigitData *rowData[NEvents];
25   AliL3DigitData *test;
26   Char_t Carry[255];
27   //Int_t event = 0;
28   Int_t DigitsTot = 0;
29   Int_t MeVDigit;
30   Int_t DigitsPerRow = 0;
31   //AliL3DigitRowData *NData;
32   Byte_t *NData;
33   Int_t tot_dig=0;
34   
35   AliL3Transform *transform = new AliL3Transform();
36   const Int_t maxdigits=50000;
37   AliL3DigitData *temp[maxdigits];
38   AliL3DigitData *temp2[maxdigits];
39   
40   
41   for(Int_t slice=0 ; slice<36 ; slice++){
42     for(Int_t event = 0; event < NEvents ; event++){
43       hand[event] = new AliL3FileHandler();
44       sprintf(Carry,"/prog/alice/data/Rawdata/1_patch/pp/event_%d/digits_%d_0.raw",event,slice);
45       if(!hand[event]->SetBinaryInput(Carry))
46         cerr<<" Error opening file :"<<Carry<<endl;
47       cout << "Open File: " << Carry;
48       hand[event]->Init(slice,0,srow);
49       data[event] = hand[event]->CompBinary2Memory(size);
50       cout<<" done"<<endl;
51       hand[event]->CloseBinaryInput();
52     }
53     DigitsTot=0;
54     Int_t sign = slice < 18 ? -1 : 1;
55     for(Int_t row = 0 ; row < 176 ; row++){
56       
57       for(Int_t event = 0; event < NEvents ; event++){
58         rowData[event] = (AliL3DigitData *)data[event]->fDigitData;
59         for(UInt_t dig=0; dig<data[event]->fNDigit; dig++)
60           {
61             Int_t time = rowData[event][dig].fTime + sign*Offset[event];
62             if(time < transform->GetNTimeBins() && time >= 0)
63               DigitsTot++;
64           }
65         //DigitsTot += data[event]->fNDigit;
66         hand[event]->UpdateRowPointer(data[event]);
67       }
68       
69     }
70     //cout << "Try to allocate : " << DigitsTot << " Digits in row : " << row << endl;
71     
72     
73     //NData = (AliL3DigitRowData*) malloc(AllDigitsSize);// Allocate !
74     
75     Int_t AllDigitsSize = sizeof(AliL3DigitData) * DigitsTot + sizeof(AliL3DigitRowData) * 176;
76     NData = new Byte_t[AllDigitsSize];
77     memset(NData,0,AllDigitsSize);
78     AliL3DigitRowData *AllRowData = (AliL3DigitRowData*)NData;
79     
80     //cout << "Allocated " << endl;
81     
82     //Reset the data pointers, because they changed when doing UpdateRowPointer.
83     for(Int_t event=0; event<NEvents; event++)
84       data[event] = (AliL3DigitRowData*)hand[event]->GetDataPointer(size);
85     
86     tot_dig=0;
87     for(Int_t row = 0 ; row < 176 ; row++){
88       DigitsPerRow = 0;
89       for(Int_t event = 0; event < NEvents ; event++){
90         rowData[event] = (AliL3DigitData *)data[event]->fDigitData;
91         for(UInt_t dig=0; dig<data[event]->fNDigit; dig++)
92           {
93             Int_t time = rowData[event][dig].fTime + sign*Offset[event];
94             if(time < transform->GetNTimeBins() && time >= 0)
95               DigitsPerRow++;
96           }
97         //DigitsPerRow += data[event]->fNDigit;// + DigitsPerRow; 
98       }
99       //tot_dig += DigitsPerRow;
100       //cout << "Found : " << DigitsPerRow << " in all Events, Row :" << row << endl;
101       
102       /*
103         AllRowData->fRow = row;
104         AllRowData->fNDigit = DigitsPerRow;
105         test = (AliL3DigitData *) AllRowData->fDigitData;
106       */
107       
108       //cout << "Got Row Pointer " << endl;
109       MeVDigit = 0;
110       for(Int_t event = 0; event < NEvents ; event++){
111         rowData[event] = (AliL3DigitData *)data[event]->fDigitData;
112         for(UInt_t digit = 0; digit < data[event]->fNDigit ; digit++){
113           
114           //Copy the digits to a temporary storage, because we have to sort them before storing them :
115           Int_t time = rowData[event][digit].fTime + sign*Offset[event];
116           if(time >= transform->GetNTimeBins() || time < 0)
117             continue;
118           
119           temp[MeVDigit] = &rowData[event][digit];
120           temp[MeVDigit]->fTime = temp[MeVDigit]->fTime + sign*Offset[event];
121           
122           //cout << MeVDigit << " : " << digit << endl;
123           MeVDigit++;
124         }
125         hand[event]->UpdateRowPointer(data[event]);
126       }
127       //Sort the digits:
128       QSort(temp,0,MeVDigit);
129       int final_count=0;
130       
131       //Now, loop and check for overlaps:
132       for(Int_t c=0; c<MeVDigit; c++)
133         {
134           
135           if(c>0 && temp[c]->fTime == temp[c-1]->fTime && temp[c]->fPad == temp[c-1]->fPad)
136             {//this one is overlapping with the previous one:
137               temp2[final_count-1]->fCharge += temp[c+1]->fCharge;
138               if(temp2[final_count-1]->fCharge>=1024)
139                 temp2[final_count-1]->fCharge=1023;//Saturation
140             }
141           else
142             {
143               temp2[final_count] = temp[c];
144               final_count++;
145             }
146         }
147       tot_dig+=final_count;
148       AllRowData->fRow = row;
149       AllRowData->fNDigit = final_count;
150       test = (AliL3DigitData *) AllRowData->fDigitData;
151       
152       for(Int_t c=0; c<final_count; c++)//and copy them to the new event:
153         {
154           test[c].fCharge = temp2[c]->fCharge;
155           test[c].fPad = temp2[c]->fPad;
156           test[c].fTime = temp2[c]->fTime;
157           if(c>0 && test[c].fTime == test[c-1].fTime && test[c].fPad == test[c-1].fPad)
158             cout<<"Overlap at row "<<row<<" pad "<<(int)test[c].fPad<<" time "<<(int)test[c].fTime<<endl;
159           //  cout<<"Copying back, pad "<<(int)test[c].fPad<<" time "<<(int)test[c].fTime<<" charge "<<(int)test[c].fCharge<<endl;
160         }
161       
162       if(MeVDigit!=DigitsPerRow)
163         cerr<<endl<<"Error: "<<MeVDigit<<" "<<DigitsPerRow<<endl;
164       if(final_count > MeVDigit)
165         cerr<<"Error; final_count "<<final_count<<" MeVDigit "<<MeVDigit<<endl;
166       Byte_t *tmp = (Byte_t *)AllRowData;
167       Int_t UpdateSize = sizeof(AliL3DigitRowData) + sizeof(AliL3DigitData)*final_count;//DigitsPerRow;
168       tmp += UpdateSize;
169       AllRowData = (AliL3DigitRowData *) tmp; 
170       
171     }
172     if(tot_dig>DigitsTot)
173       cerr<<endl<<"Mismatching digitcount "<<tot_dig<<" "<<DigitsTot<<endl;
174     
175     for(Int_t event = 0; event < NEvents ; event++)
176       {
177         data[event]=0;
178         delete hand[event];
179       }
180     
181     AliL3FileHandler *Out = new AliL3FileHandler();
182     //sprintf(Carry,"/prog/alice/data/Rawdata/PileUp/digits_1_0.raw");
183     AllRowData = (AliL3DigitRowData*)NData;
184     sprintf(Carry,"/prog/alice/data/Rawdata/1_patch/pp/pileups/event_%d/digits_%d_0.raw",triggerevent,slice);
185     Out->SetBinaryOutput(Carry);
186     cout << "Open File: " << Carry << " to wite PileUp Event" << endl;
187     Out->Memory2CompBinary(176,AllRowData);
188     Out->CloseBinaryOutput();
189     delete Out;
190     delete [] NData;
191   }
192   
193   delete transform;
194 }
195
196 void QSort(AliL3DigitData **a,Int_t first,Int_t last)
197 {
198   //General sorting routine. only sorting the pointers.
199
200   AliL3DigitData *tmp;
201   int i; 
202   int j;
203   
204   while (last - first > 1) {
205     i = first;
206     j = last;
207     for (;;) {
208       while (++i < last && CompareDigits(a[i], a[first]) < 0)
209         ;
210       while (--j > first && CompareDigits(a[j], a[first]) > 0)
211             ;
212          if (i >= j)
213             break;
214
215          tmp  = a[i];
216          a[i] = a[j];
217          a[j] = tmp;
218     }
219       if (j == first) {
220          ++first;
221          continue;
222       }
223       tmp = a[first];
224       a[first] = a[j];
225       a[j] = tmp;
226       if (j - first < last - (j + 1)) {
227          QSort(a, first, j);
228          first = j + 1;   // QSort(j + 1, last);
229       } else {
230          QSort(a, j + 1, last);
231          last = j;        // QSort(first, j);
232       }
233    }
234   
235   
236 }
237
238 Int_t CompareDigits(AliL3DigitData *a,AliL3DigitData *b)
239 {
240
241   if(a->fPad==b->fPad && a->fTime == b->fTime) return 0;
242   
243   if(a->fPad<b->fPad) return -1;
244   if(a->fPad==b->fPad && a->fTime<b->fTime) return -1;
245   
246   return 1;
247   
248 }