Cluster finding improved.
authordibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 20 Feb 2007 12:47:24 +0000 (12:47 +0000)
committerdibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 20 Feb 2007 12:47:24 +0000 (12:47 +0000)
New drawing utilities for event display.

12 files changed:
HMPID/AliHMPID.cxx
HMPID/AliHMPID.h
HMPID/AliHMPIDCluster.cxx
HMPID/AliHMPIDCluster.h
HMPID/AliHMPIDDigit.cxx
HMPID/AliHMPIDDigit.h
HMPID/AliHMPIDHit.cxx
HMPID/AliHMPIDHit.h
HMPID/AliHMPIDParam.cxx
HMPID/AliHMPIDTracker.cxx
HMPID/Hdisp.C
HMPID/Hmenu.C

index 11f25c0..56e0c32 100644 (file)
@@ -135,8 +135,10 @@ void AliHMPID::DigPrint(Int_t iEvt)const
   
   GetLoader()->TreeD()->GetEntry(0);
   DigLst()->Print();
+  Int_t totDigs=0;
+  for(Int_t i=0;i<7;i++) {totDigs+=DigLst(i)->GetEntries();}
+  Printf("totally %i Digits",totDigs);
   GetLoader()->UnloadDigits();
-  Printf("totally %i Digits",DigLst()->GetEntries());
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPID::CluPrint(Int_t iEvtN)const
index 9a71c5f..83ed19b 100644 (file)
@@ -20,6 +20,7 @@ public:
           void  BuildGeometry   (                ) {}          //from AliModule invoked from AliMC::InitGeometry() to build geometry for old event display
   virtual void  CreateMaterials (                )=0;          //from AliModule invoked from AliMC::ConstructGeometry() to define detector materials
   virtual void  CreateGeometry  (                )=0;          //from AliModule invoked from AliMC::ConstructGeometry() to build detector for simulation
+
   virtual Int_t IsVersion       (                )const=0;     //from AliModule not used        
   virtual void  Init            (                )=0;          //from AliModule invoked from AliMC::InitGeometry() after CreateGeometry() to do VolID initialization
   virtual void  DefineOpticalProperties() {}                   //from AliModule invoked from AliMC::ConstructOpGeometry() to set Cerenkov properties
index dacba95..8504030 100644 (file)
 #include "AliHMPIDCluster.h"  //class header
 #include <TMinuit.h>         //Solve()
 #include <TClonesArray.h>    //Solve()
-
+#include <TMarker.h>         //Draw()
 ClassImp(AliHMPIDCluster)
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDCluster::CoG()
 {
 // Calculates naive cluster position as a center of gravity of its digits.
 // Arguments: none 
-//   Returns: shape of the cluster i.e. the box which fully contains the cluster      
-  if(fDigs==0) return;                                  //no digits in this cluster
-  fX=fY=0;                                              //set cluster position to (0,0) to start to collect contributions
+//   Returns: none
+  
+//  if(fDigs==0) return;                                //no digits in this cluster
+  fX=fY=fQ=0;                                           //set cluster position to (0,0) to start to collect contributions
+  Int_t maxQpad=-1,maxQ=-1;                             //to calculate the pad with the highest charge
+  AliHMPIDDigit *pDig;
   for(Int_t iDig=0;iDig<fDigs->GetEntriesFast();iDig++){//digits loop
-    AliHMPIDDigit *pDig=(AliHMPIDDigit*)fDigs->At(iDig);  //get pointer to next digit
+    pDig=(AliHMPIDDigit*)fDigs->At(iDig);               //get pointer to next digit
     Float_t q=pDig->Q();                                //get QDC 
     fX += pDig->LorsX()*q;fY +=pDig->LorsY()*q;         //add digit center weighted by QDC
+    fQ+=q;                                              //increment total charge 
+    if(q>maxQ) {maxQpad = pDig->Pad();maxQ=(Int_t)q;}   // to find pad with highest charge
   }//digits loop
-  fX/=fQ;fY/=fQ;                                        //final center of gravity
-
-  CorrSin();
+  if ( fQ != 0 )   fX/=fQ;fY/=fQ;                       //final center of gravity
   
+  CorrSin();                                            //correct it by sinoid   
+  fCh=pDig->Ch();                                       //initialize chamber number
+  fMaxQpad = maxQpad; fMaxQ=maxQ;                       //store max charge pad to the field
+  fXi=fX+99; fYi=fY+99; fQi=fQ+99;                      //initial local max position is to be shifted artificially 
   fSt=kCoG;
 }//CoG()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -48,6 +56,11 @@ void AliHMPIDCluster::CorrSin()
   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;
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDCluster::Draw(Option_t*)
+{
+  TMarker *pMark=new TMarker(X(),Y(),5); pMark->SetMarkerColor(kBlue); pMark->Draw();
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDCluster::FitFunc(Int_t &iNpars, Double_t *, Double_t &chi2, Double_t *par, Int_t )
 {
 // Cluster fit function 
@@ -77,14 +90,23 @@ void AliHMPIDCluster::Print(Option_t* opt)const
 //Print current cluster  
   const char *status=0;
   switch(fSt){
-    case      kFor: status="formed"     ;break;
-    case      kUnf: status="unfolded"   ;break;
-    case      kCoG: status="coged"      ;break;
-    case      kEmp: status="empty"      ;break;
+    case        kFrm  : status="formed        "   ;break;
+    case        kUnf  : status="unfolded (fit)"   ;break;
+    case        kCoG  : status="coged         "   ;break;
+    case        kLo1  : status="locmax 1 (fit)"   ;break;
+    case        kAbn  : status="abnorm   (fit)"   ;break;
+    case        kMax  : status="exceeded (cog)"   ;break;
+    case        kNot  : status="not done (cog)"   ;break;
+    case        kEmp  : status="empty         "   ;break;
+    case        kEdg  : status="edge     (fit)"   ;break;
+    case       kSi1  : status="size 1   (cog)"   ;break;
+    case       kNoLoc: status="no LocMax(fit)"   ;break;
+    
+    default:            status="??????"          ;break;   
   }
-  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();    
+  Printf("%sCLU:(%7.3f,%7.3f) Q=%8.3f  ch=%i, FormedSize=%2i   N loc. max. %i Box %i  Chi2 %7.3f        %s",
+         opt,    X(),  Y(),   Q(),     Ch(),  Size(),           fNlocMax,       fBox,   fChi2,                 status);
+  if(fDigs) fDigs->Print();    
 }//Print()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
@@ -97,52 +119,119 @@ Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
 //Arguments: pCluLst     - cluster list pointer where to add new cluster(s)
 //           isTryUnfold - flag to switch on/off unfolding   
 //  Returns: number of local maxima of original cluster
-
+  CoG();
+  //  Printf("1 - fStatus: %d",fSt);
+  Int_t iCluCnt=pCluLst->GetEntriesFast();                                             //get current number of clusters already stored in the list by previous operations
+  if(isTryUnfold==kFALSE || Size()==1) {                                               //if cluster contains single pad there is no way to improve the knowledge 
+    (isTryUnfold)?fSt=kSi1:fSt=kNot;
+    new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);  //add this raw cluster 
+    return 1;
+  } 
+  // Printf("2 - fStatus: %d",fSt);
 //Phase 0. Initialise TMinuit  
   const Int_t kMaxLocMax=6;                                                            //max allowed number of loc max for fitting
   TMinuit *pMinuit = new TMinuit(3*kMaxLocMax);                                        //init MINUIT with this number of parameters (3 params per mathieson)
-  pMinuit->SetObjectFit((TObject*)this);  pMinuit->SetFCN(AliHMPIDCluster::FitFunc);    //set fit function
-  Double_t aArg=-1,parStart,parStep,parLow,parHigh;     Int_t iErrFlg;                 //tmp vars for TMinuit
+  pMinuit->SetObjectFit((TObject*)this);  pMinuit->SetFCN(AliHMPIDCluster::FitFunc);   //set fit function
+  Double_t aArg=-1;                                     Int_t iErrFlg;                 //tmp vars for TMinuit
   pMinuit->mnexcm("SET PRI",&aArg,1,iErrFlg);                                          //suspend all printout from TMinuit 
   pMinuit->mnexcm("SET NOW",&aArg,0,iErrFlg);                                          //suspend all warning printout from TMinuit
