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 "AliTRDrawStreamTB.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 AliTRDrawStreamTB *AliTRDtestBeam::GetTRDrawStream()
 488 {
 489   //
 490   // Get the TRD raw stream
 491   //
 492
 493   // needs AliTRDrawStreamTB
 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   //AliTRDrawStreamTB::RawBufferMissAligned(kTRUE);
 503
 504   AliTRDrawStreamTB::SetNoErrorWarning();
 505   AliTRDrawStreamTB::SetExtraWordsFix();
 506   AliTRDrawStreamTB::AllowCorruptedData();
 507
 508   AliTRDrawStreamTB *tb = new AliTRDrawStreamTB(reader);
 509   //tb->Init();
 510   return tb;
 511   /*
 512     return
 513
 514     AliEawReaderMemory *rmem = data->GetRawReader();
 515     rmem->ReadHeader();
 516
 517     AliTRDrawStreamTB tb(rmem);
 518     tb.Init();
 519     AliTRDrawStreamTB::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 //____________________________________________________________________________