]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOS.cxx
Addition of electronic noise corrected
[u/mrichter/AliRoot.git] / PHOS / AliPHOS.cxx
index 4e79940098362a9da57bbfcde48b2073d3c105ac..b72a77ab934da1cfd19501f87221971b9a376dc7 100644 (file)
@@ -30,20 +30,20 @@ class TFile;
 #include <TFolder.h> 
 #include <TTree.h>
 #include <TVirtualMC.h> 
+#include <TRandom.h> 
 
 // --- Standard library ---
 
 // --- AliRoot header files ---
 #include "AliMagF.h"
-#include "AliESDtrack.h" 
-#include "AliESD.h"
 #include "AliPHOS.h"
 #include "AliPHOSGetter.h"
-#include "AliPHOSQAChecker.h"
 #include "AliRun.h"
 #include "AliPHOSDigitizer.h"
 #include "AliPHOSSDigitizer.h"
-#include "AliPHOSReconstructioner.h"
+#include "AliPHOSDigit.h"
+#include "AliAltroBuffer.h"
+
 ClassImp(AliPHOS)
 //____________________________________________________________________________
 AliPHOS:: AliPHOS() : AliDetector()
@@ -52,6 +52,8 @@ AliPHOS:: AliPHOS() : AliDetector()
   fName="PHOS";
   fQATask = 0;
   fTreeQA = 0;
+  fDebug  = 0; 
+  fRan    = new TRandom(123456) ; 
 }
 
 //____________________________________________________________________________
@@ -61,17 +63,20 @@ AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
   
   fQATask = 0;
   fTreeQA = 0;
+  fDebug  = 0; 
+  fRan    = new TRandom(123456) ; 
 }
 
 //____________________________________________________________________________
 AliPHOS::~AliPHOS() 
 {  
-  
+  delete fRan ; 
 }
 
 //____________________________________________________________________________
 void AliPHOS::Copy(AliPHOS & phos)
 {
+  // copy method to be used byy the cpy ctor
   TObject::Copy(phos) ; 
   //  fQATask = AliPHOSQAChecker::Copy(*(phos.fQATask)) ; 
   phos.fTreeQA = fTreeQA->CloneTree() ; 
@@ -221,6 +226,7 @@ void AliPHOS::CreateMaterials()
   wGM[1] = 1. - arContent ;
   dGM    = wGM[0] * dAr + wGM[1] * dCO;
 
+  
   AliMixture(16, "ArCO2$", aGM, zGM, dGM,  2, wGM) ;
 
   // --- Stainless steel (let it be pure iron) ---
@@ -246,12 +252,14 @@ void AliPHOS::CreateMaterials()
   AliMixture(19, "Cables  $", aCA, zCA, dCA, -4, wCA) ;
 
 
-
   // --- Air ---
-  AliMaterial(99, "Air$", 14.61, 7.3, 0.001205, 30420., 67500., 0, 0) ;
-  
+  Float_t aAir[4]={12.0107,14.0067,15.9994,39.948};
+  Float_t zAir[4]={6.,7.,8.,18.};
+  Float_t wAir[4]={0.000124,0.755267,0.231781,0.012827};
+  Float_t dAir = 1.20479E-3;
  
+  AliMixture(99, "Air$", aAir, zAir, dAir, 4, wAir) ;
+
   // DEFINITION OF THE TRACKING MEDIA
 
   // for PHOS: idtmed[699->798] equivalent to fIdtmed[0->100]
@@ -375,50 +383,137 @@ void AliPHOS::CreateMaterials()
 }
 
 //____________________________________________________________________________
