]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - EMCAL/AliEMCAL.cxx
forgot to include baseLinkDef file update in yesterdays commit; sorry, and thanks...
[u/mrichter/AliRoot.git] / EMCAL / AliEMCAL.cxx
index d34a7bc66e62cfc158f0d146967fa2cdf024409a..b030d72638f9b7d3b837c07ef13c7873c1fa1ac2 100644 (file)
  **************************************************************************/
 
 /* $Id$ */
-
+/* History of cvs commits:
+ *
+ * $Log$
+
+ * Revision 1.54  2007/12/06 10:31:13  hristov
+ * Bug fix: using the mapping from CDB
+ *
+ * Revision 1.53.10.1  2007/12/06 10:29:59  hristov
+ * Bug fix: using the mapping from CDB
+ *
+ * Revision 1.53  2007/03/17 19:56:38  mvl
+ * Moved signal shape routines from AliEMCAL to separate class AliEMCALRawUtils to streamline raw data reconstruction code.
+ *
+ * Revision 1.52  2007/03/10 22:19:01  pavlinov
+ * move one varibels from AliEMCALv2 to AliEMCAL
+ *
+ * Revision 1.51  2007/02/24 20:42:35  pavlinov
+ * fixed error of Geant3 parameters initialisation
+ *
+ * Revision 1.50  2007/02/05 10:43:25  hristov
+ * Changes for correct initialization of Geant4 (Mihaela)
+ *
+ * Revision 1.49  2007/01/22 17:29:12  pavlinov
+ * EMCAL geometry can be created independently form anything now
+ *
+ * Revision 1.48  2006/12/19 02:34:13  pavlinov
+ * clean up the EMCAL name scheme : super module -> module -> tower (or cell)
+ *
+ * Revision 1.47  2006/12/05 17:12:03  gustavo
+ * Updated AliEMCAL::Digits2Raw, reads first provisional RCU mapping files to make Raw data with new AliCaloAltroMapping and AliCaloRawStream
+ *
+ *
+ */
 //_________________________________________________________________________
 // Base Class for EMCAL description:
 // This class contains material definitions    
 // --- ROOT system ---
 class TFile;
 #include <TFolder.h> 
+#include <TGeoGlobalMagField.h>
+#include <TGraph.h> 
+#include <TH1F.h> 
+#include <TRandom.h> 
 #include <TTree.h>
 #include <TVirtualMC.h> 
 
@@ -36,37 +72,140 @@ class TFile;
 // --- AliRoot header files ---
 #include "AliMagF.h"
 #include "AliEMCAL.h"
-#include "AliEMCALGeometry.h"
-#include "AliEMCALLoader.h"
 #include "AliRun.h"
+#include "AliEMCALLoader.h"
 #include "AliEMCALSDigitizer.h"
 #include "AliEMCALDigitizer.h"
+#include "AliEMCALDigit.h"
+#include "AliEMCALRawUtils.h"
+#include "AliCDBManager.h"
+#include "AliCDBEntry.h"
 
 ClassImp(AliEMCAL)
+
 //____________________________________________________________________________
-AliEMCAL::AliEMCAL():AliDetector()
+AliEMCAL::AliEMCAL()
+  : AliDetector(),
+    fBirkC0(0),
+    fBirkC1(0.),
+    fBirkC2(0.),
+    fGeometry(0)
 {
   // Default ctor 
-  fName="EMCAL";
-  //fQATask = 0;
-  fTreeQA = 0;
-  fGeom = 0 ; 
+  fName = "EMCAL" ;
+  InitConstants();
+
+  // Should call  AliEMCALGeometry::GetInstance(EMCAL->GetTitle(),"") for getting EMCAL geometry
 }
 
 //____________________________________________________________________________