-//Phase 1. Find number of local maxima. Strategy is to check if the current pad has QDC more then all neigbours   
-  Int_t iLocMaxCnt=0;
-  for(Int_t iDig1=0;iDig1<Size();iDig1++) {                                             //first digits loop
+//Phase 1. Find number of local maxima. Strategy is to check if the current pad has QDC more then all neigbours. Also find the box contaning the cluster   
+  fNlocMax=0;
+  Int_t minPadX=999,minPadY=999,maxPadX=-1,maxPadY=-1,pc=-1;                             //for box finding   
+  //Double_t  lowX,highX,lowY,highY;
+  
+  //  Printf("3 - fStatus: %d",fSt);
+ for(Int_t iDig1=0;iDig1<Size();iDig1++) {                                              //first digits loop
     AliHMPIDDigit *pDig1 = Dig(iDig1);                                                   //take next digit
-    Int_t iHowManyMoreCnt = 0;                                                          //counts how many neighbouring pads has QDC more then current one
-    for(Int_t iDig2=0;iDig2<Size();iDig2++) {                                           //loop on all digits again
-      if(iDig1==iDig2) continue;                                                        //the same digit, no need to compare 
+    pc=pDig1->Pc();                                                                      //finding the box  
+    
+    if(pDig1->PadPcX() > maxPadX) maxPadX = pDig1->PadPcX();                              
+    if(pDig1->PadPcY() > maxPadY) maxPadY = pDig1->PadPcY();
+    if(pDig1->PadPcX() < minPadX) minPadX = pDig1->PadPcX();
+    if(pDig1->PadPcY() < minPadY) minPadY = pDig1->PadPcY();
+    
+    fBox=(maxPadX-minPadX+1)*100+maxPadY-minPadY+1;
+    
+    Int_t iHowManyMoreCnt = 0;                                                           //counts how many neighbouring pads has QDC more then current one
+    for(Int_t iDig2=0;iDig2<Size();iDig2++) {                                            //loop on all digits again
+      if(iDig1==iDig2) continue;                                                         //the same digit, no need to compare 
       AliHMPIDDigit *pDig2 = Dig(iDig2);                                                 //take second digit to compare with the first one
       Int_t dist = TMath::Sign(Int_t(pDig1->PadChX()-pDig2->PadChX()),1)+TMath::Sign(Int_t(pDig1->PadChY()-pDig2->PadChY()),1);//distance between pads
-      if(dist==1)                                                                       //means dig2 is a neighbour of dig1
-         if(pDig2->Q()>=pDig1->Q()) iHowManyMoreCnt++;                                 //count number of pads with Q more then Q of current pad
+      if(dist==1)                                                                        //means dig2 is a neighbour of dig1
+         if(pDig2->Q()>=pDig1->Q()) iHowManyMoreCnt++;                                   //count number of pads with Q more then Q of current pad
     }//second digits loop
-    if(iHowManyMoreCnt==0&&iLocMaxCnt<kMaxLocMax){                                     //this pad has Q more then any neighbour so it's local maximum
-        pMinuit->mnparm(3*iLocMaxCnt  ,Form("x%i",iLocMaxCnt),parStart=pDig1->LorsX(),parStep=0.01,parLow=0,parHigh=0,iErrFlg);
-        pMinuit->mnparm(3*iLocMaxCnt+1,Form("y%i",iLocMaxCnt),parStart=pDig1->LorsY(),parStep=0.01,parLow=0,parHigh=0,iErrFlg);
-        pMinuit->mnparm(3*iLocMaxCnt+2,Form("q%i",iLocMaxCnt),parStart=pDig1->Q()    ,parStep=0.01,parLow=0,parHigh=0,iErrFlg);
-        iLocMaxCnt++;
+    if(iHowManyMoreCnt==0&&fNlocMax<kMaxLocMax){                                       //this pad has Q more then any neighbour so it's local maximum
+      
+      /*
+      lowX  = AliHMPIDDigit::LorsX(pc,minPadX) - 0.5 *AliHMPIDDigit::SizePadX();
+      highX = AliHMPIDDigit::LorsX(pc,maxPadX) + 0.5 *AliHMPIDDigit::SizePadX();
+      lowY  = AliHMPIDDigit::LorsY(pc,minPadY) - 0.5 *AliHMPIDDigit::SizePadY();
+      highY = AliHMPIDDigit::LorsY(pc,maxPadY) + 0.5 *AliHMPIDDigit::SizePadY();
+      */
+      //Double_t    lowQ=0,highQ=30000; 
+      
+      fQi=pDig1->Q();  fXi=pDig1->LorsX();  fYi=pDig1->LorsY();                          //initial position of this Mathieson is to be in the center of loc max pad                               
+      /*
+        pMinuit->mnparm(3*fNlocMax  ,Form("x%i",fNlocMax),fXi,0.01,lowX,highX,iErrFlg);
+        pMinuit->mnparm(3*fNlocMax+1,Form("y%i",fNlocMax),fYi,0.01,lowY,highY,iErrFlg);
+        pMinuit->mnparm(3*fNlocMax+2,Form("q%i",fNlocMax),fQi,0.01,lowQ,highQ,iErrFlg);
+      */
+      pMinuit->mnparm(3*fNlocMax  ,Form("x%i",fNlocMax),fXi,0.01,0,0,iErrFlg);
+      pMinuit->mnparm(3*fNlocMax+1,Form("y%i",fNlocMax),fYi,0.01,0,0,iErrFlg);
+      pMinuit->mnparm(3*fNlocMax+2,Form("q%i",fNlocMax),fQi,0.01,0,100000,iErrFlg);
+      
+      fNlocMax++;
     }//if this pad is local maximum
   }//first digits loop
+  
+ //Int_t fitChk=0;
+
 //Phase 2. Fit loc max number of Mathiesons or add this current cluster to the list
-  Int_t iCluCnt=pCluLst->GetEntriesFast();                                          //get current number of clusters already stored in the list by previous operations
-  if(isTryUnfold==kTRUE && iLocMaxCnt<kMaxLocMax){                                        //resonable number of local maxima to fit and user requested it
-    pMinuit->mnexcm("MIGRAD" ,&aArg,0,iErrFlg);                                     //start fitting
-    if (!iErrFlg) { // Only if MIGRAD converged normally
-      Double_t fitX,fitY,fitQ,d1,d2,d3; TString sName;                                //vars to get results from TMinuit
-      for(Int_t i=0;i<iLocMaxCnt;i++){//local maxima loop
-       pMinuit->mnpout(3*i   ,sName,  fitX, d1 , d2, d3, iErrFlg);
-       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,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(Ch(),X(),Y(),Q(),kCoG);  //add this raw cluster 
-  }
+// Printf("4 - fStatus: %d",fSt);
+ if ( fNlocMax == 0) { // case of no local maxima found: pads with same charge...
+   pMinuit->mnparm(3*fNlocMax  ,Form("x%i",fNlocMax),fX,0.01,0,0,iErrFlg);
+   pMinuit->mnparm(3*fNlocMax+1,Form("y%i",fNlocMax),fY,0.01,0,0,iErrFlg);
+   pMinuit->mnparm(3*fNlocMax+2,Form("q%i",fNlocMax),fQ,0.01,0,100000,iErrFlg);
+   fNlocMax = 1;
+   fSt=kNoLoc;
+ }
+ if ( fNlocMax >= kMaxLocMax)
+   {
+     fSt   = kMax;   new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);               //add this raw cluster  
+   }
+ else{                                                                                 //resonable number of local maxima to fit and user requested it
+   Double_t arglist[10];     arglist[0] = 10000;     arglist[1] = 1.;                  //number of steps and sigma on pads charges  
+   pMinuit->mnexcm("MIGRAD" ,arglist,0,iErrFlg);                                       //start fitting
+   
+   if (iErrFlg) 
+     {
+       fSt   = kAbn;                                                                   //fit fails, MINUIT returns error flag
+       new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);                             //add this raw cluster 
+     }
+   else
+     {                                                                                 //Only if MIGRAD converged normally
+       Double_t d2,d3; TString sName;                                                  //vars to get results from TMinuit
+       for(Int_t i=0;i<fNlocMax;i++){                                                  //local maxima loop
+         pMinuit->mnpout(3*i   ,sName,  fX, fXe , d2, d3, iErrFlg);
+         pMinuit->mnpout(3*i+1 ,sName,  fY, fYe , d2, d3, iErrFlg);
+         pMinuit->mnpout(3*i+2 ,sName,  fQ, fQe , d2, d3, iErrFlg);
+         pMinuit->mnstat(fChi2,d2,d2,iErrFlg,iErrFlg,iErrFlg);
+         
+         if(fNlocMax!=1)fSt=kUnf;
+         if(fNlocMax==1&&fSt!=kNoLoc) fSt=kLo1;
+         if ( !IsInPc()) fSt = kEdg;       
+         if(fSt==kNoLoc) fNlocMax=0;
+         new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);      //add new unfolded cluster
+       }
+     }
+ }
+    
+
+  
+
   delete pMinuit;
-  return iLocMaxCnt;
+  return fNlocMax;
 }//Solve()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index 60fd0ac..bfbfcb8 100644 (file)
@@ -10,39 +10,66 @@ class TClonesArray;        //Solve()
 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,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;}
+  enum EClusterStatus {kFrm,kCoG,kLo1,kUnf,kMax,kAbn,kNot,kEdg,kSi1,kNoLoc,kEmp=-1};      //status flags    
+  AliHMPIDCluster():TObject( ),fCh(-1),fSi(-1),fSt(kEmp),fBox(-1),fNlocMax(-1),fMaxQpad(-1),fMaxQ(-1),fQ(-1),fQi(-1),fQe(-1),fX(-1),fXi(-1),fXe(-1),fY(-1),fYi(-1),fYe(-1),fChi2(-1),fDigs(0) {} //empty ctor
+  
+  
+  AliHMPIDCluster           (const AliHMPIDCluster &c):TObject(c),fCh(c.fCh),fSi(c.fSi),fSt(c.fSt),fBox(c.fBox),fNlocMax(c.fNlocMax),fMaxQpad(c.fMaxQpad),fMaxQ(c.fMaxQ),
+                                                                  fQ (c.fQ ),fQi(c.fQi),fQe(c.fQe),
+                                                                  fX (c.fX ),fXi(c.fXi),fXe(c.fXe),
+                                                                  fY (c.fY ),fYi(c.fYi),fYe(c.fYe),fChi2(c.fChi2),fDigs(0)                  {}//copy ctor
+  AliHMPIDCluster &operator=(const AliHMPIDCluster &c) {if(this == &c)return *this;TObject::operator=(c);          
+                                                        fSi=c.fSi;  fSt=c.fSt; fCh=c.fCh; fBox=c.fBox;fNlocMax=c.fNlocMax;fMaxQpad=c.fMaxQpad; fMaxQ=c.fMaxQ;
+                                                        fQ=c.fQ; fQi=c.fQi;fQe=c.fQe; 
+                                                        fX=c.fX; fXi=c.fXi;fXe=c.fXe;
+                                                        fY=c.fY; fYi=c.fYi;fYe=c.fYe; fChi2=c.fChi2;fDigs=c.fDigs ? new TObjArray(*c.fDigs):0; return *this;}
     
