]> 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 ac9a7dceec373d883a00e3dd5668e0491679cb0b..b72a77ab934da1cfd19501f87221971b9a376dc7 100644 (file)
@@ -12,7 +12,6 @@
  * about the suitability of this software for any purpose. It is          *
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
-
 /* $Id$ */
 
 //_________________________________________________________________________
 //*-- Author: Laurent Aphecetche & Yves Schutz (SUBATECH) 
 //////////////////////////////////////////////////////////////////////////////
 
+
 // --- ROOT system ---
 class TFile;
-#include "TTree.h"
-#include "TTask.h"
-#include "TROOT.h"
-#include "TFolder.h"
+#include <TFolder.h> 
+#include <TTree.h>
+#include <TVirtualMC.h> 
+#include <TRandom.h> 
 
 // --- Standard library ---
 
-#include <strstream.h>
 // --- AliRoot header files ---
-
+#include "AliMagF.h"
 #include "AliPHOS.h"
-#include "AliMC.h"
+#include "AliPHOSGetter.h"
 #include "AliRun.h"
-#include "AliMagF.h"
-#include "AliPHOSGeometry.h"
+#include "AliPHOSDigitizer.h"
+#include "AliPHOSSDigitizer.h"
+#include "AliPHOSDigit.h"
+#include "AliAltroBuffer.h"
 
 ClassImp(AliPHOS)
+//____________________________________________________________________________
+AliPHOS:: AliPHOS() : AliDetector()
+{
+  // Default ctor
+  fName="PHOS";
+  fQATask = 0;
+  fTreeQA = 0;
+  fDebug  = 0; 
+  fRan    = new TRandom(123456) ; 
+}
 
+//____________________________________________________________________________
+AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title) 
+{
+  //   ctor : title is used to identify the layout
+  
+  fQATask = 0;
+  fTreeQA = 0;
+  fDebug  = 0; 
+  fRan    = new TRandom(123456) ; 
+}
+
+//____________________________________________________________________________
+AliPHOS::~AliPHOS() 
+{  
+  delete fRan ; 
+}
 
 //____________________________________________________________________________
-  AliPHOS::AliPHOS():AliDetector()
+void AliPHOS::Copy(AliPHOS & phos)
 {
-  // ctor. creates task PHOS
-  TTask * phosTask = new TTask("PHOS","");
-  TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; 
-  roottasks->Add(phosTask) ; 
+  // copy method to be used byy the cpy ctor
+  TObject::Copy(phos) ; 
+  //  fQATask = AliPHOSQAChecker::Copy(*(phos.fQATask)) ; 
+  phos.fTreeQA = fTreeQA->CloneTree() ; 
 }
 
 //____________________________________________________________________________
-  AliPHOS::AliPHOS(const char* name, const char* title=""):AliDetector(name, title)
+AliDigitizer* AliPHOS::CreateDigitizer(AliRunDigitizer* manager) const
 {
-  // ctor. creates task PHOS
-  TTask * phosTask = new TTask("PHOS","");
-  TTask * roottasks = (TTask*)gROOT->GetRootFolder()->FindObject("Tasks") ; 
-  roottasks->Add(phosTask) ; 
+  return new AliPHOSDigitizer(manager);
 }
 
 //____________________________________________________________________________
@@ -124,7 +148,7 @@ void AliPHOS::CreateMaterials()
   Float_t aTI[2] = {12.011, 1.00794} ;
   Float_t zTI[2] = {6.0, 1.0} ;
   Float_t wTI[2] = {1.0, 1.0} ;
-  Float_t dTI = 0.1 ;
+  Float_t dTI = 0.04 ;
 
   AliMixture(7, "Thermo Insul.$", aTI, zTI, dTI, -2, wTI) ;
 
@@ -180,7 +204,7 @@ void AliPHOS::CreateMaterials()
   AliMaterial(15, "Ar$", 39.948, 18.0, dAr, 14.0, 0., 0, 0) ;   
  
  // ArCo2
-  Char_t namate[21];
+  Char_t namate[21]="";
   Float_t aGM[2] ; 
   Float_t zGM[2] ; 
   Float_t wGM[2] ; 
@@ -202,15 +226,40 @@ 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) ---
   AliMaterial(17, "Steel$", 55.845, 26, 7.87, 1.76, 0., 0, 0) ;
