Bugs in clusterizing and recon fixed
authorkir <kir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 15 Dec 2006 16:49:34 +0000 (16:49 +0000)
committerkir <kir@f7af4fe6-9843-0410-8265-dc069ae4e863>
Fri, 15 Dec 2006 16:49:34 +0000 (16:49 +0000)
15 files changed:
HMPID/AliHMPID.cxx
HMPID/AliHMPID.h
HMPID/AliHMPIDCluster.cxx
HMPID/AliHMPIDCluster.h
HMPID/AliHMPIDDigit.cxx
HMPID/AliHMPIDDigit.h
HMPID/AliHMPIDDigitizer.cxx
HMPID/AliHMPIDHit.cxx
HMPID/AliHMPIDHit.h
HMPID/AliHMPIDRecon.cxx
HMPID/AliHMPIDReconstructor.cxx
HMPID/AliHMPIDReconstructor.h
HMPID/AliHMPIDSelector.C
HMPID/AliHMPIDTracker.cxx
HMPID/AliHMPIDTracker.h

index a11e9f0..11f25c0 100644 (file)
 
 #include "AliHMPID.h"
 #include "AliHMPIDHit.h"   //OccupancyPrint(), HitQa()
-#include "AliHMPIDDigit.h" //OccupancyPrint()
+#include "AliHMPIDDigit.h" //
 #include <TParticle.h>  //SummaryOfEvent(), HitQa()
 #include <TBenchmark.h>  //HitQA()
 #include <TPDGCode.h>    //HitQA()
-#include <AliStack.h>   //OccupancyPrint(), SummaryOfEvent(), HitQa()
+#include <AliStack.h>   //SummaryOfEvent(), HitQa()
 #include <AliRun.h>     //HitQa() 
 #include <AliMC.h>       //ctor
 #include <AliHeader.h>
-#include <AliGenEventHeader.h>
-#include <AliGenHijingEventHeader.h>
 #include <TH1F.h>        //HitQA()
 #include <AliLog.h>      //in many methods to print AliInfo 
 ClassImp(AliHMPID)    
@@ -141,62 +139,6 @@ void AliHMPID::DigPrint(Int_t iEvt)const
   Printf("totally %i Digits",DigLst()->GetEntries());
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPID::OccupancyPrint(Int_t iEvtNreq)
-{
-//prints occupancy for each chamber in a given event
-  Int_t iEvtNmin,iEvtNmax;
-  if(iEvtNreq==-1){
-    iEvtNmin=0;
-    iEvtNmax=gAlice->GetEventsPerRun();
-  } else { 
-    iEvtNmin=iEvtNreq;iEvtNmax=iEvtNreq+1;
-  }
-    
-  if(GetLoader()->GetRunLoader()->LoadHeader()) return;    
-  if(GetLoader()->GetRunLoader()->LoadKinematics()) return;    
-  
-//  Info("Occupancy","for event %i",iEvtN);
-  if(GetLoader()->LoadHits()) return;
-  if(GetLoader()->LoadDigits()) return;
-
-  
-  for(Int_t iEvtN=iEvtNmin;iEvtN<iEvtNmax;iEvtN++){    
-    Int_t nDigCh[7]={0,0,0,0,0,0,0};  
-    Int_t iChHits[7]={0,0,0,0,0,0,0};
-    Int_t nPrim[7]={0,0,0,0,0,0,0};
-    Int_t nSec[7]={0,0,0,0,0,0,0};
-    AliInfo(Form("events processed %i",iEvtN));
-    if(GetLoader()->GetRunLoader()->GetEvent(iEvtN)) return;    
-    AliStack *pStack = GetLoader()->GetRunLoader()->Stack();
-    for(Int_t iPrimN=0;iPrimN<GetLoader()->TreeH()->GetEntries();iPrimN++){//prims loop
-      GetLoader()->TreeH()->GetEntry(iPrimN);      
-      for(Int_t iHitN=0;iHitN<Hits()->GetEntries();iHitN++){
-        AliHMPIDHit *pHit = (AliHMPIDHit*)Hits()->At(iHitN);
-        if(pHit->E()>0){
-          iChHits[pHit->Ch()]++;
-          if(pStack->Particle(pHit->GetTrack())->Rho()<0.01) nPrim[pHit->Ch()]++;else nSec[pHit->Ch()]++;
-        }
-      }
-    }
-    
-    GetLoader()->TreeD()->GetEntry(0);
-    for(Int_t iCh=0;iCh<7;iCh++){
-      for(Int_t iDig=0;iDig<DigLst(iCh)->GetEntries();iDig++){
-        AliHMPIDDigit *pDig=(AliHMPIDDigit*)DigLst(iCh)->At(iDig);
-        nDigCh[pDig->Ch()]++;
-      }
-      Printf("Occupancy for chamber %i = %4.2f %% and charged prim tracks %i and sec. tracks %i with total %i",
-        iCh,Float_t(nDigCh[iCh])*100/AliHMPIDDigit::kPadAll,nPrim[iCh],nSec[iCh],iChHits[iCh]);
-    }      
-    
-    
-  }//events loop
-  GetLoader()->UnloadHits();
-  GetLoader()->UnloadDigits();
-  GetLoader()->GetRunLoader()->UnloadHeader();    
-  GetLoader()->GetRunLoader()->UnloadKinematics();    
-}
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPID::CluPrint(Int_t iEvtN)const
 {
 //prints a list of HMPID clusters  for a given event
@@ -213,28 +155,3 @@ void AliHMPID::CluPrint(Int_t iEvtN)const
   Printf("totally %i clusters for event %i",iCluCnt,iEvtN);
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPID::SummaryOfEvent(Int_t iEvtN) const
-{
-//prints a summary for a given event
-  AliInfo(Form("Summary of event %i",iEvtN));
-  GetLoader()->GetRunLoader()->GetEvent(iEvtN);    
-  if(GetLoader()->GetRunLoader()->LoadHeader()) return;
-  if(GetLoader()->GetRunLoader()->LoadKinematics()) return;
-  AliStack *pStack=GetLoader()->GetRunLoader()->Stack();
-  
-  AliGenEventHeader* pGenHeader =  gAlice->GetHeader()->GenEventHeader();
-  if(pGenHeader->InheritsFrom("AliGenHijingEventHeader")) {
-    AliInfo(Form(" Hijing event with impact parameter b = %.2f (fm)",((AliGenHijingEventHeader*) pGenHeader)->ImpactParameter()));
-  }
-  Int_t nChargedPrimaries=0;
-  for(Int_t i=0;i<pStack->GetNtrack();i++) {
-    TParticle *pParticle = pStack->Particle(i);
-    if(pParticle->IsPrimary()&&pParticle->GetPDG()->Charge()!=0) nChargedPrimaries++;
-    }
-  AliInfo(Form("Total number of         primaries %i",pStack->GetNprimary()));
-  AliInfo(Form("Total number of charged primaries %i",nChargedPrimaries));
-  AliInfo(Form("Total n. of tracks in stack(+sec) %i",pStack->GetNtrack()));
-  GetLoader()->GetRunLoader()->UnloadHeader();
-  GetLoader()->GetRunLoader()->UnloadKinematics();
-}
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index b7f3aef..7323f17 100644 (file)
@@ -5,7 +5,7 @@
 
 #include <AliDetector.h>  //base class
 #include <TClonesArray.h> //XxxCreate() 
-#include <TObjArray.h>    //fClu field
+#include <TObjArray.h>    //fDig,fClu field
 
 
 
@@ -40,14 +40,11 @@ public:
   void          DigReset (         )     {if(fDig)for(int i=0;i<7;i++)fDig->At(i)->Clear();                     }//clean digits list 
   void          DigPrint (Int_t evt)const;                                                                       //print digits
           
+  TObjArray*    CluLst   (         )const{return fClu;                                                          }//get clusters list for all chambers
   TClonesArray* CluLst   (Int_t c  )const{return fClu ? (TClonesArray *)fClu->At(c):0;                          }//get clusters list for chamber
   inline void   CluCreate(         )     {fClu=new TObjArray(7); for(Int_t i=0;i<7;i++)fClu->AddAt(new TClonesArray("AliHMPIDCluster"),i);}//create clusters list
          void   CluReset (         )     {if(fClu)for(int i=0;i<7;i++)fClu->At(i)->Clear();                     }//clean clusters list
          void   CluPrint (Int_t evt)const;                                                                       //print clusters list
-         
-  void          OccupancyPrint(Int_t evt=-1);                    //print chambers occupancy 
-  void          SummaryOfEvent(Int_t evt=0)const;
-
 protected:  
   TClonesArray         *fSdi;                     //! list of sdigits  
   TObjArray            *fDig;                     //! each chamber holds it's one list of digits
index 9cc00b3..dacba95 100644 (file)
@@ -43,7 +43,7 @@ void AliHMPIDCluster::CorrSin()
 // Correction of cluster x position due to sinoid, see HMPID TDR  page 30
 // Arguments: none
 //   Returns: none
-  AliHMPIDDigit dig(Ch(),100,1,fX,fY);                                               //tmp digit to get it center
+  AliHMPIDDigit dig;dig.Manual1(Ch(),fX,fY);                                               //tmp digit to get it center
   Float_t x=fX-dig.LorsX();  
   fX+=3.31267e-2*TMath::Sin(2*TMath::Pi()/0.8*x)-2.66575e-3*TMath::Sin(4*TMath::Pi()/0.8*x)+2.80553e-3*TMath::Sin(6*TMath::Pi()/0.8*x)+0.0070;
 }