-  virtual ~AliHMPIDCluster(                                     )                                                                        {DigDel();}
+  virtual ~AliHMPIDCluster(                                     )                                                                        {if(fDigs) delete fDigs; fDigs=0;}
 //framework part                   
-         void          Print  (Option_t *opt=""                                  )const;                  //overloaded TObject::Print() to print cluster info
-  static void          FitFunc(Int_t &, Double_t *, Double_t &, Double_t *, Int_t);                       //fit function to be used by MINUIT
+         void           Draw   (Option_t *opt=""                                  );                       //overloaded TObject::Print() to draw cluster in current canvas
+         void           Print  (Option_t *opt=""                                  )const;                  //overloaded TObject::Print() to print cluster info
+  static void           FitFunc(Int_t &, Double_t *, Double_t &, Double_t *, Int_t);                       //fit function to be used by MINUIT
 //private part  
-         void          CoG      (                                         );                                                      //calculates center of gravity
-         void          CorrSin  (                                         );                                                      //sinoidal correction   
-         Int_t         Ch       (                                         )const{return fCh;                                    } //chamber number
-  inline void          DigAdd   (AliHMPIDDigit *pDig                       );                                                      //add new digit ot the cluster
-         void          DigDel   (                                         )     {if(fDigs) {delete fDigs;fDigs=0;}              } //deletes the list of digits (not digits!) 
-         AliHMPIDDigit* Dig      (Int_t i                                  )const{return (AliHMPIDDigit*)fDigs->At(i);            } //pointer to i-th digit
-         TObjArray*    DigLst   (                                         )const{return fDigs;                                  } //list of digits  
-         void          Reset    (                                         )     {DigDel();fQ=fCh=-1;fX=fY=-1;fSt=kEmp;          } //cleans the cluster
-         Int_t         Solve    (TClonesArray *pCluLst,Bool_t isUnfold    );                                                      //solve cluster: MINUIT fit or CoG
-         Int_t         Size     (                                         )const{return (fDigs)?fDigs->GetEntriesFast():0;      } //number of pads in cluster
-         Int_t         Q        (                                         )const{return fQ;                                     } //cluster charge in QDC channels 
-         Float_t       X        (                                         )const{return fX;                                     } //cluster x position in LRS
-         Float_t       Y        (                                         )const{return fY;                                     } //cluster y position in LRS 
+         void           CoG      (                                         );                                                      //calculates center of gravity
+         void           CorrSin  (                                         );                                                      //sinoidal correction   
+         Int_t          Ch       (                                         )const{return fCh;                                    } //chamber number
+  inline void           DigAdd   (AliHMPIDDigit *pDig                      );                                                      //add new digit ot the cluster
+         AliHMPIDDigit* Dig      (Int_t i                                  )const{return (AliHMPIDDigit*)fDigs->At(i);           } //pointer to i-th digi 
+  inline Bool_t         IsInPc   ();                                                                                        //check if is in the current PC
+  inline void           Reset    (                                         );                                                      //cleans the cluster
+         Int_t          Size     (                                         )const{return fSi;                                    } //returns number of pads in formed cluster 
+         Int_t          Solve    (TClonesArray *pCluLst,Bool_t isUnfold    );                                                      //solve cluster: MINUIT fit or CoG
+         Int_t          Status  (                                          ) const{return fSt;}                                   //Status of cluster                                  
+         Double_t       Q        (                                         )const{return fQ;                                     } //cluster charge in QDC channels 
+         Double_t       Qi       (                                         )const{return fQi;                                    } //cluster charge in QDC channels 
+         Double_t       Qe       (                                         )const{return fQe;                                    } //cluster charge in QDC channels 
+         Double_t       X        (                                         )const{return fX;                                     } //cluster x position in LRS
+         Double_t       Xi       (                                         )const{return fXi;                                    } //cluster charge in QDC channels 
+         Double_t       Xe       (                                         )const{return fXe;                                    } //cluster charge in QDC channels 
+         Double_t       Y        (                                         )const{return fY;                                     } //cluster y position in LRS 
+         Double_t       Yi       (                                         )const{return fYi;                                    } //cluster charge in QDC channels 
+         Double_t       Ye       (                                         )const{return fYe;                                    } //cluster charge in QDC channels 
+         Double_t       Chi2     (                                         )const{return fChi2;                                  } 
 protected:
-  Int_t         fSt;          //flag to mark the quality of the cluster   
   Int_t         fCh;          //chamber number
-  Int_t         fQ;           //QDC value
-  Float_t       fX;           //local x postion, [cm] 
-  Float_t       fY;           //local y postion, [cm]  
+  Int_t         fSi;          //size of the formed cluster from which this cluster deduced  
+  Int_t         fSt;          //flag to mark the quality of the cluster   
+  Int_t         fBox;         //box contaning this cluster  
+  Int_t         fNlocMax;     //number of local maxima in formed cluster        
+  Int_t         fMaxQpad;     //abs pad number of a pad with the highest charge
+  Double_t      fMaxQ;        //that max charge value             
+  Double_t      fQ;           //QDC value
+  Double_t      fQi;          //initial Q of local maximum
+  Double_t      fQe;          //error on Q 
+  Double_t      fX;           //local x postion, [cm] 
+  Double_t      fXi;          //initial x of local maximum, [cm] 
+  Double_t      fXe;          //error on x postion, [cm] 
+  Double_t      fY;           //local y postion, [cm]  
+  Double_t      fYi;          //initial y of local maximum, [cm] 
+  Double_t      fYe;          //error on y postion, [cm] 
+  Double_t      fChi2;        //some estimatore of fit quality  
   TObjArray    *fDigs;        //! list of digits forming this cluster
-  ClassDef(AliHMPIDCluster,5)  //HMPID cluster class       
+  ClassDef(AliHMPIDCluster,5) //HMPID cluster class       
 };//class AliHMPIDCluster
 
 typedef AliHMPIDCluster AliRICHCluster; // for backward compatibility
@@ -53,11 +80,29 @@ void AliHMPIDCluster::DigAdd(AliHMPIDDigit *pDig)
 // Adds a given digit to the list of digits belonging to this cluster, cluster is not owner of digits
 // Arguments: pDig - pointer to digit to be added  
 //   Returns: none  
-  if(!fDigs) {fQ=0;fDigs = new TObjArray;}
+  if(!fDigs) {fSi=0;fDigs = new TObjArray;} //create list of digits in the first invocation
   fDigs->Add(pDig);
