]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDmcmSim.cxx
Changed AliRunLoader::GetRunLoader() into AliRunLoader::Instance()
[u/mrichter/AliRoot.git] / TRD / AliTRDmcmSim.cxx
index edcbfd8d92df8f8c1d79cd4d271d486c16c32416..534f82c23315ef0197faae8d7f4f3a9aeb5ab943 100644 (file)
@@ -101,7 +101,8 @@ The default raw version is 2.
 #include "AliTRDSimParam.h"
 #include "AliTRDgeometry.h"
 #include "AliTRDcalibDB.h"
-
+#include "AliTRDdigitsManager.h"
+#include "AliTRDarrayADC.h"
 // additional for new tail filter and/or tracklet
 #include "AliTRDtrapAlu.h"
 #include "AliTRDpadPlane.h"
@@ -608,6 +609,94 @@ Int_t AliTRDmcmSim::ProduceRawStream( UInt_t *buf, Int_t maxSize )
   if( of != 0 ) return -of; else return nw;
 }
 
+//_____________________________________________________________________________
+Int_t AliTRDmcmSim::ProduceRawStreamV2( UInt_t *buf, Int_t maxSize, UInt_t iEv )
+{
+  //
+  // 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;
+  Int_t   nActiveADC = 0;      // number of activated ADC bits in a word
+
+  if( !CheckInitialized() ) return 0;
+
+  if( fFeeParam->GetRAWstoreRaw() ) {
+    adc = fADCR;
+  } else {
+    adc = fADCF;
+  }
+
+  // Produce MCM header : xrrr mmmm eeee eeee eeee eeee eeee 1100
+  //                      x : 0 before , 1 since 10.2007
+  //                      r : Readout board position (Alice numbering)
+  //                      m : MCM posi
+  //                      e : Event counter from 1
+  //x = (1<<31) | ((fRobPos * fFeeParam->GetNmcmRob() + fMcmPos) << 24) | ((iEv % 0x100000) << 4) | 0xC;
+  x = (1<<31) | (fRobPos << 28) | (fMcmPos << 24) | ((iEv % 0x100000) << 4) | 0xC;
+  if (nw < maxSize) {
+    buf[nw++] = x;
+       //printf("\nMCM header: %X ",x);
+  }
+  else {
+    of++;
+  }
+
+  // Produce ADC mask : nncc cccm mmmm mmmm mmmm mmmm mmmm 1100
+  //                           n : unused , c : ADC count, m : selected ADCs
+  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+4) );       // last 4 digit reserved for 1100=0xc
+               nActiveADC++;           // number of 1 in mmm....m
+      }
+    }
+       x = x | (1 << 30) | ( ( 0x3FFFFFFC ) & (~(nActiveADC) << 25) ) | 0xC;   // nn = 01, ccccc are inverted, 0xc=1100
+       //printf("nActiveADC=%d=%08X, inverted=%X ",nActiveADC,nActiveADC,x );
+
+    if (nw < maxSize) {
+      buf[nw++] = x;
+         //printf("ADC mask: %X nMask=%d ADC data: ",x,nActiveADC);
+    }
+    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] != 0 ) continue; // Zero Suppression, 0 means not suppressed
+    aa = !(iAdc & 1) + 2;      // 3 for the even ADC channel , 2 for the odd ADC channel
+    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;
+       //printf("%08X ",x);
+      }
+      else {
+       of++;
+      }
+    }
+  }
+
+  if( of != 0 ) return -of; else return nw;
+}
+
 //_____________________________________________________________________________
 Int_t AliTRDmcmSim::ProduceTrackletStream( UInt_t *buf, Int_t maxSize )
 {
@@ -753,8 +842,9 @@ void AliTRDmcmSim::FilterTail()
     break;
   }
 
-  delete dtarg;
-  delete itarg;
+  delete [] dtarg;
+  delete [] itarg;
+
 }
 
 //_____________________________________________________________________________