@@ -82,8 +82,8 @@ void AliHMPIDCluster::Print(Option_t* opt)const
     case      kCoG: status="coged"      ;break;
     case      kEmp: status="empty"      ;break;
   }
-  Printf("%s cs=%2i, Size=%2i (x=%7.3f cm,y=%7.3f cm,Q=%4i qdc), %s",
-         opt,Ch(),Size(),X(),Y(),Q(),status);
+  Printf("%s ch=%i, Size=%2i        (%7.3f,%7.3f) Q=%4i          %s",
+         opt,Ch(),Size(),            X(),  Y(),   Q(),status);
   for(Int_t i=0;i<Size();i++) Dig(i)->Print();    
 }//Print()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -135,12 +135,12 @@ Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
        pMinuit->mnpout(3*i+1 ,sName,  fitY, d1 , d2, d3, iErrFlg);
        pMinuit->mnpout(3*i+2 ,sName,  fitQ, d1 , d2, d3, iErrFlg);
        if (TMath::Abs(fitQ)>2147483647.0) fitQ = TMath::Sign((Double_t)2147483647,fitQ);//???????????????
-       new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(Ch(),fitX,fitY,(Int_t)fitQ);            //add new unfolded clusters
+       new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(Ch(),fitX,fitY,(Int_t)fitQ,kUnf);       //add new unfolded clusters
       }//local maxima loop
     }
   }else{//do not unfold since number of loc max is unresonably high or user's baned unfolding 
     CoG();
-    new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);  //add this raw cluster 
+    new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(Ch(),X(),Y(),Q(),kCoG);  //add this raw cluster 
   }
   delete pMinuit;
   return iLocMaxCnt;