-void AliPHOS::FillESD(AliESD* esd) const 
+void AliPHOS::Digits2Raw()
 {
-  // Called by AliReconstruct after Reconstruct() and global tracking and vertxing
-
-  //Creates the tracksegments and Recparticles
-  AliPHOSReconstructioner * rec = new AliPHOSReconstructioner((fLoader->GetRunLoader()->GetFileName()).Data()) ; 
-  TList * taskslist = rec->GetListOfTasks() ; 
-  Int_t index ;
-  TTask * task ; 
-  TString name ; 
-  // set clusterizer task inactive 
-  for (index = 0; index < taskslist->GetSize(); index++) {
-    task = dynamic_cast<TTask *>(taskslist->At(index)) ; 
-    name = task->GetName() ; 
-    if ( name.Contains(AliConfig::Instance()->GetReconstructionerTaskName()))
-      task->SetActive(kFALSE) ; 
-    if ( name.Contains(AliConfig::Instance()->GetTrackerTaskName()))
-      (dynamic_cast<AliPHOSTrackSegmentMaker *> (task))->SetESD(esd) ; 
+// convert digits of the current event to raw data
+
+  // get the digits
+  AliPHOSGetter * gime = AliPHOSGetter::Instance(AliRunLoader::GetGAliceName()) ; 
+  if (!gime) {
+    Error("Digits2Raw", "PHOS Getter not instantiated") ;
+    return ; 
+  }
+  gime->Event(gime->EventNumber(), "D") ; 
+  TClonesArray* digits = gime->Digits() ;
+
+  if (!digits) {
+    Error("Digits2Raw", "no digits found !");
+    return;
+  }
+
+  // get the geometry
+  AliPHOSGeometry* geom = gime->PHOSGeometry();
+  if (!geom) {
+    Error("Digits2Raw", "no geometry found !");
+    return;
   }
 
-  AliPHOSGetter *gime = AliPHOSGetter::Instance( (fLoader->GetRunLoader()->GetFileName()).Data() ) ;
-  Int_t eventNumber = gime->EventNumber();
-  rec->SetEventRange(eventNumber, eventNumber) ; // do all the events
-  rec->ExecuteTask("deb") ; 
-
-  // Creates AliESDtrack from AliPHOSRecParticles 
-  gime->Event(eventNumber, "P") ; 
-  TClonesArray *recParticles = gime->RecParticles();
-  Int_t nOfRecParticles = recParticles->GetEntries();
-  for (Int_t recpart = 0 ; recpart < nOfRecParticles ; recpart++) {
-    AliPHOSRecParticle * rp = dynamic_cast<AliPHOSRecParticle*>(recParticles->At(recpart));
-    rp->Print();
-    AliESDtrack * et = new AliESDtrack() ; 
-    // fills the ESDtrack
-    Double_t xyz[3];
-    for (Int_t ixyz=0; ixyz<3; ixyz++) xyz[ixyz] = rp->GetPos()[ixyz];
-    et->SetPHOSposition(xyz) ; 
-    et->SetPHOSsignal  (rp->Energy()) ; 
-    et->SetPHOSpid     (rp->GetPID()) ; 
-    // add the track to the esd object
-    esd->AddTrack(et);
-    delete et;
+  // some digitization constants
+  const Int_t    kDDLOffset = 0x600;
+  const Double_t kTimeMax = 1.28E-5;
+  const Int_t    kTimeBins = 256;
+  const Double_t kTimePeak = 2.0E-6;
+  const Double_t kTimeRes = 1.5E-6;
+  const Int_t    kThreshold = 3;
+  const Int_t    kHighGainFactor = 40;
+  const Int_t    kHighGainOffset = 0x200;
+
+  AliAltroBuffer* buffer = NULL;
+  Int_t prevDDL = -1;
+  Int_t adcValuesLow[kTimeBins];
+  Int_t adcValuesHigh[kTimeBins];
+
+  // loop over digits (assume ordered digits)
+  for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
+    AliPHOSDigit* digit = gime->Digit(iDigit);
+    if (digit->GetAmp() < kThreshold) 
+      continue;
+    Int_t relId[4];
+    geom->AbsToRelNumbering(digit->GetId(), relId);
+    Int_t module = relId[0];
+   // Begin FIXME 
+    if (relId[1] != 0) 
+      continue;    // ignore digits from CPV
+   // End FIXME 
+
+    // PHOS EMCA has 4 DDL per module. Splitting is done based on the row number
+    Int_t iDDL = 4 * (module - 1) + (4 * (relId[2] - 1)) / geom->GetNPhi();
+
+    // new DDL
+    if (iDDL != prevDDL) {
+      // write real header and close previous file
+      if (buffer) {
+       buffer->Flush();
+       buffer->WriteDataHeader(kFALSE, kFALSE);
+       delete buffer;
+      }
+
+      // open new file and write dummy header
+      TString fileName("PHOS_") ;
+      fileName += (iDDL + kDDLOffset) ; 
+      fileName += ".ddl" ; 
+      buffer = new AliAltroBuffer(fileName.Data(), 1);
+      buffer->WriteDataHeader(kTRUE, kFALSE);  //Dummy;
+
+      prevDDL = iDDL;
+    }
+
+    // out of time range signal (?)
+    if (digit->GetTime() > kTimeMax) {
+      buffer->FillBuffer(digit->GetAmp());
+      buffer->FillBuffer(kTimeBins);  // time bin
+      buffer->FillBuffer(3);          // bunch length
+      buffer->WriteTrailer(3, relId[3], relId[2], module);  // trailer
+      
+    // simulate linear rise and gaussian decay of signal
+    } else {
+      Bool_t highGain = kFALSE;
+
+      // fill time bin values :
+      // 1. the signal starts at the time given by the digit
+      // 2. the rise is linear and the maximum is reached kTimePeak after start
+      // 3. the decay is gaussian with a sigma of kTimeRes
+      // 4. the signal is binned into kTimeBins bins 
+      for (Int_t iTime = 0; iTime < kTimeBins; iTime++) {
+       Double_t time = iTime * kTimeMax/kTimeBins;
+       Int_t signal = 0;
+       if (time < digit->GetTime() + kTimePeak) {      // signal is rising
+         signal = static_cast<Int_t>(fRan->Rndm() + digit->GetAmp() * 
+                        (time - digit->GetTime()) / kTimePeak);
+       } else {                                        // signal is decaying
+         signal = static_cast<Int_t>(fRan->Rndm() + digit->GetAmp() * 
+                 TMath::Gaus(time, digit->GetTime() + kTimePeak, kTimeRes));
+       }
+       if (signal < 0) 
+         signal = 0;
+       adcValuesLow[iTime] = signal;
+       if (signal > 0x3FF) // larger than 10 bits 
+         adcValuesLow[iTime] = 0x3FF;
+       adcValuesHigh[iTime] = signal / kHighGainFactor;
+       if (adcValuesHigh[iTime] > 0) 
+         highGain = kTRUE;
+      }
+
+      // write low and eventually high gain channel
+      buffer->WriteChannel(relId[3], relId[2], module, 
+                          kTimeBins, adcValuesLow, kThreshold);
+      if (highGain) {
+       buffer->WriteChannel(relId[3], relId[2], module + kHighGainOffset, 
+                            kTimeBins, adcValuesHigh, 1);
+      }
+    }
+  }
+
+  // write real header and close last file
+  if (buffer) {
+    buffer->Flush();
+    buffer->WriteDataHeader(kFALSE, kFALSE);
+    delete buffer;
   }
-}       
+
+  gime->PhosLoader()->UnloadDigits();
+}
 
 //____________________________________________________________________________
 void AliPHOS::Hits2SDigits()  
