]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PHOS/AliPHOS.cxx
Charged jets (pPb): Improved trackcut analysis
[u/mrichter/AliRoot.git] / PHOS / AliPHOS.cxx
index 69ebedf1913e7686665b51e97eaa01656e91e8a0..8c0e9e8fddc460ccdd7e0dbe4b28c7607744d97a 100644 (file)
 
 // --- ROOT system ---
 class TFile;
+#include <TF1.h> 
 #include <TFolder.h> 
-#include <TTree.h>
-#include <TVirtualMC.h> 
+#include <TGeoGlobalMagField.h>
 #include <TH1F.h> 
-#include <TF1.h> 
 #include <TRandom.h> 
+#include <TTree.h>
+#include <TVirtualMC.h> 
 
 // --- Standard library ---
 
@@ -80,6 +81,7 @@ class TFile;
 #include "AliPHOS.h"
 #include "AliPHOSLoader.h"
 #include "AliRun.h"
+#include "AliRawReader.h"
 #include "AliPHOSDigitizer.h"
 #include "AliPHOSSDigitizer.h"
 #include "AliPHOSDigit.h"
@@ -93,7 +95,8 @@ class TFile;
 #include "AliPHOSCalibData.h"
 #include "AliPHOSPulseGenerator.h"
 #include "AliDAQ.h"
-#include "AliPHOSRawDecoder.h"
+#include "AliPHOSRawFitterv0.h"
+#include "AliPHOSCalibData.h"
 #include "AliPHOSRawDigiProducer.h"
 #include "AliPHOSQAChecker.h"
 #include "AliPHOSRecoParam.h"
@@ -102,7 +105,7 @@ class TFile;
 ClassImp(AliPHOS)
 
 //____________________________________________________________________________
-  AliPHOS:: AliPHOS() : AliDetector()
+  AliPHOS:: AliPHOS() : AliDetector(),fgCalibData(0)
 {
   // Default ctor
   fName   = "PHOS" ;
@@ -110,7 +113,8 @@ ClassImp(AliPHOS)
 }
 
 //____________________________________________________________________________
-AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
+AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title),
+fgCalibData(0)
 {
   //   ctor : title is used to identify the layout
 }
@@ -118,12 +122,13 @@ AliPHOS::AliPHOS(const char* name, const char* title): AliDetector(name, title)
 //____________________________________________________________________________
 AliPHOS::~AliPHOS() 
 {  
+  if(fgCalibData) delete fgCalibData ;
 }
 
 //____________________________________________________________________________
-AliDigitizer* AliPHOS::CreateDigitizer(AliRunDigitizer* manager) const
+AliDigitizer* AliPHOS::CreateDigitizer(AliDigitizationInput* digInput) const
 {
-  return new AliPHOSDigitizer(manager);
+  return new AliPHOSDigitizer(digInput);
 }
 
 //____________________________________________________________________________
@@ -288,8 +293,8 @@ void AliPHOS::CreateMaterials()
   // DEFINITION OF THE TRACKING MEDIA
 
   // for PHOS: idtmed[699->798] equivalent to fIdtmed[0->100]