index 77e9f95..60fd0ac 100644 (file)
@@ -11,9 +11,9 @@ class AliHMPIDCluster :public TObject
 {
 public:
   enum EClusterStatus {kFor,kCoG,kUnf,kEmp=-1};      //status flags    
-  AliHMPIDCluster(                                   ):TObject( ),fSt(kEmp ),fCh(-1   ),fQ(-1  ),fX(-1  ),fY(-1  ),fDigs(0                                  ) {} 
-  AliHMPIDCluster(Int_t c,Float_t x,Float_t y,Int_t q):TObject( ),fSt(kUnf ),fCh(c    ),fQ(q   ),fX(x   ),fY(y   ),fDigs(0                                  ) {} 
-  AliHMPIDCluster(const AliHMPIDCluster &c            ):TObject(c),fSt(c.fSt),fCh(c.fCh),fQ(c.fQ),fX(c.fX),fY(c.fY),fDigs(c.fDigs ? new TObjArray(*c.fDigs):0) {}
+  AliHMPIDCluster(                                           ):TObject( ),fSt(kEmp ),fCh(-1   ),fQ(-1  ),fX(-1  ),fY(-1  ),fDigs(0                                  ) {} 
+  AliHMPIDCluster(Int_t c,Float_t x,Float_t y,Int_t q,Int_t s):TObject( ),fSt(s    ),fCh(c    ),fQ(q   ),fX(x   ),fY(y   ),fDigs(0                                  ) {} 
+  AliHMPIDCluster(const AliHMPIDCluster &c                   ):TObject(c),fSt(c.fSt),fCh(c.fCh),fQ(c.fQ),fX(c.fX),fY(c.fY),fDigs(c.fDigs ? new TObjArray(*c.fDigs):0) {}
   AliHMPIDCluster &operator=(const AliHMPIDCluster &c) {
     if(this == &c)return *this;TObject::operator=(c);            fSt=c.fSt; fCh=c.fCh; fQ=c.fQ; fX=c.fX; fY=c.fY; fDigs=c.fDigs ? new TObjArray(*c.fDigs):0; return *this;}
     
index 171336b..18ddfe7 100644 (file)
@@ -72,28 +72,19 @@ HMPID raw word is 32 bits with the structure:
  5 bits zero  5 bits row number (1..24)  4 bits DILOGIC chip number (1..10) 6 bits DILOGIC address (0..47)  12 bits QDC value (0..4095)
 */
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Float_t AliHMPIDDigit::Hit2Sdi(AliHMPIDHit *pHit,TClonesArray *pSdiLst)
+void AliHMPIDDigit::Hit2Sdi(AliHMPIDHit *pHit,TClonesArray *pSdiLst)
 {
 // Creates a list of sdigits out of provided hit
 // Arguments: pHit- hit
 //   Returns: none
   
-  Float_t  x=(pHit->LorsX() > SizePcX())? pHit->LorsX()-SizePcX()-SizeDead():pHit->LorsX();          //sagita is for PC (0-64) and not for chamber   
-  Float_t  qdcEle=34.06311+0.2337070*x+5.807476e-3*x*x-2.956471e-04*x*x*x+2.310001e-06*x*x*x*x;      //reparametrised from DiMauro
-  
-  Int_t iNele=Int_t(pHit->E()/26e-9);  if(iNele<1) iNele=1;                                          //number of electrons created by hit
-  Float_t qdcTot=0; for(Int_t i=1;i<=iNele;i++) qdcTot-=qdcEle*TMath::Log(gRandom->Rndm()+1e-6);     //total qdc fro hit, 1e-6 is a protection against 0 from rndm
-  
-  AliHMPIDDigit dd(1,1,1,pHit->LorsX(),pHit->LorsY());                                                //tmp digit to shift hit y to the nearest anod wire  
-  Float_t y=  (pHit->LorsY() > dd.LorsY()) ? dd.LorsY()+0.21 :  dd.LorsY()-0.21;
-  
-  Int_t iSdiCnt=pSdiLst->GetEntries();
+  Int_t iSdiCnt=pSdiLst->GetEntries(); //list of sdigits contains sdigits from previous ivocations of Hit2Sdi, do not override them
+  AliHMPIDDigit dig;
   for(Int_t i=0;i<9;i++){                                      //affected pads loop
-    AliHMPIDDigit dig(pHit->Ch(),qdcTot,pHit->Tid(),pHit->LorsX(),y,i); //c,q,tid,x,y   create tmp sdigit for pad i around hit position
+    dig.Set(pHit,i); //c,q,tid,x,y   create tmp sdigit for pad i around hit position
     if(dig.PadPcX()==-1) continue;
     new((*pSdiLst)[iSdiCnt++]) AliHMPIDDigit(dig);
   }
-  return qdcTot;
 }//Hit2Sdi()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Print(Option_t*)const
@@ -102,8 +93,8 @@ void AliHMPIDDigit::Print(Option_t*)const
 // Arguments: option string not used
 //   Returns: none    
   UInt_t w32; Raw(w32);
-  Printf("(ch=%1i,pc=%1i,x=%2i,y=%2i),pos=(%7.2f,%7.2f) Q=%8.3f TID=(%5i,%5i,%5i) ddl=%i raw=0x%x (r=%2i,d=%2i,a=%2i)",
-               A2C(fPad),A2P(fPad),A2X(fPad),A2Y(fPad),LorsX(),LorsY(), Q(),  fTracks[0],fTracks[1],fTracks[2],DdlIdx(),w32,Row(),Dilogic(),Addr());
+  Printf("(ch=%1i,pc=%1i,x=%2i,y=%2i) (%7.3f,%7.3f) Q=%8.3f TID=(%5i,%5i,%5i) ddl=%i raw=0x%x (r=%2i,d=%2i,a=%2i) %s",
+               A2C(fPad),A2P(fPad),A2X(fPad),A2Y(fPad),LorsX(),LorsY(), Q(),  fTracks[0],fTracks[1],fTracks[2],DdlIdx(),w32,Row(),Dilogic(),Addr(), (IsOverTh(Q()))?"":"below thr");
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::PrintSize()
@@ -140,7 +131,7 @@ void AliHMPIDDigit::DrawZoom()
   pBar->SetParts(5);
   Float_t x=gPad->AbsPixeltoX(gPad->GetEventX());
   Float_t y=gPad->AbsPixeltoY(gPad->GetEventY());
-  AliHMPIDDigit dig(1,0,1,x,y); UInt_t w32=0; 
+  AliHMPIDDigit dig;dig.Manual1(1,x,y); UInt_t w32=0; 
   if(IsInDead(x,y))
     pBar->SetText("Out of sensitive area",4);    
   else{
@@ -196,11 +187,11 @@ void AliHMPIDDigit::DrawPc(Bool_t isFill)
     if(iPc==0) pc=new TPolyLine(5,xL,yD); if(iPc==1) pc=new TPolyLine(5,xR,yD);
     (iPc%2)? pc->SetFillColor(iColLeft): pc->SetFillColor(iColRight);
     if(isFill) pc->Draw("f"); else pc->Draw();
-    if(iPc%2) {dig.Set(0,iPc,79,25); txt.DrawText(dig.LorsX()+2,dig.LorsY(),Form("PC%i",dig.Pc()));}//print PC#    
+    if(iPc%2) {dig.Manual2(0,iPc,79,25); txt.DrawText(dig.LorsX()+2,dig.LorsY(),Form("PC%i",dig.Pc()));}//print PC#    
     
     txt.SetTextAlign(32);
     for(Int_t iRow=0;iRow<8 ;iRow++){//draw row lines (horizontal)
-      dig.Set(0,iPc,0,iRow*6);   //set digit to the left-down pad of this row
+      dig.Manual2(0,iPc,0,iRow*6);   //set digit to the left-down pad of this row
       if(iPc%2) txt.DrawText(dig.LorsX()-1           ,dig.LorsY(),Form("%i",dig.PadPcY())); //print PadY#    
                 txt.DrawText(dig.LorsX()-1+(iPc%2)*67,dig.LorsY()+2,Form("r%i",dig.Row())); //print Row#    
       pL=new TLine(dig.LorsX()-0.5*SizePadX(),dig.LorsY()-0.5*SizePadY(),dig.LorsX()+SizePcX()-0.5*SizePadX(),dig.LorsY()-0.5*SizePadY()); 
@@ -209,7 +200,7 @@ void AliHMPIDDigit::DrawPc(Bool_t isFill)
     
     txt.SetTextAlign(13);
     for(Int_t iDil=0;iDil<10;iDil++){//draw dilogic lines (vertical)
-      dig.Set(0,iPc,iDil*8,0);       //set this digit to the left-down pad of this dilogic        
+      dig.Manual2(0,iPc,iDil*8,0);       //set this digit to the left-down pad of this dilogic        
                            txt.DrawText(dig.LorsX()  ,dig.LorsY()-1,Form("%i",dig.PadPcX()));   //print PadX# 
       if(iPc==4 || iPc==5) txt.DrawText(dig.LorsX()+2,dig.LorsY()+42,Form("d%i",dig.Dilogic())); //print Dilogic#    
       pL=new TLine(dig.LorsX()-0.5*SizePadX(),dig.LorsY()-0.5*SizePadY(),dig.LorsX()-0.5*SizePadX(),dig.LorsY()+SizePcY()-0.5*SizePadY()); 
@@ -226,7 +217,7 @@ void AliHMPIDDigit::Test()
     Int_t pc=Int_t(gRandom->Rndm()*6);
     Int_t px=Int_t(gRandom->Rndm()*80);
     Int_t py=Int_t(gRandom->Rndm()*48);
-    d1.Set(ch,pc,px,py);                
+    d1.Manual2(ch,pc,px,py);                
     ddl=d1.Raw(w32);    d2.Raw(ddl,w32);
     if(d1.Compare(&d2)) Printf("Problem!!!");
   }
index 514d345..3f55a75 100644 (file)
@@ -7,10 +7,10 @@
 #include <TMath.h>         //Mathieson()
 #include <TRandom.h>       //IsOverTh()  
 #include <AliBitPacking.h> //Raw()
+#include "AliHMPIDHit.h"   //Hit2Sdi(), ctor()
 
-class AliHMPIDHit;          //Hit2Sdi()
 class TClonesArray;        //Hit2Sdi()
+  
 class AliHMPIDDigit :public AliDigit //TObject-AliDigit-AliHMPIDDigit
 {
 public:
@@ -18,10 +18,9 @@ public:
   enum ERawData{kNddls=14};                                                      //RAW data structure
   enum EPadData{kPcX=2,kPcY=3,kPadPcX=80,kPadPcY=48,kPadAllX=kPadPcX*kPcX,kPadAllY=kPadPcY*kPcY,kPcAll=kPcX*kPcY,kPadAll=kPadAllX*kPadAllY};   //Segmentation structure 
 //ctor&dtor    
-           AliHMPIDDigit(                                                     ):AliDigit( ),fPad(Abs(-1,-1,-1,-1)),fQ(-1)  {}                //default ctor
-           AliHMPIDDigit(Int_t pad,Int_t q,Int_t *t                           ):AliDigit(t),fPad(pad             ),fQ(q )  {}                //digit ctor 
-  inline   AliHMPIDDigit(Int_t c,Float_t q,Int_t t,Float_t x,Float_t y,Int_t f=0);                                                               //sdigit ctor 
-  virtual ~AliHMPIDDigit()                                {} //dtor
+           AliHMPIDDigit(                                      ):AliDigit( ),fPad(Abs(-1,-1,-1,-1)),fQ(-1)  {}                                      //default ctor
+           AliHMPIDDigit(Int_t pad,Int_t q,Int_t *t            ):AliDigit(t),fPad(pad             ),fQ(q )  {}                                      //ctor used in digitizer
+  virtual ~AliHMPIDDigit(                                      )                                            {}                                      //dtor
 //framework part    
          Bool_t  IsSortable  (                               )const{return kTRUE;}                                                     //provision to use TObject::Sort() 
   inline Int_t   Compare     (const TObject *pObj            )const;                                                                   //provision to use TObject::Sort()
@@ -41,17 +40,19 @@ public:
          void    DrawZoom    (                               ); 
          Int_t   DdlIdx      (                               )const{return 2*Ch()+Pc()%2;                           }                  //DDL# 0..13
          Int_t   DdlId       (                               )const{return (6<<8)+DdlIdx();                         }                  //DDL ID 0x600..0x60d
-  static Float_t Hit2Sdi     (AliHMPIDHit *pHit,TClonesArray*);                                                                        //hit -> 9 sdigits, returns total QDC   
+  static void    Hit2Sdi     (AliHMPIDHit *pHit,TClonesArray*);                                                                        //hit -> 9 sdigits  
   static Bool_t  IsOverTh    (Float_t q                      )     {return q > 6;                                   }                  //is digit over threshold????
   static Bool_t  IsInside    (Float_t x,Float_t y            )     {return x>0&&y>0&&x<SizeAllX()&&y<SizeAllY();    }                  //is point inside pc boundary?
   inline static Bool_t IsInDead  (Float_t x,Float_t y        );                                                                        //is point in dead area?
          Float_t LorsX       (                               )const{return (PadPcX()+0.5)*SizePadX()+(Pc()%2)*(SizePcX()+SizeDead());} //center of the pad x, [cm]
          Float_t LorsY       (                               )const{return (PadPcY()+0.5)*SizePadY()+(Pc()/2)*(SizePcY()+SizeDead());} //center of the pad y, [cm]
+         void    Manual1     (Int_t c,Float_t x,Float_t y,Int_t q=33){AliHMPIDHit h(c,q,x,y); Set(&h,0);}                              //manual creation
+         void    Manual2     (Int_t c,Int_t p,Int_t x,Int_t y)     {fPad=Abs(c,p,x,y);}                                                //manual creation 
   inline Float_t Mathieson   (Float_t x,Float_t y            )const;                                                                   //Mathieson distribution 
-         Int_t   PadPcX      (                               )const{return A2X(fPad);}                                                 //pad x within PC 0..79
-         Int_t   PadPcY      (                               )const{return A2Y(fPad);}                                                 //pad y within PC 0..47
-         Int_t   PadChX      (                               )const{return A2X(fPad);}                                                 //pad x within chamber
-         Int_t   PadChY      (                               )const{return A2Y(fPad);}                                                 //pad y within chamber
+         Int_t   PadPcX      (                               )const{return A2X(fPad);}                                                 //pad pc x # 0..79
+         Int_t   PadPcY      (                               )const{return A2Y(fPad);}                                                 //pad pc y # 0..47
+         Int_t   PadChX      (                               )const{return (Pc()%2)*kPadPcX+PadPcX();}                                 //pad ch x # 0..159
+         Int_t   PadChY      (                               )const{return (Pc()/2)*kPadPcY+PadPcY();}                                 //pad ch y # 0..143
          Int_t   Pad         (                               )const{return fPad;}                                                      //absolute id of this pad
          Int_t   Pc          (                               )const{return A2P(fPad);}                                                 //PC position number
   static void    PrintSize   (                               );                                                                        //print all segmentation sizes      
@@ -63,7 +64,7 @@ public:
   static Int_t   Raw2X       (         UInt_t d,UInt_t a     )     {                                              return (d-1)*8+a/6;} //padx=f(d,a)
   static Int_t   Raw2Y       (UInt_t l,UInt_t r,UInt_t a     )     {Int_t a2y[6]={3,2,4,1,5,0};r=(l%2)?(24-r):r-1;return 6*(r%8)+a2y[a%6];}//pady=f(ddl,r,a)
          Int_t   Row         (                               )const{Int_t r=1+Pc()/2*8+PadPcY()/6; return (Pc()%2)? 25-r:r;}           //row r=1..24
-         void    Set         (Int_t c,Int_t s,Int_t x,Int_t y)     {fPad=Abs(c,s,x,y);fQ=0xa3;}                                        //set new digit
+  inline Bool_t  Set         (AliHMPIDHit *pHit,Int_t pad    );                                                                        //sdigit from hit in given pad
   
   static Float_t SizeAllX    (                               )     {return SizePadX()*kPadAllX+SizeDead();}                            //all PCs size x, [cm]        
   static Float_t SizeAllY    (                               )     {return SizePadY()*kPadAllY+2*SizeDead();}                          //all PCs size y, [cm]    
@@ -78,7 +79,7 @@ public:
   static Float_t SizeRad     (                               )     {return 1.5;}                                                       //Rad width   
   static void    Test        (                               );                                                                        //Test conversions
 protected:                  //AliDigit has fTracks[3]
-  Int_t   fPad;             //absolute pad number is chamber*kCham
+  Int_t    fPad;            //absolute pad number
   Float_t  fQ;              //QDC value, fractions are permitted for summable procedure  
   ClassDef(AliHMPIDDigit,4) //HMPID digit class       
 };//class AliHMPIDDigitN
@@ -86,37 +87,39 @@ protected:                  //AliDigit has fTracks[3]
 typedef AliHMPIDDigit AliRICHDigit; // for backward compatibility
 
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-AliHMPIDDigit::AliHMPIDDigit(Int_t c,Float_t q,Int_t t,Float_t x,Float_t y,Int_t flag):AliDigit(),fPad(Abs(-1,-1,-1,-1)),fQ(-1)  
+Bool_t AliHMPIDDigit::Set(AliHMPIDHit *pHit,Int_t pad) 
 {
-// Creation of sdigit  
-// Arguments: c- chamber 
-//            q- total QDC
-//            t -TID
-//            x,y - hit position, LORS      
-//            flag- which pad to try     
+// Creates digit
+// Arguments: pHit- pointer to the hit
+//            pad - for which pad to create
 //   Returns: none    
-  Int_t pc,padx,pady;
-  if     (x>=          0          && x<=  SizePcX()            ) {pc=0; padx=Int_t( x                           / SizePadX());}//PC 0 or 2 or 4
-  else if(x>=SizePcX()+SizeDead() && x<=  SizeAllX()           ) {pc=1; padx=Int_t((x-  SizePcX()-  SizeDead()) / SizePadX());}//PC 2 or 4 or 6
-  else return;
-  if     (y>=          0          && y<=  SizePcY()            ) {      pady=Int_t( y                           / SizePadY());}//PC 0 or 1
-  else if(y>=SizePcY()+SizeDead() && y<=2*SizePcY()+SizeDead() ) {pc+=2;pady=Int_t((y-  SizePcY()-  SizeDead()) / SizePadY());}//PC 2 or 3
-  else if(y>=SizeAllY()-SizePcY() && y<=  SizeAllY()           ) {pc+=4;pady=Int_t((y-2*SizePcY()-2*SizeDead()) / SizePadY());}//PC 4 or 5
-  else return;
   
-  switch(flag){
-    case 8:padx--;pady++;break;    case 1:pady++;break;    case 2:padx++; pady++;break;
+  fPad=Abs(-1,-1,-1,-1); fQ=-1; //reset
+  Int_t pc,px,py;
+  Float_t x=pHit->LorsX(),y=pHit->LorsY();
+  
+  if     (x>=          0          && x<=  SizePcX()            ) {pc=0; px=Int_t( x                           / SizePadX());}//PC 0 or 2 or 4
+  else if(x>=SizePcX()+SizeDead() && x<=  SizeAllX()           ) {pc=1; px=Int_t((x-  SizePcX()-  SizeDead()) / SizePadX());}//PC 2 or 4 or 6
+  else return kFALSE;
+  if     (y>=          0          && y<=  SizePcY()            ) {      py=Int_t( y                           / SizePadY());}//PC 0 or 1
+  else if(y>=SizePcY()+SizeDead() && y<=2*SizePcY()+SizeDead() ) {pc+=2;py=Int_t((y-  SizePcY()-  SizeDead()) / SizePadY());}//PC 2 or 3
+  else if(y>=SizeAllY()-SizePcY() && y<=  SizeAllY()           ) {pc+=4;py=Int_t((y-2*SizePcY()-2*SizeDead()) / SizePadY());}//PC 4 or 5
+  else return kFALSE;
+  
+  switch(pad){
+    case 8: px--;py++;break;    case 1:py++;break;    case 2:px++; py++;break;
                                               
-    case 7: padx--;      break;    case 0:       break;    case 3:padx++;        break;
+    case 7: px--;     break;    case 0:     break;    case 3:px++;      break;
                                                  
-    case 6:padx--;pady--;break;    case 5:pady--;break;    case 4:padx++; pady--;break;                                            
+    case 6: px--;py--;break;    case 5:py--;break;    case 4:px++; py--;break;                                            
   }
-  if(padx<0 || padx>=kPadPcX) return;
-  if(pady<0 || pady>=kPadPcY) return;
-  fPad=Abs(c,pc,padx,pady);
-  fQ=q*Mathieson(x,y);
-  fTracks[0]=t; 
-}    
+  if(px<0 || px>=kPadPcX) return kFALSE;
+  if(py<0 || py>=kPadPcY) return kFALSE;
+  fPad=Abs(pHit->Ch(),pc,px,py);
+  fQ=pHit->Q()*Mathieson(x,y);
+  fTracks[0]=pHit->Tid(); 
+  return kTRUE;
+}//Set()    
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 Int_t AliHMPIDDigit::Compare(const TObject *pObj) const
 {
index 6b41fb2..2667f14 100644 (file)
@@ -78,7 +78,7 @@ void AliHMPIDDigitizer::Sdi2Dig(TClonesArray *pSdiLst,TObjArray *pDigLst)
   
   for(Int_t i=0;i<7;i++){
     pLst[i]=(TClonesArray*)(*pDigLst)[i];
-    iCnt[i]=pLst[i]->GetEntries();         //in principle those lists should be empty                                                                       
+    iCnt[i]=0; if(pLst[i]->GetEntries()!=0) AliErrorClass("Some of digits lists is not empty");         //in principle those lists should be empty                                                                       
   }
   
   pSdiLst->Sort();  
@@ -98,7 +98,7 @@ void AliHMPIDDigitizer::Sdi2Dig(TClonesArray *pSdiLst,TObjArray *pDigLst)
     q=pSdig->Q();    
   }//sdigits loop (sorted)
   
-  if(AliHMPIDDigit::IsOverTh(q))  new((*pLst[iCh])[iCnt[iCh]++]) AliHMPIDDigit(iPad,(Int_t)q,aTids);             //add the last one, in case of empty sdigits list q=-1
+  if(AliHMPIDDigit::IsOverTh(q))  new((*pLst[iCh])[iCnt[iCh]++]) AliHMPIDDigit(iPad,(Int_t)q,aTids);           //add the last one, in case of empty sdigits list q=-1
                                                                                                                //so digit is not created    
 }//Sdi2Dig()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index 6390076..68ad3ff 100644 (file)
@@ -13,9 +13,9 @@
 //  * provided "as is" without express or implied warranty.                  *
 //  **************************************************************************
 
-#include "AliHMPIDHit.h" //class header
-#include <TPDGCode.h>   //Print() 
-#include <TString.h>
+#include "AliHMPIDHit.h"  //class header
+#include "AliHMPIDDigit.h"
+#include <TPDGCode.h>    
  
 ClassImp(AliHMPIDHit)
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -38,6 +38,8 @@ void AliHMPIDHit::Print(Option_t*)const
     case 50000051:     sPart="feed";break;
   }
 
-  Printf("%s Ch:%i TID:%6i,E:%9.3f eV, LORS:(%7.2f,%7.2f) MARS:(%7.2f,%7.2f,%7.2f)cm",sPart,Ch(),Tid(),E()*1e9,LorsX(),LorsY(),X(),Y(),Z());
+  Printf(" ch=%i                 (%7.3f,%7.3f) Q=%8.3f TID= %5i, MARS=(%7.2f,%7.2f,%7.2f) %s  %s",
+           Ch(),                LorsX(),LorsY(), Q(),  Tid(),         X(),  Y(),  Z(),   sPart, 
+                        (AliHMPIDDigit::IsInDead(LorsX(),LorsY()))? "IN DEAD ZONE":"");
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index dbf6bdd..c768945 100644 (file)
@@ -5,33 +5,50 @@
 
 #include <AliHit.h>           //base class
 #include <TVector3.h>         //ctor
+#include <TRandom.h>          //QdcTot()
 
 class AliHMPIDHit : public AliHit //   TObject-AliHit-AliHMPIDHit
 {
 public:
-  AliHMPIDHit(                                                                             ):AliHit(     ),fCh(-1),fPid(-1 ),fE(-1),fLorsX(-1),fLorsY(-1) {} //default ctor
-  AliHMPIDHit(Int_t c,Float_t e,Int_t pid,Int_t tid,Float_t xl,Float_t yl,const TVector3 &p):AliHit(0,tid),fCh(c ),fPid(pid),fE(e ),fLorsX(xl),fLorsY(yl) {fX=p.X();fY=p.Y();fZ=p.Z();}           
-  AliHMPIDHit(Int_t c,Float_t e,Int_t pid,Int_t tid,Float_t xl,Float_t yl                  ):              fCh(c ),fPid(pid),fE(e ),fLorsX(xl),fLorsY(yl) {fTrack=tid;}           
+  AliHMPIDHit(                                                                             ):AliHit(     ),fCh(-1),fPid(-1  ),fQ(-1),fLorsX(-1),fLorsY(-1) {} //default ctor
+  AliHMPIDHit(Int_t c,Float_t e,Int_t pid,Int_t tid,Float_t xl,Float_t yl,const TVector3 &p):AliHit(0,tid),fCh(c ),fPid(pid ),fQ(-1),fLorsX(xl),fLorsY(yl) {QdcTot(e);fX=p.X();fY=p.Y();fZ=p.Z();}           
+  AliHMPIDHit(Int_t c,Int_t q  ,                    Float_t xl,Float_t yl                  ):              fCh(c ),fPid(2212),fQ(q ),fLorsX(xl),fLorsY(yl) {fTrack=Int_t(1000*gRandom->Rndm());}//manual ctor           
   virtual ~AliHMPIDHit()                                                                                                {}
 //framework part
   void     Print(Option_t *option="")const;                                    //from TObject to print current status
 //private part  
-  Int_t   Ch    ()const{return fCh;                                    }       //Chamber
-  Float_t E     ()const{return fE;                                     }       //Eloss for MIP hit or Etot for photon hit, [GeV]
-  Float_t LorsX ()const{return fLorsX;                                 }       //hit X position in LORS, [cm]
-  Float_t LorsY ()const{return fLorsY;                                 }       //hit Y position in LORS, [cm]
-  Int_t   Pid   ()const{return fPid;                                   }       //PID
-  Int_t   Tid   ()const{return fTrack;                                 }       //TID
+         Int_t   Ch    (         )const{return fCh;                                    }       //Chamber
+         Float_t LorsX (         )const{return fLorsX;                                 }       //hit X position in LORS, [cm]
+         Float_t LorsY (         )const{return fLorsY;                                 }       //hit Y position in LORS, [cm]
+         Int_t   Pid   (         )const{return fPid;                                   }       //PID
+         Float_t Q     (         )const{return fQ;                                     }       //Eloss for MIP hit or Etot for photon hit, [GeV]
+  inline Float_t QdcTot(Float_t e);                                                            //calculate total charge of the hit          
+         Int_t   Tid   (         )const{return fTrack;                                 }       //TID
          
 protected:                                                                     //AliHit has fTid,fX,fY,fZ 
   Int_t    fCh;                                                                //Chamber
   Int_t    fPid;                                                               //PID
-  Float_t  fE;                                                                 //Eloss for MIP or Etot for photon [GeV]
+  Float_t  fQ;                                                                 //total charge [QDC]
   Float_t  fLorsX;                                                             //hit X position in chamber LORS, [cm]
   Float_t  fLorsY;                                                             //hit Y position in chamber LORS, [cm]
-  ClassDef(AliHMPIDHit,4)                                                       //HMPID hit class 
+  ClassDef(AliHMPIDHit,5)                                                      //HMPID hit class 
 };//class AliHMPIDhit
-
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++  
+Float_t AliHMPIDHit::QdcTot(Float_t e)
+{
+// Samples total charge of the hit
+// Arguments: e- hit energy [GeV]; for mip Eloss for photon Etot   
+//   Returns: total QDC
+  Float_t  x=(LorsX() > 66.6)? LorsX()-66.6:LorsX();                                                 //sagita is for PC (0-64) and not for chamber   
+  Float_t  qdcEle=34.06311+0.2337070*x+5.807476e-3*x*x-2.956471e-04*x*x*x+2.310001e-06*x*x*x*x;      //reparametrised from DiMauro
+  
+  Int_t iNele=Int_t(e/26e-9);  if(iNele<1) iNele=1;                                                  //number of electrons created by hit
+  for(Int_t i=1;i<=iNele;i++) fQ-=qdcEle*TMath::Log(gRandom->Rndm()+1e-6);                           //1e-6 is a protection against 0 from rndm  
+  return fQ;
+}  
+  
+  
+  
 typedef AliHMPIDHit AliRICHHit; // for backward compatibility
     
 #endif
index 2059c27..f44e55a 100644 (file)
@@ -150,13 +150,11 @@ Double_t AliHMPIDRecon::FindRingCkov(Int_t)
   
   for(Int_t i=0;i<fPhotCnt;i++){//candidates loop
     if(fPhotFlag[i] == 2){
+      if(fPhotCkov[i]<=0) continue;//?????????????????Flag photos = 2 may imply CkovEta = 0?????????????? 
       if(fPhotCkov[i]<ckovMin) ckovMin=fPhotCkov[i];  //find max and min Theta ckov from all candidates within probable window
       if(fPhotCkov[i]>ckovMax) ckovMax=fPhotCkov[i]; 
       weightThetaCerenkov += fPhotCkov[i]*fPhotWei[i];   wei += fPhotWei[i];                 //collect weight as sum of all candidate weghts   
       
-     //Double_t phiref=(GetPhiPoint()-GetTrackPhi());
-       if(fPhotCkov[i]<=0) continue;//?????????????????Flag photos = 2 may imply CkovEta = 0?????????????? 
-                                     
       sigma2 += 1./Sigma2(fPhotCkov[i],fPhotPhi[i]);
     }
   }//candidates loop
@@ -164,7 +162,6 @@ Double_t AliHMPIDRecon::FindRingCkov(Int_t)
   if(sigma2>0) fCkovSigma2=1./sigma2;
   else         fCkovSigma2=1e10;  
   
-
   if(wei != 0.) weightThetaCerenkov /= wei; else weightThetaCerenkov = 0.;  
   return weightThetaCerenkov;
 }//FindCkovRing()