-AliEMCAL::AliEMCAL(const char* name, const char* title): AliDetector(name,title)
+AliEMCAL::AliEMCAL(const char* name, const char* title)
+  : AliDetector(name,title),
+    fBirkC0(0),
+    fBirkC1(0.),
+    fBirkC2(0.),
+    fGeometry(0)
 {
   //   ctor : title is used to identify the layout
-  
-  //fQATask = 0;
-  fTreeQA = 0;
-  fGeom = 0;
+  InitConstants();
 }
 
 //____________________________________________________________________________
 AliEMCAL::~AliEMCAL()
 {
+  //dtor
+}
+
+//____________________________________________________________________________
+void AliEMCAL::InitConstants()
+{
+  //initialize EMCAL values
+  fBirkC0 = 1;
+  fBirkC1 = 0.013/1.032;
+  fBirkC2 = 9.6e-6/(1.032 * 1.032);
+  }
+
+//____________________________________________________________________________
+void AliEMCAL::DefineMediumParameters()
+{
+  //
+  // EMCAL cuts (Geant3) 
+  // 
+  Int_t * idtmed = fIdtmed->GetArray() - 1599 ; 
+// --- Set decent energy thresholds for gamma and electron tracking
+
+  // Tracking threshold for photons and electrons in Lead 
+  Float_t cutgam=10.e-5; // 100 kev;
+  Float_t cutele=10.e-5; // 100 kev;
+  TString ntmp(GetTitle()); 
+  ntmp.ToUpper();
+  if(ntmp.Contains("10KEV")) {
+    cutele = cutgam = 1.e-5;
+  } else if(ntmp.Contains("50KEV")) {
+    cutele = cutgam = 5.e-5;
+  } else if(ntmp.Contains("100KEV")) {
+    cutele = cutgam = 1.e-4;
+  } else if(ntmp.Contains("200KEV")) {
+    cutele = cutgam = 2.e-4;
+  } else if(ntmp.Contains("500KEV")) {
+    cutele = cutgam = 5.e-4;
+  }
+
+  gMC->Gstpar(idtmed[1600],"CUTGAM", cutgam);
+  gMC->Gstpar(idtmed[1600],"CUTELE", cutele); // 1MEV -> 0.1MEV; 15-aug-05
+  gMC->Gstpar(idtmed[1600],"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1600],"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  // --- Generate explicitly delta rays in Lead ---
+  gMC->Gstpar(idtmed[1600], "LOSS", 3) ;
+  gMC->Gstpar(idtmed[1600], "DRAY", 1) ;
+  gMC->Gstpar(idtmed[1600], "DCUTE", cutele) ;
+  gMC->Gstpar(idtmed[1600], "DCUTM", cutele) ;
 
+// --- in aluminium parts ---
+  gMC->Gstpar(idtmed[1602],"CUTGAM", cutgam) ;
+  gMC->Gstpar(idtmed[1602],"CUTELE", cutele) ;
+  gMC->Gstpar(idtmed[1602],"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1602],"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1602], "LOSS",3.) ;
+  gMC->Gstpar(idtmed[1602], "DRAY",1.) ;
+  gMC->Gstpar(idtmed[1602], "DCUTE", cutele) ;
+  gMC->Gstpar(idtmed[1602], "DCUTM", cutele) ;
+
+// --- and finally thresholds for photons and electrons in the scintillator ---
+  gMC->Gstpar(idtmed[1601],"CUTGAM", cutgam) ;
+  gMC->Gstpar(idtmed[1601],"CUTELE", cutele) ;// 1MEV -> 0.1MEV; 15-aug-05
+  gMC->Gstpar(idtmed[1601],"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1601],"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1601], "LOSS",3) ; // generate delta rays 
+  gMC->Gstpar(idtmed[1601], "DRAY",1) ;
+  gMC->Gstpar(idtmed[1601], "DCUTE", cutele) ;
+  gMC->Gstpar(idtmed[1601], "DCUTM", cutele) ;
+
+  // S steel - 
+  gMC->Gstpar(idtmed[1603],"CUTGAM", cutgam);
+  gMC->Gstpar(idtmed[1603],"CUTELE", cutele);
+  gMC->Gstpar(idtmed[1603],"BCUTE",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  gMC->Gstpar(idtmed[1603],"BCUTM",  cutgam);  // BCUTE and BCUTM start from GUTGUM
+  // --- Generate explicitly delta rays 
+  gMC->Gstpar(idtmed[1603], "LOSS",3);
+  gMC->Gstpar(idtmed[1603], "DRAY",1);
+  gMC->Gstpar(idtmed[1603], "DCUTE", cutele) ;
+  gMC->Gstpar(idtmed[1603], "DCUTM", cutele) ;
+
+  AliEMCALGeometry* geom = GetGeometry();
+  if(geom->GetILOSS()>=0) {
+    for(int i=1600; i<=1603; i++) gMC->Gstpar(idtmed[i], "LOSS", geom->GetILOSS()) ; 
+  } 
+  if(geom->GetIHADR()>=0) {
+    for(int i=1600; i<=1603; i++) gMC->Gstpar(idtmed[i], "HADR", geom->GetIHADR()) ; 
+  }
+}
+
+//____________________________________________________________________________
+AliDigitizer* AliEMCAL::CreateDigitizer(AliRunDigitizer* manager) const
+{
+  //create and return the digitizer
+  return new AliEMCALDigitizer(manager);
 }
 
 //____________________________________________________________________________
