From Andrea Dainese: protection in AliAODTrack::PropagateToDCA for usage in case...

[u/mrichter/AliRoot.git] / EMCAL / AliEMCALReconstructor.cxx

diff --git a/EMCAL/AliEMCALReconstructor.cxx b/EMCAL/AliEMCALReconstructor.cxx
index 0a1e6660a49cd3b3ab08bbc6db47ae316765cbf7..55eb4a40a211157f4bdb79fb6fa5b0ac42f9274a 100644 (file)
--- a/EMCAL/AliEMCALReconstructor.cxx
+++ b/EMCAL/AliEMCALReconstructor.cxx
@@ -135,8 +135,8 @@ AliEMCALReconstructor::AliEMCALReconstructor()
   //Init temporary list of digits
   fgDigitsArr     = new TClonesArray("AliEMCALDigit",1000);
   fgClustersArr   = new TObjArray(1000);
-  fgTriggerDigits = new TClonesArray("AliEMCALTriggerRawDigit",1000);  
-
+  fgTriggerDigits = new TClonesArray("AliEMCALTriggerRawDigit", 32 * 96);      
+       
   //Track matching
   fMatches = new TList();
   fMatches->SetOwner(kTRUE);
@@ -354,18 +354,15 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
   
   if (esdV0) 
     {
-      for (Int_t i = 0; i < 32; i++)
-       {
-         v0M[0] += (Int_t)esdV0->GetAdcV0C(i);
-         v0M[1] += (Int_t)esdV0->GetAdcV0A(i);
-       }
+               v0M[0] = esdV0->GetTriggerChargeC();
+               v0M[1] = esdV0->GetTriggerChargeA();
     }
   else
     {
-      AliWarning("Cannot retrieve V0 ESD! Run w/ null V0 charges");
+      AliWarning("No V0 ESD! Run trigger processor w/ null V0 charges");
     }
   
-  if (fgTriggerDigits) fgTriggerDigits->Clear();
+  if (fgTriggerDigits && fgTriggerDigits->GetEntriesFast()) fgTriggerDigits->Delete();
   
   TBranch *branchtrg = digitsTree->GetBranch("EMTRG");
   
@@ -379,7 +376,7 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
   branchtrg->GetEntry(0);
   
   // Note: fgTriggerProcessor reset done at the end of this method
-//   fgTriggerProcessor->Digits2Trigger(fgTriggerDigits, v0M, fTriggerData);
+  fgTriggerProcessor->Digits2Trigger(fgTriggerDigits, v0M, fTriggerData);
   
   // Fill ESD
   AliESDCaloTrigger* trgESD = esd->GetCaloTrigger("EMCAL");
@@ -391,7 +388,8 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
       for (Int_t i = 0; i < fgTriggerDigits->GetEntriesFast(); i++)
        {         
          AliEMCALTriggerRawDigit* rdig = (AliEMCALTriggerRawDigit*)fgTriggerDigits->At(i);
-         
+         if (AliDebugLevel() > 999) rdig->Print("");
+               
          Int_t px, py;
          if (fGeom->GetPositionInEMCALFromAbsFastORIndex(rdig->GetId(), px, py))
            {
@@ -399,14 +397,15 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
              
              rdig->GetMaximum(a, t);
              rdig->GetL0Times(times);
-             
+                       
              trgESD->Add(px, py, a, t, times, rdig->GetNL0Times(), rdig->GetL1TimeSum(), rdig->GetTriggerBits());
            }
        }
       
-      trgESD->SetL1Threshold(0, fTriggerData->GetL1GammaThreshold());
-      
-      trgESD->SetL1Threshold(1, fTriggerData->GetL1JetThreshold()  );
+               for (int i = 0; i < 2; i++) {
+                       trgESD->SetL1Threshold(2 * i    , fTriggerData->GetL1JetThreshold(  i));
+                       trgESD->SetL1Threshold(2 * i + 1, fTriggerData->GetL1GammaThreshold(i));
+               }
       
       Int_t v0[2];
       fTriggerData->GetL1V0(v0);