index e338f37..385ecc8 100644 (file)
 #include <AliRunLoader.h>         //Reconstruct() for simulated digits
 #include <AliRawReader.h>         //Reconstruct() for raw digits
 #include <AliRun.h>               //Reconstruct()
-#include <TH1F.h>                 //CluQA() 
-#include <TH2F.h>                 //CluQA() 
-#include <TCanvas.h>              //CluQA() 
-#include <TNtupleD.h>             //CheckPR()
 ClassImp(AliHMPIDReconstructor)
 
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDReconstructor::CluQA(AliRunLoader *pAL)
-{
-// Quality assesment plots for clusters. 
-// This methode takes list of digits and form list of clusters again in order to 
-// calculate cluster shape and cluster particle mixture    
-  AliLoader *pRL=pAL->GetDetectorLoader("HMPID");  AliHMPID *pRich=(AliHMPID*)pAL->GetAliRun()->GetDetector("HMPID");//get pointers for HMPID and HMPID loader
-  Int_t iNevt=pAL->GetNumberOfEvents();  if(iNevt==0)             {AliInfoClass("No events");return;}   
-                                         if(pRL->LoadDigits())    {AliInfoClass("No digits file");return;}
-                                            pAL->LoadHeader();
-                                            pAL->LoadKinematics(); 
-//                                            AliStack *pStack=pAL->Stack();
-  TH1::AddDirectory(kFALSE);
-  
-        
-  TH1F*    pQ=new TH1F("RiAllQ"  ,"Charge All"           ,4000 ,0  ,4000);// Q hists
-  TH1F* pCerQ=new TH1F("RiCerQ"  ,"Charge Ckov"          ,4000 ,0  ,4000);
-  TH1F* pMipQ=new TH1F("RiMipQ"  ,"Charge MIP"           ,4000 ,0  ,4000);
-  
-  TH1F*    pS=new TH1F("RichCluSize"    ,"Cluster size;size"         ,100  ,0  ,100 );// size hists
-  TH1F* pCerS=new TH1F("RichCluCerSize" ,"Ckov size;size"            ,100  ,0  ,100 );
-  TH1F* pMipS=new TH1F("RichCluMipSize" ,"MIP size;size"             ,100  ,0  ,100 );
-  
-  TH2F*    pM=new TH2F("RichCluMap"     ,"Cluster map;x [cm];y [cm]" ,1000 ,0  ,AliHMPIDDigit::SizePcX(),1000,0,AliHMPIDDigit::SizePcY()); // maps
-  TH2F* pMipM=new TH2F("RichCluMipMap"  ,"MIP map;x [cm];y [cm]"     ,1000 ,0  ,AliHMPIDDigit::SizePcX(),1000,0,AliHMPIDDigit::SizePcY());
-  TH2F* pCerM=new TH2F("RichCluCerMap"  ,"Ckov map;x [cm];y [cm]"    ,1000 ,0  ,AliHMPIDDigit::SizePcX(),1000,0,AliHMPIDDigit::SizePcY());
-  
-  
-  for(Int_t iEvt=0;iEvt<iNevt;iEvt++){
-    pAL->GetEvent(iEvt);               
-    pRL->TreeD()->GetEntry(0); 
-    TClonesArray *pCluLst=new TClonesArray("AliHMPIDCluster");//tmp list of clusters for this event
-    
-    for(Int_t iCh=0;iCh<7;iCh++) Dig2Clu(pRich->DigLst(iCh),pCluLst,kFALSE);//cluster finder for all chamber if any digits present
-    
-    for(Int_t iClu=0;iClu<pCluLst->GetEntriesFast();iClu++){
-      AliHMPIDCluster *pClu = (AliHMPIDCluster*)pCluLst->At(iClu);
-      Int_t cfm=0; for(Int_t iDig=0;iDig<pClu->Size();iDig++)  cfm+=pClu->Dig(iDig)->Ch(); //collect ckov-fee-mip structure of current cluster ?????
-      Int_t iNckov=cfm/1000000;      Int_t iNfee =cfm%1000000/1000;      Int_t iNmip =cfm%1000000%1000; 
-
-                                             pQ   ->Fill(pClu->Q()) ; pS   ->Fill(pClu->Size()) ; pM    ->Fill(pClu->X(),pClu->Y()); //all clusters                                      
-      if(iNckov!=0 && iNfee==0 && iNmip==0) {pCerQ->Fill(pClu->Q()) ; pCerS->Fill(pClu->Size()) ; pCerM ->Fill(pClu->X(),pClu->Y());}//ckov only cluster
-      if(iNckov==0 && iNfee==0 && iNmip!=0) {pMipQ->Fill(pClu->Q()) ; pMipS->Fill(pClu->Size()) ; pMipM ->Fill(pClu->X(),pClu->Y());}//mip only cluster
-                                       
-    }//clusters loop   
-    pCluLst->Clear();delete pCluLst;
-  }//events loop
-  
-  pRL->UnloadDigits(); pAL->UnloadKinematics(); pAL->UnloadHeader();
-  TCanvas *pC=new TCanvas("RichCluQA",Form("QA for cluster from %i events",iNevt),1000,900); pC->Divide(3,3);
-  pC->cd(1);    pM->Draw();          pC->cd(2);    pQ->Draw();       pC->cd(3);    pS->Draw();        
-  pC->cd(4); pMipM->Draw();          pC->cd(5); pMipQ->Draw();       pC->cd(6); pMipS->Draw();        
-  pC->cd(7); pCerM->Draw();          pC->cd(8); pCerQ->Draw();       pC->cd(9); pCerS->Draw();        
-}//CluQA()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDReconstructor::Dig2Clu(TClonesArray *pDigLst,TClonesArray *pCluLst,Bool_t isTryUnfold)
+void AliHMPIDReconstructor::Dig2Clu(TObjArray *pDigAll,TObjArray *pCluAll,Bool_t isTryUnfold)
 {
 // Finds all clusters for a given digits list provided not empty. Currently digits list is a list of all digits for a single chamber.
 // Puts all found clusters in separate lists, one per clusters. 
-// Arguments: pDigLst     - list of digits provided not empty
-//            pCluLst     - list of clusters, provided empty     
+// Arguments: pDigAll     - list of digits for all chambers 
+//            pCluAll     - list of clusters for all chambers
 //            isTryUnfold - flag to choose between CoG and Mathieson fitting  
 //  Returns: none    
-  TMatrixF digMap(AliHMPIDDigit::kPadAllX,AliHMPIDDigit::kPadAllY);  digMap=(Float_t)-1;   //digit map for single chamber reseted to -1
-  for(Int_t iDig=0;iDig<pDigLst->GetEntriesFast();iDig++){                                 //digits loop to fill digits map
-    AliHMPIDDigit *pDig= (AliHMPIDDigit*)pDigLst->At(iDig);                                //get current digit
-    digMap( pDig->PadChX(), pDig->PadChY() )=iDig;                                             //fill the map, (padx,pady) cell takes digit index
-  }                                                                                        //digits loop to fill digits map 
+  TMatrixF padMap(AliHMPIDDigit::kPadAllX,AliHMPIDDigit::kPadAllY);                   //pads map for single chamber 0..159 x 0..143 
+  
+  for(Int_t iCh=0;iCh<7;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
   
-  AliHMPIDCluster clu;                                                                      //tmp cluster to be used as current
+    padMap=(Float_t)-1;                                                                   //reset map to -1 (means no digit for this pad)  
+    TClonesArray *pCluCur=(TClonesArray*)pCluAll->At(iCh);                                //get list of clusters for current chamber
+    
+    for(Int_t iDig=0;iDig<pDigCur->GetEntriesFast();iDig++){                              //digits loop to fill pads map
+      AliHMPIDDigit *pDig= (AliHMPIDDigit*)pDigCur->At(iDig);                             //get current digit
+      padMap( pDig->PadChX(), pDig->PadChY() )=iDig;                                      //fill the map, (padx,pady) cell takes digit index
+    }//digits loop to fill digits map 
+    
+    AliHMPIDCluster clu;                                                                  //tmp cluster to be used as current
   
-  for(Int_t iDig=0;iDig<pDigLst->GetEntriesFast();iDig++){                                 //digits loop to form clusters list
-    AliHMPIDDigit *pDig=(AliHMPIDDigit*)pDigLst->At(iDig);                                 //take current digit
-    if(!(pDig=UseDig(pDig->PadChX(),pDig->PadChY(),pDigLst,&digMap))) continue;            //this digit is already taken in FormClu(), go after next digit
-    FormClu(&clu,pDig,pDigLst,&digMap);                                                    //form cluster starting from this digit by recursion
-    clu.Solve(pCluLst,isTryUnfold);                                                        //solve this cluster and add all unfolded clusters to provided list  
-    clu.Reset();                                                                           //empty current cluster
-  }                                                                                        //digits loop to form clusters list
+    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
+      FormClu(&clu,pDig,pDigCur,&padMap);                                                 //form cluster starting from this digit by recursion
+      
+      clu.Solve(pCluCur,isTryUnfold);                                                     //solve this cluster and add all unfolded clusters to provided list  
+      clu.Reset();                                                                        //empty current cluster
+    }//digits loop for current chamber
+  }//chambers loop
 }//Dig2Clu()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void  AliHMPIDReconstructor::FormClu(AliHMPIDCluster *pClu,AliHMPIDDigit *pDig,TClonesArray *pDigLst,TMatrixF *pDigMap)
