541223d59fbfd1e5406e50282a30fbcb134ebfb9
[u/mrichter/AliRoot.git] / PHOS / AliPHOSRaw2Digits.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 // Short description  
20 //
21 /*-- Author: Maxim Volkov (RRC KI)
22              Dmitri Peressounko (RRC KI & SUBATECH)
23              Yuri Kharlov (IHEP & SUBATECH)     */
24
25 //////////////////////////////////////////////////////////////////////////////
26
27 // --- ROOT system ---
28 #include "TClonesArray.h"
29 #include "TFile.h"
30 #include "TTree.h"
31
32 // --- Standard library ---
33 //#include <sys/types.h>
34 #include <sys/stat.h>
35 #include <fcntl.h>
36 #include <unistd.h>
37 #include <ctype.h>
38 #include <stdio.h>
39 #include <stdlib.h>
40 #include <unistd.h>
41 #include <string.h>
42 #include <netinet/in.h>
43
44 // --- AliRoot header files ---
45 #include "AliPHOSDigit.h"
46 #include "AliPHOSConTableDB.h"
47 #include "AliPHOSBeamTestEvent.h"
48 #include "AliPHOSRaw2Digits.h"
49 #include "AliPHOSv1.h"
50 #include "../EVGEN/AliGenBox.h"
51 #include "AliRun.h"
52
53 ClassImp(AliPHOSRaw2Digits)
54   
55   
56 //____________________________________________________________________________ 
57   AliPHOSRaw2Digits::AliPHOSRaw2Digits():TTask() 
58 {
59   fInName="";  
60   fMK1 = 0x0123CDEF ;
61   fMK2 = 0x80708070 ;
62   fMK3 = 0x4321ABCD ;
63   fMK4 = 0x80618061 ;
64   fCKW = 0x4640E400 ;
65   fDebug = kFALSE;             //  Debug flag
66   fIsInitialized = kFALSE ;
67   fTarget[0] = 0 ;
68   fTarget[1] = 0 ;
69   fTarget[2] = 0 ;
70   fDigits = 0 ;
71   fPHOSHeader =0 ;
72   fEvent = 0 ;
73   fctdb = 0;
74 }
75 //____________________________________________________________________________ 
76   AliPHOSRaw2Digits::AliPHOSRaw2Digits(const char * filename):TTask("Default","") 
77 {
78   fInName=filename;
79   TString outname(fInName) ;
80   outname.ToLower() ;
81   outname.ReplaceAll(".fz",".root") ;
82   outname.ReplaceAll(".gz","") ;
83   SetTitle(outname) ;
84
85   fMK1 = 0x0123CDEF ;
86   fMK2 = 0x80708070 ;
87   fMK3 = 0x4321ABCD ;
88   fMK4 = 0x80618061 ;
89   fCKW = 0x4640E400 ;
90   fDebug = kFALSE;             //  Debug flag
91   fIsInitialized = kFALSE ;
92   fTarget[0] = 0 ;
93   fTarget[1] = 0 ;
94   fTarget[2] = 0 ;
95   fDigits = 0 ;
96   fPHOSHeader =0 ;
97   fEvent = 0 ;
98   fctdb = 0;
99 }
100
101 //____________________________________________________________________________ 
102 AliPHOSRaw2Digits::~AliPHOSRaw2Digits()
103 {
104   if(fPHOSHeader)
105     fPHOSHeader->Delete() ;
106   if(fDigits){
107     fDigits->Delete() ;
108     delete fDigits ;
109   }
110   
111 }
112 //____________________________________________________________________________ 
113 void AliPHOSRaw2Digits::Exec(Option_t * option){
114   //This is steering method performing all the conversion
115
116   if(!fIsInitialized) //need initialization
117     if(!Init())       //failed to initialize
118       return ;
119
120   ProcessRawFile() ;
121
122   FinishRun() ;
123
124 //____________________________________________________________________________ 
125 Bool_t AliPHOSRaw2Digits::Init(void){
126   //Create PHOS geometry, sets magnetic field to zero, 
127   //create Generator - to store target position, 
128   //opens out file, creates TreeE and make initialization of contaniers
129
130
131   if(fIsInitialized)
132     return kTRUE;
133
134   //Create PHOS
135   new AliPHOSv1("PHOS","GPS2") ;
136
137   //Set Magnetic field
138   gAlice->SetField(0,2);  
139
140   //Set positin of the virtex
141   AliGenBox *gener = new AliGenBox(1);
142   Float_t ox = fTarget[1]; 
143   Float_t oy = fTarget[2]-460.; 
144   Float_t oz = fTarget[0];
145   gener->SetOrigin(ox, oy, oz);
146
147   //  Create the output file
148   TString outname("") ;
149   if(strstr(GetTitle(),"root")){
150     outname=GetTitle();
151   }
152   else{
153     outname = fInName ;
154     outname.ToLower() ;
155     outname.ReplaceAll(".fz",".root") ;
156   }
157   TFile *rootfile = new TFile(outname,"recreate");
158   rootfile->SetCompressionLevel(2);
159
160   // Create the Root Trees
161   gAlice->MakeTree("E") ;
162
163   //Make container for digits
164   fDigits = new TClonesArray("AliPHOSDigit",1000) ;
165
166   //Fill now TreeE
167   fPHOSHeader = new  AliPHOSBeamTestEvent() ;
168   Int_t splitlevel = 0 ;
169   Int_t bufferSize = 32000 ;    
170   TBranch * headerBranch = gAlice->TreeE()->Branch("AliPHOSBeamTestEvent", 
171                                                    "AliPHOSBeamTestEvent", 
172                                                    &fPHOSHeader,bufferSize,splitlevel);
173   headerBranch->SetName("AliPHOSBeamTestEvent") ;
174
175   fIsInitialized = kTRUE ;
176   return kTRUE ;
177 }
178 //____________________________________________________________________________ 
179 Bool_t AliPHOSRaw2Digits::ProcessRawFile(){
180
181   //Method opens zebra file and reads successively bytes from it,
182   //filling corresponding fields in header and digits.
183
184
185   fStatus= -3 ;
186   //First of all, open file and check if it is a zebra file
187
188   char command[256];
189   sprintf(command,"zcat %s",fInName.Data());
190   FILE *dataFile = popen(command, "r");
191   if (dataFile == NULL) {
192     Warning("ProcessRawFile", " Cannot open file %s\n", fInName.Data() ) ;
193     perror(fInName.Data()) ;
194     fStatus = -1 ;
195     return kFALSE ;
196   }
197   printf("Open pipe: %s\n",command);
198
199   // Check if byte ordering is little-endian 
200   UInt_t w = 0x12345678;
201   Int_t swapo = memcmp(&w, "\x78\x56\x34\x12", sizeof(UInt_t)) == 0;
202   if(fDebug)
203     Info("ProcessRawFile", "swapo=%f\n", swapo ) ;
204   
205   
206   UInt_t recBuf[300] ;
207
208   // Read physical record control words 
209   UInt_t nb = 8*sizeof(UInt_t);
210   Int_t n = fread(recBuf, nb, 1, dataFile);
211   if(static_cast<UInt_t>(n) != 1) {
212     if (n < 0 )
213       perror(fInName.Data());
214     else
215       Error("ProcessRawFile", "Could not read physical record control words" ) ;
216     fStatus = -2 ;
217     return kFALSE;
218   }
219   
220   if(fDebug)
221     Info("ProcessRawFile", "recbuf[0] = %d\n", recBuf[0] );
222   
223   // Check if it is a ZEBRA file and if the words are swapped 
224   UInt_t swapi = 0 ;
225   if (recBuf[0] != fMK1) {
226     Swab4(recBuf, &w, 1);
227     if (w != fMK1) {
228       Error("ProcessRawFile", "Does not look like a ZEBRA file\n" ) ;
229       pclose(dataFile) ;
230       fStatus = -2 ;
231       return kFALSE;
232     }
233     swapi=1 ;
234   }
235   
236   if(fDebug){
237     TString message ; 
238     message  = "        w = %f\n" ; 
239     message += "    swapi = %f\n" ; 
240     Info("ProcessRawFile", message.Data(), w, swapi ) ; 
241   }
242   
243   // Get number of words in physical record 
244   UInt_t  nwphr ;
245   if (swapi)
246     Swab4(&recBuf[4],&nwphr,1);
247   else 
248     nwphr = recBuf[4];
249   nwphr*=2; // 1998 -- Now we have 2 records 150 words each 
250
251
252   //start loop over data  
253   // Read rest of record 
254   nb = (nwphr-8)*sizeof(UInt_t);
255   n = fread(&recBuf[8], nb, 1, dataFile) ;
256   if (static_cast<UInt_t>(n) != 1) {
257     if (n < 0 ){
258       perror(fInName.Data());
259       fStatus = -2 ;
260       return kFALSE;
261     }
262   }
263   nb = nwphr *sizeof(UInt_t);
264
265   UInt_t userVector[16] ;
266   UInt_t zheader[12];    
267   UShort_t pattern ;
268   UShort_t scanning[32] ;
269   UShort_t charge[12];
270   UInt_t scaler[12]; 
271   UShort_t tdc2228[32];
272   
273   //read untill the end of file
274   fEvent=0 ;
275   while(1){
276
277     //    StartNewEvent() ;
278     fDigits->Delete() ;
279     gAlice->SetEvent(fEvent) ;
280           
281     Int_t i ;
282     for(i=0;i<16;i++)
283       userVector[i]=*(recBuf+21+i);
284     if(!swapi)
285       Swab4(userVector, userVector, 16);     
286     fPHOSHeader->SetUserVector(userVector) ;
287     
288     
289     // ZEBRA event header
290     for(i=0;i<12;i++)
291       zheader[i]=*(recBuf+47+i);
292     if(swapi)
293       Swab4(zheader, zheader, 12);
294     fPHOSHeader->SetHeader(zheader) ;
295     
296     // Swap input 
297     if (swapi)
298       Swab4(recBuf, recBuf, nwphr);
299     
300     /* Physical record control words */
301     UInt_t * recptr = recBuf;  //Pointer to current position
302
303     if(recptr[7] != 1) {
304       Error("ProcessRawFile", "Cannot handle fast blocks" ) ; 
305       fStatus = -2 ;
306       return kFALSE;
307     }    
308     recptr += 8;
309     
310     // Logical record control words   
311     UInt_t lrtyp = recptr[1];
312     if (lrtyp != 3) {
313       Error("ProcessRawFile", "Can not handle logical record type %d", lrtyp ) ;
314       fStatus = -2 ;
315       return kFALSE;
316     }
317     
318     recptr += 2;
319     if (recptr[0] != fCKW) {
320       Error("ProcessRawFile", "Bad check word" ) ;
321       fStatus = -2 ;
322       return kFALSE;
323     }
324     
325     UInt_t  nwuh = recptr[9];
326     recptr += 10+nwuh;
327     
328     // Bank system words 
329     UInt_t nd = recptr[8];            /* Number of data words */
330     recptr += 10;
331     
332     // Data words 
333     UInt_t evtno = recptr[2];                   /* Event number */
334     
335     if(fDebug)
336        Info("ProcessRawFile", "evtno= %d", evtno);
337     
338     UInt_t nh = recptr[4];                   /* Number of header words in data bank */
339     recptr += nh;
340     
341     // Unswap data from VME 
342     if (swapi)
343       Swab4(recptr, recptr, nd-nh-3);
344     
345     // Give buffer to monitor program 
346     //  UInt_t esize = nd-nh-3;
347     //    if (swapo)
348     //       Swab2(recptr, recptr, esize); 
349     // Two byte data are correct after this operation. 
350     //But we're in trouble if the data array contains 4 byte data!
351     
352     // From now on deal with VME data (MSB first, or network byte order).
353     
354     
355     // Fill the event with data from ADCs
356     UChar_t *byteptr=(UChar_t*)recptr;
357     
358     // Trigger bit register  
359     pattern=ntohs(*(UShort_t*)byteptr);
360     fPHOSHeader->SetPattern(pattern) ;
361     byteptr+=sizeof(UShort_t);
362     
363     // Either peak ADCs, 10 modulesX8=80 channels, 
364     //or Kurchatov 64+2 channel ADC 
365     //(the rest of the channels padded with 0xffff) 
366     for(i=0;i<80;i++){
367       Int_t peak = static_cast<Int_t>(ntohs(*(UShort_t*)byteptr));
368       //make digit
369       Int_t absID = fctdb->Raw2AbsId(i) ;
370       if(absID > 0)
371         new((*fDigits)[i])AliPHOSDigit(-1,absID,peak,0.,i) ;
372       if(fDebug){
373         if(peak>(UShort_t)1000)
374           Info("ProcessRawFile", "event= %d peak[%d] = %f", fEvent, i, peak);
375       }
376       byteptr+=sizeof(UShort_t);
377     }
378     
379     // Scanning ADCs, 4 modulesX8=32 channels
380     for(i=0;i<32;i++){
381       scanning[i]=ntohs(*(UShort_t*)byteptr);
382       byteptr+=sizeof(UShort_t);
383     }
384     fPHOSHeader->SetScanning(scanning) ;
385     
386     // Charge ADCs, 1 moduleX12=12 channels
387     for(i=0;i<12;i++){
388       charge[i]=ntohs(*(UShort_t*)byteptr);
389       byteptr+=sizeof(UShort_t);
390     }
391     fPHOSHeader->SetCharge(charge) ;
392     
393     // Scalers, 1 moduleX12=12 (4 byte) channels
394     for(i=0;i<12;i++){
395       scaler[i]=ntohl(*(UInt_t*)byteptr);
396       byteptr+=sizeof(UInt_t);
397     }
398     fPHOSHeader->SetScaler(scaler) ;
399     
400     // LeCroy TDC 2228A, 4 moduleX8=32 channels
401     for(i=0;i<8;i++){
402       tdc2228[i]=ntohs(*(UShort_t*)byteptr);
403       byteptr+=sizeof(UShort_t);
404     }
405     fPHOSHeader->SetTDC(tdc2228) ;
406
407     WriteDigits() ;
408     if(fDebug)
409       Info("ProcessRawFile", "event= %d written", fEvent) ;
410  
411     // Read next record 
412     UInt_t nb = nwphr *sizeof(UInt_t);
413     n = fread( recBuf, nb,1,dataFile);
414     if (n < 0 ){
415       perror(fInName);
416       fStatus = -2 ;
417       return kFALSE;
418     }
419     if (static_cast<UInt_t>(n) != 1) {
420       pclose(dataFile) ;
421       fStatus = 1 ;
422       return kTRUE ; //all read
423     }
424     fEvent++ ;
425   }
426   
427   fStatus = 1 ;  
428   return kTRUE ;  
429 }
430 //____________________________________________________________________________ 
431 void AliPHOSRaw2Digits::Swab4(void *from, void *to, size_t nwords){
432   // The function swaps 4 bytes: byte#3<-->byte#0, byte#2<-->byte#1 
433   register char *pf=static_cast<char*>(from) ;
434   register char *pt=static_cast<char*>(to) ;
435   register char c;
436   while (nwords-- > 0 ) {
437     c = pf[0];
438     pt[0] = pf[3];
439     pt[3] = c;
440     c = pf[1];
441     pt[1] = pf[2];
442     pt[2] = c;
443     pf += 4;
444     pt += 4;
445   }
446 }
447
448 //____________________________________________________________________________ 
449 void AliPHOSRaw2Digits::Swab2(void *from, void *to, size_t nwords)
450 { //The function swaps 2x2 bytes: byte#0<-->byte#1, byte#2<-->byte#3 
451   register char *pf=static_cast<char*>(from) ;
452   register char *pt=static_cast<char*>(to);
453   register char c;   
454   while (nwords-- > 0 ) {
455     c = pf[0];
456     pt[0] = pf[1];
457     pt[1] = c;
458     c = pf[2];
459     pt[2] = pf[3];
460     pt[3] = c;
461     pf += 4;
462     pt += 4;
463   }
464 }
465
466 //____________________________________________________________________________ 
467 void AliPHOSRaw2Digits::FinishRun(){
468   //Write geometry and header tree
469   gAlice->Write(0,TObject::kOverwrite);
470   gAlice->TreeE()->Write(0,TObject::kOverwrite);
471   
472 }
473 //____________________________________________________________________________ 
474 void AliPHOSRaw2Digits::WriteDigits(void){
475   //In this method we create TreeD, write digits and Raw2Digits to it
476   // and write Header to TreeE. Finally we write TreeD to root file 
477   
478   //Start from Digits
479   fDigits->Sort() ;
480   fDigits->Expand(fDigits->GetEntriesFast()) ;
481   for(Int_t i=0;i<fDigits->GetEntriesFast(); i++)
482     static_cast<AliPHOSDigit*>(fDigits->At(i))->SetIndexInList(i) ;
483
484   char hname[30];
485   sprintf(hname,"TreeD%d",fEvent);
486   TTree * treeD = new TTree(hname,"Digits");
487   //treeD->Write(0,TObject::kOverwrite);
488   
489   // -- create Digits branch
490   Int_t bufferSize = 32000 ;    
491   TBranch * digitsBranch = treeD->Branch("PHOS",&fDigits,bufferSize);
492   digitsBranch->SetTitle("Default");
493   
494   // -- Create Digitizer branch
495   Int_t splitlevel = 0 ;
496   const AliPHOSRaw2Digits * d = this ;
497   TBranch * digitizerBranch = treeD->Branch("AliPHOSRaw2Digits", 
498                                             "AliPHOSRaw2Digits", &d,bufferSize,splitlevel); 
499   digitizerBranch->SetTitle("Default");
500   
501   digitsBranch->Fill() ;
502   digitizerBranch->Fill() ; 
503   treeD->Write(0,TObject::kOverwrite);
504  
505   delete treeD ;
506
507   //Write header
508   gAlice->TreeE()->Fill();
509 }
510 //____________________________________________________________________________ 
511 void AliPHOSRaw2Digits::Print(Option_t * option)const{
512   
513   TString message ;
514   message  = "----------AliPHOSRaw2Digits---------- \n" ;
515   message += "Input stream \n" ;
516   message += "Current input  File: %s\n" ; 
517   message += "Current output File: %s\n" ; 
518   message += "Events processes in the last file %d\n" ; 
519   message += "Input file status\n" ;
520   switch (fStatus){
521   case 0: message += "`Have not processed yet'\n" ;
522     break ;
523   case 1: message += "`Processed normally'\n" ;
524     break ;
525   case -1: message += "`File not found'\n" ;
526     break ;
527   case -2: message += "`Error in reading'\n" ;
528     break ;
529   case -3: message += "'Interupted'\n" ;
530   default: ;
531   }
532   message += "Connection table: " ;
533   if(fctdb)
534     message += "%s %s \n" ; 
535   else
536     message += " no DB \n" ;
537
538   Info("Print", message.Data(),  
539        fInName.Data(), 
540        GetTitle(), 
541        fEvent, 
542        fctdb->GetName(), fctdb->GetTitle() ) ; 
543 }