-  fQ+=(Int_t)pDig->Q(); 
-  fCh=pDig->Ch();
-  fSt=kFor;
+  fSt=kFrm;
+  fSi++;
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDCluster::Reset()
+{
+  if(fDigs) delete fDigs; 
+  fDigs=0; 
+  fQ=fCh=fSi=-1;fX=fY=-1;fSt=kEmp;      
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Bool_t AliHMPIDCluster::IsInPc()
+{
+  Int_t pc = ((AliHMPIDDigit*)fDigs->At(0))->Pc();
+  
+  if ( fX < AliHMPIDDigit::fMinPcX[pc] || fX > AliHMPIDDigit::fMaxPcX[pc] || 
+       fY < AliHMPIDDigit::fMinPcY[pc] || fY > AliHMPIDDigit::fMaxPcY[pc] ) return kFALSE;
+  
+  return kTRUE;
+  
+}
+
 #endif
index 18ddfe7..3b14f0f 100644 (file)
 //  **************************************************************************
 
 #include "AliHMPIDDigit.h" //class header
-#include "AliHMPIDHit.h"   //Hit2Sdi() 
-#include <TClonesArray.h> //Hit2Sdi() 
-#include <TCanvas.h>      //TestSeg()
-#include <TLine.h>        //TestSeg()
-#include <TLatex.h>       //TestSeg()
-#include <TPolyLine.h>    //DrawPc() 
-#include <TGStatusBar.h>  //Zoom()
-#include <TRootCanvas.h>  //Zoom() 
+#include <TClonesArray.h>  //Hit2Sdi() 
+#include <TMarker.h>       //Draw() 
 ClassImp(AliHMPIDDigit)
 
+const Float_t AliHMPIDDigit::fMinPcX[]={0,SizePcX() + SizeDead(),0,SizePcX() + SizeDead(),       0,SizePcX() + SizeDead()};
+const Float_t AliHMPIDDigit::fMinPcY[]={0,                     0,SizePcY()+SizeDead(),SizePcY() + SizeDead(),2*(SizePcY()+SizeDead()),2*(SizePcY()+SizeDead())};
+
+const Float_t AliHMPIDDigit::fMaxPcX[]={SizePcX(),SizeAllX(),SizePcX(),SizeAllX(),SizePcX(),SizeAllX()};
+const Float_t AliHMPIDDigit::fMaxPcY[]={SizePcY(),SizePcY(),SizeAllY() - SizePcY(),SizeAllY() - SizePcY(),SizeAllY(),SizeAllY()};
+
+
 /*
 Preface: all geometrical information (like left-right sides) is reported as seen from electronic side.
 
@@ -72,6 +73,11 @@ 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)
 */
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDDigit::Draw(Option_t*)
+{
+  TMarker *pMark=new TMarker(LorsX(),LorsY(),25); pMark->SetMarkerColor(kGreen); pMark->Draw();
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Hit2Sdi(AliHMPIDHit *pHit,TClonesArray *pSdiLst)
 {
 // Creates a list of sdigits out of provided hit
@@ -93,8 +99,10 @@ 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) (%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");
+  Printf("DIG:(%7.3f,%7.3f) Q=%8.3f (ch=%1i,pc=%1i,x=%2i,y=%2i)  TID=(%5i,%5i,%5i) ddl=%i raw=0x%x (r=%2i,d=%2i,a=%2i) %s",
+              LorsX(),LorsY(),Q(), A2C(fPad),A2P(fPad),A2X(fPad),A2Y(fPad),   
+              fTracks[0],fTracks[1],fTracks[2],DdlIdx(),w32,Row(),Dilogic(),Addr(), 
+              (IsOverTh(Q()))?"":"!!!");
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::PrintSize()
@@ -110,105 +118,6 @@ void AliHMPIDDigit::PrintSize()
                SizeAllX(),SizeAllY(),kPadAllX,kPadAllY);
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDDigit::DrawSeg()
-{
-// Draws the picture of segmentation   
-// Arguments: none
-//   Returns: none     
-  TCanvas *pC=new TCanvas("pads","View from electronics side, IP is behind the picture.",1000,900);pC->ToggleEventStatus();
-  gPad->AddExec("test","AliHMPIDDigit::DrawZoom()");
-  DrawPc();  
-}//TestSeg()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDDigit::DrawZoom()
-{
-// Show info about current cursur position in status bar of the canvas
-// Arguments: none
-//   Returns: none     
-  TCanvas *pC=(TCanvas*)gPad; 
-  TRootCanvas *pRC= (TRootCanvas*)pC->GetCanvasImp();
-  TGStatusBar *pBar=pRC->GetStatusBar();
-  pBar->SetParts(5);
-  Float_t x=gPad->AbsPixeltoX(gPad->GetEventX());
-  Float_t y=gPad->AbsPixeltoY(gPad->GetEventY());
-  AliHMPIDDigit dig;dig.Manual1(1,x,y); UInt_t w32=0; 
-  if(IsInDead(x,y))
-    pBar->SetText("Out of sensitive area",4);    
-  else{
-    Int_t ddl=dig.Raw(w32);
-    pBar->SetText(Form("(p%i,x%i,y%i) ddl=%i 0x%x (r%i,d%i,a%i) (%.2f,%.2f)",
-        dig.Pc(),dig.PadPcX(),dig.PadPcY(),
-        ddl,w32,
-        dig.Row(),dig.Dilogic(),dig.Addr(),
-        dig.LorsX(),dig.LorsY()                            ),4);
-  }    
-  if(gPad->GetEvent()==1){
-    new TCanvas("zoom",Form("Row %i DILOGIC %i",dig.Row(),dig.Dilogic()));  
-  }
-}//Zoom()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void AliHMPIDDigit::DrawPc(Bool_t isFill) 
-{ 
-// Utility methode draws HMPID chamber PCs on event display.
-// Arguments: none
-//   Returns: none      
-//  y6  ----------  ----------
-//      |        |  |        |
-//      |    4   |  |    5   |
-//  y5  ----------  ----------
-//
-//  y4  ----------  ----------
-//      |        |  |        |
-//      |    2   |  |    3   |   view from electronics side
-//  y3  ----------  ----------
-//          
-//  y2  ----------  ----------
-//      |        |  |        |
-//      |    0   |  |    1   |
-//  y1  ----------  ----------
-//      x1      x2  x3       x4
-  gPad->Range(-5,-5,SizeAllX()+5,SizeAllY()+5); 
-  Float_t x1=0,x2=SizePcX(),x3=SizePcX()+SizeDead(),                                                  x4=SizeAllX();
-  Float_t y1=0,y2=SizePcY(),y3=SizePcY()+SizeDead(),y4=2*SizePcY()+SizeDead(),y5=SizeAllY()-SizePcY(),y6=SizeAllY();
-
-  Float_t xL[5]={x1,x1,x2,x2,x1}; //clockwise
-  Float_t xR[5]={x3,x3,x4,x4,x3};  
-  Float_t yD[5]={y1,y2,y2,y1,y1};
-  Float_t yC[5]={y3,y4,y4,y3,y3};  
-  Float_t yU[5]={y5,y6,y6,y5,y5};
-    
-  TLatex txt; txt.SetTextSize(0.01);
-  Int_t iColLeft=29,iColRight=41;
-  TPolyLine *pc=0;  TLine *pL;
-  AliHMPIDDigit dig;
-  for(Int_t iPc=0;iPc<kPcAll;iPc++){
-    if(iPc==4) pc=new TPolyLine(5,xL,yU); if(iPc==5) pc=new TPolyLine(5,xR,yU); //draw PCs
-    if(iPc==2) pc=new TPolyLine(5,xL,yC); if(iPc==3) pc=new TPolyLine(5,xR,yC);
-    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.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.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()); 
-      if(iRow!=0) pL->Draw(); 
-    }//row loop  
-    
-    txt.SetTextAlign(13);
-    for(Int_t iDil=0;iDil<10;iDil++){//draw dilogic lines (vertical)
-      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()); 
-      if(iDil!=0)pL->Draw();
-    }//dilogic loop        
-  }//PC loop      
-}//DrawPc()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Test()
 {
   AliHMPIDDigit d1,d2; Int_t ddl;UInt_t w32;
index 3f55a75..df44e6e 100644 (file)
@@ -7,7 +7,7 @@
 #include <TMath.h>         //Mathieson()
 #include <TRandom.h>       //IsOverTh()  
 #include <AliBitPacking.h> //Raw()
-#include "AliHMPIDHit.h"   //Hit2Sdi(), ctor()
+#include "AliHMPIDHit.h"   //Set()
 
 class TClonesArray;        //Hit2Sdi()
   
@@ -18,12 +18,13 @@ 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 )  {}                                      //ctor used in digitizer
-  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()
+         void    Draw        (Option_t *opt=""               );                                                                        //TObject::Draw() overloaded
          void    Print       (Option_t *opt=""               )const;                                                                   //TObject::Print() overloaded
 //private part  
   static Int_t   Abs         (Int_t c,Int_t s,Int_t x,Int_t y)     {return c*kChAbs+s*kPcAbs+x*kPadAbsX+y*kPadAbsY; }                  //(ch,pc,padx,pady)-> abs pad
@@ -35,19 +36,18 @@ public:
          void    AddTidOffset(Int_t offset                   )     {for (Int_t i=0; i<3; i++) if (fTracks[i]>0) fTracks[i]+=offset;};  //needed for merging
          Int_t   Ch          (                               )const{return A2C(fPad);                               }                  //chamber number
          Int_t   Dilogic     (                               )const{return 1+PadPcX()/8;                            }                  //DILOGIC# 1..10
