gcc 4.3 warnings fixed
[u/mrichter/AliRoot.git] / TRD / AliTRDtestBeam.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id$ */
17
18 ////////////////////////////////////////////////////////////////////////////
19 //                                                                        //
20 // Class to handle the test beam data of 2007                             //
21 //                                                                        //
22 // Authors:                                                               //
23 //   Sylwester Radomski (radomski@physi.uni-heidelberg.de)                //
24 //   Anton Andronic (A.Andronic@gsi.de)                                   //
25 //                                                                        //
26 ////////////////////////////////////////////////////////////////////////////
27
28 //#include <iostream>
29 //#include <fstream>
30 //#include <sys/types.h>
31 //#include <sys/stat.h>
32 //#include <fcntl.h>
33 //include <unistd.h>
34
35 #include "AliTRDrawStream.h"
36 #include "AliRawReaderMemory.h"
37 #include "AliTRDtestBeam.h"
38
39 ClassImp(AliTRDtestBeam)
40
41 const Long_t AliTRDtestBeam::fgkFileHeadSize = 544; // ?
42 const Long_t AliTRDtestBeam::fgkEventHeadSize = 68; //?
43 const Long_t AliTRDtestBeam::fgkLdcHeadSize = 68; //?
44 const Long_t AliTRDtestBeam::fgkEquipHeadSize = 28; //
45 const Int_t AliTRDtestBeam::fgkVmeIn =1; //VME event in
46 const Int_t AliTRDtestBeam::fgkSimIn =1; //Si-strips in
47
48 //typedef char byte;
49
50 //offsets in bytes
51 const Int_t AliTRDtestBeam::fgkPosRun = 20; //run nr. (in file and event header)
52 const Int_t AliTRDtestBeam::fgkPosLength = 0; //event/equip. length
53 const Int_t AliTRDtestBeam::fgkEqId = 8;  //equipment id.
54 const Int_t AliTRDtestBeam::fgkPosSiOff = 12;  //Si data size offset (3 extra words!!!)      
55      
56 using namespace std;
57
58 //____________________________________________________________________________ 
59 AliTRDtestBeam::AliTRDtestBeam() :
60   fDataStream(0),
61   fHeaderIsRead(0),
62   fEventCount(0),
63   fLimit(4), 
64   fCurrent(0),
65   fDdlOff(0),
66   fSiOff(0),
67   fQdcOff(0),
68   fDdlSize(0),
69   fFileHeader(0),
70   fEventHeader(0),
71   fEventData(0),
72   fNSi1(0),
73   fNSi2(0),
74   fCher(0),
75   fPb(0)
76 {
77   //
78   // Standard construction
79   //
80
81 }
82 //____________________________________________________________________________ 
83 AliTRDtestBeam::AliTRDtestBeam(const char *filename) :
84   fDataStream(0),
85   fHeaderIsRead(0),
86   fEventCount(0),
87   fLimit(4), 
88   fCurrent(0),
89   fDdlOff(0),
90   fSiOff(0),
91   fQdcOff(0),
92   fDdlSize(0),
93   fFileHeader(0),
94   fEventHeader(0),
95   fEventData(0),
96   fNSi1(0),
97   fNSi2(0),
98   fCher(0),
99   fPb(0)
100 {
101   //
102   // AliTRDtestBeam constructor
103   //
104
105   fDataStream = new ifstream(filename, ifstream::in | ifstream::binary );
106   cout << fDataStream->is_open() << endl;
107   //fHeaderIsRead = kTRUE;
108   fHeaderIsRead = kTRUE;
109
110   fFileHeader = new Char_t[fgkFileHeadSize];
111   fEventHeader = new Char_t[fgkEventHeadSize];
112   fEventData = new Char_t[fLimit];
113
114 }
115
116 //____________________________________________________________________________
117 AliTRDtestBeam::AliTRDtestBeam(const AliTRDtestBeam &tb)
118  :TObject(tb),
119   fDataStream(0),
120   fHeaderIsRead(0),
121   fEventCount(0),
122   fLimit(4),
123   fCurrent(0),
124   fDdlOff(0),
125   fSiOff(0),
126   fQdcOff(0),
127   fDdlSize(0),
128   fFileHeader(0),
129   fEventHeader(0),
130   fEventData(0),
131   fNSi1(0),
132   fNSi2(0),
133   fCher(0),
134   fPb(0)
135 {
136   //
137   // Copy constructor
138   //
139
140 }
141
142 //____________________________________________________________________________ 
143 AliTRDtestBeam::~AliTRDtestBeam() 
144 {
145   //
146   // Destructor
147   //
148
149   if (fDataStream)  delete fDataStream;
150   if (fEventHeader) delete fEventHeader;
151   if (fFileHeader)  delete fFileHeader;
152   if (fEventData)   delete fEventData;
153
154 }
155
156 //____________________________________________________________________________ 
157 Int_t AliTRDtestBeam::NextEvent() 
158 {
159   //
160   // Read the next event
161   //
162   
163   Long_t dataSize=0,ldcOff; //,ldc_id,ldc2_id;
164   Long_t ldcSize,eqId; //,ev_l2;
165   Long_t eventNr,evL1;
166   Long_t word;
167   
168   if ( !fHeaderIsRead ) {
169     fDataStream->read(fFileHeader, fgkFileHeadSize);
170     if(fDataStream->fail()) {
171       cerr << "Error reading file header! " << endl;    
172       return false;
173     }
174     cout  << " Run nr.  " << Int(fgkPosRun, fFileHeader) << endl;    
175     fHeaderIsRead=kTRUE;
176   }
177
178   fDataStream->read(fEventHeader, fgkEventHeadSize);
179   if(fDataStream->fail()) {
180     cerr << "End of file, Event " << fEventCount  << endl;      
181     return false;
182   }
183   
184   dataSize = Int(fgkPosLength, fEventHeader)-fgkEventHeadSize; //?
185   eventNr = Int((4+fgkPosRun), fEventHeader); //ev.nr.
186   //cout << " Event " << eventNr <<" size "<< dataSize <<endl;
187   
188     if (eventNr <= fEventCount-1) { //watch-out ...event counter starts at 1?
189       cout << fEventCount << " End of file?, Event " << fEventCount << endl;    
190       return false;
191     }
192     //cout <<  "Run " << Int(fgkPosRun, header)<< " , Event " <<eventNr <<endl;
193     
194     // enough space for data?
195     if (fLimit < dataSize) {
196       delete[] fEventData;
197       fEventData = new Char_t[dataSize];
198       fLimit = dataSize;
199     }
200     
201     fDataStream->read(fEventData, dataSize);
202     
203     if(fDataStream->fail()) {
204       cerr << "End of file, Event " << fEventCount; // << endl; 
205         return false;
206     }
207     
208     //cout  << " ...IDs (size) : ";
209     
210     ldcOff=0; // size of data from one DDL link
211     
212     for ( size_t k = 0; k < 2; k++ ) { // 2 LDCs (DDL & VME)
213       
214       ldcSize = Int(ldcOff+fgkPosLength, fEventData); //
215       //ldcSize1=(ldcSize-fgkLdcHeadSize);
216       eqId = Int(ldcOff+fgkLdcHeadSize+fgkEqId, fEventData);        
217       //cout  << eqId <<" ("<<ldcSize<<") ";        
218       
219       evL1 = Int((4+ldcOff+fgkPosRun), fEventData); //ev.nr.
220       if ( evL1 != eventNr ){
221         //cerr << "eqId " <<eqId<<" event nr. mismatch? " << eventNr <<" / "<< evL1 <<" ...LDC data size (header:68) " <<ldcSize<<endl;
222       }
223       
224       if (eqId == 1024) {  //DDL data
225         fDdlOff = ldcOff; //+fgkLdcHeadSize+fgkEquipHeadSize + 32;
226         fDdlSize = ldcSize;
227       }
228       
229       if (eqId == 550) {  //Si-strip data (+QDC)
230         //cout << "550" << endl;
231         fSiOff=ldcOff+fgkLdcHeadSize+fgkEquipHeadSize+fgkPosSiOff;
232         word = Int(fSiOff, fEventData);
233         Short_t lenSi1 = (word >> 16) & 0xffff;
234         Short_t lenSi2 = word & 0xffff;
235         fQdcOff=fSiOff+4*(lenSi1+lenSi2+1)+fgkEquipHeadSize+4; 
236       } 
237       else if (eqId == 1182) {  //QDC first...
238         //cout << "1182" << endl;
239         fQdcOff=ldcOff+fgkLdcHeadSize+fgkEquipHeadSize+fgkPosSiOff;
240         fSiOff=fQdcOff+fgkEquipHeadSize+4;
241       }
242       
243       ldcOff=ldcSize;
244       
245     }
246     //cout << endl;
247
248     //cout << "DDL = " << fDdlOff << endl;
249     // cout << "Si  = " << fSiOff << endl;
250     //cout << "QDC = " << fQdcOff << endl;
251     
252     DecodeSi();
253
254     fEventCount++; //event counter
255     return true;
256 }
257
258 //____________________________________________________________________________
259 Int_t AliTRDtestBeam::DecodeSi() 
260 {
261   //
262   // Decode the silicon detector
263   //
264   
265   if (fSiOff < 0) return 0;
266   
267   // cout << "decoding Si data" << endl;
268
269   Long_t word;
270   
271   word=Int(fSiOff, fEventData);
272   fNSi1 = (word >> 16) & 0xffff;
273   fNSi2 = word & 0xffff;
274   
275   Int_t cSi=fSiOff; //   
276   for (int i = 0; i < fNSi1; i++) {
277     fSi1Address[i] =  ( Int(cSi, fEventData) >> 12 ) & 0x7ff;
278     fSi1Charge[i] = Int(cSi, fEventData)  & 0xfff;
279     cSi+=4;
280   }
281     
282   for (int i = 0; i < fNSi2; i++) {  //1,for Date!
283     fSi2Address[i] =  ( Int(cSi, fEventData) >> 12 ) & 0x7ff;
284     fSi2Charge[i] = Int(cSi, fEventData)  & 0xfff;
285     cSi+=4;
286   }  
287   
288   // reconstruction
289
290   int aLenSiX = 640;
291
292   int amaxX=0;
293   int amaxY=0;
294
295   Int_t q, a;  
296   Int_t nst1=0,nst2=0;
297   Int_t qclX=0,qclY=0, nclX=0,nclY=0, nstX=0,nstY=0;
298   const Int_t kThr = 20;
299
300   nst1=0;
301   nstX=0;
302   nstY=0;
303   nclX=0;
304   nclY=0;
305   qclX=0;
306   qclY=0;
307   
308   for( int i = 0; i < GetNSi1(); i++ ) {
309  
310     if (fSi1Address[i] == 0) continue; // noise
311
312     q = fSi1Charge[i];
313     a = fSi1Address[i];
314
315     if ( q > kThr ) 
316     {
317         if ( i > 0 && i < (GetNSi1()-1) ) {
318
319             if ( (a-fSi1Address[i+1]) == -1 &&
320                  (a-fSi1Address[i-1]) == 1) 
321             {  
322                 nst1++;   
323                 if (a < aLenSiX) {
324                     qclX = q+fSi1Charge[i+1]+fSi1Charge[i-1];
325                     nclX++;
326                     nstX+=3;
327                     amaxX = a;
328                 }
329                 else {
330                     qclY = q+fSi1Charge[i+1]+fSi1Charge[i-1];
331                     nclY++;
332                     nstY+=3;
333                     amaxY = a;
334                 }
335                 i+=1;
336             }
337             else if ( (a-fSi1Address[i-1]) == 1)
338             {  
339                 nst1++;   
340                 if (a < aLenSiX) {
341                     qclX = q+fSi1Charge[i-1];
342                     nclX++;
343                     nstX+=2;
344                     amaxX = a;
345                 }
346                 else {
347                     qclY = q+fSi1Charge[i-1];
348                     nclY++;
349                     nstY+=2;
350                     amaxY = a;
351                 }
352             }
353             else if ( (a-fSi1Address[i+1]) == -1)
354             {  
355                 nst1++;   
356                 if (a < aLenSiX) {
357                     qclX = q+fSi1Charge[i+1];
358                     nclX++;
359                     nstX+=2;
360                     amaxX = a;
361                 }
362                 else {
363                     qclY = q+fSi1Charge[i+1];
364                     nclY++;
365                     nstY+=2;
366                     amaxY = a;
367                 }
368                 i+=1;
369             }
370         }
371     }
372   }
373   if (nst1==2 && nstX<4 && nstY<4 ) {
374       fX[0] = (float)(amaxX*0.05);  // [mm]
375       fY[0] = (float)((amaxY-aLenSiX)*0.05);
376       fQx[0] = (float)qclX;
377       fQy[0] = (float)qclY;
378   }
379   else {
380       fX[0] = -1.;
381       fY[0] = -1.;
382       fQx[0] = 0.;
383       fQy[0] = 0.;
384   }
385   
386   // ...and Si2
387
388   nst2=0;
389   nstX=0;
390   nstY=0;
391   nclX=0;
392   nclY=0;
393   qclX=0;
394   qclY=0;
395
396   for( int i = 0; i < GetNSi2(); i++ ) {
397     
398     if (fSi2Address[i] == 1279) continue; // noise
399     if (fSi2Address[i] == 0) continue;    // noise
400     
401     q = fSi2Charge[i];
402     a = fSi2Address[i];
403
404     if ( q > kThr/2 ) //...as Si2 has 1/2 gain! 
405     {
406         if ( i > 0 && i < (GetNSi2()-1) ) {
407
408             if ( (a-fSi2Address[i+1]) == -1 &&
409                  (a-fSi2Address[i-1]) == 1) 
410             {  
411                 nst2++;   
412                 if (a < aLenSiX) {
413                     qclX = q+fSi2Charge[i+1]+fSi2Charge[i-1];
414                     nclX++;
415                     nstX+=3;
416                     amaxX = a;
417                 }
418                 else {
419                     qclY = q+fSi2Charge[i+1]+fSi2Charge[i-1];
420                     nclY++;
421                     nstY+=3;
422                     amaxY = a;
423                 }
424                 i+=1;
425             }
426             else if ( (a-fSi2Address[i-1]) == 1)
427             {  
428                 nst2++;   
429                 if (a < aLenSiX) {
430                     qclX = q+fSi2Charge[i-1];
431                     nclX++;
432                     nstX+=2;
433                     amaxX = a;
434                 }
435                 else {
436                     qclY = q+fSi2Charge[i-1];
437                     nclY++;
438                     nstY+=2;
439                     amaxY = a;
440                 }
441             }
442             else if ( (a-fSi2Address[i+1]) == -1)
443             {  
444                 nst2++;   
445                 if (a < aLenSiX) {
446                     qclX = q+fSi2Charge[i+1];
447                     nclX++;
448                     nstX+=2;
449                     amaxX = a;
450                 }
451                 else {
452                     qclY = q+fSi2Charge[i+1];
453                     nclY++;
454                     nstY+=2;
455                     amaxY = a;
456                 }
457                 i+=1;
458             }
459         }
460     }
461   }
462   
463   if (nst2==2 && nstX<4 && nstY<4 ) {
464       fX[1] = (float)(amaxX*0.05);  // [mm]
465       fY[1] = (float)((amaxY-aLenSiX)*0.05);
466       fQx[1] = (float)qclX;
467       fQy[1] = (float)qclY;
468   }
469   else {
470       fX[1] = -1.;
471       fY[1] = -1.;
472       fQx[1] = 0.;
473       fQy[1] = 0.;
474   }
475   
476   if (fQdcOff < 0) return 0;
477  
478   word=Int(fQdcOff, fEventData);
479   fPb   = (Double_t)((word >> 16) & 0xFFF);
480   fCher = (Double_t)((word ) & 0xFFF);
481
482   //cout << fCher << " " << fPb << endl;
483   return 1;
484
485 }
486 //____________________________________________________________________________ 
487 AliTRDrawStream *AliTRDtestBeam::GetTRDrawStream() 
488 {
489   //
490   // Get the TRD raw stream
491   //
492   
493   // needs AliTRDrawStream  
494   //cout << "Chamber reader:" << (Int_t)(fEventData+fDdlOff) << " " << fDdlSize << endl;
495   //int ifout = open("dump.dat", O_WRONLY | O_TRUNC | O_CREAT);
496   //write(ifout, (void*)(fEventData+fDdlOff+16), fDdlSize);
497   //close(ifout);
498
499   AliRawReaderMemory *reader = new AliRawReaderMemory((UChar_t*)(fEventData+fDdlOff), (UInt_t)fDdlSize);
500   reader->SetEquipmentID(1024);
501   reader->ReadHeader();
502   //AliTRDrawStream::RawBufferMissAligned(kTRUE);
503
504   AliTRDrawStream::SetNoErrorWarning();
505   AliTRDrawStream::SetExtraWordsFix();
506   AliTRDrawStream::AllowCorruptedData();
507   
508   AliTRDrawStream *tb = new AliTRDrawStream(reader); 
509   //tb->Init();
510   return tb;
511   /*
512     return 
513
514     AliEawReaderMemory *rmem = data->GetRawReader();
515     rmem->ReadHeader();
516     
517     AliTRDrawStream tb(rmem);
518     tb.Init();
519     AliTRDrawStream::SupressWarnings(kTRUE);
520     
521   */
522 }
523
524 //____________________________________________________________________________ 
525 Int_t AliTRDtestBeam::Int(Int_t i, Char_t *start) const
526 {
527   //
528   // ?????
529   //
530   
531   bool swap = kFALSE;
532
533   if(swap) {
534     char *q=(char*)(start+i); 
535     char p[] = {q[3], q[2], q[1], q[0]};
536     return *((int*) p);
537   } else return *((int*)(start+i));
538
539 }
540
541 //____________________________________________________________________________