1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //_________________________________________________________________________
19 // Class designed to convert raw data to aliroot format. We assume, that
20 // prototype is situated in the center of 3 PHOS module and convert prototype
21 // outpur to AliPHOSDigits. In addition, we fill branch of TreeE with
22 // AliPHOSBeamTestEvent, contaning description of event(triggers etc).
23 // Note, that one byte per channel in raw data is transvormed to class
24 // AliPHOSDigit, so finale zise increase ~ 100 times. So output can be split
25 // into peases of reasonable size: each file can not contain more than
26 // fMaxPerFile: if there are still events in raw file, then new directory
27 // is created and header+digits files are written to it.
30 // AliPHOSRaw2Digits * r = new AliPHOSRaw2Digits("path/input.file") ;
31 // //note, that it can be gzipped file!
32 // //Set position of the target in the given run.
33 // //Z axis along beam direction, from target to prototype (0-surface of prototype)
34 // //X axis along columns of prototype (0-center of prototype)
35 // //Y axis along raws of prototype (0-center of prototype)
36 // Double_t pos[3]={0,0,-452.} ;
37 // r->SetTargetPosition(pos) ;
38 // //Read/create connection Table:
39 // TFile f("ConTableDB.root") ;
40 // AliPHOSConTableDB * cdb = f.Get("AliPHOSConTableDB") ;
42 // r->SetConTableDB(cdb) ;
45 // As a result files galice.root and PHOS.Digits.root should be produced in
46 // current dir, and, possibly, dirs 1,2,3... each with galice.root and PHOS.Digits.root,
47 // where the rest of data are written.
49 /*-- Author: Maxim Volkov (RRC KI)
50 Dmitri Peressounko (RRC KI & SUBATECH)
51 Yuri Kharlov (IHEP & SUBATECH) */
53 //////////////////////////////////////////////////////////////////////////////
55 // --- ROOT system ---
56 #include "TClonesArray.h"
62 // --- Standard library ---
71 #include <netinet/in.h>
73 // --- AliRoot header files ---
74 #include "AliPHOSDigit.h"
75 #include "AliPHOSConTableDB.h"
76 #include "AliPHOSBeamTestEvent.h"
77 #include "AliPHOSRaw2Digits.h"
78 #include "AliPHOSv1.h"
79 #include "../EVGEN/AliGenBox.h"
82 ClassImp(AliPHOSRaw2Digits)
85 //____________________________________________________________________________
86 AliPHOSRaw2Digits::AliPHOSRaw2Digits():TTask()
94 fDebug = kFALSE; // Debug flag
96 fIsInitialized = kFALSE ;
104 fMaxPerFile = 20000 ;
108 //____________________________________________________________________________
109 AliPHOSRaw2Digits::AliPHOSRaw2Digits(const char * filename,Bool_t toSplit):TTask("Default","")
113 TString outname("") ;
115 outname = "galice.root" ;
119 outname.ReplaceAll(".fz",".root") ;
120 outname.ReplaceAll(".gz","") ;
129 fDebug = kFALSE; // Debug flag
130 fIsInitialized = kFALSE ;
138 fMaxPerFile = 20000 ;
143 //____________________________________________________________________________
144 AliPHOSRaw2Digits::~AliPHOSRaw2Digits()
147 fPHOSHeader->Delete() ;
154 //____________________________________________________________________________
155 void AliPHOSRaw2Digits::Exec(Option_t * option){
156 //This is steering method performing all the conversion
158 if(!fIsInitialized) //need initialization
159 if(!Init()) //failed to initialize
166 //____________________________________________________________________________
167 Bool_t AliPHOSRaw2Digits::Init(void){
168 //Makes initialization of contaniers
173 //Make container for digits
174 fDigits = new TClonesArray("AliPHOSDigit",1000) ;
175 fPHOSHeader = new AliPHOSBeamTestEvent() ;
176 fIsInitialized = kTRUE ;
177 return StartRootFiles() ;
180 //____________________________________________________________________________
181 Bool_t AliPHOSRaw2Digits::StartRootFiles(void ){
182 //Create PHOS geometry, sets magnetic field to zero,
183 //create Generator - to store target position,
184 //opens out file, creates TreeE
186 //create gAlice if nececcary
188 new AliRun("gAlice","The ALICE Off-line Simulation Framework") ;
191 if(!gAlice->GetModule("PHOS"))
192 new AliPHOSv1("PHOS","GPS2") ;
195 gAlice->SetField(0,2);
197 //Set positin of the virtex
198 AliGenerator * gener = gAlice->Generator() ;
200 gener = new AliGenBox(1);
201 Float_t ox = fTarget[1];
202 Float_t oy = fTarget[2]-460.;
203 Float_t oz = fTarget[0];
204 gener->SetOrigin(ox, oy, oz);
207 Int_t nRootFile = (fEvent+1)/fMaxPerFile ;
210 sprintf(dname,"%d",nRootFile) ;
211 if(gSystem->AccessPathName(dname)) //strange return: 0 if exists
212 if(gSystem->MakeDirectory(dname)!=0)
213 Fatal("StartRootFiles","Can not make directory %s \n",dname) ;
215 if(!gSystem->ChangeDirectory(dname))
216 Fatal("StartRootFiles","Can not cd to %s\n",dname) ;
219 // Create the output file
220 TString outname("") ;
221 if(strstr(GetTitle(),"root")){
227 outname.ReplaceAll(".fz",".root") ;
230 fHeaderFile = new TFile(outname,"recreate");
231 fHeaderFile->SetCompressionLevel(2);
233 // Create the Root Trees
234 gAlice->MakeTree("E") ;
237 Int_t splitlevel = 0 ;
238 Int_t bufferSize = 32000 ;
239 TBranch * headerBranch = gAlice->TreeE()->Branch("AliPHOSBeamTestEvent",
240 "AliPHOSBeamTestEvent",
241 &fPHOSHeader,bufferSize,splitlevel);
242 headerBranch->SetName("AliPHOSBeamTestEvent") ;
245 fDigitsFile = new TFile("PHOS.Digits.root","recreate") ;
246 fDigitsFile->SetCompressionLevel(2) ;
250 //____________________________________________________________________________
251 Bool_t AliPHOSRaw2Digits::CloseRootFiles(void ){
252 //cleans everething to start next root file
257 fHeaderFile->Close() ;
263 fDigitsFile->Close() ;
268 Int_t nRootFile = (fEvent-1)/fMaxPerFile ;
270 if(!gSystem->ChangeDirectory("../")){
271 Fatal("CloseRootFile","Can not return to initial dir \n") ;
277 //____________________________________________________________________________
278 Bool_t AliPHOSRaw2Digits::ProcessRawFile(){
280 //Method opens zebra file and reads successively bytes from it,
281 //filling corresponding fields in header and digits.
285 //First of all, open file and check if it is a zebra file
288 sprintf(command,"zcat %s",fInName.Data());
289 FILE *dataFile = popen(command, "r");
290 if (dataFile == NULL) {
291 Warning("ProcessRawFile", " Cannot open file %s\n", fInName.Data() ) ;
292 perror(fInName.Data()) ;
297 // Check if byte ordering is little-endian
298 UInt_t w = 0x12345678;
299 Int_t swapo = memcmp(&w, "\x78\x56\x34\x12", sizeof(UInt_t)) == 0;
301 Info("ProcessRawFile", "swapo=%f\n", swapo ) ;
306 // Read physical record control words
307 UInt_t nb = 8*sizeof(UInt_t);
308 Int_t n = fread(recBuf, nb, 1, dataFile);
309 if(static_cast<UInt_t>(n) != 1) {
311 perror(fInName.Data());
313 Error("ProcessRawFile", "Could not read physical record control words" ) ;
319 Info("ProcessRawFile", "recbuf[0] = %d\n", recBuf[0] );
321 // Check if it is a ZEBRA file and if the words are swapped
323 if (recBuf[0] != fMK1) {
324 Swab4(recBuf, &w, 1);
326 Error("ProcessRawFile", "Does not look like a ZEBRA file\n" ) ;
336 message = " w = %f\n" ;
337 message += " swapi = %f\n" ;
338 Info("ProcessRawFile", message.Data(), w, swapi ) ;
341 // Get number of words in physical record
344 Swab4(&recBuf[4],&nwphr,1);
347 nwphr*=2; // 1998 -- Now we have 2 records 150 words each
350 //start loop over data
351 // Read rest of record
352 nb = (nwphr-8)*sizeof(UInt_t);
353 n = fread(&recBuf[8], nb, 1, dataFile) ;
354 if (static_cast<UInt_t>(n) != 1) {
356 perror(fInName.Data());
361 nb = nwphr *sizeof(UInt_t);
363 UInt_t userVector[16] ;
366 UShort_t scanning[32] ;
369 UShort_t tdc2228[32];
371 //read untill the end of file
377 if((fEvent%fMaxPerFile == 0) && fEvent ){
381 gAlice->SetEvent(fEvent%fMaxPerFile) ;
384 fPHOSHeader->SetBeamEnergy(fBeamEnergy) ;
388 userVector[i]=*(recBuf+21+i);
390 Swab4(userVector, userVector, 16);
391 fPHOSHeader->SetUserVector(userVector) ;
394 // ZEBRA event header
396 zheader[i]=*(recBuf+47+i);
398 Swab4(zheader, zheader, 12);
399 fPHOSHeader->SetHeader(zheader) ;
403 Swab4(recBuf, recBuf, nwphr);
405 /* Physical record control words */
406 UInt_t * recptr = recBuf; //Pointer to current position
409 Error("ProcessRawFile", "Cannot handle fast blocks" ) ;
415 // Logical record control words
416 UInt_t lrtyp = recptr[1];
418 Error("ProcessRawFile", "Can not handle logical record type %d", lrtyp ) ;
424 if (recptr[0] != fCKW) {
425 Error("ProcessRawFile", "Bad check word" ) ;
430 UInt_t nwuh = recptr[9];
434 UInt_t nd = recptr[8]; /* Number of data words */
438 UInt_t evtno = recptr[2]; /* Event number */
441 Info("ProcessRawFile", "evtno= %d", evtno);
443 UInt_t nh = recptr[4]; /* Number of header words in data bank */
446 // Unswap data from VME
448 Swab4(recptr, recptr, nd-nh-3);
450 // Give buffer to monitor program
451 // UInt_t esize = nd-nh-3;
453 // Swab2(recptr, recptr, esize);
454 // Two byte data are correct after this operation.
455 //But we're in trouble if the data array contains 4 byte data!
457 // From now on deal with VME data (MSB first, or network byte order).
460 // Fill the event with data from ADCs
461 UChar_t *byteptr=(UChar_t*)recptr;
463 // Trigger bit register
464 pattern=net2host(*(UShort_t*)byteptr);
465 fPHOSHeader->SetPattern(pattern) ;
466 byteptr+=sizeof(UShort_t);
468 // Either peak ADCs, 10 modulesX8=80 channels,
469 //or Kurchatov 64+2 channel ADC
470 //(the rest of the channels padded with 0xffff)
472 Int_t peak = static_cast<Int_t>(net2host(*(UShort_t*)byteptr));
474 Int_t absID = fctdb->Raw2AbsId(i) ;
476 new((*fDigits)[i])AliPHOSDigit(-1,absID,peak,0.,i) ;
478 if(peak>(UShort_t)1000)
479 Info("ProcessRawFile", "event= %d peak[%d] = %f", fEvent, i, peak);
481 byteptr+=sizeof(UShort_t);
484 // Scanning ADCs, 4 modulesX8=32 channels
486 scanning[i]=net2host(*(UShort_t*)byteptr);
487 byteptr+=sizeof(UShort_t);
489 fPHOSHeader->SetScanning(scanning) ;
491 // Charge ADCs, 1 moduleX12=12 channels
493 charge[i]=net2host(*(UShort_t*)byteptr);
494 byteptr+=sizeof(UShort_t);
496 fPHOSHeader->SetCharge(charge) ;
498 // Scalers, 1 moduleX12=12 (4 byte) channels
500 scaler[i]=net2host(*(UInt_t*)byteptr);
501 byteptr+=sizeof(UInt_t);
503 fPHOSHeader->SetScaler(scaler) ;
505 // LeCroy TDC 2228A, 4 moduleX8=32 channels
507 tdc2228[i]=net2host(*(UShort_t*)byteptr);
508 byteptr+=sizeof(UShort_t);
510 fPHOSHeader->SetTDC(tdc2228) ;
514 Info("ProcessRawFile", "event= %d written", fEvent) ;
517 UInt_t nb = nwphr *sizeof(UInt_t);
518 n = fread( recBuf, nb,1,dataFile);
524 if (static_cast<UInt_t>(n) != 1) {
527 return kTRUE ; //all read
535 //____________________________________________________________________________
536 void AliPHOSRaw2Digits::Swab4(void *from, void *to, size_t nwords){
537 // The function swaps 4 bytes: byte#3<-->byte#0, byte#2<-->byte#1
538 register char *pf=static_cast<char*>(from) ;
539 register char *pt=static_cast<char*>(to) ;
541 while (nwords-- > 0 ) {
553 //____________________________________________________________________________
554 void AliPHOSRaw2Digits::Swab2(void *from, void *to, size_t nwords)
555 { //The function swaps 2x2 bytes: byte#0<-->byte#1, byte#2<-->byte#3
556 register char *pf=static_cast<char*>(from) ;
557 register char *pt=static_cast<char*>(to);
559 while (nwords-- > 0 ) {
571 //____________________________________________________________________________
572 void AliPHOSRaw2Digits::FinishRun(){
573 //Write geometry and header tree
574 gAlice->Write(0,TObject::kOverwrite);
575 gAlice->TreeE()->Write(0,TObject::kOverwrite);
578 //____________________________________________________________________________
579 void AliPHOSRaw2Digits::WriteDigits(void){
580 //In this method we create TreeD, write digits and Raw2Digits to it
581 // and write Header to TreeE. Finally we write TreeD to root file
585 fDigits->Expand(fDigits->GetEntriesFast()) ;
586 for(Int_t i=0;i<fDigits->GetEntriesFast(); i++)
587 static_cast<AliPHOSDigit*>(fDigits->At(i))->SetIndexInList(i) ;
590 sprintf(hname,"TreeD%d",fEvent%fMaxPerFile);
591 TTree * treeD = new TTree(hname,"Digits");
592 //treeD->Write(0,TObject::kOverwrite);
594 // -- create Digits branch
595 Int_t bufferSize = 32000 ;
596 TBranch * digitsBranch = treeD->Branch("PHOS",&fDigits,bufferSize);
597 digitsBranch->SetTitle("Default");
599 // -- Create Digitizer branch
600 Int_t splitlevel = 0 ;
601 const AliPHOSRaw2Digits * d = this ;
602 TBranch * digitizerBranch = treeD->Branch("AliPHOSRaw2Digits",
603 "AliPHOSRaw2Digits", &d,bufferSize,splitlevel);
604 digitizerBranch->SetTitle("Default");
608 digitsBranch->Fill() ;
609 digitizerBranch->Fill() ;
610 treeD->Write(0,TObject::kOverwrite);
616 gAlice->TreeE()->Fill();
618 //____________________________________________________________________________
619 void AliPHOSRaw2Digits::Print(Option_t * option)const{
622 printf("----------AliPHOSRaw2Digits---------- \n") ;
623 printf("Current input File: %s\n",fInName.Data()) ;
624 printf("Current output File: %s\n", GetTitle());
625 printf("Events processes in the last file %d\n",fEvent) ;
626 printf("Input file status\n") ;
628 case 0: printf("`Have not processed yet'\n") ;
630 case 1: printf("`Processed normally'\n") ;
632 case -1: printf("`File not found'\n") ;
634 case -2: printf("`Error in reading'\n") ;
636 case -3: printf("'Interupted'\n") ;
639 printf("Connection table: " );
641 printf("%s %s \n",fctdb->GetName(), fctdb->GetTitle() ) ;
643 printf(" no DB \n" );