-  static void    DrawPc      (Bool_t isFill=kTRUE            );                                                                        //draw PCs
-  static void    DrawSeg     (                               );                                                                        //draw segmentation
-         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 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?
+  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 chamber 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]
+  static Float_t LorsX       (Int_t pc,Int_t padx            )     {return (padx    +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 
+  static Float_t LorsY       (Int_t pc,Int_t pady            )     {return (pady    +0.5)*SizePadY()+(pc  /2)*(SizePcY()+SizeDead());} //center of the pad y, [cm]
+         void    Manual1     (Int_t c,Float_t x,Float_t y    )     {AliHMPIDHit h(c,200e-9,2212,3,x,y); Set(&h,0);}                    //manual from hit
+         void    Manual2     (Int_t c,Int_t p,Int_t x,Int_t y,Float_t q=0)     {fPad=Abs(c,p,x,y);fQ=q;}                               //manual creation 
   inline Float_t Mathieson   (Float_t x,Float_t y            )const;                                                                   //Mathieson distribution 
          Int_t   PadPcX      (                               )const{return A2X(fPad);}                                                 //pad pc x # 0..79
          Int_t   PadPcY      (                               )const{return A2Y(fPad);}                                                 //pad pc y # 0..47
@@ -78,6 +78,12 @@ public:
   static Float_t SizeWin     (                               )     {return 0.5;}                                                       //Quartz window width
   static Float_t SizeRad     (                               )     {return 1.5;}                                                       //Rad width   
   static void    Test        (                               );                                                                        //Test conversions
+  static const Float_t fMinPcX[6];
+  static const Float_t fMinPcY[6];
+  static const Float_t fMaxPcX[6];
+  static const Float_t fMaxPcY[6];
+  
+  
 protected:                  //AliDigit has fTracks[3]
   Int_t    fPad;            //absolute pad number
   Float_t  fQ;              //QDC value, fractions are permitted for summable procedure  
index 68ad3ff..6b289fa 100644 (file)
 #include "AliHMPIDHit.h"  //class header
 #include "AliHMPIDDigit.h"
 #include <TPDGCode.h>    
+#include <TMarker.h>
  
 ClassImp(AliHMPIDHit)
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDHit::Draw(Option_t*)
+{
+  Int_t iMark;
+  switch(Pid()){
+    case 50000050:   iMark=4;  break;
+    case 50000051:   iMark=27; break;
+    default:         iMark=26; break;
+  }    
+  TMarker *pMark=new TMarker(LorsX(),LorsY(),iMark); pMark->SetMarkerColor(kRed); pMark->Draw();
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDHit::Print(Option_t*)const
 {
 //Print hit
@@ -38,8 +50,8 @@ void AliHMPIDHit::Print(Option_t*)const
     case 50000051:     sPart="feed";break;
   }
 
-  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, 
+  Printf("HIT:(%7.3f,%7.3f) Q=%8.3f ch=%i                   TID= %5i, MARS=(%7.2f,%7.2f,%7.2f) %s  %s",
+              LorsX(),LorsY(), Q(), Ch(),                  Tid(),         X(),  Y(),  Z(),   sPart, 
                         (AliHMPIDDigit::IsInDead(LorsX(),LorsY()))? "IN DEAD ZONE":"");
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index dd30365..645249b 100644 (file)
 class AliHMPIDHit : public AliHit //   TObject-AliHit-AliHMPIDHit
 {
 public:
-  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()                                                                                                {}
+  AliHMPIDHit(                                                                         ):AliHit(    ),fCh(-1),fPid(-1),fQ(-1),fLx(0),fLy(0){} //default ctor
+  AliHMPIDHit(Int_t c,Float_t e,Int_t id,Int_t tn,Float_t x,Float_t y,const TVector3 &p):AliHit(0,tn),fCh(c ),fPid(id),fQ(-1),fLx(x),fLy(y){QdcTot(e);fX=p.X();fY=p.Y();fZ=p.Z();}
+  AliHMPIDHit(Int_t c,Float_t e,Int_t id,Int_t tn,Float_t x,Float_t y                  ):             fCh(c ),fPid(id),fQ(-1),fLx(x),fLy(y){QdcTot(e);fTrack=tn;}//manual ctor 
+  virtual ~AliHMPIDHit()                                                                                                                                 {}
 //framework part
-  void     Print(Option_t *option="")const;                                    //from TObject to print current status
+         void    Print(Option_t *opt="")const;                                                 //from TObject to print current status
+         void    Draw (Option_t *opt="");                                                      //from TObject to Draw this hit in current canvas
 //private part  
          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]
+         //clm: hit LorsX and Y were changed for simulation hits
+         Float_t LorsX (         )const{return fLx;                                    }       //hit X position in LORS, [cm]
+         Float_t LorsY (         )const{return fLy;                                    }       //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          
@@ -30,8 +32,8 @@ protected:                                                                     /
   Int_t    fCh;                                                                //Chamber
   Int_t    fPid;                                                               //PID
   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]
+  Float_t  fLx;                                                                //hit X position in chamber LORS, [cm]
+  Float_t  fLy;                                                                //hit Y position in chamber LORS, [cm]
   ClassDef(AliHMPIDHit,5)                                                      //HMPID hit class 
 };//class AliHMPIDhit
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++  
@@ -47,9 +49,8 @@ Float_t AliHMPIDHit::QdcTot(Float_t e)
   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 c1ab264..1fa574c 100644 (file)
@@ -48,7 +48,7 @@ void AliHMPIDParam::Print(Option_t* opt) const
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 Int_t AliHMPIDParam::Stack(Int_t evt,Int_t tid)
 {
-// Prints some usefull info from stack
+// Prints some useful info from stack
 // Arguments: evt - event number. if not -1 print info only for that event
 //            tid - track id. if not -1 then print it and all it's mothers if any   
 //   Returns: mother tid of the given tid if any
index 4a14ac1..064eebc 100644 (file)
@@ -88,13 +88,13 @@ Int_t AliHMPIDTracker::Recon(AliESD *pESD,TObjArray *pCluAll)
     
     AliHMPIDCluster *pMipClu=(AliHMPIDCluster*)pCluLst->At(iMip);                            //take mip cluster 
     
-                   pTrk->SetHMPIDmip      (pMipClu->X(),pMipClu->Y(),pMipClu->Q());          //store mip info 
+                   pTrk->SetHMPIDmip      ((Float_t)pMipClu->X(),(Float_t)pMipClu->Y(),(Int_t)pMipClu->Q());          //store mip info 
     if(dMin>1)    {pTrk->SetHMPIDsignal   (kMipDistCut); continue;}                          //closest cluster with enough charge is still too far 
                    pTrk->SetHMPIDcluIdx   (iCh,iMip);                                        //set mip cluster index
   recon.SetTrack(th,ph,xRa,yRa); Int_t iNphot=0;                                            //initialize track parameters  
                    pTrk->SetHMPIDsignal   (recon.CkovAngle(pCluLst,iNphot));                 //search for Cerenkov angle for this track
                    pTrk->SetHMPIDchi2     (recon.CkovSigma2());                              //error squared 
-                   pTrk->SetHMPIDmip      (pMipClu->X(),pMipClu->Y(),pMipClu->Q(),iMip);     //info on mip cluster + n. phot.
+                   pTrk->SetHMPIDmip      ((Float_t)pMipClu->X(),(Float_t)pMipClu->Y(),(Int_t)pMipClu->Q(),iMip);     //info on mip cluster + n. phot.
  }//ESD tracks loop
   AliDebugClass(1,"Stop pattern recognition");
   return 0; // error code: 0=no error;
index 3470a99..03a7394 100644 (file)
@@ -1,4 +1,4 @@
-red()
+Hdisp()
 {
   gSystem->Load("libMinuit.so");
   gSystem->Load("libVMC.so");
index 1713c2e..9e6c56c 100644 (file)
@@ -10,7 +10,7 @@ void GetParam()
   if(isGeomType) g=TGeoManager::Import("geometry.root");
   else           g=TGeoManager::Import("misaligned_geometry.root");
   rp=AliHMPIDParam::Instance();
-}
+}//GetParam()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void Hmenu()
 {   
@@ -30,45 +30,44 @@ void Hmenu()
   
   TControlBar *pMenu = new TControlBar("horizontal",status.Data(),0,0);
     pMenu->AddButton("                     ","","");
-    pMenu->AddButton(" General  "           ,"General()"  ,"general items which do not depend on any files");
+    pMenu->AddButton("       General       ","General()"  ,"general items which do not depend on any files");
     pMenu->AddButton("                     ",""           ,"");
-    pMenu->AddButton(" Sim data "           ,"SimData()"  ,"items which expect to have simulated files"    );
+    pMenu->AddButton("       Sim data      ","SimData()"  ,"items which expect to have simulated files"    );
     pMenu->AddButton("                     ",""           ,"");
-    pMenu->AddButton(" Raw data "           ,"RawData()"  ,"items which expect to have raw files"          );
+    pMenu->AddButton("       Raw data      ","RawData()"  ,"items which expect to have raw files"          );
     pMenu->AddButton("                     ","print()"    ,"");
-    pMenu->AddButton("Test"                 ,"Test()"     ,"all test utilities");
+    pMenu->AddButton("         Test        ","Test()"     ,"all test utilities");
     pMenu->AddButton("                     ","GetParam()" ,"");
-    pMenu->AddButton("Quit"                 ,".q"         ,"close session"                                 );
+    pMenu->AddButton("         Quit        ",".q"         ,"close session"                                 );
   pMenu->Show();
-}
+}//Menu()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void General()
 {         
-  TControlBar *pMenu = new TControlBar("vertical","Sim data",60,50);  
-    pMenu->AddButton("Debug ON"           ,"don();"                    ,"Switch debug on-off"                        );   
-    pMenu->AddButton("Debug OFF"          ,"doff();"                   ,"Switch debug on-off"                        );   
-    pMenu->AddButton("Geo GUI"            ,"geo();"                    ,"Shows geometry"                             ); 
-    pMenu->AddButton("Browser"            ,"new TBrowser;"             ,"Start ROOT TBrowser"                        );
+  TControlBar *pMenu = new TControlBar("vertical","General purpose",100,50);  
+    pMenu->AddButton("       Debug ON      ","don();"                    ,"Switch debug on-off"                        );   
+    pMenu->AddButton("       Debug OFF     ","doff();"                   ,"Switch debug on-off"                        );   
+    pMenu->AddButton("        Geo GUI      ","geo();"                    ,"Shows geometry"                             ); 
+    pMenu->AddButton("        Browser      ","new TBrowser;"             ,"Start ROOT TBrowser"                        );
   pMenu->Show();  
-}//menu()
+}//General()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void SimData()
 {
-  TControlBar *pMenu = new TControlBar("vertical","Sim",190,50);  
-    pMenu->AddButton("Display ALL chambers"            ,"ed();"     , "Display Fast");
-    pMenu->AddButton("HITS QA"                         ,"hqa()"     ,"QA plots for hits: hqa()");
-    
-    pMenu->AddButton("Print hits"                      ,"hp();"      ,"To print hits:     hp(evt)");
-    pMenu->AddButton("Print sdigits"                   ,"sp();"      ,"To print sdigits:  sp(evt)");
-    pMenu->AddButton("Print digits"                    ,"dp();"      ,"To print digits:   dp(evt)");
-    pMenu->AddButton("Print clusters"                  ,"cp();"      ,"To print clusters: cp(evt)");
-    
+  TControlBar *pMenu = new TControlBar("vertical","Sim data",340,50);  
+    pMenu->AddButton("Display ALL chambers ","ed();"      ,"Display Fast");
+    pMenu->AddButton("      HITS QA        ","hqa()"      ,"QA plots for hits: hqa()");
+    pMenu->AddButton("     Print Stack     ","stack();"   ,"To print hits:     hp(evt)");
+    pMenu->AddButton("     Print hits      ","hp();"      ,"To print hits:     hp(evt)");
+    pMenu->AddButton("    Print sdigits    ","sp();"      ,"To print sdigits:  sp(evt)");
+    pMenu->AddButton("    Print digits     ","dp();"      ,"To print digits:   dp(evt)");
+    pMenu->AddButton("   Print clusters    ","cp();"      ,"To print clusters: cp(evt)");
   pMenu->Show();         
-}
+}//SimData()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void RawData()
 {
-  TControlBar *pMenu = new TControlBar("vertical","Raw",350,50);  
+  TControlBar *pMenu = new TControlBar("vertical","Raw data",580,50);  
     pMenu->AddButton("ESD print"                       ,"ep();"                  ,"To print ESD info: ep()"         );  
     pMenu->AddButton("ESD QA"                          ,"eq();"                  ,"To draw ESD hists: eq()"         );  
     pMenu->AddButton("Clusters print"                  ,"cp();"                  ,"To print clusters: cp()"         );  
@@ -77,11 +76,11 @@ void RawData()
     pMenu->AddButton("Print occupancy"                 ,"r->OccupancyPrint(-1);" ,"To print occupancy"              );  
     pMenu->AddButton("Print event summary  "           ,"r->SummaryOfEvent();"   ,"To print a summary of the event" );  
   pMenu->Show();         
-}//RawData
+}//RawData()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void Test()
 {         
-  TControlBar *pMenu = new TControlBar("vertical","Test",400,50);  
+  TControlBar *pMenu = new TControlBar("vertical","Test",820,50);  
     pMenu->AddButton("Hits->Digits"       ,"thd();"                    ,"test hits->sdigits->digits"                 );   
     pMenu->AddButton("Segmentation"       ,"ts()"                      ,"test segmentation methods"                  );
     pMenu->AddButton("Test response"      ,"AliHMPIDParam::TestResp();","Test AliHMPIDParam response methods"         );
@@ -113,59 +112,83 @@ void eq  (                       ) {AliHMPIDTracker::EsdQA();             }
 void mp  (Double_t probCut=0.7   ) {AliHMPIDTracker::MatrixPrint(probCut);}                   
 
 
-void t   (Int_t evt=0          )   {AliHMPIDParam::Stack(evt);}    
-void tid (Int_t tid,Int_t evt=0)   {AliHMPIDParam::Stack(evt,tid);} 
+void stack(                     )   {AliHMPIDParam::Stack();}    
+void tid  (Int_t tid,Int_t evt=0)   {AliHMPIDParam::Stack(evt,tid);} 
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void tst()
 {
-//create all lists  
-  TClonesArray hit("AliHMPIDHit"),sdi("AliHMPIDDigit");
-  TObjArray    dig,clu,cog;      for(Int_t i=0;i<7;i++) {dig.AddAt(new TClonesArray("AliHMPIDDigit"),i); 
-                                                         clu.AddAt(new TClonesArray("AliHMPIDCluster"),i);
-                                                         cog.AddAt(new TClonesArray("AliHMPIDCluster"),i);}
-//simulate track  
-  TLorentzVector mom; mom.SetPtEtaPhiM(3,0,30*TMath::DegToRad(),0.938);                                                                    
-  TLorentzVector vtx; vtx.SetXYZT(0,0,0,0);                                                                 
-
-  AliESD *pESD=new AliESD;  pESD->SetMagneticField(0.2);
-  pESD->AddTrack(new AliESDtrack(new TParticle(kProton,0,-1,-1,-1,-1,mom,vtx)));
-  pESD->Print();
-  
-  AliTracker::SetFieldMap(new AliMagF,1);
-  AliHMPIDTracker *pTracker=new AliHMPIDTracker;
+  TCanvas *pC=new TCanvas("pads","View from electronics side, IP is behind the picture.",1000,900);  pC->ToggleEventStatus();
+  SimEvt();
+}//tst()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void SimEvt()
+{
+//sample track         
+  Int_t ch=1;
+  Float_t th=8*TMath::DegToRad();                     
+  Float_t ph=gRandom->Rndm()*TMath::TwoPi(); 
+  Float_t radx=gRandom->Rndm()*AliHMPIDDigit::SizeAllX(); 
+  Float_t rady=gRandom->Rndm()*AliHMPIDDigit::SizeAllY();
+  Float_t simCkov=0.63,simErr=0; Int_t simN=8,recN=0;
+  
+  AliHMPIDRecon rec; rec.SetTrack(th,ph,radx,rady);                                                                
+//sample list of hits  
+  TClonesArray hl("AliHMPIDHit"); Int_t hc=0;  
+  TVector2 pc;//tmp position
+  
+  Int_t kCerenkov=50000050,kFeedback=50000051;
+  Int_t tn=0;
+                                                                                   new(hl[hc++]) AliHMPIDHit(ch,200e-9,kProton  ,tn++,radx  ,rady    );   //mip hit
+  for(int i=0;i<simN;i++){rec.TracePhot(simCkov,gRandom->Rndm()*TMath::TwoPi(),pc);new(hl[hc++]) AliHMPIDHit(ch,7.5e-9,kCerenkov,tn++,pc.X(),pc.Y());}  //photon hits
+  for(int i=0;i<10;i++)  { Float_t x=gRandom->Rndm()*130,y=gRandom->Rndm()*126;    new(hl[hc++]) AliHMPIDHit(ch,7.5e-9,kFeedback,tn++,   x  ,   y  );}  //bkg hits
+//do reconstruction
+                                                        TClonesArray sl("AliHMPIDDigit");              AliHMPIDv1::Hit2Sdi(&hl,&sl);                               
+  TObjArray dl(7);  for(Int_t i=0;i<7;i++) dl.AddAt(new TClonesArray("AliHMPIDDigit"),i);       AliHMPIDDigitizer::Sdi2Dig(&sl,&dl);     
+  TObjArray cl(7);  for(Int_t i=0;i<7;i++) cl.AddAt(new TClonesArray("AliHMPIDCluster"),i); AliHMPIDReconstructor::Dig2Clu(&dl,&cl);
+  
+  
+                             TClonesArray* pDigLst=(TClonesArray*)dl.At(ch);
+                             TClonesArray* pCluLst=(TClonesArray*)cl.At(ch);
+                    
+          Int_t recN=0;  Float_t recCkov=rec.CkovAngle(pCluLst,recN);    //reconstruct the ring
   
-        
-  return;                                                                      
-  Double_t th=8*TMath::DegToRad();                              //gRandom->Rndm()*TMath::PiOver4();
-  Double_t ph=gRandom->Rndm()*TMath::TwoPi(); 
-  Double_t radx=gRandom->Rndm()*AliHMPIDDigit::SizeAllX(); 
-  Double_t rady=gRandom->Rndm()*AliHMPIDDigit::SizeAllY();
+  AliESDtrack trk; trk.SetHMPIDtrk(radx,rady,th,ph);  trk.SetHMPIDmip(radx,rady,340,recN); trk.SetHMPIDsignal(recCkov); trk.SetHMPIDchi2(rec.CkovSigma2());
   
-  Int_t iHitCnt=0;  
+  Printf("Start of EVENT\n");                 
+  Printf("SDI------SDI---------SDI--------SDI------SDI------SDI");sl.Print();Printf("");
+  Printf("DIG------DIG---------DIG--------DIG------DIG------DIG");dl.Print();Printf("");                   
+  Printf("HIT------HIT---------HIT--------HIT------HIT------HIT");hl.Print();Printf("");
+  Printf("CLU------CLU---------CLU--------CLU------CLU------CLU");cl.Print();Printf("");                     
+  Printf("End of EVENT\n");                 
   
   
+  DrawEvt(&hl,pDigLst,pCluLst,&trk);  
+}//SimEvt()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void DrawEvt(TClonesArray *pHitLst,TClonesArray *pDigLst,TClonesArray *pCluLst,AliESDtrack *pTrk)
+{
+//visualization  
+  TPolyLine    *pRing  =new TPolyLine;     pRing->SetLineColor(kMagenta);     
   
+  Float_t recCkov=pTrk->GetHMPIDsignal();  Float_t recErr=TMath::Sqrt(pTrk->GetHMPIDchi2());
   
-      
-                                                                      
-             AliHMPIDv1::Hit2Sdi(&hitLst,&sdiLst);                               
-      AliHMPIDDigitizer::Sdi2Dig(&sdiLst,&digLst);     
-  AliHMPIDReconstructor::Dig2Clu(&digLst,&cluLst);   AliHMPIDReconstructor::Dig2Clu(&digLst,&cogLst,0);  
-  
-  
-  Int_t iDigN=0,iCluN=0,iCogN=0; for(Int_t i=0;i<7;i++){ iDigN+=((TClonesArray*)digLst.At(i))->GetEntries();                  
-                                                         iCluN+=((TClonesArray*)cluLst.At(i))->GetEntries();
-                                                         iCogN+=((TClonesArray*)cogLst.At(i))->GetEntries(); }                   
-  
-  Printf("SDI------SDI---------SDI--------SDI------SDI------SDI #%i",sdiLst.GetEntries());sdiLst.Print();Printf("");
-  Printf("DIG------DIG---------DIG--------DIG------DIG------DIG #%i",iDigN              );digLst.Print();Printf("");                   
-  Printf("HIT------HIT---------HIT--------HIT------HIT------HIT #%i",hitLst.GetEntries());hitLst.Print();Printf("");
-  Printf("CLU------CLU---------CLU--------CLU------CLU------CLU #%i",iCluN              );cluLst.Print();Printf("");                     
-  Printf("COG------COG---------COG--------COG------COG------COG #%i",iCogN              );cogLst.Print();Printf("");                     
-}
-
-
+  Float_t trkTh,trkPh,trkX,trkY,mipX,mipY; pTrk->GetHMPIDtrk(trkX,trkY,trkTh,trkPh);  Int_t mipQ,recN; pTrk->GetHMPIDmip(mipX,mipY,mipQ,recN);
+  AliHMPIDRecon rec; rec.SetTrack(trkTh,trkPh,trkX,trkY);
 
+  TVector2 pos;    for(int j=0;j<100;j++){rec.TracePhot(recCkov,j*0.0628,pos); pRing->SetNextPoint(pos.X(),pos.Y());}  
+  
+  gPad->Clear(); TButton *pBtn=new TButton("Next","SimEvt()",0,0,0.07,0.05);   pBtn->Draw(); DrawPc(0);
+  
+  pRing->Draw();
+  pDigLst->Draw();
+  pCluLst->Draw();
+  pHitLst->Draw();  
+  TLatex txt; txt.SetTextSize(0.02);
+//  txt.DrawLatex(20,-5,Form("#theta=%.4f#pm%.5f with %2i #check{C}"          ,simCkov,simErr,simN));
+  txt.DrawLatex(25,-5,Form("#theta=%.4f#pm%.5f with %2i #check{C}"          ,recCkov,recErr,recN));
+//  txt.DrawLatex(0 ,127,Form("#theta=%.2f#circ   #phi=%.2f#circ @(%.2f,%.2f) ",th*TMath::RadToDeg(),ph*TMath::RadToDeg(),radx,rady));
+}//DrawEvt()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void PrintMap()
 {
  
@@ -224,74 +247,7 @@ void PrintMap()
   
   
 }//PrintMap()