@@ -416,7 +415,7 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
       
       if (!saveOnce && fTriggerData->GetL1DataDecoded()) 
        {
-         int type[8] = {0};
+         int type[15] = {0};
          fTriggerData->GetL1TriggerType(type);
          
          esd->SetCaloTriggerType(type);
@@ -437,23 +436,51 @@ void AliEMCALReconstructor::FillESD(TTree* digitsTree, TTree* clustersTree,
   
   Int_t nDigits = fgDigitsArr->GetEntries(), idignew = 0 ;
   AliDebug(1,Form("%d digits",nDigits));
+  
   AliESDCaloCells &emcCells = *(esd->GetEMCALCells());
   emcCells.CreateContainer(nDigits);
   emcCells.SetType(AliVCaloCells::kEMCALCell);
+  
   Float_t energy = 0;
   Float_t time   = 0;
-  for (Int_t idig = 0 ; idig < nDigits ; idig++) {
+  for (Int_t idig = 0 ; idig < nDigits ; idig++) 
+  {
     const AliEMCALDigit * dig = (const AliEMCALDigit*)fgDigitsArr->At(idig);
     time   = dig->GetTime();      // Time already calibrated in clusterizer
     energy = dig->GetAmplitude(); // energy calibrated in clusterizer
-    if(energy > 0 ){
+    
+    if(energy > 0 )
+    {
       fgClusterizer->Calibrate(energy,time,dig->GetId()); //Digits already calibrated in clusterizers
-      if(energy > 0){ //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them
-        emcCells.SetCell(idignew,dig->GetId(),energy, time);   
+      
+      if(energy > 0) //Digits tagged as bad (dead, hot, not alive) are set to 0 in calibrate, remove them
+      { 
+        // Only for MC
+        // Get the label of the primary particle that generated the cell
+        // Assign the particle that deposited more energy
+        Int_t   nprimaries = dig->GetNprimary() ;
+        Int_t   digLabel   =-1 ;
+        Float_t edep       =-1.;
+        if ( nprimaries > 0 )
+        {
+          Int_t jndex ;
+          for ( jndex = 0 ; jndex < nprimaries ; jndex++ ) { // all primaries in digit
+ 
+            if(edep < dig->GetDEParent(jndex+1))
+            {
+              digLabel = dig->GetIparent (jndex+1);
+              edep     = dig->GetDEParent(jndex+1);
+            }
+                   
+          } // all primaries in digit      
+        } // select primary label
+        
+        emcCells.SetCell(idignew,dig->GetId(),energy, time,digLabel);   
         idignew++;
       }
     }
   }
+  
   emcCells.SetNumberOfCells(idignew);
   emcCells.Sort();
   
@@ -647,7 +674,8 @@ void AliEMCALReconstructor::FillMisalMatrixes(AliESDEvent* esd)const{
   TGeoHMatrix * m = 0x0;
   for(Int_t sm = 0; sm < fGeom->GetNumberOfSuperModules(); sm++){
     snprintf(path,bufsize,"/ALIC_1/XEN1_1/SMOD_%d",sm+1) ; //In Geometry modules numbered 1,2,.,5
-    if(sm >= 10) snprintf(path,bufsize,"/ALIC_1/XEN1_1/SM10_%d",sm-10+1) ;
+    if(sm >= 10 && !((fGeom->GetEMCGeometry()->GetGeoName()).Contains("12SMV1"))) snprintf(path,bufsize,"/ALIC_1/XEN1_1/SM10_%d",sm-10+1) ;
+    if(sm >= 10 &&  ((fGeom->GetEMCGeometry()->GetGeoName()).Contains("12SMV1"))) snprintf(path,bufsize,"/ALIC_1/XEN1_1/SM3rd_%d",sm-10+1) ;
     
     if (gGeoManager->CheckPath(path)){
       gGeoManager->cd(path);
@@ -701,18 +729,28 @@ Bool_t AliEMCALReconstructor::CalculateResidual(AliESDtrack *track, AliESDCaloCl
   // Otherwise use the TPCInner point
 
   dEta = -999, dPhi = -999;
+  Bool_t ITSTrackSA = 0;
 
   AliExternalTrackParam *trkParam = 0;
+  
   const AliESDfriendTrack*  friendTrack = track->GetFriendTrack();
   if(friendTrack && friendTrack->GetTPCOut())
     trkParam = const_cast<AliExternalTrackParam*>(friendTrack->GetTPCOut());
-  else
+  else if(track->GetInnerParam())
     trkParam = const_cast<AliExternalTrackParam*>(track->GetInnerParam());
+  else{
+    trkParam = new AliExternalTrackParam(*track); //If there is ITSSa track 
+    ITSTrackSA = 1;    
+  }
   if(!trkParam) return kFALSE;
-
+  
   AliExternalTrackParam trkParamTmp (*trkParam);
-  if(!AliEMCALRecoUtils::ExtrapolateTrackToCluster(&trkParamTmp, cluster, track->GetMass(), GetRecParam()->GetExtrapolateStep(), dEta, dPhi)) return kFALSE;
+  if(!AliEMCALRecoUtils::ExtrapolateTrackToCluster(&trkParamTmp, cluster, track->GetMass(kTRUE), GetRecParam()->GetExtrapolateStep(), dEta, dPhi)){
+       if(ITSTrackSA) delete trkParam;
+       return kFALSE;
+  }
 
+  if(ITSTrackSA) delete trkParam;
   return kTRUE;
 }