@@ -113,19 +64,18 @@ void  AliHMPIDReconstructor::FormClu(AliHMPIDCluster *pClu,AliHMPIDDigit *pDig,T
 // Forms the initial cluster as a combination of all adjascent digits. Starts from the given digit
 // then calls itself recursevly  for all possible neighbours.
 // Arguments: pClu - pointer to cluster being formed
-//  Returns: none   ???????????????????
+//  Returns: none   
   pClu->DigAdd(pDig);//take this digit in cluster
 
-  Int_t x[4],y[4];
+  Int_t cnt=0,cx[4],cy[4];
   
-  Int_t cnt=0;  Int_t iPadX=pDig->PadPcX(); Int_t iPadY=pDig->PadPcY();
-  if(iPadX != 0)                        {x[cnt]=iPadX-1; y[cnt]=iPadY;   cnt++;}       //left
-  if(iPadX != AliHMPIDDigit::kPadPcX-1) {x[cnt]=iPadX+1; y[cnt]=iPadY;   cnt++;}       //right
-  if(iPadY != 0)                        {x[cnt]=iPadX;   y[cnt]=iPadY-1; cnt++;}       //down
-  if(iPadY != AliHMPIDDigit::kPadPcY-1) {x[cnt]=iPadX;   y[cnt]=iPadY+1; cnt++;}       //up
+  if(pDig->PadPcX() != 0                       ){cx[cnt]=pDig->PadChX()-1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //left
+  if(pDig->PadPcX() != AliHMPIDDigit::kPadPcX-1){cx[cnt]=pDig->PadChX()+1; cy[cnt]=pDig->PadChY()  ;cnt++;}       //right
+  if(pDig->PadPcY() != 0                       ){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()-1;cnt++;}       //down
+  if(pDig->PadPcY() != AliHMPIDDigit::kPadPcY-1){cx[cnt]=pDig->PadChX()  ; cy[cnt]=pDig->PadChY()+1;cnt++;}       //up
   
   for (Int_t i=0;i<cnt;i++)
-    if((pDig=UseDig(x[i],y[i],pDigLst,pDigMap))) FormClu(pClu,pDig,pDigLst,pDigMap);   //check if this neighbour pad fired and mark it as taken  
+    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  
   
 }//FormClu()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -144,7 +94,7 @@ void AliHMPIDReconstructor::Reconstruct(AliRunLoader *pAL)const
     pAL->GetEvent(iEvtN); AliDebug(1,Form("Processing event %i...",iEvtN)); //switch current directory to next event    
     pRL->TreeD()->GetEntry(0);  pRL->MakeTree("R");  pRich->MakeBranch("R");  //load digits to memory  and create branches for clusters              
     
-    for(Int_t iCh=0;iCh<7;iCh++) Dig2Clu(pRich->DigLst(iCh),pRich->CluLst(iCh));//cluster finder 
+    Dig2Clu(pRich->DigLst(),pRich->CluLst());//cluster finder 
       
     pRL->TreeR()->Fill();            //fill tree for current event
     pRL->WriteRecPoints("OVERWRITE");//write out clusters for current event
@@ -167,25 +117,28 @@ void AliHMPIDReconstructor::Reconstruct(AliRunLoader *pAL,AliRawReader* pRR)cons
   
   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 all chambers
+  
   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
-      TClonesArray *pDigLst=new TClonesArray("AliHMPIDDigit"); Int_t iDigCnt=0; //tmp list of digits for single chamber
       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(ddl,w32);  
         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((*pDigLst)[iDigCnt++]) AliHMPIDDigit(dig); //add this digit to the tmp list
+        new((*((TClonesArray*)digLst.At(iCh)))[iDigCnt[iCh]++]) AliHMPIDDigit(dig); //add this digit to the tmp list
       }//raw records loop