@@ -74,10 +213,12 @@ void AliEMCAL::CreateMaterials()
 {
   // Definitions of materials to build EMCAL and associated tracking media.
   // media number in idtmed are 1599 to 1698.
-
-  // --- Air ---
-  AliMaterial(0, "Air$", 14.61, 7.3, 0.001205, 30420., 67500., 0, 0) ;
-
+  // --- Air ---               
+  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(0, "Air$", aAir, zAir, dAir, 4, wAir) ;
 
   // --- Lead ---                                                                     
   AliMaterial(1, "Pb$", 207.2, 82, 11.35, 0.56, 0., 0, 0) ;
@@ -95,167 +236,73 @@ void AliEMCAL::CreateMaterials()
   AliMaterial(3, "Al$", 26.98, 13., 2.7, 8.9, 999., 0, 0) ;
   // ---         Absorption length is ignored ^
 
-  // --- Copper ---
-  AliMaterial(4, "Cu$", 63.546, 29, 8.96, 1.43, 14.8, 0, 0) ; 
-  // ---         Absorption length is ignored ^
-
+  // 25-aug-04 by PAI - see  PMD/AliPMDv0.cxx for STEEL definition
+  Float_t asteel[4] = { 55.847,51.9961,58.6934,28.0855 };
+  Float_t zsteel[4] = { 26.,24.,28.,14. };
+  Float_t wsteel[4] = { .715,.18,.1,.005 };
+  AliMixture(4, "STAINLESS STEEL$", asteel, zsteel, 7.88, 4, wsteel);
 
   // DEFINITION OF THE TRACKING MEDIA
 
   // for EMCAL: idtmed[1599->1698] equivalent to fIdtmed[0->100]
-  Int_t * idtmed = fIdtmed->GetArray() - 1599 ; 
-  Int_t   isxfld = gAlice->Field()->Integ() ;
-  Float_t sxmgmx = gAlice->Field()->Max() ;
-
+  Int_t   isxfld = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Integ() ;
+  Float_t sxmgmx = ((AliMagF*)TGeoGlobalMagField::Instance()->GetField())->Max() ;
 
-   // Air                                                                           -> idtmed[1599] 
 AliMedium(0, "Air          $", 0, 0,
+  // Air                                                                         -> idtmed[1599]
AliMedium(0, "Air$", 0, 0,
             isxfld, sxmgmx, 10.0, 1.0, 0.1, 0.1, 10.0, 0, 0) ;
 
   // The Lead                                                                      -> idtmed[1600]
  
-  AliMedium(1, "Lead      $", 1, 0,
+  AliMedium(1, "Lead$", 1, 0,
             isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.1, 0.1, 0, 0) ;
 
  // The scintillator of the CPV made of Polystyrene scintillator                   -> idtmed[1601]
-  AliMedium(2, "CPV scint.   $", 2, 1,
-            isxfld, sxmgmx, 10.0, 0.001, 0.1, 0.001, 0.001, 0, 0) ;
+  float deemax = 0.1; // maximum fractional energy loss in one step (0 < DEEMAX < deemax )
+  AliMedium(2, "Scintillator$", 2, 1,
+            isxfld, sxmgmx, 10.0, 0.001, deemax, 0.001, 0.001, 0, 0) ;
 
   // Various Aluminium parts made of Al                                            -> idtmed[1602]
-  AliMedium(3, "Al parts     $", 3, 0,
+  AliMedium(3, "Al$", 3, 0,
              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
 
-  // Copper for HCal (post shower)                                                 -> idtmed[1603]
-  AliMedium(4, "Copper       $", 4, 0,
+  // 25-aug-04 by PAI : see  PMD/AliPMDv0.cxx for STEEL definition                 -> idtmed[1603]
+  AliMedium(4, "S steel$", 4, 0, 
              isxfld, sxmgmx, 10.0, 0.1, 0.1, 0.001, 0.001, 0, 0) ;
 
 
+  //set constants for Birk's Law implentation
+  fBirkC0 =  1;
+  fBirkC1 =  0.013/dP;
+  fBirkC2 =  9.6e-6/(dP * dP);
 
-// --- Set decent energy thresholds for gamma and electron tracking
-
-  // Tracking threshold for photons and electrons in Lead 
-  gMC->Gstpar(idtmed[1600],"CUTGAM",0.00008) ;
-  gMC->Gstpar(idtmed[1600],"CUTELE",0.001) ;
-  gMC->Gstpar(idtmed[1600],"BCUTE",0.0001) ;
-
-  // --- Generate explicitly delta rays in Lead ---
-  gMC->Gstpar(idtmed[1600], "LOSS",3.) ;
-  gMC->Gstpar(idtmed[1600], "DRAY",1.) ;
-  gMC->Gstpar(idtmed[1600], "DCUTE",0.00001) ;
-  gMC->Gstpar(idtmed[1600], "DCUTM",0.00001) ;
-
-// --- in aluminium parts ---
-  gMC->Gstpar(idtmed[1602], "LOSS",3.) ;
-  gMC->Gstpar(idtmed[1602], "DRAY",1.) ;
-  gMC->Gstpar(idtmed[1602], "DCUTE",0.00001) ;
-  gMC->Gstpar(idtmed[1602], "DCUTM",0.00001) ;
-
-// --- in copper parts ---
-  gMC->Gstpar(idtmed[1603], "LOSS",3.) ;
-  gMC->Gstpar(idtmed[1603], "DRAY",1.) ;
-  gMC->Gstpar(idtmed[1603], "DCUTE",0.00001) ;
-  gMC->Gstpar(idtmed[1603], "DCUTM",0.00001) ;
-
-
-
-// --- and finally thresholds for photons and electrons in the scintillator ---
-  gMC->Gstpar(idtmed[1601],"CUTGAM",0.00008) ;
-  gMC->Gstpar(idtmed[1601],"CUTELE",0.001) ;
-  gMC->Gstpar(idtmed[1601],"BCUTE",0.0001) ;
-
-
+  // Call just in case of Geant3; What to do in case of Geant4 ?
+  if(gMC->InheritsFrom("TGeant3")) DefineMediumParameters(); // Feb 20, 2007
 }
-
 //____________________________________________________________________________
-AliEMCALGeometry * AliEMCAL::GetGeometry() const 
-{  
-  // gets the pointer to the AliEMCALGeometry unique instance 
-
-  if (fGeom) 
-    return fGeom ; 
-  else 
-    return AliEMCALGeometry::GetInstance(GetTitle(),"") ;  
-}
+void AliEMCAL::Digits2Raw() {
 
-//____________________________________________________________________________
-void AliEMCAL::SetTreeAddress()
-{ 
-  // Linking Hits in Tree to Hits array
-  TBranch *branch;
-  char branchname[20];
-  sprintf(branchname,"%s",GetName());
-  
-  // Branch address for hit tree
-  TTree *treeH = TreeH();
-  if (treeH) {
-    branch = treeH->GetBranch(branchname);
-    if (branch) 
-      { 
-       if (fHits == 0x0) 
-         fHits= new TClonesArray("AliEMCALHit",1000);
-       //Info("SetTreeAddress","<%s> Setting Hits Address",GetName());
-       branch->SetAddress(&fHits);
-      }
-    else
-      {
-       Warning("SetTreeAddress","<%s> Failed",GetName());
-      }
-  }
-}
-//____________________________________________________________________________
-void AliEMCAL::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 emcalqafn(AliConfig::Instance()->GetQAFolderName()+"/"); //get name of QAaut folder relative to top event; skowron
-  emcalqafn+=GetName(); //hard wired string!!! add the detector name to the pathname; skowron 
-  TFolder * alarmsF = (TFolder*)topfold->FindObjectAny(emcalqafn); //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
-}
+  static AliEMCALRawUtils rawUtils;
+  rawUtils.Digits2Raw();
 
-//____________________________________________________________________________
-AliLoader* AliEMCAL::MakeLoader(const char* topfoldername)
-{
-//different behaviour than standard (singleton getter)
-// --> to be discussed and made eventually coherent
- fLoader = new AliEMCALLoader(GetName(),topfoldername);
- return fLoader;
 }
-
 //____________________________________________________________________________
 void AliEMCAL::Hits2SDigits()  
 { 
 // create summable digits
 
-  AliEMCALSDigitizer* emcalDigitizer = 
-    new AliEMCALSDigitizer(fLoader->GetRunLoader()->GetFileName().Data());
-  emcalDigitizer->ExecuteTask();
+  GetGeometry();
+  AliEMCALSDigitizer emcalDigitizer(fLoader->GetRunLoader()->GetFileName().Data()) ;
+  emcalDigitizer.SetEventRange(0, -1) ; // do all the events
+  emcalDigitizer.ExecuteTask() ;
 }
 
 //____________________________________________________________________________
-AliDigitizer* AliEMCAL::CreateDigitizer(AliRunDigitizer* manager)
+
+AliLoader* AliEMCAL::MakeLoader(const char* topfoldername)
 {
-  return new AliEMCALDigitizer(manager);
+//different behaviour than standard (singleton getter)
+// --> to be discussed and made eventually coherent
+ fLoader = new AliEMCALLoader(GetName(),topfoldername);
+ return fLoader;
 }