-  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() ;
 
   // The scintillator of the calorimeter made of PBW04                              -> idtmed[699]
   AliMedium(0, "PHOS Xtal    $", 0, 1,
@@ -375,42 +380,6 @@ void AliPHOS::CreateMaterials()
 //_____________________________________________________________________________
 void AliPHOS::Init()
 {
-  //
-  // Initialises cuts for PHOS
-  //
-  // --- Set decent energy thresholds for gamma and electron tracking
-  Int_t * idtmed = fIdtmed->GetArray() - 699 ; 
-
-  // Tracking threshold for photons and electrons in the scintillator crystal 
-  gMC->Gstpar(idtmed[699], "CUTGAM",0.5E-4) ; 
-  gMC->Gstpar(idtmed[699], "CUTELE",1.0E-4) ;
-  // --- Generate explicitly delta rays in the titan cover ---
-  gMC->Gstpar(idtmed[704], "LOSS",3.) ;
-  gMC->Gstpar(idtmed[704], "DRAY",1.) ;
-  // --- and in aluminium parts ---
-  gMC->Gstpar(idtmed[701], "LOSS",3.) ;
-  gMC->Gstpar(idtmed[701], "DRAY",1.) ;
-  // --- and in PIN diode
-  gMC->Gstpar(idtmed[705], "LOSS",3) ;
-  gMC->Gstpar(idtmed[705], "DRAY",1) ;
-  // --- and in the passive convertor
-  gMC->Gstpar(idtmed[712], "LOSS",3) ;
-  gMC->Gstpar(idtmed[712], "DRAY",1) ;
-  // Tracking threshold for photons and electrons in the gas ArC02 
-  gMC->Gstpar(idtmed[715], "CUTGAM",1.E-5) ; 
-  gMC->Gstpar(idtmed[715], "CUTELE",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "CUTNEU",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "CUTHAD",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "CUTMUO",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "BCUTE",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "BCUTM",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "DCUTE",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "DCUTM",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "PPCUTM",1.E-5) ;
-  gMC->Gstpar(idtmed[715], "LOSS",2.) ;
-  gMC->Gstpar(idtmed[715], "DRAY",0.) ;
-  gMC->Gstpar(idtmed[715], "STRA",2.) ;
 }
 
 //____________________________________________________________________________
@@ -422,7 +391,7 @@ void AliPHOS::Digits2Raw()
     AliError("Energy digitization should be OFF if use Digits2Raw") ;
   }
 
-  AliPHOSLoader * loader = dynamic_cast<AliPHOSLoader*>(fLoader) ; 
+  AliPHOSLoader * loader = static_cast<AliPHOSLoader*>(fLoader) ; 
 
   // get the digits
   loader->LoadDigits();
@@ -445,13 +414,19 @@ void AliPHOS::Digits2Raw()
   if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
 
   // some digitization constants
-  const Float_t    kThreshold = 0.001; // skip digits below 1 MeV
+  const Float_t    kThreshold = 1.; // skip digits below 1 ADC channel
   const Int_t      kAdcThreshold = 1;  // Lower ADC threshold to write to raw data
 
   Int_t prevDDL = -1;
 
+  if(fgCalibData==0)
+    fgCalibData= new AliPHOSCalibData(-1) ;
+
   // Create a shaper pulse object
   AliPHOSPulseGenerator *pulse = new AliPHOSPulseGenerator();
+
+  //Set Time step of sample
+  pulse->SetTimeStep(TMath::Abs(fgCalibData->GetSampleTimeStep())) ;
   
   Int_t *adcValuesLow = new Int_t[pulse->GetRawFormatTimeBins()];
   Int_t *adcValuesHigh= new Int_t[pulse->GetRawFormatTimeBins()];