-      if(iDigCnt) Dig2Clu(pDigLst,pRich->CluLst(iCh));//cluster finder for the current chamber if any digits present
       pRR->Reset();        
-      pDigLst->Delete();  iDigCnt=0;//clean up list of digits for the current chamber
     }//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();
index b453c44..33c6395 100644 (file)
@@ -29,10 +29,9 @@ public:
    using AliReconstructor::Reconstruct;                                                                                 //to get rid of virtual hidden warning 
 
   //private part  
-  static        void          Dig2Clu (TClonesArray*pDigLst,TClonesArray *pCluLst,Bool_t isTryUnfold=kTRUE            );//digits list -> clusters list
-  static        void          CluQA   (AliRunLoader* pAL                                                              );//QA for clusters
-  static        void          FormClu (AliHMPIDCluster *pClu,AliHMPIDDigit *pDig,TClonesArray *pDigLst,TMatrixF *pDigMap);//cluster formation recursive algorithm
-  static inline AliHMPIDDigit* UseDig  (Int_t padX,Int_t padY,TClonesArray *pDigList,TMatrixF *pDigMap                 );//use this pad's digit to form a cluster
+  static        void           Dig2Clu (TObjArray *pDigLst,TObjArray *pCluLst,Bool_t isUnfold=kTRUE                      );//digits->clusters
+  static        void           FormClu (AliHMPIDCluster *pClu,AliHMPIDDigit *pDig,TClonesArray *pDigLst,TMatrixF *pPadMap);//cluster formation recursive algorithm
+  static inline AliHMPIDDigit* UseDig  (Int_t padX,Int_t padY,                    TClonesArray *pDigLst,TMatrixF *pDigMap);//use this pad's digit to form a cluster
 
   protected:
   ClassDef(AliHMPIDReconstructor, 0)   //class for the HMPID reconstruction