@@ -1525,7 +1615,7 @@ void AliTRDmcmSim::Tracklet(){
     Int_t j = 1; // selection number
     while(i<fNADC-2 && j<=3){
       i = i + 1;
-      if((mark>>(i-1)) & 1 == 1) {
+      if( ((mark>>(i-1)) & 1) == 1) {
        selection[j].iadc  = fNADC-1-i;
        selection[j].value = qsum[fNADC-1-i]>>6;   // for hit-selection only the first 8 out of the 14 Bits are used for comparison
        
@@ -1561,7 +1651,7 @@ void AliTRDmcmSim::Tracklet(){
     // read three from left side
     Int_t k = fNADC-2;
     while(k>i && j<=6) {
-      if((mark>>(k-1)) & 1 == 1) {
+      if( ((mark>>(k-1)) & 1) == 1) {
        selection[j].iadc  = fNADC-1-k;
        selection[j].value = qsum[fNADC-1-k]>>6;
        
@@ -1768,7 +1858,7 @@ void AliTRDmcmSim::Tracklet(){
   // insert marked channels into list and sort according to hit-sum
   while(adcL < fNADC-3 && selNr < fNADC-3){
      
-    if((mPair>>((fNADC-4)-(adcL))) & 1 == 1) {
+    if( ((mPair>>((fNADC-4)-(adcL))) & 1) == 1) {
       selectPair[selNr].iadc  = adcL;
       selectPair[selNr].value = hitSum[adcL];   
       
@@ -2355,7 +2445,7 @@ void AliTRDmcmSim::Tracklet(){
        shift2 = 1;
        for(Int_t iBit = 1; iBit < 7; iBit++) {
            shift2  = shift2<<1;
-           shift2 |= (1-((shift)>>(6-iBit))&1);
+           shift2 |= (1- (((shift)>>(6-iBit))&1) );
            //printf("%d",(1-((-mSlope[i])>>(6-iBit))&1));
        }
        shift2 = shift2 + 1;
@@ -2384,7 +2474,7 @@ void AliTRDmcmSim::Tracklet(){
        shift2 = 1;
        for(Int_t iBit = 1; iBit < 13; iBit++) {
            shift2  = shift2<<1;
-           shift2 |= (1-((shift)>>(12-iBit))&1);
+           shift2 |= (1-(((shift)>>(12-iBit))&1));
            //printf("%d",(1-((-mOffset[i])>>(12-iBit))&1));
        }
        shift2 = shift2 + 1;
@@ -2436,25 +2526,26 @@ void AliTRDmcmSim::Tracklet(){
 
 //if you want to activate the MC tracklet output, set fgkMCTrackletOutput=kTRUE in AliTRDfeeParam
        
-  if (!fFeeParam->GetMCTrackletOutput()) return;
+  if (!fFeeParam->GetMCTrackletOutput()) 
+      return;
  
-  
   AliLog::SetClassDebugLevel("AliTRDmcmSim", 10);
   AliLog::SetFileOutput("../log/tracklet.log");
   
-   UInt_t* trackletWord;
-  Int_t*  adcChannel;
-
-  Int_t u = 0;
-
   // testing for wordnr in order to speed up the simulation
   if (wordnr == 0) 
     return;
    
-  //Int_t mcmNr = fRobPos * (fGeo->MCMmax()) + fMcmPos;
-  
-  trackletWord = new UInt_t[fMaxTracklets];
-  adcChannel   = new Int_t[fMaxTracklets];
+  UInt_t       *trackletWord = new UInt_t[fMaxTracklets];
+  Int_t        *adcChannel   = new Int_t[fMaxTracklets];
+  Int_t        *trackRef     = new Int_t[fMaxTracklets];
+
+  Int_t u = 0;
+
+  AliTRDdigitsManager *digman = new AliTRDdigitsManager();
+  digman->ReadDigits(AliRunLoader::Instance()->GetLoader("TRDLoader")->TreeD());
+  digman->SetUseDictionaries(kTRUE);
+  AliTRDfeeParam *feeParam = AliTRDfeeParam::Instance();
 
   for (Int_t j = 0; j < fMaxTracklets; j++) {
       Int_t i = order[j];
@@ -2463,13 +2554,41 @@ void AliTRDmcmSim::Tracklet(){
       if (bitWord[j]!=0) {
          trackletWord[u] = bitWord[j];
          adcChannel[u]   = mADC[i];   // mapping onto the original adc-array to be in line with the digits-adc-ordering (21 channels in total on 1 mcm, 18 belonging to pads); mADC[i] should be >-1 in case bitWord[i]>0
-         
-         //fMCMT[u] = bitWord[j];
+
+// Finding label of MC track
+         TH1F *hTrkRef = new TH1F("trackref", "trackref", 100000, 0, 100000);
+         Int_t track[3];
+         Int_t padcol = feeParam->GetPadColFromADC(fRobPos, fMcmPos, adcChannel[u]);
+         Int_t padcol_ngb = feeParam->GetPadColFromADC(fRobPos, fMcmPos, adcChannel[u] - 1);
+         Int_t padrow = 4 * (fRobPos / 2) + fMcmPos / 4;
+         Int_t det = 30 * fSector + 6 * fStack + fLayer;
+         for(Int_t iTimebin = feeParam->GetLinearFitStart(); iTimebin < feeParam->GetLinearFitEnd(); iTimebin++) {
+             track[0] = digman->GetTrack(0, padrow, padcol, iTimebin, det);
+             track[1] = digman->GetTrack(1, padrow, padcol, iTimebin, det);
+             track[2] = digman->GetTrack(2, padrow, padcol, iTimebin, det);
+             hTrkRef->Fill(track[0]);
+             if (track[1] != track[0] && track[1] != -1)
+                 hTrkRef->Fill(track[1]);
+             if (track[2] != track[0] && track[2] != track[1] && track[2] != -1)
+                 hTrkRef->Fill(track[2]);
+             if (padcol_ngb >= 0) {
+                 track[0] = digman->GetTrack(0, padrow, padcol, iTimebin, det);
+                 track[1] = digman->GetTrack(1, padrow, padcol, iTimebin, det);
+                 track[2] = digman->GetTrack(2, padrow, padcol, iTimebin, det);
+                 hTrkRef->Fill(track[0]);
+                 if (track[1] != track[0] && track[1] != -1)
+                     hTrkRef->Fill(track[1]);
+                 if (track[2] != track[0] && track[2] != track[1] && track[2] != -1)
+                     hTrkRef->Fill(track[2]);
+             }
+         }
+         trackRef[u] = hTrkRef->GetMaximumBin() - 1;
+         delete hTrkRef;
          u = u + 1;
       }
   }
 
AliDataLoader *dl = gAlice->GetRunLoader()->GetLoader("TRDLoader")->GetDataLoader("tracklets");
 AliDataLoader *dl = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets");
   if (!dl) {
     AliError("Could not get the tracklets data loader!");
   }
@@ -2492,22 +2611,154 @@ void AliTRDmcmSim::Tracklet(){
        trkl->SetDetector(30*fSector + 6*fStack + fLayer);
        trkl->SetROB(fRobPos);
        trkl->SetMCM(fMcmPos);
+       trkl->SetLabel(trackRef[iTracklet]);
        trackletTree->Fill();
-//     AliInfo(Form("Filling tracklet tree with trkl: %i", iTracklet));
     }
     delete trkl;
     dl->WriteData("OVERWRITE");
   }
 
+  delete [] trackletWord;
+  delete [] adcChannel; 
+  delete [] trackRef;
+  delete digman;
 
   // to be done:
   // error measure for quality of fit (not necessarily needed for the trigger)
   // cluster quality threshold (not yet set)
   // electron probability
 }
+//_____________________________________________________________________________________
+void AliTRDmcmSim::GeneratefZSM1Dim()
+{
+  //
+  // Generate the array fZSM1Dim necessary
+  // for the method ProduceRawStream
+  //
+
+  // Fill the mapping
+  // Supressed zeros indicated by -1 in digits array
+  for( Int_t iadc = 1 ; iadc < fNADC-1; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         
+         if(fADCF[iadc][it]==-1)  // If is a supressed value
+           {
+             fZSM[iadc][it]=1;
+           }
+         else                    // Not suppressed
+           {
+             fZSM[iadc][it]=0;
+           }
+       }
+    }
+
+  // Make the 1 dim projection
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         fZSM1Dim[iadc] &= fZSM[iadc][it];
+       }
+    }
+}
+//_______________________________________________________________________________________
+void AliTRDmcmSim::CopyArrays()
+{
+  //
+  // Initialize filtered data array with raw data
+  // Method added for internal consistency
+  //
+
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         fADCF[iadc][it] = fADCR[iadc][it]; 
+       }
+    }
+}
+//_______________________________________________________________________________________
+void AliTRDmcmSim::StartfastZS(Int_t pads, Int_t timebins)
+{
+  //
+  // Initialize just the necessary elements to perform
+  // the zero suppression in the digitizer
+  //
+   
+  fFeeParam  = AliTRDfeeParam::Instance();
+  fSimParam  = AliTRDSimParam::Instance();
+  fNADC      = pads;      
+  fNTimeBin  = timebins; 
 
+  if( fADCR == NULL ) 
+    {
+      fADCR    = new Int_t *[fNADC];
+      fADCF    = new Int_t *[fNADC];
+      fADCT    = 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];
+       fADCT[iadc] = new Int_t[fNTimeBin]; 
+       fZSM [iadc] = new Int_t[fNTimeBin];
+      }
+    }
 
+  for( Int_t iadc = 0 ; iadc < fNADC; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         fADCR[iadc][it] =  0;
+         fADCF[iadc][it] =  0;
+         fADCT[iadc][it] = -1;  
+         fZSM [iadc][it] =  1;   
+       }
+      fZSM1Dim[iadc] = 1;      
+    }
+  
+  fInitialized = kTRUE;
+}
+//_______________________________________________________________________________________
+void AliTRDmcmSim::FlagDigitsArray(AliTRDarrayADC *tempdigs, Int_t valrow)
+{
+  //
+  // Modify the digits array to flag suppressed values
+  //
 
+  for( Int_t iadc = 1 ; iadc < fNADC-1; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         if(fZSM[iadc][it]==1)
+           {
+             tempdigs->SetData(valrow,iadc,it,-1);
+           }
+       }
+    }
+}
+//_______________________________________________________________________________________
+void AliTRDmcmSim::RestoreZeros()
+{
+  //
+  // Restore the zero-suppressed values (set as -1) to the value 0
+  //
 
+  for( Int_t iadc = 1 ; iadc < fNADC-1; iadc++ ) 
+    {
+      for( Int_t it = 0 ; it < fNTimeBin ; it++ ) 
+       {
+         
+         if(fADCF[iadc][it]==-1)  //if is a supressed zero, reset to zero
+           {
+             fADCF[iadc][it]=0;
+             fADCR[iadc][it]=0;
+           }     
+       }
+    }
 
+}