]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - HMPID/AliHMPIDReconstructor.cxx
Bug in area evaluation for reconstruction of ring. Important for PbPb events (high...
[u/mrichter/AliRoot.git] / HMPID / AliHMPIDReconstructor.cxx
index 0c1c106e275bd7dfa1a590a085d3350998d5df18..4a4f076b1bb8502718e7e2e6a23211f728b51549 100644 (file)
 #include "AliHMPID.h"              //Reconstruct() 
 #include "AliHMPIDCluster.h"       //Dig2Clu()
 #include "AliHMPIDParam.h"         //FillEsd() 
-#include <AliESD.h>               //FillEsd()
-#include <AliRawReader.h>         //Reconstruct() for raw digits
+#include <AliCDBEntry.h>           //ctor
+#include <AliCDBManager.h>         //ctor
+#include <AliESDEvent.h>           //FillEsd()
+#include <AliRawReader.h>          //Reconstruct() for raw digits
+#include "AliHMPIDRawStream.h"     //ConvertDigits()
 ClassImp(AliHMPIDReconstructor)
 
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-AliHMPIDReconstructor::AliHMPIDReconstructor():AliReconstructor()
+AliHMPIDReconstructor::AliHMPIDReconstructor():AliReconstructor(),fUserCut(0),fDaqSig(0),fDig(0),fClu(0)
 {
 //
 //ctor
 //
-  fClu=new TObjArray(AliHMPIDDigit::kMaxCh+1);
-  fDig=new TObjArray(AliHMPIDDigit::kMaxCh+1);
-  fClu->SetOwner(kTRUE);
-  for(int i=AliHMPIDDigit::kMinCh;i<=AliHMPIDDigit::kMaxCh;i++){ 
+  AliHMPIDParam::Instance();                                                        //geometry loaded for reconstruction
+  fUserCut = new Int_t[7];
+  fClu=new TObjArray(AliHMPIDParam::kMaxCh+1); fClu->SetOwner(kTRUE);
+  fDig=new TObjArray(AliHMPIDParam::kMaxCh+1); fDig->SetOwner(kTRUE);
+  
+  for(int i=AliHMPIDParam::kMinCh;i<=AliHMPIDParam::kMaxCh;i++){ 
     fDig->AddAt(new TClonesArray("AliHMPIDDigit"),i);
-    fClu->AddAt(new TClonesArray("AliHMPIDCluster"),i);
+    TClonesArray *pClus = new TClonesArray("AliHMPIDCluster");
+    pClus->SetUniqueID(i);
+    fClu->AddAt(pClus,i);
+  }
+  
+  AliCDBEntry *pUserCutEnt =AliCDBManager::Instance()->Get("HMPID/Calib/UserCut");    //contains TObjArray of 14 TObject with n. of sigmas to cut charge 
+  if(pUserCutEnt) {
+    TObjArray *pUserCut = (TObjArray*)pUserCutEnt->GetObject(); 
+    for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++){                  //chambers loop 
+      fUserCut[iCh] = pUserCut->At(iCh)->GetUniqueID();
+      Printf("HMPID: UserCut successfully loaded for chamber %i -> %i ",iCh,fUserCut[iCh]);
+    }
+  }
+
+  
+  AliCDBEntry *pDaqSigEnt =AliCDBManager::Instance()->Get("HMPID/Calib/DaqSig");  //contains TObjArray of TObjArray 14 TMatrixF sigmas values for pads 
+  if(!pDaqSigEnt) AliFatal("No pedestals from DAQ!");
+  fDaqSig = (TObjArray*)pDaqSigEnt->GetObject();
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++){                  //chambers loop 
+    Printf(" HMPID: DaqSigCut successfully loaded for chamber %i -> %i ",iCh,(Int_t)fDaqSig->At(iCh)->GetUniqueID());
   }
-  fDig->SetOwner(kTRUE);
 }//AliHMPIDReconstructor
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDReconstructor::Dig2Clu(TObjArray *pDigAll,TObjArray *pCluAll,Bool_t isTryUnfold)
@@ -49,9 +72,9 @@ void AliHMPIDReconstructor::Dig2Clu(TObjArray *pDigAll,TObjArray *pCluAll,Bool_t
 //            pCluAll     - list of clusters for all chambers
 //            isTryUnfold - flag to choose between CoG and Mathieson fitting  
 //  Returns: none    
-  TMatrixF padMap(AliHMPIDDigit::kMinPx,AliHMPIDDigit::kMaxPcx,AliHMPIDDigit::kMinPy,AliHMPIDDigit::kMaxPcy);  //pads map for single chamber 0..159 x 0..143 
+  TMatrixF padMap(AliHMPIDParam::kMinPx,AliHMPIDParam::kMaxPcx,AliHMPIDParam::kMinPy,AliHMPIDParam::kMaxPcy);  //pads map for single chamber 0..159 x 0..143 
   
-  for(Int_t iCh=AliHMPIDDigit::kMinCh;iCh<=AliHMPIDDigit::kMaxCh;iCh++){                  //chambers loop 
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++){                  //chambers loop 
     TClonesArray *pDigCur=(TClonesArray*)pDigAll->At(iCh);                                //get list of digits for current chamber
     if(pDigCur->GetEntriesFast()==0) continue;                                            //no digits for current chamber
   
@@ -64,7 +87,6 @@ void AliHMPIDReconstructor::Dig2Clu(TObjArray *pDigAll,TObjArray *pCluAll,Bool_t
     }//digits loop to fill digits map 
     
     AliHMPIDCluster clu;                                                                  //tmp cluster to be used as current
-  
     for(Int_t iDig=0;iDig<pDigCur->GetEntriesFast();iDig++){                              //digits loop for current chamber
       AliHMPIDDigit *pDig=(AliHMPIDDigit*)pDigCur->At(iDig);                              //take current digit
       if(!(pDig=UseDig(pDig->PadChX(),pDig->PadChY(),pDigCur,&padMap))) continue;         //this digit is already taken in FormClu(), go after next digit
@@ -85,10 +107,10 @@ void  AliHMPIDReconstructor::FormClu(AliHMPIDCluster *pClu,AliHMPIDDigit *pDig,T
 
   Int_t cnt=0,cx[4],cy[4];
   
-  if(pDig->PadPcX() != AliHMPIDDigit::kMinPx){cx[cnt]=pDig->PadChX()-1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //left
-  if(pDig->PadPcX() != AliHMPIDDigit::kMaxPx){cx[cnt]=pDig->PadChX()+1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //right
-  if(pDig->PadPcY() != AliHMPIDDigit::kMinPy){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()-1;cnt++;}       //down
-  if(pDig->PadPcY() != AliHMPIDDigit::kMaxPy){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()+1;cnt++;}       //up
+  if(pDig->PadPcX() != AliHMPIDParam::kMinPx){cx[cnt]=pDig->PadChX()-1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //left
+  if(pDig->PadPcX() != AliHMPIDParam::kMaxPx){cx[cnt]=pDig->PadChX()+1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //right
+  if(pDig->PadPcY() != AliHMPIDParam::kMinPy){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()-1;cnt++;}       //down
+  if(pDig->PadPcY() != AliHMPIDParam::kMaxPy){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()+1;cnt++;}       //up
   
   for (Int_t i=0;i<cnt;i++){//neighbours loop
     if((pDig=UseDig(cx[i],cy[i],pDigLst,pDigMap))) FormClu(pClu,pDig,pDigLst,pDigMap);   //check if this neighbour pad fired and mark it as taken  
@@ -102,7 +124,7 @@ void AliHMPIDReconstructor::Reconstruct(TTree *pDigTree,TTree *pCluTree)const
 //           pCluTree - poitner to Cluster tree
 //  Returns: none    
   AliDebug(1,"Start.");
-  for(Int_t iCh=AliHMPIDDigit::kMinCh;iCh<=AliHMPIDDigit::kMaxCh;iCh++) {
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++) {
     pCluTree->Branch(Form("HMPID%d",iCh),&((*fClu)[iCh]),4000,0);
     pDigTree->SetBranchAddress(Form("HMPID%d",iCh),&((*fDig)[iCh]));
   }   
@@ -110,7 +132,7 @@ void AliHMPIDReconstructor::Reconstruct(TTree *pDigTree,TTree *pCluTree)const
   Dig2Clu(fDig,fClu);              //cluster finder 
   pCluTree->Fill();                //fill tree for current event
   
-  for(Int_t iCh=AliHMPIDDigit::kMinCh;iCh<=AliHMPIDDigit::kMaxCh;iCh++){
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++){
     fDig->At(iCh)->Clear();
     fClu->At(iCh)->Clear();
   }
@@ -118,45 +140,6 @@ void AliHMPIDReconstructor::Reconstruct(TTree *pDigTree,TTree *pCluTree)const
   AliDebug(1,"Stop.");      
 }//Reconstruct(for simulated digits)
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDReconstructor::Reconstruct(AliRunLoader *pAL,AliRawReader* pRR)const
-{
-//Invoked  by AliReconstruction to convert raw digits from DDL files to clusters
-//Arguments: pAL - ALICE run loader pointer
-//           pRR - ALICE raw reader pointer
-//  Returns: none
-  AliLoader *pRL=pAL->GetDetectorLoader("HMPID");  AliHMPID *pRich=(AliHMPID*)pAL->GetAliRun()->GetDetector("HMPID");//get pointers for HMPID and HMPID loader
-
-  AliHMPIDDigit dig; //tmp digit, raw digit will be converted to it
-
-  TObjArray digLst; Int_t iDigCnt[7]; for(Int_t i=0;i<7;i++){digLst.AddAt(new TClonesArray("AliHMPIDDigit"),i); iDigCnt[i]=0;} //tmp list of digits for allchambers
-
-  Int_t iEvtN=0;
-  while(pRR->NextEvent()){//events loop
-    pAL->GetEvent(iEvtN++);
-    pRL->MakeTree("R");  pRich->MakeBranch("R");
-
-    for(Int_t iCh=0;iCh<7;iCh++){//chambers loop
-      pRR->Select("HMPID",2*iCh,2*iCh+1);//select only DDL files for the current chamber
-      UInt_t w32=0;
-      while(pRR->ReadNextInt(w32)){//raw records loop (in selected DDL files)
-        UInt_t ddl=pRR->GetDDLID(); //returns 0,1,2 ... 13
-        dig.Raw(w32,ddl);
-        AliDebug(1,Form("Ch=%i DDL=%i raw=0x%x digit=(%3i,%3i,%3i,%3i) Q=%5.2f",iCh,ddl,w32,dig.Ch(),dig.Pc(),dig.PadPcX(),dig.PadPcY(),dig.Q()));
-        new((*((TClonesArray*)digLst.At(iCh)))[iDigCnt[iCh]++]) AliHMPIDDigit(dig); //add this digit to the tmp list
-      }//raw records loop
-      pRR->Reset();
-    }//chambers loop
-
-    Dig2Clu(&digLst,pRich->CluLst());//cluster finder for all chambers
-    for(Int_t i=0;i<7;i++){digLst.At(i)->Delete(); iDigCnt[i]=0;}                    //clean up list of digits for all chambers
-
-    pRL->TreeR()->Fill();            //fill tree for current event
-    pRL->WriteRecPoints("OVERWRITE");//write out clusters for current event
-    pRich->CluReset();
-  }//events loop
-  pRL->UnloadRecPoints();
-}//Reconstruct raw data
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDReconstructor::ConvertDigits(AliRawReader *pRR,TTree *pDigTree)const
 {
 //Invoked  by AliReconstruction to convert raw digits from DDL files to digits
@@ -164,26 +147,41 @@ void AliHMPIDReconstructor::ConvertDigits(AliRawReader *pRR,TTree *pDigTree)cons
 //           pDigTree - pointer to Digit tree
 //  Returns: none    
   AliDebug(1,"Start.");
-  AliHMPIDDigit dig; //tmp digit, raw digit will be converted to it
-  for(Int_t iCh=AliHMPIDDigit::kMinCh;iCh<=AliHMPIDDigit::kMaxCh;iCh++) {
-    pDigTree->Branch(Form("HMPID%d",iCh),&((*fDig)[iCh]),4000,0);
-    pRR->Select("HMPID",2*iCh,2*iCh+1);//select only DDL files for the current chamber      
+//  Int_t digcnt=0;
+  
+
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++) {
+    pDigTree->Branch(Form("HMPID%d",iCh),&((*fDig)[iCh]),4000,0); 
+    
     Int_t iDigCnt=0;
-    UInt_t w32=0;
-                while(pRR->ReadNextInt(w32)){//raw records loop (in selected DDL files)
-      UInt_t ddl=pRR->GetDDLID(); //returns 0,1,2 ... 13
-      dig.Raw(w32,ddl);  
-      AliDebug(1,Form("Ch=%i DDL=%i raw=0x%x digit=(%3i,%3i,%3i,%3i) Q=%5.2f",iCh,ddl,w32,dig.Ch(),dig.Pc(),dig.PadPcX(),dig.PadPcY(),dig.Q()));
-      new((*((TClonesArray*)fDig->At(iCh)))[iDigCnt++]) AliHMPIDDigit(dig); //add this digit to the tmp list
-    }//raw records loop
-    pRR->Reset();        
-  }//chambers loop
+    AliHMPIDRawStream stream(pRR);    
+    while(stream.Next())
+    {
+      
+      UInt_t ddl=stream.GetDDLNumber(); //returns 0,1,2 ... 13
+      if((UInt_t)(2*iCh)==ddl || (UInt_t)(2*iCh+1)==ddl) {
+       for(Int_t row = 1; row <= AliHMPIDRawStream::kNRows; row++){
+        for(Int_t dil = 1; dil <= AliHMPIDRawStream::kNDILOGICAdd; dil++){
+          for(Int_t pad = 0; pad < AliHMPIDRawStream::kNPadAdd; pad++){
+            if(stream.GetCharge(ddl,row,dil,pad) < 1) continue; 
+              AliHMPIDDigit dig(stream.GetPad(ddl,row,dil,pad),stream.GetCharge(ddl,row,dil,pad));
+              if(!IsDigSurvive(&dig)) continue; 
+              new((*((TClonesArray*)fDig->At(iCh)))[iDigCnt++]) AliHMPIDDigit(dig); //add this digit to the tmp list 
+            }//pad
+          }//dil
+        }//row
+      }//while
+    }
+    stream.Delete();
+  }
   pDigTree->Fill();
-  for(Int_t iCh=AliHMPIDDigit::kMinCh;iCh<=AliHMPIDDigit::kMaxCh;iCh++)fDig->At(iCh)->Clear();
+  
+  for(Int_t iCh=AliHMPIDParam::kMinCh;iCh<=AliHMPIDParam::kMaxCh;iCh++)fDig->At(iCh)->Clear();
+  
   AliDebug(1,"Stop.");
 }//Reconstruct digits from raw digits
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDReconstructor::FillESD(AliRunLoader *, AliESD *pESD) const
+void AliHMPIDReconstructor::FillESD(TTree */*digitsTree*/, TTree */*clustersTree*/, AliESDEvent *pESD) const
 {
 // Calculates probability to be a electron-muon-pion-kaon-proton
 // from the given Cerenkov angle and momentum assuming no initial particle composition
@@ -221,4 +219,3 @@ void AliHMPIDReconstructor::FillESD(AliRunLoader *, AliESD *pESD) const
     pTrk->SetHMPIDpid(pid);
   }//ESD tracks loop
 }//FillESD()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++