@@ -474,18 +449,15 @@ void AliPHOS::Digits2Raw()
 
   // loop over digits (assume ordered digits)
   for (Int_t iDigit = 0; iDigit < digits->GetEntries(); iDigit++) {
-    AliPHOSDigit* digit = dynamic_cast<AliPHOSDigit *>(digits->At(iDigit)) ;
+    AliPHOSDigit* digit = static_cast<AliPHOSDigit *>(digits->At(iDigit)) ;
 
     // Skip small energy below treshold
     if (digit->GetEnergy() < kThreshold) 
       continue;
-    // Skip CPV digits
-    if (digit->GetId() > geom->GetNModules() * geom->GetNCristalsInModule()) 
-      continue;
 
     Int_t relId[4];
     geom->AbsToRelNumbering(digit->GetId(), relId);
-    Int_t module = relId[0];
+    Int_t module = 5-relId[0];
  
     // Begin FIXME 
     if (relId[1] != 0) 
@@ -494,25 +466,25 @@ void AliPHOS::Digits2Raw()
 
     Int_t row = relId[2]-1;
     Int_t col = relId[3]-1;
-
+    
     Int_t iRCU = -111;
-
+    
     //RCU0
-    if(0<=row&&row<32 && 0<=col&&col<28) iRCU=0;
-
+    if(0<=row&&row<16 && 0<=col&&col<56) iRCU=0;
+    
     //RCU1
-    if(0<=row&&row<32 && 28<=col&&col<56) iRCU=1;
+    if(16<=row&&row<32 && 0<=col&&col<56) iRCU=1;
 
     //RCU2
-    if(32<=row&&row<64 && 0<=col&&col<28) iRCU=2;
+    if(32<=row&&row<48 && 0<=col&&col<56) iRCU=2;
 
     //RCU3
-    if(32<=row&&row<64 && 28<=col&&col<56) iRCU=3;
-
-
+    if(48<=row&&row<64 && 0<=col&&col<56) iRCU=3;
+    
+    
     // PHOS EMCA has 4 DDL per module. Splitting is based on the (row,column) numbers.
-    // PHOS internal convention: 1<module<5.
-    Int_t iDDL = 4 * (module - 1) + iRCU;
+    // here module already in PHOS online convention: 0<module<4 and inverse order.
+    Int_t iDDL = 4 * module  + iRCU;
 
     // new DDL
     if (iDDL != prevDDL) {
@@ -520,12 +492,7 @@ void AliPHOS::Digits2Raw()
        // open new file and write dummy header
        TString fileName = AliDAQ::DdlFileName("PHOS",iDDL);
 
-       TString path = gSystem->Getenv("ALICE_ROOT");
-       path += "/PHOS/mapping/RCU";
-       path += iRCU;
-       path += ".data";
-
-       mapping[iDDL] = (AliAltroMapping*)maps->At(iRCU);
+       mapping[iDDL] = (AliAltroMapping*)maps->At(iDDL);
        buffer[iDDL]  = new AliAltroBuffer(fileName.Data(),mapping[iDDL]);
        buffer[iDDL]->WriteDataHeader(kTRUE, kFALSE);  //Dummy;
       }
@@ -546,25 +513,29 @@ void AliPHOS::Digits2Raw()
     // calculate the time response function
     } else {
       Double_t energy = 0 ;
-      module = relId[0];
       if (digit->GetId() <= geom->GetNModules() * geom->GetNCristalsInModule()) {
        energy=digit->GetEnergy();
-       if(energy>eMax) {eMax=energy; modMax=module; colMax=col; rowMax=row;}
+       if(energy>eMax) {eMax=energy; modMax=relId[0]; colMax=col; rowMax=row;}
       }
       else {
        energy = 0; // CPV raw data format is now know yet
       }
       pulse->SetAmplitude(energy);
       pulse->SetTZero(digit->GetTimeR());
+      Double_t r =fgCalibData->GetHighLowRatioEmc(relId[0],relId[3],relId[2]) ;
+      pulse->SetHG2LGRatio(r) ;
       pulse->MakeSamples();
       pulse->GetSamples(adcValuesHigh, adcValuesLow) ; 
+
       buffer[iDDL]->WriteChannel(relId[3]-1, relId[2]-1, 0, 
                           pulse->GetRawFormatTimeBins(), adcValuesLow , kAdcThreshold);
       buffer[iDDL]->WriteChannel(relId[3]-1, relId[2]-1, 1, 
                           pulse->GetRawFormatTimeBins(), adcValuesHigh, kAdcThreshold);
     }
   }