+
+
+  // --- Fiberglass ---
+  Float_t aFG[4] = {16.0, 28.09, 12.011, 1.00794} ;
+  Float_t zFG[4] = {8.0, 14.0, 6.0, 1.0} ;
+  Float_t wFG[4] = {292.0, 68.0, 462.0, 736.0} ;
+  Float_t dFG    = 1.9 ;
+
+  AliMixture(18, "Fibergla$", aFG, zFG, dFG, -4, wFG) ;
+
+  // --- Cables in Air box  ---
+  // SERVICES
+
+  Float_t aCA[4] = { 1.,12.,55.8,63.5 };
+  Float_t zCA[4] = { 1.,6.,26.,29. }; 
+  Float_t wCA[4] = { .014,.086,.42,.48 };
+  Float_t dCA    = 0.8 ;  //this density is raw estimation, if you know better - correct
+
+  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]
@@ -286,6 +335,14 @@ void AliPHOS::CreateMaterials()
   AliMedium(17, "Steel     $", 17, 0,
             isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.0001, 0, 0) ;
 
+  // Fibergalss                                                                     -> idtmed[717]
+  AliMedium(18, "Fiberglass$", 18, 0,
+            isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
+
+  // Cables in air                                                                  -> idtmed[718]
+  AliMedium(19, "Cables    $", 19, 0,
+            isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
+
   // Air                                                                            -> idtmed[798] 
   AliMedium(99, "Air          $", 99, 0,
             isxfld, sxmgmx, 10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
@@ -324,23 +381,209 @@ void AliPHOS::CreateMaterials()
   gMC->Gstpar(idtmed[715], "STRA",2.) ;
 
 }
+
 //____________________________________________________________________________
-void AliPHOS::SetTreeAddress()
+void AliPHOS::Digits2Raw()
+{
+// 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;
+  }
+
+  // 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()  
 { 
-  // called by AliRun
-  // TBranch *branch;
-  //  AliDetector::SetTreeAddress();
+// create summable digits
 
+  AliPHOSSDigitizer* phosDigitizer = 
+    new AliPHOSSDigitizer(fLoader->GetRunLoader()->GetFileName().Data()) ;
+  phosDigitizer->SetEventRange(0, -1) ; // do all the events
+  phosDigitizer->ExecuteTask("all") ;
+}
+
+//____________________________________________________________________________
+AliLoader* AliPHOS::MakeLoader(const char* topfoldername)
+{
+//different behaviour than standard (singleton getter)
+// --> to be discussed and made eventually coherent
+ fLoader = new AliPHOSLoader(GetName(),topfoldername);
+ return fLoader;
+}
+
+//____________________________________________________________________________
+void AliPHOS::SetTreeAddress()
+{ 
+  // Links Hits in the Tree to Hits array
   TBranch *branch;
   char branchname[20];
   sprintf(branchname,"%s",GetName());
   
   // Branch address for hit tree
-  TTree *treeH = gAlice->TreeH();
-  if (treeH && fHits) {
+  TTree *treeH = TreeH();
+  if (treeH) {
     branch = treeH->GetBranch(branchname);
-    if (branch) branch->SetAddress(&fHits);
+    if (branch) 
+     { 
+       if (fHits == 0x0) fHits= new TClonesArray("AliPHOSHit",1000);
+       //Info("SetTreeAddress","<%s> Setting Hits Address",GetName());
+       branch->SetAddress(&fHits);
+     }
   }
 }
 
+//____________________________________________________________________________
+void AliPHOS::WriteQA()
+{
+
+  // Make TreeQA in the output file. 
+
+  if(fTreeQA == 0)
+    fTreeQA = new TTree("TreeQA", "QA Alarms") ;    
+  // Create Alarms branches
+  Int_t bufferSize = 32000 ;    
+  Int_t splitlevel = 0 ; 
+
+  TFolder* topfold = GetLoader()->GetTopFolder(); //get top aliroot folder; skowron
+  TString phosqafn(AliConfig::Instance()->GetQAFolderName()+"/"); //get name of QAaut folder relative to top event; skowron
+  phosqafn+=GetName(); //hard wired string!!! add the detector name to the pathname; skowron 
+  TFolder * alarmsF = (TFolder*)topfold->FindObjectAny(phosqafn); //get the folder
+  if (alarmsF == 0x0)
+   {
+     Error("WriteQA","Can not find folder with qa alarms");
+     return;
+   }
+  TString branchName(alarmsF->GetName());
+  TBranch * alarmsBranch = fTreeQA->Branch(branchName,"TFolder", &alarmsF, bufferSize, splitlevel);
+  TString branchTitle = branchName + " QA alarms" ; 
+  alarmsBranch->SetTitle(branchTitle);
+  alarmsBranch->Fill() ; 
+
+  //fTreeQA->Fill() ; 
+}