@@ -41,7 +40,7 @@ public:
 typedef AliHMPIDReconstructor AliRICHReconstructor; // for backward compatibility
 
 //__________________________________________________________________________________________________
-AliHMPIDDigit* AliHMPIDReconstructor::UseDig(Int_t padX,Int_t padY,TClonesArray *pDigLst,TMatrixF *pDigMap)
+AliHMPIDDigit* AliHMPIDReconstructor::UseDig(Int_t padX,Int_t padY,TClonesArray *pDigLst,TMatrixF *pPadMap)
 {
 //Digit map contains a matrix if digit numbers.
 //Main operation in forming initial cluster is done here. Requested digit pointer is returned and this digit marked as taken.
@@ -49,8 +48,8 @@ AliHMPIDDigit* AliHMPIDReconstructor::UseDig(Int_t padX,Int_t padY,TClonesArray
 //           pDigLst   - list of digits for one sector
 //           pDigMap   - map of those digits
 //  Returns: pointer to digit if not yet used or 0 if used
-  Int_t iDig=(Int_t)(*pDigMap)(padX,padY);(*pDigMap)(padX,padY)=-1;//take digit number from the map and reset this map cell to -1
-  if(iDig!=-1)    return (AliHMPIDDigit*)pDigLst->At(iDig);         //digit pointer
+  Int_t iDig=(Int_t)(*pPadMap)(padX,padY);(*pPadMap)(padX,padY)=-1;//take digit number from the map and reset this map cell to -1
+  if(iDig!=-1)    return (AliHMPIDDigit*)pDigLst->At(iDig);        //digit pointer
   else            return 0;
 }
 
index 083892c..b0d9b49 100644 (file)
@@ -235,11 +235,14 @@ void HmpidHeader(ifstream *pFile,Bool_t isPrint=kFALSE)
   for(Int_t i=13;i<=15;i++){ pFile->read((char*)&w32,4); Printf("Word #%2i=%12i empty",i,w32);}
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void HmpidPayload(ifstream *pFile,Int_t iDdl,TClonesArray *pDigLst)
+void HmpidPayload(ifstream *pFile,Int_t iDdl,TObjArray *pDigAll)
 {
 // payload is 8 sequences with structure: WC36A8 then WC number of w32
+  
+  TClonesArray *pDig1=(TClonesArray*)pDigAll->At(iDdl/2); //get list of digits for requested chamber
+  
   UInt_t w32=0;
-  Int_t iDigCnt=pDigLst->GetEntriesFast();
+  Int_t iDigCnt=pDig1->GetEntriesFast();
   for(Int_t row=1;row<=8;row++){  
     pFile->read((char*)&w32,4);  Int_t wc=AliBitPacking::UnpackWord(w32,16,31); Int_t ma=AliBitPacking::UnpackWord(w32, 0,15);    
     if(ma!=0x36a8) Printf("ERROR ERROR ERROR WRONG Marker=0x%x ERROR ERROR ERROR",ma);    
@@ -248,12 +251,12 @@ void HmpidPayload(ifstream *pFile,Int_t iDdl,TClonesArray *pDigLst)
       if(w32&0x08000000) continue; //it's DILOGIC CW
       AliHMPIDDigit *pDig=new AliHMPIDDigit;
       pDig->Raw(iDdl,w32);   
-      new ((*pDigLst)[iDigCnt++]) AliHMPIDDigit(*pDig);
+      new ((*pDig1)[iDigCnt++]) AliHMPIDDigit(*pDig);
     }//words loop 
   }//rows loop
 }//HmpidPayload()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void docosmic(const char* name)
+void docosmic(const char* name,Int_t iMaxEvt=9999999)
 {
   gBenchmark->Start("Cosmic converter");
   ifstream in(name);
@@ -267,70 +270,85 @@ void docosmic(const char* name)
   TFile *pOut=new TFile(rooName,"recreate");
   TTree *pTr=new TTree("cosmic","some for time being");
   
-  TClonesArray *pDigLst=new TClonesArray("AliHMPIDDigit");   pTr->Branch("Digs",&pDigLst);
-  TClonesArray *pCluLst=new TClonesArray("AliHMPIDCluster"); pTr->Branch("Clus",&pCluLst);
+  
+  TObjArray *pDigAll=new TObjArray(7); pDigAll->SetOwner(); for(Int_t i=0;i<7;i++) pDigAll->AddAt(new TClonesArray("AliHMPIDDigit")  ,i); pTr->Branch("Digs",&pDigAll,64000,0); 
+  TObjArray *pCluAll=new TObjArray(7); pCluAll->SetOwner(); for(Int_t i=0;i<7;i++) pCluAll->AddAt(new TClonesArray("AliHMPIDCluster"),i); pTr->Branch("Clus",&pCluAll,64000,0);
 
+  Int_t iEvt=0;
   while(1){      
     Int_t iSize=DateHeader(&in,      isPrint);  if(iSize==68) continue;  //DATE header 
     if(in.eof()) break;
     HmpidHeader (&in,      isPrint);    //HMPID header 
-    HmpidPayload(&in,ddl+1,pDigLst);    //HMPID payload
+    HmpidPayload(&in,ddl+1,pDigAll);    //HMPID payload
     HmpidHeader (&in,      isPrint);    //HMPID header 
-    HmpidPayload(&in,ddl  ,pDigLst);    //HMPID payload
+    HmpidPayload(&in,ddl  ,pDigAll);    //HMPID payload
     
-    AliHMPIDReconstructor::Dig2Clu(pDigLst,pCluLst);
+    AliHMPIDReconstructor::Dig2Clu(pDigAll,pCluAll);
     pTr->Fill();
-    pDigLst->Clear(); pCluLst->Clear();
-  }
+    for(Int_t i=0;i<7;i++){
+      ((TClonesArray*)pDigAll->At(i))->Clear(); 
+      ((TClonesArray*)pCluAll->At(i))->Clear(); 
+    }
+    
+    if(!(iEvt%200)) Printf("Event %i processed",iEvt);
+    iEvt++;
+    if(iEvt==iMaxEvt) break;
+  }//events loop
   
-  pTr->Write();
-  pOut->Close();
-  delete pDigLst;
+  pTr->Write();  pOut->Close();
+  pDigAll->Delete(); pCluAll->Delete();
   in.close();