-  
+  delete [] adcValuesLow;
+  delete [] adcValuesHigh;
+
   // write real header and close last file
   for (Int_t iDDL=0; iDDL<maxDDL; iDDL++) {
     if (buffer[iDDL]) {
@@ -590,7 +561,7 @@ void AliPHOS::Hits2SDigits()
   AliPHOSSDigitizer phosDigitizer(fLoader->GetRunLoader()->GetFileName().Data()) ;
   phosDigitizer.SetEventRange(0, -1) ; // do all the events
  
-  phosDigitizer.ExecuteTask("all") ; 
+  phosDigitizer.Digitize("all") ; 
 }
 
 
@@ -609,7 +580,7 @@ void AliPHOS::SetTreeAddress()
   // Links Hits in the Tree to Hits array
   TBranch *branch;
   char branchname[20];
-  sprintf(branchname,"%s",GetName());
+  snprintf(branchname,20,"%s",GetName());
   // Branch address for hit tree
     TTree *treeH = fLoader->TreeH();
   if (treeH) {
@@ -623,36 +594,52 @@ void AliPHOS::SetTreeAddress()
   }
 }
 
-//____________________________________________________________________________
-Bool_t AliPHOS::Raw2SDigits(AliRawReader* rawReader)
-{
+//____________________________________________________________________________          
+Bool_t AliPHOS::Raw2SDigits(AliRawReader* rawReader)    
+{       
+                
+  AliPHOSLoader * loader = static_cast<AliPHOSLoader*>(fLoader) ;       
+                
+  TTree * tree = 0 ;    
+  tree = loader->TreeS() ;      
+  if ( !tree ) {        
+    loader->MakeTree("S");      
+    tree = loader->TreeS() ;    
+  }     
+                
+  TClonesArray * sdigits = loader->SDigits() ;          
+  if(!sdigits) {        
+    loader->MakeSDigitsArray();         
+    sdigits = loader->SDigits();        
+  }     
+  sdigits->Clear();     
+                
+  rawReader->Reset() ;
 
-  AliPHOSLoader * loader = dynamic_cast<AliPHOSLoader*>(fLoader) ;
-
-  TTree * tree = 0 ;
-  tree = loader->TreeS() ;
-  if ( !tree ) {
-    loader->MakeTree("S");
-    tree = loader->TreeS() ;
-  }
+  const TObjArray* maps = AliPHOSRecoParam::GetMappings();
+  if(!maps) AliFatal("Cannot retrieve ALTRO mappings!!");
 
-  TClonesArray * sdigits = loader->SDigits() ;
-  if(!sdigits) { 
-    loader->MakeSDigitsArray();
-    sdigits = loader->SDigits();
+  AliAltroMapping *mapping[20];
+  for(Int_t i = 0; i < 20; i++) {
+    mapping[i] = (AliAltroMapping*)maps->At(i);
   }
-  sdigits->Clear();
 
-  AliPHOSRawDecoder dc(rawReader);
-  AliPHOSRawDigiProducer pr;
-  pr.MakeDigits(sdigits,&dc);
-
-  Int_t bufferSize = 32000 ;
-  // TBranch * sdigitsBranch = tree->Branch("PHOS",&sdigits,bufferSize);
-  tree->Branch("PHOS",&sdigits,bufferSize);
-  tree->Fill();
-
-  fLoader->WriteSDigits("OVERWRITE");
-  return kTRUE;
-    
+  AliPHOSRawFitterv0 fitter;
+
+  fitter.SubtractPedestals(AliPHOSSimParam::GetInstance()->EMCSubtractPedestals());
+  fitter.SetAmpOffset(AliPHOSSimParam::GetInstance()->GetGlobalAltroOffset());
+  fitter.SetAmpThreshold(AliPHOSSimParam::GetInstance()->GetGlobalAltroThreshold());
+
+  AliPHOSRawDigiProducer pr(rawReader,mapping);
+  pr.SetSampleQualityCut(AliPHOSSimParam::GetInstance()->GetEMCSampleQualityCut());     
+  pr.MakeDigits(sdigits,&fitter);
+                
+  Int_t bufferSize = 32000 ;    
+  // TBranch * sdigitsBranch = tree->Branch("PHOS",&sdigits,bufferSize);        
+  tree->Branch("PHOS",&sdigits,bufferSize);     
+  tree->Fill();         
+                
+  fLoader->WriteSDigits("OVERWRITE");   
+  return kTRUE;         
+  
 }