-
-
-void rec()
-{
-  
-  
-  Double_t ckovMax=0.75,ckovSim;
-  Int_t nSim=0;
-  while(nSim<3){
-    ckovSim=gRandom->Rndm()*ckovMax;//0.6468;
-    nSim=20*TMath::Power(TMath::Sin(ckovSim)/TMath::Sin(ckovMax),2); //scale number of photons 
-  }
-  
-  
-  TClonesArray *pCluLst=new TClonesArray("AliHMPIDCluster");
-  TPolyMarker  *pMipMap=new TPolyMarker();   pMipMap->SetMarkerStyle(8);  pMipMap->SetMarkerColor(kRed); 
-  TPolyMarker  *pPhoMap=new TPolyMarker();   pPhoMap->SetMarkerStyle(4);  pPhoMap->SetMarkerColor(kRed);
-  TPolyMarker  *pBkgMap=new TPolyMarker();   pBkgMap->SetMarkerStyle(25); pBkgMap->SetMarkerColor(kRed);
-  TPolyLine    *pRing  =new TPolyLine;                                    pRing->SetLineColor(kGreen);     
-  TPolyLine    *pOld   =new TPolyLine;                                    pOld->SetLineColor(kBlue);     
-  
-  Int_t iCluCnt=0; pMipMap->SetPoint(iCluCnt,x,y); new((*pCluLst)[iCluCnt++])   AliHMPIDCluster(1,x,y,200); //add mip cluster
-  
-//  for(int j=0;j<30;j++){                                                                                   //add bkg photons  
-//    Float_t bkgX=gRandom->Rndm()*AliHMPIDDigit::SizeAllX();
-//    Float_t bkgY=gRandom->Rndm()*AliHMPIDDigit::SizeAllY();
-//    pBkgMap->SetPoint(iCluCnt,bkgX,bkgY); new((*pCluLst)[iCluCnt++]) AliHMPIDCluster(1,bkgX,bkgY,35);
-//  }   
-
-  
-  
-  
-  AliHMPIDRecon    rec; rec.SetTrack(th,ph,x,y);                                                                
-  
-  TVector2 pos;
-  for(int i=0;i<nSim;i++){
-    rec.TracePhot(ckovSim,gRandom->Rndm()*2*TMath::Pi(),pos);                                   //add photons 
-    if(AliHMPIDDigit::IsInDead(pos.X(),pos.Y())) continue; 
-    pPhoMap->SetPoint(iCluCnt,pos.X(),pos.Y());    new((*pCluLst)[iCluCnt++]) AliHMPIDCluster(1,pos.X(),pos.Y(),35); 
-  }  
-
-  
-  Int_t nRec=0,nOld=0;
-  Double_t ckovRec=rec.CkovAngle(pCluLst,nRec); Double_t err=TMath::Sqrt(rec.CkovSigma2());   
-  Double_t ckovOld=old.CkovAngle(pCluLst,nOld);
-  
-  Printf("---------------- Now reconstructed --------------------");
-  
-  
-  for(int j=0;j<100;j++){rec.TracePhot(ckovRec,j*0.0628,pos); pRing->SetPoint(j,pos.X(),pos.Y());}  
-  for(int j=0;j<100;j++){rec.TracePhot(ckovOld,j*0.0628,pos); pOld->SetPoint(j,pos.X(),pos.Y());}  
-    
-  new TCanvas;  AliHMPIDDigit::DrawPc();  pMipMap->Draw(); pPhoMap->Draw(); pBkgMap->Draw(); pRing->Draw();  pOld->Draw(); 
-  
-  TLatex txt; txt.SetTextSize(0.03);
-  txt.DrawLatex(65,127,Form("#theta=%.4f#pm%.5f with %2i #check{C}"                             ,ckovSim, 0.,nSim             ));
-  txt.DrawLatex(65,122,Form("#theta=%.4f#pm%.5f with %2i #check{C} Old=%.4f with %i #check{C}"  ,ckovRec,err,nRec,ckovOld,nOld));
-  txt.DrawLatex(0 ,127,Form("#theta=%.2f#circ   #phi=%.2f#circ @(%.2f,%.2f) ",th*TMath::RadToDeg(),ph*TMath::RadToDeg(),x,y));
-                   
-//  for(int i=0;i<35;i++){
-//    Double_t ckov=0.1+i*0.02;
-//    Printf("Ckov=%.2f Old=%.3f New=%.3f",ckov,old.FindRingArea(ckov),rec.FindRingArea(ckov));
-//  }
-  
-  
-}//rec()
-
-
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void HitQA(Double_t cut,Double_t cutele,Double_t cutR)
 {
 // Provides a set of control plots intended primarily for charged particle flux analisys
@@ -313,8 +269,8 @@ void HitQA(Double_t cut,Double_t cutele,Double_t cutR)
   Double_t cutPproton        =cut;
                        
 
-  TH2F *pEleHitRZ    =new TH2F("EleHitRZ"    ,Form("e^{+} e^{-} hit %s;z[cm];R[cm]" ,GetName())     , 400,-300,300 ,400,-500,500);   //R-z plot 0cm<R<550cm -300cm<z<300cm  
-  TH2F *pEleHitRP    =new TH2F("EleHitRP"    ,Form("e^{+} e^{-} hit %s;p[GeV];R[cm]",GetName())     ,1000,-1  ,1   ,400,   0,550);   //R-p plot 0cm<R<550cm -1GeV<p<1GeV 
+  TH2F *pEleHitRZ    =new TH2F("EleHitRZ"    ,Form("e^{+} e^{-} hit %s;z[cm];R[cm]" ,GetName())     , 400,-300,300 ,400,-500,500);   //R-z
+  TH2F *pEleHitRP    =new TH2F("EleHitRP"    ,Form("e^{+} e^{-} hit %s;p[GeV];R[cm]",GetName())     ,1000,-1  ,1   ,400,   0,550);   //R-p
   TH1F *pEleAllP     =new TH1F("EleAllP"     ,     "e^{+} e^{-} all;p[GeV]"                         ,1000,-1  ,1                );  
   TH1F *pEleHitP     =new TH1F("EleHitP"     ,Form("e^{+} e^{-} hit %s;p[GeV]"      ,GetName())     ,1000,-1  ,1                );   
   TH1F *pMuoHitP     =new TH1F("MuoHitP"     ,Form("#mu^{-} #mu^{+} hit %s;p[GeV]"  ,GetName())     ,1000,-4  ,4                ); 
@@ -447,8 +403,7 @@ void HitQA(Double_t cut,Double_t cutele,Double_t cutR)
   
   gBenchmark->Show("HitsPlots");
 }//HitQA()