+  
+  Printf("Total %i events processed",iEvt);
   gBenchmark->Show("Cosmic converter");
 }//docosmic()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void cosmic()
 {
-  TChain *pCh=new TChain("cosmic");
-  pCh->Add("test.root");
+  TChain *pCh=new TChain("cosmic");  pCh->Add("test.root");
 
   
   TH1F *pDigQ=new TH1F("digQ","Digits QDC",500,0,4100);
   TH1F *pDigO=new TH1F("digO","Digits # per event",500,0,8000);
   TH2F *pDigM=new TH2F("digM","Digits map",500,0,131,500,0,127);
     
-  TH1F *pCluQ=new TH1F("cluQ","Clusters QDC",500,0,12100);
+  TH1F *pCluQ   =new TH1F("cluQ","Clusters QDC",500,0,12100);
+  TH1F *pCluQmax=new TH1F("cluQmax","Clusters Max QDC",500,0,12100); pCluQmax->SetLineColor(kRed);
   TH1F *pCluO=new TH1F("cluO","Clusters # per event",500,0,5000);
   TH2F *pCluM=new TH2F("cluM","Clusters map",500,0,131,500,0,127);
   
-  TClonesArray *pDigLst=new TClonesArray("AliHMPIDDigit");   pCh->SetBranchAddress("Digs",&pDigLst);
-  TClonesArray *pCluLst=new TClonesArray("AliHMPIDCluster"); pCh->SetBranchAddress("Clus",&pCluLst);
-  
+  TObjArray *pDigAll=0; pCh->SetBranchAddress("Digs",&pDigAll);
+  TObjArray *pCluAll=0; pCh->SetBranchAddress("Clus",&pCluAll);
+
+    
   for(Int_t iEvt=0;iEvt<pCh->GetEntries();iEvt++){
     pCh->GetEntry(iEvt);
     
-    pDigO->Fill(pDigLst->GetEntriesFast());
-    pCluO->Fill(pCluLst->GetEntriesFast());
+    TClonesArray *pDigCh=(TClonesArray*)pDigAll->At(0);
+    TClonesArray *pCluCh=(TClonesArray*)pCluAll->At(0);
+    pDigO->Fill(pDigCh->GetEntriesFast());
+    pCluO->Fill(pCluCh->GetEntriesFast());
     
-    for(Int_t iDig=0;iDig<pDigLst->GetEntriesFast();iDig++){//digits loop
-      AliHMPIDDigit *pDig=(AliHMPIDDigit*)pDigLst->UncheckedAt(iDig);
+    for(Int_t iDig=0;iDig<pDigCh->GetEntriesFast();iDig++){//digits loop
+      AliHMPIDDigit *pDig=(AliHMPIDDigit*)pDigCh->UncheckedAt(iDig);
       pDigQ->Fill(pDig->Q());
       pDigM->Fill(pDig->LorsX(),pDig->LorsY());
     }//digits loop
-    
-    for(Int_t iClu=0;iClu<pCluLst->GetEntriesFast();iClu++){//clusters loop
-      AliHMPIDCluster *pClu=(AliHMPIDCluster*)pCluLst->UncheckedAt(iClu);
+    Int_t qmax=0;    
+    for(Int_t iClu=0;iClu<pCluCh->GetEntriesFast();iClu++){//clusters loop
+      AliHMPIDCluster *pClu=(AliHMPIDCluster*)pCluCh->UncheckedAt(iClu);
       pCluQ->Fill(pClu->Q());
+      if(pClu->Q()>qmax) qmax=pClu->Q();
       pCluM->Fill(pClu->X(),pClu->Y());
     }//digits loop
+    pCluQmax->Fill(qmax);
   }//entries loop
   
   TCanvas *pC=new TCanvas("comic","cosmic"); pC->Divide(2,3);
   
   pC->cd(1); pDigM->Draw(); pC->cd(2); pCluM->Draw();
-  pC->cd(3); pDigQ->Draw(); pC->cd(4); pCluQ->Draw();
-  pC->cd(5); pDigO->Draw(); pC->cd(6); pCluO->Draw();
+  pC->cd(3); gPad->SetLogy(); pDigQ->Draw(); pC->cd(4); gPad->SetLogy(); pCluQ->Draw(); pCluQmax->Draw("same");
+  pC->cd(5); gPad->SetLogy(); pDigO->Draw(); pC->cd(6); gPad->SetLogy(); pCluO->Draw();
 }//cosmic()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
index 41ef32f..4a14ac1 100644 (file)
@@ -1,20 +1,12 @@
-#include "AliHMPIDTracker.h" //class header
-#include "AliHMPID.h"             //GetTrackPoint(),PropagateBack()   
-#include "AliHMPIDCluster.h"      //GetTrackPoint(),PropagateBack() 
-#include "AliHMPIDParam.h"        //GetTrackPoint(),PropagateBack()
-#include "AliHMPIDRecon.h"        //PropagateBack()
+#include "AliHMPIDTracker.h"     //class header
+#include "AliHMPIDCluster.h"     //GetTrackPoint(),PropagateBack() 
+#include "AliHMPIDParam.h"       //GetTrackPoint(),PropagateBack()
+#include "AliHMPIDRecon.h"       //PropagateBack()
 #include <AliESD.h>              //PropagateBack()  
 #include <AliRun.h>              //GetTrackPoint(),PropagateBack()  
 #include <AliTrackPointArray.h>  //GetTrackPoint()
 #include <AliAlignObj.h>         //GetTrackPoint()
 ClassImp(AliHMPIDTracker)
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-AliHMPIDTracker::AliHMPIDTracker():AliTracker()
-{
-// AliHMPIDTracker is created from AliReconstraction::Run() which invokes AliReconstraction::CreateTrackers() 
-// which in turn invokes AliHMPIDReconstructor::CreateTracker(). 
-// Note that this is done just once per session before AliReconstruction::Run() goes to events loop.
-}
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++    
 Bool_t AliHMPIDTracker::GetTrackPoint(Int_t idx, AliTrackPoint& point) const
 {
@@ -43,13 +35,12 @@ Int_t AliHMPIDTracker::LoadClusters(TTree *pCluTree)
   AliDebug(1,"Start.");  pCluTree->GetEntry(0);  AliDebug(1,"Stop."); return 0;
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Int_t AliHMPIDTracker::PropagateBack(AliESD *pESD)
+Int_t AliHMPIDTracker::Recon(AliESD *pESD,TObjArray *pCluAll)
 {
 // Interface callback methode invoked by AliRecontruction::RunTracking() during tracking after TOF. It's done just once per event
 // Arguments: pESD - pointer to Event Summary Data class instance which contains a list of tracks
 //   Returns: error code, 0 if no errors   
-  Int_t iNtracks=pESD->GetNumberOfTracks();  AliDebug(1,Form("Start with %i tracks",iNtracks));
-  AliHMPID *pRich=((AliHMPID*)gAlice->GetDetector("HMPID"));  
+  Int_t iNtracks=pESD->GetNumberOfTracks();  AliDebugClass(1,Form("Start with %i tracks",iNtracks));
   AliHMPIDRecon recon;                                                        //instance of reconstruction class, nothing important in ctor
    
   AliHMPIDParam *pParam=AliHMPIDParam::Instance();
@@ -77,7 +68,7 @@ Int_t AliHMPIDTracker::PropagateBack(AliESD *pESD)
     
     if(iCh==-1) continue;                                                                  //no intersection at all, go after next track
     
-    TClonesArray *pCluLst=pRich->CluLst(iCh);                                              //get clusters list for intersected chamber
+    TClonesArray *pCluLst=(TClonesArray *)pCluAll->At(iCh);                                //get clusters list for intersected chamber
     
     Double_t    dMin=999;                                                                  //distance between track-PC intersection point and current cluster
     Int_t   iMip=-1;                                                                       //index of cluster nearest to intersection point
@@ -105,7 +96,7 @@ Int_t AliHMPIDTracker::PropagateBack(AliESD *pESD)
                    pTrk->SetHMPIDchi2     (recon.CkovSigma2());                              //error squared 
                    pTrk->SetHMPIDmip      (pMipClu->X(),pMipClu->Y(),pMipClu->Q(),iMip);     //info on mip cluster + n. phot.
  }//ESD tracks loop
-  AliDebug(1,"Stop pattern recognition");
+  AliDebugClass(1,"Stop pattern recognition");
   return 0; // error code: 0=no error;
 }//PropagateBack()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index 7849eb0..e1f836f 100644 (file)
@@ -2,28 +2,41 @@
 #define AliHMPIDTracker_h
 
 #include <AliTracker.h> //base class
+#include "AliHMPID.h"   //Recon()
+#include <AliRun.h>     //Recon()
 
-class TNtupleD;         //RecWithStack()   
-class AliESD;           //Clusters2Tracks(), RefitInward(), PropagateBack(), RecWithESD()
+class AliESD;           
 
 class AliHMPIDTracker : public AliTracker
 {
 public:
-           AliHMPIDTracker(); 
-  virtual ~AliHMPIDTracker()                    {}
+           AliHMPIDTracker():AliTracker()                               {} 
+  virtual ~AliHMPIDTracker()                                            {}
 //framework part  
-  AliCluster *GetCluster     (Int_t                      )const  {return 0;} //pure virtual from AliTracker 
-  Bool_t      GetTrackPoint  (Int_t idx,AliTrackPoint &pt)const;             //             from AliTracker  
-  Int_t       Clusters2Tracks(AliESD *                   )       {return 0;} //pure virtual from AliTracker 
-  Int_t       LoadClusters   (TTree *pCluTr              );                  //pure virtual from AliTracker   
-  Int_t       PropagateBack  (AliESD *                   );                  //pure virtual from AliTracker invoked from AliReconstruction::RunTracking()
-  Int_t       RefitInward    (AliESD *                   )       {return 0;} //pure virtual from AliTracker 
-  void        UnloadClusters (                           )       {         } //pure virtual from AliTracker 
+         AliCluster *GetCluster     (Int_t                      )const  {return 0;} //pure virtual from AliTracker 
+         Bool_t      GetTrackPoint  (Int_t idx,AliTrackPoint &pt)const;             //             from AliTracker  
+         Int_t       Clusters2Tracks(AliESD *                   )       {return 0;} //pure virtual from AliTracker 
+         Int_t       LoadClusters   (TTree *pCluTr              );                  //pure virtual from AliTracker   
+  inline Int_t       PropagateBack  (AliESD *                   );                  //pure virtual from AliTracker   
+         Int_t       RefitInward    (AliESD *                   )       {return 0;} //pure virtual from AliTracker 
+         void        UnloadClusters (                           )       {         } //pure virtual from AliTracker 
 //private part  
   enum ETrackingFlags {kMipDistCut=-9,kMipQdcCut=-5};
+  static Int_t       Recon(AliESD *pEsd,TObjArray *pCluAll);                        //do actual job 
 protected:
   ClassDef(AliHMPIDTracker,0)
 };//class AliHMPIDTracker
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Int_t AliHMPIDTracker::PropagateBack(AliESD *pEsd)
+{
+// This method defined as pure virtual in AliTracker. It is invoked from AliReconstruction::RunTracking() after invocation of AliTracker::LoadClusters()
+// Agruments: pEsd - pointer to ESD
+//   Returns: error code    
+  AliHMPID *pHmpid=((AliHMPID*)gAlice->GetDetector("HMPID"));  
+  return Recon(pEsd,pHmpid->CluLst());  
+}
+
+
 
 typedef AliHMPIDTracker AliRICHTracker; // for backward compatibility