@@ -461,31 +556,6 @@ void AliPHOS::SetTreeAddress()
   }
 }
 
-//____________________________________________________________________________
-void AliPHOS::Reconstruct() const 
-{ 
-  // method called by AliReconstruction; it should not be called otherwise but the AliPHOSReconstructionner 
-  // should be called for stand alone reconstruction.
-  // Only the clusterization is performed,; the rest of the reconstruction is done in FillESD because the track
-  // segment maker needs access to the AliESD object to retrieve the tracks reconstructed by 
-  // the global tracking.
-  AliPHOSReconstructioner * rec = new AliPHOSReconstructioner((fLoader->GetRunLoader()->GetFileName()).Data()) ; 
-  TList * taskslist = rec->GetListOfTasks() ; 
-  Int_t index ;
-  TTask * task ; 
-  TString name ; 
-  // set all tasks inactive except clusterizer
-  for (index = 0; index < taskslist->GetSize(); index++) {
-    task = dynamic_cast<TTask *>(taskslist->At(index)) ; 
-    name = task->GetName() ; 
-    if ( !name.Contains(AliConfig::Instance()->GetReconstructionerTaskName()))
-      task->SetActive(kFALSE) ; 
-  }
-  rec->SetEventRange(0, -1) ; // do all the events
-  rec->ExecuteTask("deb all") ; 
-}
-
 //____________________________________________________________________________
 void AliPHOS::WriteQA()
 {
@@ -516,3 +586,4 @@ void AliPHOS::WriteQA()
 
   //fTreeQA->Fill() ; 
 }
+