-
-
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void RecWithStack()
 {
   al->LoadHeader();al->LoadKinematics();
@@ -468,10 +423,9 @@ void RecWithStack()
   t.LoadClusters(rl->TreeR()); 
   t.PropagateBack(pEsd);
   rl->UnloadRecPoints();
-}
-
-
-void AliHMPIDReconstructor::CluQA(AliRunLoader *pAL)
+}//RecWithStack()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void CluQA(AliRunLoader *pAL)
 {
 // Quality assesment plots for clusters. 
 // This methode takes list of digits and form list of clusters again in order to 
@@ -485,17 +439,17 @@ void AliHMPIDReconstructor::CluQA(AliRunLoader *pAL)
   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*    pQ=new TH1F("HmpAllQ"  ,"Charge All"           ,4000 ,0  ,4000);// Q hists
+  TH1F* pCerQ=new TH1F("HmpCerQ"  ,"Charge Ckov"          ,4000 ,0  ,4000);
+  TH1F* pMipQ=new TH1F("HmpMipQ"  ,"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 );
+  TH1F*    pS=new TH1F("HmpCluSize"    ,"Cluster size;size"         ,100  ,0  ,100 );// size hists
+  TH1F* pCerS=new TH1F("HmpCluCerSize" ,"Ckov size;size"            ,100  ,0  ,100 );
+  TH1F* pMipS=new TH1F("HmpCluMipSize" ,"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());
+  TH2F*    pM=new TH2F("HmpCluMap"     ,"Cluster map;x [cm];y [cm]" ,1000 ,0  ,AliHMPIDDigit::SizeAllX(),1000,0,AliHMPIDDigit::SizeAllY()); // maps
+  TH2F* pMipM=new TH2F("HmpCluMipMap"  ,"MIP map;x [cm];y [cm]"     ,1000 ,0  ,AliHMPIDDigit::SizeAllX(),1000,0,AliHMPIDDigit::SizeAllY());
+  TH2F* pCerM=new TH2F("HmpCluCerMap"  ,"Ckov map;x [cm];y [cm]"    ,1000 ,0  ,AliHMPIDDigit::SizeAllX(),1000,0,AliHMPIDDigit::SizeAllY());
  
   
   
@@ -504,14 +458,14 @@ void AliHMPIDReconstructor::CluQA(AliRunLoader *pAL)
     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 iCh=0;iCh<7;iCh++) AliHMPIDReconstructor::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                                      
+                                             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
                                        
@@ -525,10 +479,8 @@ void AliHMPIDReconstructor::CluQA(AliRunLoader *pAL)
   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 AliHMPID::OccupancyPrint(Int_t iEvtNreq)
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void OccupancyPrint(Int_t iEvtNreq)
 {
 //prints occupancy for each chamber in a given event
   Int_t iEvtNmin,iEvtNmax;
@@ -582,9 +534,7 @@ void AliHMPID::OccupancyPrint(Int_t iEvtNreq)
   GetLoader()->UnloadDigits();
   GetLoader()->GetRunLoader()->UnloadHeader();    
   GetLoader()->GetRunLoader()->UnloadKinematics();    
-}
-
-
+}//OccupancyPrint()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPID::SummaryOfEvent(Int_t iEvtN) const
 {
@@ -609,5 +559,170 @@ void AliHMPID::SummaryOfEvent(Int_t iEvtN) const
   AliInfo(Form("Total n. of tracks in stack(+sec) %i",pStack->GetNtrack()));
   GetLoader()->GetRunLoader()->UnloadHeader();
   GetLoader()->GetRunLoader()->UnloadKinematics();
-}
+}//SummaryOfEvent()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void DrawZoom()
+{
+// Show info about current cursur position in status bar of the canvas
+// Arguments: none
+//   Returns: none     
+  TCanvas *pC=(TCanvas*)gPad; 
+  TRootCanvas *pRC= (TRootCanvas*)pC->GetCanvasImp();
+  TGStatusBar *pBar=pRC->GetStatusBar();
+  pBar->SetParts(5);
+  Float_t x=gPad->AbsPixeltoX(gPad->GetEventX());
+  Float_t y=gPad->AbsPixeltoY(gPad->GetEventY());
+  AliHMPIDDigit dig;dig.Manual1(1,x,y); UInt_t w32=0; 
+  if(IsInDead(x,y))
+    pBar->SetText("Out of sensitive area",4);    
+  else{
+    Int_t ddl=dig.Raw(w32);
+    pBar->SetText(Form("(p%i,x%i,y%i) ddl=%i 0x%x (r%i,d%i,a%i) (%.2f,%.2f)",
+        dig.Pc(),dig.PadPcX(),dig.PadPcY(),
+        ddl,w32,
+        dig.Row(),dig.Dilogic(),dig.Addr(),
+        dig.LorsX(),dig.LorsY()                            ),4);
+  }    
+  if(gPad->GetEvent()==1){
+    new TCanvas("zoom",Form("Row %i DILOGIC %i",dig.Row(),dig.Dilogic()));  
+  }
+}//DrawZoom()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void DrawPc(Bool_t isFill) 
+{ 
+// Utility methode draws HMPID chamber PCs on event display.
+// Arguments: none
+//   Returns: none      
+//  y6  ----------  ----------
+//      |        |  |        |
+//      |    4   |  |    5   |
+//  y5  ----------  ----------
+//
+//  y4  ----------  ----------
+//      |        |  |        |
+//      |    2   |  |    3   |   view from electronics side
+//  y3  ----------  ----------
+//          
+//  y2  ----------  ----------
+//      |        |  |        |
+//      |    0   |  |    1   |
+//  y1  ----------  ----------
+//      x1      x2  x3       x4
+  gPad->Range(-10,-10,AliHMPIDDigit::SizeAllX()+5,AliHMPIDDigit::SizeAllY()+5); 
+  Float_t x1=0,
+          x2=AliHMPIDDigit::SizePcX(),
+          x3=AliHMPIDDigit::SizePcX()+AliHMPIDDigit::SizeDead(),   
+          x4=AliHMPIDDigit::SizeAllX();
+  Float_t y1=0,
+          y2=  AliHMPIDDigit::SizePcY(),
+          y3=  AliHMPIDDigit::SizePcY()+AliHMPIDDigit::SizeDead(),
+          y4=2*AliHMPIDDigit::SizePcY()+AliHMPIDDigit::SizeDead(),
+          y5=  AliHMPIDDigit::SizeAllY()-AliHMPIDDigit::SizePcY(),
+          y6=  AliHMPIDDigit::SizeAllY();
+
+  Float_t xL[5]={x1,x1,x2,x2,x1}; //clockwise
+  Float_t xR[5]={x3,x3,x4,x4,x3};  
+  Float_t yD[5]={y1,y2,y2,y1,y1};
+  Float_t yC[5]={y3,y4,y4,y3,y3};  
+  Float_t yU[5]={y5,y6,y6,y5,y5};
+    
+  Float_t dX2=0.5*AliHMPIDDigit::SizePadX(),
+          dY2=0.5*AliHMPIDDigit::SizePadY() ;
+  
+  TLatex txt; txt.SetTextSize(0.01);
+  Int_t iColLeft=29,iColRight=41;
+  TPolyLine *pc=0;  TLine *pL;
+  AliHMPIDDigit dig;
+  for(Int_t iPc=0;iPc<AliHMPIDDigit::kPcAll;iPc++){
+    if(iPc==4) pc=new TPolyLine(5,xL,yU); if(iPc==5) pc=new TPolyLine(5,xR,yU); //draw PCs
+    if(iPc==2) pc=new TPolyLine(5,xL,yC); if(iPc==3) pc=new TPolyLine(5,xR,yC);
+    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();continue;}
+    
+    pc->Draw("f");
+    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.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()-dX2,dig.LorsY()-dY2,dig.LorsX()+AliHMPIDDigit::SizePcX()-dX2,dig.LorsY()-dY2);//draw horizontal line 
+      if(iRow!=0) pL->Draw(); 
+    }//row loop  
+    
+    txt.SetTextAlign(13);
+    for(Int_t iDil=0;iDil<10;iDil++){//draw dilogic lines (vertical)
+      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()-dX2,dig.LorsY()-dY2,dig.LorsX()-dX2,dig.LorsY()+AliHMPIDDigit::SizePcY()-dY2); //draw vertical line
+      if(iDil!=0)pL->Draw();
+    }//dilogic loop        
+  }//PC loop      
+}//DrawPc()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void ReadEvt() 
+{
+// Display digits, reconstructed tracks intersections and HMPID rings if available 
+// Arguments: none
+//   Returns: none    
+  TFile *pEsdFl=TFile::Open("AliESDs.root");     if(!pEsdFl || !pEsdFl->IsOpen()) return;//open AliESDs.root
+  TTree *pEsdTr=(TTree*) pEsdFl->Get("esdTree"); if(!pEsdTr)                      return;//get ESD tree
+                                                                 
+  AliESD *pEsd=0;  pEsdTr->SetBranchAddress("ESD", &pEsd);
+  
+  for(Int_t iEvt=0;iEvt<pEsdTr->GetEntries();iEvt++) {                //events loop
+    
+    
+    
+    pEsdTr->GetEntry(iEvt);                                       //get ESD for this event   
+    for(Int_t iTrk=0;iTrk<pEsd->GetNumberOfTracks();iTrk++){//ESD tracks loop
+      AliESDtrack *pTrk = pEsd->GetTrack(iTrk);             //
+    }//ESD tracks loop
+    
+    al->GetEvent(iEvt);   rl->TreeD()->GetEntry(0); //get digits list
+    
+    
+  }//events loop
+  delete pEsd;  pEsdFl->Close();//close AliESDs.root
+//  rl->UnloadDigits();
+}//ReadEvt()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void t1(Int_t case=1)
+{
+  AliHMPIDDigit *d[10]; for(Int_t i=0;i<10;i++) d[i]=new AliHMPIDDigit;
+  
+  
+  Int_t iNdig;
+  
+  if(case==1){
+    iNdig=9;  
+  
+                                                              d[0]->Manual2(1,2,67,26, 33); 
+                                d[1]->Manual2(1,2,66,25,431); d[2]->Manual2(1,2,67,25, 21);
+  d[3]->Manual2(1,2,65,24,127); d[4]->Manual2(1,2,66,24, 54); d[5]->Manual2(1,2,67,24,  5);
+  d[6]->Manual2(1,2,65,23, 20); d[7]->Manual2(1,2,66,23,  5); d[8]->Manual2(1,2,67,23,  6);
+  }else if(case==2){
+    iNdig=3;
+    d[0]->Manual2(0,0,36,14,  8); 
+    d[1]->Manual2(0,0,36,13, 33); d[2]->Manual2(0,0,37,13, 22);
+  }
+  
+  AliHMPIDCluster c;
+  for(int i=0;i<iNdig;i++) c.DigAdd(d[i]);  c.Print();
+  
+  
+  TClonesArray *cl=new TClonesArray("AliHMPIDCluster");
+  
+  c.Solve(cl,kTRUE);
+  Printf("");
+  
+  cl->Print();  
+}//t1()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+
+
+
+