Ssemi-final version of TRD raw data simulation
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 9 Aug 2007 10:34:04 +0000 (10:34 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 9 Aug 2007 10:34:04 +0000 (10:34 +0000)
code with zero suppression similar to TRD FEE.
It contains the following classes:

AliTRDfeeParam
AliTRDmcmSim
AliTRDrawData

AliTRDfeeParam has been modified to have more parameters
like raw data production version and so on.

AliTRDmcmSim is new class and this is the core of MCM
(PASA+TRAP) simulator. It has still very simple function
and it will be a half year project to improve this to make
it closer to the reall FEE.

AliTRDrawData has been modified to use new class AliTRDmcmSim.

These modifications were tested on Aug. 02 HEAD
version that code itself compiles. I'm sure there must
be still bugs and we need testing by as many as possible
persons now. Especially it seems HLT part is impacted by
problems because some parameters were moved from
AliTRDrawData to AliTRDfeeParam (like fRawVersion
disappeared from AliTRDrawData).

In TRD definition, we have now 4 raw data versions.

0 very old offline version (by Bogdan)
1 test version (no zero suppression)
2 final version (no zero suppression)
3 test version (with zero suppression)

The default is still set to 2 in AliTRDfeeParam::fgkRAWversion
and it uses previously existing codes. If you set this to 3,
AliTRDrawData changes behavior to use my new class with zero
suppression, but might crash yet.

Plan is after we make sure it works stably, we delete AliTRDmcm
which is obsolete. However it still take time because tarcklet
part is not yet touched.

The default raw version is 2.

Ken Oyama

TRD/AliTRDfeeParam.cxx
TRD/AliTRDfeeParam.h
TRD/AliTRDmcmSim.cxx [new file with mode: 0644]
TRD/AliTRDmcmSim.h [new file with mode: 0644]
TRD/AliTRDrawData.cxx
TRD/AliTRDrawData.h
TRD/TRDbaseLinkDef.h
TRD/libTRDbase.pkg

index 2da453c..f4b540a 100644 (file)
@@ -71,27 +71,27 @@ AliTRDfeeParam::AliTRDfeeParam()
   :TObject()
   //  ,fGeo(0)
   ,fCP(0)
-  ,fTFaR1(0)
-  ,fTFaR2(0)
-  ,fTFaC1(0)
-  ,fTFaC2(0)
+  ,fTFr1(0)
+  ,fTFr2(0)
+  ,fTFc1(0)
+  ,fTFc2(0)
 {
   //
   // Default constructor
   //
   
   // PASA V.4
-  if      (fgkTFaNExp == 1) {
-    fTFaR1 = 1.1563;
-    fTFaR2 = 0.1299;
-    fTFaC1 = 0.0657;
-    fTFaC2 = 0.0000;
+  if      (fgkTFnExp == 1) {
+    fTFr1 = 1.1563;
+    fTFr2 = 0.1299;
+    fTFc1 = 0.0657;
+    fTFc2 = 0.0000;
   }
-  else if (fgkTFaNExp == 2) {
-    fTFaR1 = 1.1563;
-    fTFaR2 = 0.1299;
-    fTFaC1 = 0.1141;
-    fTFaC2 = 0.6241;
+  else if (fgkTFnExp == 2) {
+    fTFr1 = 1.1563;
+    fTFr2 = 0.1299;
+    fTFc1 = 0.1141;
+    fTFc2 = 0.6241;
   }
 
   //  fGeo = AliTRDgeometry::Instance();
@@ -104,10 +104,10 @@ AliTRDfeeParam::AliTRDfeeParam(const AliTRDfeeParam &p)
   :TObject(p)
   //  ,fGeo(p.fGeo)
   ,fCP(p.fCP)
-  ,fTFaR1(p.fTFaR1)
-  ,fTFaR2(p.fTFaR2)
-  ,fTFaC1(p.fTFaC1)
-  ,fTFaC2(p.fTFaC2)
+  ,fTFr1(p.fTFr1)
+  ,fTFr2(p.fTFr2)
+  ,fTFc1(p.fTFc1)
+  ,fTFc2(p.fTFc2)
 {
   //
   // AliTRDfeeParam copy constructor
@@ -144,10 +144,10 @@ void AliTRDfeeParam::Copy(TObject &p) const
 
   //  ((AliTRDfeeParam &) p).fGeo     = fGeo;
   ((AliTRDfeeParam &) p).fCP      = fCP;
-  ((AliTRDfeeParam &) p).fTFaR1   = fTFaR1;
-  ((AliTRDfeeParam &) p).fTFaR2   = fTFaR2;
-  ((AliTRDfeeParam &) p).fTFaC1   = fTFaC1;
-  ((AliTRDfeeParam &) p).fTFaC2   = fTFaC2;
+  ((AliTRDfeeParam &) p).fTFr1   = fTFr1;
+  ((AliTRDfeeParam &) p).fTFr2   = fTFr2;
+  ((AliTRDfeeParam &) p).fTFc1   = fTFc1;
+  ((AliTRDfeeParam &) p).fTFc2   = fTFc2;
   
   TObject::Copy(p);
 }
index 96b632b..92c116e 100644 (file)
@@ -54,6 +54,17 @@ class AliTRDfeeParam : public TObject
   static  Int_t    GetNrowC0()            { return fgkNrowC0; }
   static  Int_t    GetNrowC1()            { return fgkNrowC1; }
 
+  static  Int_t    GetADCpedestal()       { return fgkADCpedestal; }
+  static  Int_t    GetADCnoise()          { return fgkADCnoise; }
+  static  Int_t    GetADCDAC()            { return fgkADCDAC; }
+
+  static  Bool_t   isPFon()               { return fgkPFon; }
+  static  Bool_t   isGFon()               { return fgkGFon; }
+  static  Bool_t   isTFon()               { return fgkTFon; }
+
+  static  Int_t    GetPFtimeConstant()    {  return fgkPFtimeConstant; }
+  static  Int_t    GetPFeffectPedestal()  {  return fgkPFeffectPedestal; }
+
   //          Float_t  GetClusThr()           { return fClusThr; };
   //        Float_t  GetPadThr() const { return fPadThr; };
   //        Int_t    GetTailCancelation() const { return fTCOn; };
@@ -61,8 +72,23 @@ class AliTRDfeeParam : public TObject
   //virtual void     GetFilterParam(Float_t &r1, Float_t &r2, Float_t &c1, Float_t &c2, Float_t &ped) const;
   //        Int_t    GetFilterType() const { return fFilterType; };
 
-  static  Float_t  GetTFattenuationParam() { return ((Float_t)fgkTFattenuationParameter1) / ((Float_t)fgkTFattenuationParameter2) ; }
-  static  Float_t  GetTFf0()               { return 1 + fgkTFon*(-1+GetTFattenuationParam()); }   // 1 if TC off
+  static  Int_t    GetTFtype()            { return fgkTFtype; }
+  static  Int_t    GetTFnExp()            { return fgkTFnExp; }
+          Float_t  GetTFr1()              { return fTFr1; }
+          Float_t  GetTFr2()              { return fTFr2; }
+          Float_t  GetTFc1()              { return fTFc1; }
+          Float_t  GetTFc2()              { return fTFc2; }
+
+  static  Float_t  GetTFattPar()          { return ((Float_t)fgkTFattPar1) / ((Float_t)fgkTFattPar2) ; }
+  static  Float_t  GetTFf0()              { return 1 + fgkTFon*(-1+GetTFattPar()); }   // 1 if TC off
+  
+  static  Int_t    GetEBsglIndThr()       { return fgkEBsglIndThr; }
+  static  Int_t    GetEBsumIndThr()       { return fgkEBsumIndThr; }
+  static  Int_t    GetEBindLUT()          { return fgkEBindLUT; }
+  static  Int_t    GetEBignoreNeighbour() { return fgkEBignoreNeighbour; }
+
+  static  Int_t    GetRAWversion()        { return fgkRAWversion; }
+  static  Bool_t   GetRAWstoreRaw()       { return fgkRAWstoreRaw; }
 
  protected:
 
@@ -89,7 +115,7 @@ class AliTRDfeeParam : public TObject
   static const Int_t    fgkNrowC1            = 16;        // Number of Rows per C1 chamber  (old fgkRowmaxC1)
 
   // ADC intrinsic parameters
-  static const Int_t    fgkADCpedestal       = 100000;    // ADC baseline * 100 (old name fPedestal)
+  static const Int_t    fgkADCpedestal       = 0;         // ADC baseline * 100 (old name fPedestal) it is not PFeffectPedestal
   static const Int_t    fgkADCnoise          = 10;        // ADC noise    * 100 (not contained in the digitizer) [in ADC] 
   static const Int_t    fgkADCDAC            = 0;         // 5 bit ADC gain parameter
 
@@ -106,40 +132,44 @@ class AliTRDfeeParam : public TObject
   static const Int_t    fgkGFnoise           =  0;        // Noise level increased by gain filter x 100 [in ADC] (to be measured)
 
   // TF setup
-  static const Int_t    fgkTFtype                 = 2;    // TC type (0=analog, 1=digital, 2=TRAPsim) (old name fFilterType)
-  static const Int_t    fgkTFlongDecayWeight      = 270;  // 0 to 1024 corresponds to 0 to 0.5
-  static const Int_t    fgkTFlongDecayParameter   = 348;  // 0 to 511 corresponds to 0.75 to 1
-  static const Int_t    fgkTFshortDecayParameter  = 449;  // 0 to 511 correponds to 0.25 to 0.5
-  static const Int_t    fgkTFattenuationParameter1= 45;   // attenuationParameter = fgkTFattenuationParameter1/fgkTFattenuationParameter2
-  static const Int_t    fgkTFattenuationParameter2= 14;   //                      = -alphaL/ln(lambdaL)-(1-alphaL)/ln(lambdaS)
-
-  // OLD TF setup (calculated from above)  (valid only for fgkTFtype = 0 or 1)
-  static const Int_t    fgkTFaNExp                 = 1;    // Number of exponential
-               Float_t  fTFaR1;                        // Time constant [microseconds] long (old name fR1)
-               Float_t  fTFaR2;                        // Time constant [microseconds] short(old name fR2)
-               Float_t  fTFaC1;                        // Weight long  (old name fC1)
-               Float_t  fTFaC2;                        // Weight short (old name fC2)
-
-  // Zero suppression parameters
-  static const Int_t    fgkEBsingleIndicatorThreshold = 3;    // used in EBIS, in ADC units above the pedestal
-  static const Int_t    fgkEBsumIndicatorThreshold    = 4;    // used in EBIT, in ADC units above the pedestal
-  static const Int_t    fgkEBindicatorLookupTable     = 0xF0; // see the TRAP user manual, used in EBIL
-  static const Int_t    fgkEBmarkIgnoreNeighbour      = 1;    // used in EBIN
+  static const Int_t    fgkTFtype            = 1;    // TC type (0=analog, 1=digital, 2=MI) (old name fFilterType)
+
+  // OLD TF setup (calculated from above)  (valid only for fgkTFsimType = 0 or 1)
+  static const Int_t    fgkTFnExp          = 1;    // Number of exponential for simType 0 and 1
+
+  // following need Instance because initialized in constructor
+               Float_t  fTFr1;                        // Time constant [us] long (old name fR1)
+               Float_t  fTFr2;                        // Time constant [us] short(old name fR2)
+               Float_t  fTFc1;                        // Weight long  (old name fC1)
+               Float_t  fTFc2;                        // Weight short (old name fC2)
+
+  // here is for TRAP simulation (not yet used)
+  static const Int_t    fgkTFdecayWeightL     = 270;  // 0 to 1024 corresponds to 0 to 0.5
+  static const Int_t    fgkTFdecayParL        = 348;  // 0 to 511 corresponds to 0.75 to 1
+  static const Int_t    fgkTFdecayParS        = 449;  // 0 to 511 correponds to 0.25 to 0.5
+  static const Int_t    fgkTFattPar1          = 45;   // attenuationParameter = fgkTFattenuationParameter1/fgkTFattenuationParameter2
+  static const Int_t    fgkTFattPar2          = 14;   //                      = -alphaL/ln(lambdaL)-(1-alphaL)/ln(lambdaS)
+
+  // ZS parameters
+  static const Int_t    fgkEBsglIndThr        = 3;    // EBIS in ADC units
+  static const Int_t    fgkEBsumIndThr        = 4;    // EBIT in ADC units
+  static const Int_t    fgkEBindLUT           = 0xF0; // EBIL lookup table
+  static const Int_t    fgkEBignoreNeighbour  = 1;    // EBIN 0:include neighbor
 
   // Charge accumulators
-  static const Int_t    fgkPREPqAcc0Start             =  0;   // Preprocessor Charge Accumulator 0 Start
-  static const Int_t    fgkPREPqAcc0End               = 10;   // Preprocessor Charge Accumulator 0 End
-  static const Int_t    fgkPREPqAcc1Start             = 11;   // Preprocessor Charge Accumulator 1 Start
-  static const Int_t    fgkPREPqAcc1End               = 20;   // Preprocessor Charge Accumulator 1 End
-  static const Int_t    fgkMinClusterCharge           = 20;   // Hit detection [in ADC units]
+  static const Int_t    fgkPREPqAcc0Start      =  0;   // Preprocessor Charge Accumulator 0 Start
+  static const Int_t    fgkPREPqAcc0End        = 10;   // Preprocessor Charge Accumulator 0 End
+  static const Int_t    fgkPREPqAcc1Start      = 11;   // Preprocessor Charge Accumulator 1 Start
+  static const Int_t    fgkPREPqAcc1End        = 20;   // Preprocessor Charge Accumulator 1 End
+  static const Int_t    fgkMinClusterCharge    = 20;   // Hit detection [in ADC units]
 
   // OLD TRAP processing parameters calculated from above
   //static const Float_t  fClusThr;                     // Cluster threshold
   //static const Float_t  fPadThr;                      // Pad threshold
 
   // For raw production
-  static const Int_t    fgkRAWversion            = 1;         // Raw data production version
-  static const Bool_t   fgkRAWstoreRaw           = kTRUE;     // Store unfiltered data for raw data stream
+  static const Int_t    fgkRAWversion          = 2;     // Raw data production version
+  static const Bool_t   fgkRAWstoreRaw         = kTRUE; // Store unfiltered data for raw data stream
 
  private:
 
diff --git a/TRD/AliTRDmcmSim.cxx b/TRD/AliTRDmcmSim.cxx
new file mode 100644 (file)
index 0000000..4c48a53
--- /dev/null
@@ -0,0 +1,694 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/*
+$Log$
+*/
+
+///////////////////////////////////////////////////////////////////////////////
+//                                                                           //
+//  TRD MCM (Multi Chip Module) simulator                                    //
+//                                                                           //
+///////////////////////////////////////////////////////////////////////////////
+
+#include <TMath.h>
+
+#include "AliLog.h"
+#include "AliTRDmcmSim.h"
+#include "AliTRDfeeParam.h"
+#include "AliTRDgeometry.h"
+
+ClassImp(AliTRDmcmSim)
+
+//_____________________________________________________________________________
+AliTRDmcmSim::AliTRDmcmSim() :TObject()
+  ,fInitialized(kFALSE)
+  ,fFeeParam(NULL)
+  ,fGeo(NULL)
+  ,fChaId(-1)
+  ,fSector(-1)
+  ,fStack(-1)
+  ,fLayer(-1)
+  ,fRobPos(-1)
+  ,fMcmPos(-1)
+  ,fNADC(-1)
+  ,fNTimeBin(-1)
+  ,fRow(-1)
+  ,fColOfADCbeg(-1)
+  ,fColOfADCend(-1)
+  ,fADCR(NULL)
+  ,fADCF(NULL)
+  ,fZSM(NULL)
+  ,fZSM1Dim(NULL)
+{
+  //
+  // AliTRDmcmSim default constructor
+  //
+
+  // By default, nothing is initialized.
+  // It is necessary to issue Init before use.
+}
+
+//_____________________________________________________________________________
+AliTRDmcmSim::~AliTRDmcmSim() 
+{
+  //
+  // AliTRDmcmSim destructor
+  //
+  if( fADCR != NULL ) {
+    for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+      delete fADCR[iadc];
+      delete fADCF[iadc];
+      delete fZSM [iadc];
+    }
+    delete fADCR;
+    delete fADCF;
+    delete fZSM;
+    delete fZSM1Dim;
+  }
+  delete fGeo;
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::Init( Int_t cha_id, Int_t rob_pos, Int_t mcm_pos )
+{
+  // Initialize the class with new geometry information
+  // fADC array will be reused with filled by zero
+
+  fFeeParam      = AliTRDfeeParam::Instance();
+  fGeo           = new AliTRDgeometry();
+  fChaId         = cha_id;
+  fSector        = fGeo->GetSector( fChaId );
+  fStack         = fGeo->GetChamber( fChaId );
+  fLayer         = fGeo->GetPlane( fChaId );
+  fRobPos        = rob_pos;
+  fMcmPos        = mcm_pos;
+  fNADC          = fFeeParam->GetNadcMcm();
+  fNTimeBin      = fFeeParam->GetNtimebin();
+  fRow           = fFeeParam->GetPadRowFromMCM( fRobPos, fMcmPos );
+  fColOfADCbeg   = fFeeParam->GetPadColFromADC( fRobPos, fMcmPos, 0 );
+  fColOfADCend   = fFeeParam->GetPadColFromADC( fRobPos, fMcmPos, fNADC-1 );
+
+  // Allocate ADC data memory if not yet done
+  if( fADCR == NULL ) {
+    fADCR    = new Int_t *[fNADC];
+    fADCF    = new Int_t *[fNADC];
+    fZSM     = new Int_t *[fNADC];
+    fZSM1Dim = new Int_t  [fNADC];
+    for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+      fADCR[iadc] = new Int_t[fNTimeBin];
+      fADCF[iadc] = new Int_t[fNTimeBin];
+      fZSM [iadc] = new Int_t[fNTimeBin];
+    }
+  }
+
+  // Initialize ADC data
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+    for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+      fADCR[iadc][it] = 0;
+      fADCF[iadc][it] = 0;
+      fZSM [iadc][it] = 1;   // Default unread = 1
+    }
+    fZSM1Dim[iadc] = 1;      // Default unread = 1
+  }
+
+  fInitialized = kTRUE;
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::SetData( Int_t iadc, Int_t *adc )
+{
+  // Store ADC data into array of raw data
+
+  if( ! fInitialized ) {
+    //Log (Form ("Error: AliTRDmcmSim is not initialized but setData is called."));
+    return;
+  }
+
+  if( iadc < 0 || iadc >= fNADC ) {
+    //Log (Form ("Error: iadc is out of range (should be 0 to %d).", fNADC-1));
+    return;
+  }
+
+  for( int it = 0 ;  it < fNTimeBin ; it++ ) {
+    fADCR[iadc][it] = (Int_t)(adc[it]);
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::SetData( Int_t iadc, Int_t it, Int_t adc )
+{
+  // Store ADC data into array of raw data
+
+  if( ! fInitialized ) {
+    //Log (Form ("Error: AliTRDmcmSim is not initialized but setData is called."));
+    return;
+  }
+
+  if( iadc < 0 || iadc >= fNADC ) {
+    //Log (Form ("Error: iadc is out of range (should be 0 to %d).", fNADC-1));
+    return;
+  }
+
+  fADCR[iadc][it] = adc;
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::SetDataPedestal( Int_t iadc )
+{
+  // Store ADC data into array of raw data
+
+  if( ! fInitialized ) {
+    //Log (Form ("Error: AliTRDmcmSim is not initialized but setData is called."));
+    return;
+  }
+
+  if( iadc < 0 || iadc >= fNADC ) {
+    //Log (Form ("Error: iadc is out of range (should be 0 to %d).", fNADC-1));
+    return;
+  }
+
+  for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+    fADCR[iadc][it] = fFeeParam->GetADCpedestal();
+  }
+}
+
+//_____________________________________________________________________________
+Int_t AliTRDmcmSim::GetCol( Int_t iadc )
+{
+  // Return column id of the pad for the given ADC channel
+
+  return (fColOfADCbeg - iadc);
+}
+
+
+
+
+//_____________________________________________________________________________
+Int_t AliTRDmcmSim::ProduceRawStream( UInt_t *buf, Int_t maxSize )
+{
+  // Produce raw data stream from this MCM and put in buf
+  // Returns number of words filled, or negative value with -1 * number of overflowed words
+
+  UInt_t  x;
+  UInt_t  iEv = 0;
+  Int_t   nw  = 0;  // Number of written words
+  Int_t   of  = 0;  // Number of overflowed words
+  Int_t   rawVer   = fFeeParam->GetRAWversion();
+  Int_t **adc;
+
+  if( fFeeParam->GetRAWstoreRaw() ) {
+    adc = fADCR;
+  } else {
+    adc = fADCF;
+  }
+
+  // Produce MCM header
+  x = ((fRobPos * fFeeParam->GetNmcmRob() + fMcmPos) << 24) | ((iEv % 0x100000) << 4) | 0xC;
+  if (nw < maxSize) {
+    buf[nw++] = x;
+  }
+  else {
+    of++;
+  }
+
+  // Produce ADC mask
+  if( rawVer >= 3 ) {
+    x = 0;
+    for( Int_t iAdc = 0 ; iAdc < fNADC ; iAdc++ ) {
+      if( fZSM1Dim[iAdc] == 0 ) { //  0 means not suppressed
+       x = x | (1 << iAdc);
+      }
+    }
+    if (nw < maxSize) {
+      buf[nw++] = x;
+    }
+    else {
+      of++;
+    }
+  }
+
+  // Produce ADC data. 3 timebins are packed into one 32 bits word
+  // In this version, different ADC channel will NOT share the same word
+
+  UInt_t aa=0, a1=0, a2=0, a3=0;
+
+  for (Int_t iAdc = 0; iAdc < 21; iAdc++ ) {
+    if( rawVer>= 3 && fZSM1Dim[iAdc] == 1 ) continue; // suppressed
+    aa = !(iAdc & 1) + 2;
+    for (Int_t iT = 0; iT < fNTimeBin; iT+=3 ) {
+      a1 = ((iT    ) < fNTimeBin ) ? adc[iAdc][iT  ] : 0;
+      a2 = ((iT + 1) < fNTimeBin ) ? adc[iAdc][iT+1] : 0;
+      a3 = ((iT + 2) < fNTimeBin ) ? adc[iAdc][iT+2] : 0;
+    }
+    x = (a3 << 22) | (a2 << 12) | (a1 << 2) | aa;
+    if (nw < maxSize) {
+      buf[nw++] = x;
+    }
+    else {
+      of++;
+    }
+  }
+
+  if( of != 0 ) return -of; else return nw;
+}
+
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::Filter()
+{
+  // Apply digital filter
+
+  if( ! fInitialized ) {
+    // Log (Form ("Error: AliTRDmcmSim is not initialized but setData is called."));
+    return;
+  }
+
+  // Initialize filtered data array with raw data
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+    for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+      fADCF[iadc][it] = fADCR[iadc][it]; 
+    }
+  }
+
+  // Then apply fileters one by one to filtered data array
+  if( fFeeParam->isPFon() ) FilterPedestal();
+  if( fFeeParam->isGFon() ) FilterGain();
+  if( fFeeParam->isTFon() ) FilterTail();
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterPedestal()
+{
+  // Apply pedestal filter
+
+  Int_t ap = fFeeParam->GetADCpedestal();      // ADC instrinsic pedestal
+  Int_t ep = fFeeParam->GetPFeffectPedestal(); // effective pedestal
+  Int_t tc = fFeeParam->GetPFtimeConstant();   // this makes no sense yet
+
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+    for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+      fADCF[iadc][it] = fADCF[iadc][it] - ap + ep;
+    }
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterGain()
+{
+  // Apply gain filter (not implemented)
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterTail()
+{
+  // Apply exponential tail filter (Bogdan's version)
+
+  Double_t *dtarg  = new Double_t[fNTimeBin];
+  Int_t    *itarg  = new Int_t[fNTimeBin];
+  Int_t     nexp   = fFeeParam->GetTFnExp();
+  Int_t     tftype = fFeeParam->GetTFtype();
+
+  switch( tftype ) {
+    
+  case 0: // Exponential Filter Analog Bogdan
+    for (Int_t iCol = 0; iCol < fNADC; iCol++) {
+      FilterSimDeConvExpA( fADCF[iCol], dtarg, fNTimeBin, nexp);
+      for (Int_t iTime = 0; iTime < fNTimeBin; iTime++) {
+       fADCF[iCol][iTime] = (Int_t) TMath::Max(0.0,dtarg[iTime]);
+      }
+    }
+    break;
+
+  case 1: // Exponential filter digital Bogdan
+    for (Int_t iCol = 0; iCol < fNADC; iCol++) {
+      FilterSimDeConvExpD( fADCF[iCol], itarg, fNTimeBin, nexp);
+      for (Int_t iTime = 0; iTime < fNTimeBin; iTime++) {
+       fADCF[iCol][iTime] = itarg[iTime];
+      }
+    }
+    break;
+    
+  case 2: // Exponential filter Marian special
+    for (Int_t iCol = 0; iCol < fNADC; iCol++) {
+      FilterSimDeConvExpMI( fADCF[iCol], dtarg, fNTimeBin);
+      for (Int_t iTime = 0; iTime < fNTimeBin; iTime++) {
+       fADCF[iCol][iTime] = (Int_t) TMath::Max(0.0,dtarg[iTime]);
+      }
+    }
+    break;
+    
+  default:
+    AliError(Form("Invalid filter type %d ! \n", tftype ));
+    break;
+  }
+
+  delete dtarg;
+  delete itarg;
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::ZSMapping()
+{
+  // Zero Suppression Mapping implemented in TRAP chip
+  //
+  // See detail TRAP manual "Data Indication" section:
+  // http://www.kip.uni-heidelberg.de/ti/TRD/doc/trap/TRAP-UserManual.pdf
+
+  Int_t EBIS = fFeeParam->GetEBsglIndThr();       // TRAP default = 0x4  (Tis=4)
+  Int_t EBIT = fFeeParam->GetEBsumIndThr();       // TRAP default = 0x28 (Tit=40)
+  Int_t EBIL = fFeeParam->GetEBindLUT();          // TRAP default = 0xf0 (lookup table accept (I2,I1,I0)=(111) or (110) or (101) or (100))
+  Int_t EBIN = fFeeParam->GetEBignoreNeighbour(); // TRAP default = 1 (no neighbor sensitivity)
+
+  for( Int_t iadc = 1 ; iadc < fNADC-1; iadc++ ) {
+    for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+
+      // evaluate three conditions
+      Int_t I0 = ( fADCF[iadc][it] >=  fADCF[iadc-1][it] && fADCF[iadc][it] >=  fADCF[iadc+1][it] ) ? 0 : 1; // peak center detection
+      Int_t I1 = ( (fADCF[iadc-1][it] + fADCF[iadc][it] + fADCF[iadc+1][it]) > EBIT )               ? 0 : 1; // cluster
+      Int_t I2 = ( fADCF[iadc][it] > EBIS )                                                         ? 0 : 1; // absolute large peak
+
+      Int_t I = I2 * 4 + I1 * 2 + I0; // Bit position in lookup table
+      Int_t D = (EBIL >> I) & 1;      // Looking up  (here D=0 means true and D=1 means false according to TRAP manual)
+
+      fZSM[iadc][it] &= D;
+      if( EBIN == 0 ) {  // turn on neighboring ADCs
+       fZSM[iadc-1][it] &= D;
+       fZSM[iadc+1][it] &= D;
+      }
+    }
+  }
+
+  // do 1 dim projection
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) {
+    for( Int_t it = 0 ; it < fNTimeBin ; it++ ) {
+      fZSM1Dim[iadc] &= fZSM[iadc+1][it];
+    }
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterSimDeConvExpA(Int_t *source, Double_t *target, Int_t n, Int_t nexp) 
+{
+  //
+  // Exponential filter "analog"
+  // source will not be changed
+  //
+
+  Int_t    i = 0;
+  Int_t    k = 0;
+  Double_t reminder[2];
+  Double_t correction;
+  Double_t result;
+  Double_t rates[2];
+  Double_t coefficients[2];
+
+  // Initialize (coefficient = alpha, rates = lambda)
+  // FilterOpt.C (aliroot@pel:/homel/aliroot/root/work/beamt/CERN02)
+
+  Double_t r1 = (Double_t)fFeeParam->GetTFr1();
+  Double_t r2 = (Double_t)fFeeParam->GetTFr2();
+  Double_t c1 = (Double_t)fFeeParam->GetTFc1();
+  Double_t c2 = (Double_t)fFeeParam->GetTFc2();
+  
+  coefficients[0] = c1;
+  coefficients[1] = c2;
+
+  Double_t dt = 0.1;
+  rates[0] = TMath::Exp(-dt/(r1));
+  rates[1] = TMath::Exp(-dt/(r2));
+
+  // Attention: computation order is important
+  correction = 0.0;
+  for (k = 0; k < nexp; k++) {
+    reminder[k] = 0.0;
+  }
+    
+  for (i = 0; i < n; i++) {
+
+    result    = ((Double_t)source[i] - correction);    // no rescaling
+    target[i] = result;
+    
+    for (k = 0; k < nexp; k++) {
+      reminder[k] = rates[k] * (reminder[k] + coefficients[k] * result);
+    }
+      
+    correction = 0.0;
+    for (k = 0; k < nexp; k++) {
+      correction += reminder[k];
+    }
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterSimDeConvExpD(Int_t *source, Int_t *target, Int_t n, Int_t nexp) 
+{
+  //
+  // Exponential filter "digital"
+  // source will not be changed
+  //
+
+  Int_t i = 0;
+
+  Int_t fAlphaL  = 0;
+  Int_t fAlphaS  = 0;
+  Int_t fLambdaL = 0;
+  Int_t fLambdaS = 0;
+  Int_t fTailPed = 0;
+
+  Int_t iAlphaL  = 0;
+  Int_t iAlphaS  = 0;
+  Int_t iLambdaL = 0;
+  Int_t iLambdaS = 0;
+
+  // FilterOpt.C (aliroot@pel:/homel/aliroot/root/work/beamt/CERN02)
+  // initialize (coefficient = alpha, rates = lambda)
+
+  Double_t dt = 0.1;
+  Double_t r1 = (Double_t)fFeeParam->GetTFr1();
+  Double_t r2 = (Double_t)fFeeParam->GetTFr2();
+  Double_t c1 = (Double_t)fFeeParam->GetTFc1();
+  Double_t c2 = (Double_t)fFeeParam->GetTFc2();
+
+  fLambdaL = (Int_t)((TMath::Exp(-dt/r1) - 0.75) * 2048.0);
+  fLambdaS = (Int_t)((TMath::Exp(-dt/r2) - 0.25) * 2048.0);
+  iLambdaL = fLambdaL & 0x01FF; iLambdaL |= 0x0600;    //  9 bit paramter + fixed bits
+  iLambdaS = fLambdaS & 0x01FF; iLambdaS |= 0x0200;    //  9 bit paramter + fixed bits
+
+  if (nexp == 1) {
+    fAlphaL = (Int_t) (c1 * 2048.0);
+    iAlphaL = fAlphaL & 0x03FF;                                // 10 bit paramter
+  }
+  if (nexp == 2) {
+    fAlphaL = (Int_t) (c1 * 2048.0);
+    fAlphaS = (Int_t) ((c2 - 0.5) * 2048.0);
+    iAlphaL = fAlphaL & 0x03FF;                                // 10 bit paramter
+    iAlphaS = fAlphaS & 0x03FF; iAlphaS |= 0x0400;             // 10 bit paramter + fixed bits
+  }
+  
+  Double_t iAl = iAlphaL  / 2048.0;           // alpha  L: correspondence to floating point numbers
+  Double_t iAs = iAlphaS  / 2048.0;           // alpha  S: correspondence to floating point numbers
+  Double_t iLl = iLambdaL / 2048.0;           // lambda L: correspondence to floating point numbers
+  Double_t iLs = iLambdaS / 2048.0;           // lambda S: correspondence to floating point numbers
+
+  Int_t h1;
+  Int_t h2;
+  Int_t rem1;
+  Int_t rem2;
+  Int_t correction;
+  Int_t result;
+  Int_t iFactor = ((Int_t) fFeeParam->GetPFeffectPedestal() ) << 2;
+
+  Double_t xi = 1 - (iLl*iAs + iLs*iAl);            // Calculation of equilibrium values of the
+  rem1 = (Int_t) ((iFactor/xi) * ((1-iLs)*iLl*iAl)); // Internal registers to prevent switch on effects.
+  rem2 = (Int_t) ((iFactor/xi) * ((1-iLl)*iLs*iAs));
+  
+  // further initialization
+  if ((rem1 + rem2) > 0x0FFF) {
+    correction = 0x0FFF;
+  } 
+  else {
+    correction = (rem1 + rem2) & 0x0FFF;
+  }
+
+  fTailPed = iFactor - correction;
+
+  for (i = 0; i < n; i++) {
+
+    result = (source[i] - correction);
+    if (result < 0) {                          
+      result = 0;
+    }
+
+    target[i] = result;
+                                                        
+    h1 = (rem1 + ((iAlphaL * result) >> 11));
+    if (h1 > 0x0FFF) {
+      h1 = 0x0FFF;
+    } 
+    else {
+      h1 &= 0x0FFF;
+    }
+
+    h2 = (rem2 + ((iAlphaS * result) >> 11));
+    if (h2 > 0x0FFF) {
+      h2 = 0x0FFF;
+    } 
+    else {
+      h2 &= 0x0FFF;
+    }
+  
+    rem1 = (iLambdaL * h1 ) >> 11;
+    rem2 = (iLambdaS * h2 ) >> 11;
+    
+    if ((rem1 + rem2) > 0x0FFF) {
+      correction = 0x0FFF;
+    } 
+    else {
+      correction = (rem1 + rem2) & 0x0FFF;
+    }
+
+  }
+}
+
+//_____________________________________________________________________________
+void AliTRDmcmSim::FilterSimDeConvExpMI(Int_t *source, Double_t *target, Int_t n) 
+{
+  //
+  // Exponential filter (M. Ivanov)
+  // source will not be changed
+  //
+
+  Int_t i = 0;
+  Double_t sig1[100];
+  Double_t sig2[100];
+  Double_t sig3[100];
+
+  for (i = 0; i < n; i++) {
+    sig1[i] = (Double_t)source[i];
+  }
+
+  Float_t dt      = 0.1;
+  Float_t lambda0 = (1.0 / fFeeParam->GetTFr2()) * dt;
+  Float_t lambda1 = (1.0 / fFeeParam->GetTFr1()) * dt;
+
+  FilterSimTailMakerSpline( sig1, sig2, lambda0, n);
+  FilterSimTailCancelationMI( sig2, sig3, 0.7, lambda1, n);
+
+  for (i = 0; i < n; i++) {
+    target[i] = sig3[i];
+  }
+
+}
+
+//______________________________________________________________________________
+void AliTRDmcmSim::FilterSimTailMakerSpline(Double_t *ampin, Double_t *ampout, Double_t lambda, Int_t n) 
+{
+  //
+  // Special filter (M. Ivanov)
+  //
+
+  Int_t    i = 0;
+  Double_t l = TMath::Exp(-lambda*0.5);
+  Double_t in[1000];
+  Double_t out[1000];
+
+  // Initialize in[] and out[] goes 0 ... 2*n+19
+  for (i = 0; i < n*2+20; i++) {
+    in[i] = out[i] = 0;
+  }
+
+  // in[] goes 0, 1
+  in[0] = ampin[0];
+  in[1] = (ampin[0] + ampin[1]) * 0.5;
+   
+  // Add charge to the end
+  for (i = 0; i < 22; i++) {
+    in[2*(n-1)+i] = ampin[n-1]; // in[] goes 2*n-2, 2*n-1, ... , 2*n+19 
+  }
+
+  // Use arithmetic mean
+  for (i = 1; i < n-1; i++) {
+    in[2*i]   = ampin[i];    // in[] goes 2, 3, ... , 2*n-4, 2*n-3
+    in[2*i+1] = ((ampin[i]+ampin[i+1]))/2.;
+  }
+
+  Double_t temp;
+  out[2*n]    = in[2*n];
+  temp        = 0;
+  for (i = 2*n; i >= 0; i--) {
+    out[i]    = in[i] + temp;
+    temp      = l*(temp+in[i]);
+  }
+
+  for (i = 0; i < n; i++){
+    //ampout[i] = out[2*i+1];  // org
+    ampout[i] = out[2*i];
+  }
+
+}
+
+//______________________________________________________________________________
+void AliTRDmcmSim::FilterSimTailCancelationMI(Double_t *ampin, Double_t *ampout, Double_t norm, Double_t lambda, Int_t n) 
+{
+  //
+  // Special filter (M. Ivanov)
+  //
+
+  Int_t    i = 0;
+
+  Double_t l = TMath::Exp(-lambda*0.5);
+  Double_t k = l*(1.0 - norm*lambda*0.5);
+  Double_t in[1000];
+  Double_t out[1000];
+
+  // Initialize in[] and out[] goes 0 ... 2*n+19
+  for (i = 0; i < n*2+20; i++) {
+    in[i] = out[i] = 0;
+  }
+
+  // in[] goes 0, 1
+  in[0] = ampin[0];
+  in[1] = (ampin[0]+ampin[1])*0.5;
+
+  // Add charge to the end
+  for (i =-2; i < 22; i++) {
+    // in[] goes 2*n-4, 2*n-3, ... , 2*n+19 
+    in[2*(n-1)+i] = ampin[n-1];
+  }
+
+  for (i = 1; i < n-2; i++) {
+    // in[] goes 2, 3, ... , 2*n-6, 2*n-5
+    in[2*i]    = ampin[i];
+    in[2*i+1]  = (9.0 * (ampin[i]+ampin[i+1]) - (ampin[i-1]+ampin[i+2])) / 16.0;
+    //in[2*i+1]  = ((ampin[i]+ampin[i+1]))/2.0;
+  }
+
+  Double_t temp;
+  out[0] = in[0];
+  temp   = in[0];
+  for (i = 1; i <= 2*n; i++) {
+    out[i] = in[i] + (k-l)*temp;
+    temp   = in[i] +  k   *temp;
+  }
+
+  for (i = 0; i < n; i++) {
+    //ampout[i] = out[2*i+1];  // org
+    //ampout[i] = TMath::Max(out[2*i+1],0.0);  // org
+    ampout[i] = TMath::Max(out[2*i],0.0);
+  }
+}
+
+// EOF
diff --git a/TRD/AliTRDmcmSim.h b/TRD/AliTRDmcmSim.h
new file mode 100644 (file)
index 0000000..0f4690d
--- /dev/null
@@ -0,0 +1,81 @@
+#ifndef ALITRDMCMSIM_H
+#define ALITRDMCMSIM_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id$ */
+
+///////////////////////////////////////////////////////
+//                                                   //
+//  Multi Chip Module Simulation Class               //
+//                                                   //
+///////////////////////////////////////////////////////
+
+#include <TObject.h>
+
+class AliTRDfeeParam;
+class AliTRDgeometry;
+
+class AliTRDmcmSim : public TObject {
+
+ public:
+
+  AliTRDmcmSim();
+  virtual         ~AliTRDmcmSim();
+
+  //PH  virtual void      Copy(TObject &m) const;
+
+          void      Init( Int_t cha_id, Int_t rob_pos, Int_t mcm_pos ); // Initialize MCM by the position parameters
+          void      SetData(Int_t iadc, Int_t *adc);           // Set ADC data with array 
+          void      SetData(Int_t iadc, Int_t it, Int_t adc ); // Set ADC data one by one
+          void      SetDataPedestal(Int_t iadc );              // Fill ADC data with pedestal values
+
+          Int_t     GetChaId()  { return fChaId;  };    // Returns Chamber ID (0-539)
+          Int_t     GetRobPos() { return fRobPos; };    // Returns ROB position (0-7)
+          Int_t     GetMcmPos() { return fMcmPos; };    // Returns MCM position (0-17) (16,17 are mergers)
+          Int_t     GetRow()    { return fRow; };       // Returns Row number on chamber where it is sitting
+         Int_t     GetCol( Int_t iadc );               // Get corresponding column (0-143) from adc ID (= 0-20) 
+                    // for the ADC/Col mapping, see: http://wiki.kip.uni-heidelberg.de/ti/TRD/index.php/Image:ROB_MCM_numbering.pdf
+
+         Int_t     ProduceRawStream( UInt_t *buf, Int_t bufsize );    // Produce raw data stream from this MCM
+         void      Filter();                           //  Apply digital filters
+         void      ZSMapping();                        //  Do ZS mapping
+
+ protected:
+
+          Bool_t    fInitialized;                       // Status
+
+         AliTRDfeeParam *fFeeParam;
+         AliTRDgeometry *fGeo;
+
+         Int_t     fChaId;                             // Chamber ID
+         Int_t     fSector;                            // Sector number of the supermodule
+         Int_t     fStack;                             // Chamber stack ID
+         Int_t     fLayer;                             // Chamber layer ID
+          Int_t     fRobPos;                            // ROB Position
+          Int_t     fMcmPos;                            // MCM Position
+         Int_t     fNADC;                              // Number of ADC (usually 21)
+         Int_t     fNTimeBin;                          // Number of Timebin (valiable)
+          Int_t     fRow;                               // Pad row number (0-11 or 0-15)
+          Int_t     fColOfADCbeg;                       // Column corresponds to ADC channel 0
+          Int_t     fColOfADCend;                       // Column corresponds to ADC channel [fNADC-1] (usually 20)
+          Int_t   **fADCR;                              // Array with MCM ADC values (Raw)
+          Int_t   **fADCF;                              // Array with MCM ADC values (Filtered)
+          Int_t   **fZSM;                               // Zero Suppression Map
+          Int_t    *fZSM1Dim;                           // Zero Suppression (1 dimensional projection)
+
+         void      FilterPedestal();                   // Apply pedestal filter
+         void      FilterGain();                       // Apply gain filter
+         void      FilterTail();                       // Apply tail filter
+
+         void      FilterSimDeConvExpA(Int_t *source, Double_t *target, Int_t n, Int_t nexp);
+         void      FilterSimDeConvExpD(Int_t *source, Int_t *target, Int_t n, Int_t nexp);
+         void      FilterSimDeConvExpMI(Int_t *source, Double_t *target, Int_t n);
+         void      FilterSimTailMakerSpline(Double_t *ampin, Double_t *ampout, Double_t lambda, Int_t n);
+         void      FilterSimTailCancelationMI(Double_t *ampin, Double_t *ampout, Double_t norm, Double_t lambda, Int_t n);
+
+  ClassDef(AliTRDmcmSim,3)
+
+};
+
+#endif
index 53f6a3b..aed37f9 100644 (file)
 #include "AliFstream.h"
 
 #include "AliTRDSignalIndex.h"
+
+#include "AliTRDfeeParam.h"
+#include "AliTRDmcmSim.h"
+
 ClassImp(AliTRDrawData)
 
 //_____________________________________________________________________________
 AliTRDrawData::AliTRDrawData()
   :TObject()
-  ,fRawVersion(2)    // Default Raw Data version set here
-  ,fGeo(0)
+  ,fGeo(NULL)
+  ,fFee(NULL)
   ,fNumberOfDDLs(0)
 {
   //
   // Default constructor
   //
-
+  fFee = AliTRDfeeParam::Instance();
+  fNumberOfDDLs = AliDAQ::NumberOfDdls("TRD");
 }
 
 //_____________________________________________________________________________
 AliTRDrawData::AliTRDrawData(const AliTRDrawData &r)
   :TObject(r)
-  ,fRawVersion(2)    // Default Raw Data version set here
-  ,fGeo(0)
+  ,fGeo(NULL)
+  ,fFee(NULL)
   ,fNumberOfDDLs(0)
 {
   //
   // Copy constructor
   //
-
+  fFee = AliTRDfeeParam::Instance();
+  fNumberOfDDLs = AliDAQ::NumberOfDdls("TRD");
 }
 
 //_____________________________________________________________________________
@@ -73,23 +79,6 @@ AliTRDrawData::~AliTRDrawData()
   //
   // Destructor
   //
-
-}
-
-//_____________________________________________________________________________
-Bool_t AliTRDrawData::SetRawVersion(Int_t v)
-{
-  //
-  // Set the raw data version (Currently only version 0, 1 and 2 are available)
-  //
-
-  if ( (v >= 0) && (v <= 2) ) {
-    fRawVersion = v;
-    return kTRUE;
-  }
-
-  return kFALSE;
-
 }
 
 //_____________________________________________________________________________
@@ -103,8 +92,6 @@ Bool_t AliTRDrawData::Digits2Raw(TTree *digitsTree, TTree *tracks )
   // will be supported in higher version simulator.
   //
 
-  fNumberOfDDLs = AliDAQ::NumberOfDdls("TRD");
-
   AliTRDdigitsManager* digitsManager = new AliTRDdigitsManager();
 
   if (!digitsManager->ReadDigits(digitsTree)) {
@@ -128,12 +115,13 @@ Bool_t AliTRDrawData::Digits2Raw(TTree *digitsTree, TTree *tracks )
   }
 
   Int_t retval = kTRUE;
+  Int_t rv     = fFee->GetRAWversion();
 
   // Call appropriate Raw Simulator
-  if ( fRawVersion > 0 && fRawVersion <= 2 ) retval = Digits2Raw(digitsManager); 
+  if ( rv > 0 && rv <= 3 ) retval = Digits2Raw(digitsManager); 
   else {
     retval = kFALSE;
-    AliWarning(Form("Unsupported raw version (fRawVersion=%d).",fRawVersion));
+    AliWarning(Form("Unsupported raw version (%d).", rv));
   }
 
   // Cleanup
@@ -198,17 +186,21 @@ Bool_t AliTRDrawData::Digits2Raw(AliTRDdigitsManager *digitsManager)
         digits->Expand();
 
         Int_t hcwords = 0;
+       Int_t rv = fFee->GetRAWversion();
 
         // Process A side of the chamber
-       if ( fRawVersion >= 1 && fRawVersion <= 2 ) hcwords = ProduceHcDataV1andV2(digits,0,iDet,hc_buffer,kMaxHcWords);
+       if ( rv >= 1 && rv <= 2 ) hcwords = ProduceHcDataV1andV2(digits,0,iDet,hc_buffer,kMaxHcWords);
+       if ( rv == 3 )            hcwords = ProduceHcDataV3     (digits,0,iDet,hc_buffer,kMaxHcWords);
+
         of->WriteBuffer((char *) hc_buffer, hcwords*4);
         npayloadbyte += hcwords*4;
 
         // Process B side of the chamber
-       if ( fRawVersion >= 1 && fRawVersion <= 2 ) hcwords = ProduceHcDataV1andV2(digits,1,iDet,hc_buffer,kMaxHcWords);
+       if ( rv >= 1 && rv <= 2 ) hcwords = ProduceHcDataV1andV2(digits,1,iDet,hc_buffer,kMaxHcWords);
+       if ( rv >= 3 )            hcwords = ProduceHcDataV3     (digits,1,iDet,hc_buffer,kMaxHcWords);
+
         of->WriteBuffer((char *) hc_buffer, hcwords*4);
         npayloadbyte += hcwords*4;
-
       }
     }
 
@@ -230,8 +222,8 @@ Int_t AliTRDrawData::ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side
                                         , Int_t det, UInt_t *buf, Int_t maxSize)
 {
   //
-  // This function simulates: 1) SM-I commissiong data Oct. 06 (fRawVersion == 1).
-  //                          2) Full Raw Production Version   (fRawVersion == 2)
+  // This function simulates: 1) SM-I commissiong data Oct. 06 (Raw Version == 1).
+  //                          2) Full Raw Production Version   (Raw Version == 2)
   //
   // Produce half chamber data (= an ORI data) for the given chamber (det) and side (side)
   // where
@@ -259,6 +251,7 @@ Int_t AliTRDrawData::ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side
   Int_t      kCtype = 0;                       // Chamber type (0:C0, 1:C1)
   Int_t         iEv = 0xA;                     // Event ID. Now fixed to 10, how do I get event id?
   UInt_t          x = 0;                       // General used number
+  Int_t          rv = fFee->GetRAWversion();
 
   // Check the nCol and nRow.
   if ((nCol == 144) && 
@@ -285,7 +278,7 @@ Int_t AliTRDrawData::ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side
   }
 
   // Half Chamber header
-  if      ( fRawVersion == 1 ) {
+  if      ( rv == 1 ) {
     // Now it is the same version as used in SM-I commissioning.
     Int_t  dcs = det+100;      // DCS Serial (in simulation, it is meaningless
     x = (dcs<<20) | (sect<<15) | (plan<<12) | (cham<<9) | (side<<8) | 1;
@@ -296,11 +289,11 @@ Int_t AliTRDrawData::ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side
       of++;
     }
   } 
-  else if ( fRawVersion == 2 ) {
+  else if ( rv == 2 ) {
     // h[0] (there are 3 HC header)
     Int_t minorv = 0;      // The minor version number
     Int_t add    = 1;      // The number of additional header words to follow
-    x = (1<<31) | (fRawVersion<<24) | (minorv<<17) | (add<<14) | (sect<<9) | (plan<<6) | (cham<<3) | (side<<2) | 1;
+    x = (1<<31) | (rv<<24) | (minorv<<17) | (add<<14) | (sect<<9) | (plan<<6) | (cham<<3) | (side<<2) | 1;
     if (nw < maxSize) {
       buf[nw++] = x; 
     }
@@ -409,6 +402,143 @@ Int_t AliTRDrawData::ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side
 
 }
 
+
+//_____________________________________________________________________________
+Int_t AliTRDrawData::ProduceHcDataV3(AliTRDdataArrayI *digits, Int_t side
+                                    , Int_t det, UInt_t *buf, Int_t maxSize)
+{
+  //
+  // This function simulates: Raw Version == 3 (Zero Suppression Prototype)
+  //
+
+  Int_t          nw = 0;                       // Number of written    words
+  Int_t          of = 0;                       // Number of overflowed words
+  Int_t        plan = fGeo->GetPlane( det );   // Plane
+  Int_t        cham = fGeo->GetChamber( det ); // Chamber
+  Int_t        sect = fGeo->GetSector( det );  // Sector (=iDDL)
+  Int_t        nRow = fGeo->GetRowMax( plan, cham, sect );
+  Int_t        nCol = fGeo->GetColMax( plan );
+  const Int_t nTBin = AliTRDcalibDB::Instance()->GetNumberOfTimeBins();
+  Int_t      kCtype = 0;                       // Chamber type (0:C0, 1:C1)
+  Int_t         iEv = 0xA;                     // Event ID. Now fixed to 10, how do I get event id?
+  UInt_t          x = 0;                       // General used number
+  Int_t          rv = fFee->GetRAWversion();
+
+  // Check the nCol and nRow.
+  if ((nCol == 144) && 
+      (nRow == 16 || nRow == 12)) {
+    kCtype = (nRow-12) / 4;
+  } 
+  else {
+    AliError(Form("This type of chamber is not supported (nRow=%d, nCol=%d)."
+                 ,nRow,nCol));
+    return 0;
+  }
+
+  AliDebug(1,Form("Producing raw data for sect=%d plan=%d cham=%d side=%d"
+                 ,sect,plan,cham,side));
+
+  // Tracklet should be processed here but not implemented yet
+
+  // Write end of tracklet marker
+  if (nw < maxSize) {
+    buf[nw++] = kEndoftrackletmarker;
+  } 
+  else {
+    of++;
+  }
+
+  // Half Chamber header
+  // h[0] (there are 3 HC header)
+  Int_t minorv = 0;    // The minor version number
+  Int_t add    = 1;    // The number of additional header words to follow
+  x = (1<<31) | (rv<<24) | (minorv<<17) | (add<<14) | (sect<<9) | (plan<<6) | (cham<<3) | (side<<2) | 1;
+  if (nw < maxSize) {
+    buf[nw++] = x; 
+  }
+  else {
+    of++;
+  }
+  // h[1]
+  Int_t bc_ctr   = 99; // bunch crossing counter. Here it is set to 99 always for no reason
+  Int_t pt_ctr   = 15; // pretrigger counter. Here it is set to 15 always for no reason
+  Int_t pt_phase = 11; // pretrigger phase. Here it is set to 11 always for no reason
+  x = (bc_ctr<<16) | (pt_ctr<<12) | (pt_phase<<8) | ((nTBin-1)<<2) | 1;
+  if (nw < maxSize) {
+    buf[nw++] = x; 
+  }
+  else {
+    of++;
+  }
+  // h[2]
+  Int_t ped_setup       = 1;  // Pedestal filter setup (0:1). Here it is always 1 for no reason
+  Int_t gain_setup      = 1;  // Gain filter setup (0:1). Here it is always 1 for no reason
+  Int_t tail_setup      = 1;  // Tail filter setup (0:1). Here it is always 1 for no reason
+  Int_t xt_setup        = 0;  // Cross talk filter setup (0:1). Here it is always 0 for no reason
+  Int_t nonlin_setup    = 0;  // Nonlinearity filter setup (0:1). Here it is always 0 for no reason
+  Int_t bypass_setup    = 0;  // Filter bypass (for raw data) setup (0:1). Here it is always 0 for no reason
+  Int_t common_additive = 10; // Digital filter common additive (0:63). Here it is always 10 for no reason
+  x = (ped_setup<<31) | (gain_setup<<30) | (tail_setup<<29) | (xt_setup<<28) | (nonlin_setup<<27)
+    | (bypass_setup<<26) | (common_additive<<20) | 1;
+  if (nw < maxSize) {
+    buf[nw++] = x; 
+  }
+  else {
+    of++;
+  }
+
+  AliTRDmcmSim *mcm = new AliTRDmcmSim();
+
+  // Scan for ROB and MCM
+  for (Int_t iRobRow = 0; iRobRow < (kCtype + 3); iRobRow++ ) {
+    Int_t iRob = iRobRow * 2 + side;
+    for (Int_t iMcm = 0; iMcm < fGeo->MCMmax(); iMcm++ ) {
+
+      mcm->Init( det, iRob, iMcm );
+      Int_t padrow = mcm->GetRow();
+
+      // Copy ADC data to MCM simulator
+      for (Int_t iAdc = 0; iAdc < 21; iAdc++ ) {
+       Int_t padcol = mcm->GetCol( iAdc );
+       if ((padcol >=    0) && (padcol <  nCol)) {
+         for (Int_t iT = 0; iT < nTBin; iT++) { 
+           mcm->SetData( iAdc, iT, digits->GetDataUnchecked( padrow, padcol, iT) );
+         } 
+       } else {  // this means it is out of chamber, and masked ADC
+         mcm->SetDataPedestal( iAdc );
+       }
+      }
+      // Simulate process in MCM
+      mcm->Filter();     // Apply filter
+      mcm->ZSMapping();  // Calculate zero suppression mapping
+
+      // Write MCM data to buffer
+      Int_t temp_nw =  mcm->ProduceRawStream( &buf[nw], maxSize - nw );
+      if( temp_nw < 0 ) {
+       of += temp_nw;
+       nw += maxSize - nw;
+       AliError(Form("Buffer overflow detected. Please increase the buffer size and recompile."));
+      } else {
+       nw += temp_nw;
+      }
+    }
+  }
+
+  // Write end of raw data marker
+  if (nw < maxSize) {
+    buf[nw++] = kEndofrawdatamarker; 
+  }
+  else {
+    of++;
+  }
+  if (of != 0) {
+    AliWarning("Buffer overflow. Data is truncated. Please increase buffer size and recompile.");
+  }
+
+  return nw;
+
+}
+
 //_____________________________________________________________________________
 AliTRDdigitsManager *AliTRDrawData::Raw2Digits(AliRawReader *rawReader)
 {
@@ -427,7 +557,7 @@ AliTRDdigitsManager *AliTRDrawData::Raw2Digits(AliRawReader *rawReader)
   digitsManager->CreateArrays();
 
   AliTRDRawStream input(rawReader);
-  input.SetRawVersion( fRawVersion );
+  input.SetRawVersion( fFee->GetRAWversion() );
   input.Init();
 
   // Loop through the digits
index 772b595..6dc74b4 100644 (file)
@@ -19,6 +19,7 @@ class AliRawReader;
 
 class AliTRDdigitsManager;
 class AliTRDgeometry;
+class AliTRDfeeParam;
 class AliTRDdataArrayI;
 
 class AliTRDrawData : public TObject {
@@ -32,21 +33,20 @@ class AliTRDrawData : public TObject {
   AliTRDrawData &operator=(const AliTRDrawData &/*r*/) { return *this; }
 
   virtual Bool_t       Digits2Raw(TTree *digits, TTree *tracks = NULL);
-  virtual Bool_t       SetRawVersion(Int_t v);
 
   virtual AliTRDdigitsManager* Raw2Digits(AliRawReader *rawReader);
 
  protected:
 
   virtual Bool_t       Digits2Raw(AliTRDdigitsManager* digitsManager); // for fRawVersion > 0
-
   virtual Int_t        ProduceHcDataV1andV2(AliTRDdataArrayI *digits, Int_t side, Int_t det, UInt_t *buf, Int_t maxSize);
+  virtual Int_t        ProduceHcDataV3(AliTRDdataArrayI *digits, Int_t side, Int_t det, UInt_t *buf, Int_t maxSize);
   
-  Int_t                fRawVersion;     //  Which version of raw simulator is used
   AliTRDgeometry      *fGeo;            //! Geometry
+  AliTRDfeeParam      *fFee;            //! Fee Parameters
   Int_t                fNumberOfDDLs;   //  Number of DDLs
 
-  ClassDef(AliTRDrawData,4)             //  TRD raw data class
+  ClassDef(AliTRDrawData,5)             //  TRD raw data class
 
 };
 #endif
index 58222ed..9680688 100644 (file)
@@ -66,6 +66,7 @@
 #pragma link C++ class  AliTRDtrigParam+;
 #pragma link C++ class  AliTRDmcmTracklet+;
 #pragma link C++ class  AliTRDmcm+;
+#pragma link C++ class  AliTRDmcmSim+;
 #pragma link C++ class  AliTRDltuTracklet+;
 #pragma link C++ class  AliTRDgtuTrack+;
 #pragma link C++ class  AliTRDmodule+;
index 2c6b5c7..5c4293b 100644 (file)
@@ -40,6 +40,7 @@ SRCS= AliTRDarrayI.cxx \
       AliTRDTriggerL1.cxx \
       AliTRDmcmTracklet.cxx \
       AliTRDmcm.cxx \
+      AliTRDmcmSim.cxx \
       AliTRDtrigParam.cxx \
       AliTRDtrigger.cxx \
       AliTRDltuTracklet.cxx \