]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
Cluster improvements + drawing of digits like TBox.
authordibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 26 Feb 2007 11:23:11 +0000 (11:23 +0000)
committerdibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 26 Feb 2007 11:23:11 +0000 (11:23 +0000)
HMPID/AliHMPIDCluster.cxx
HMPID/AliHMPIDCluster.h
HMPID/AliHMPIDDigit.cxx
HMPID/AliHMPIDDigit.h
HMPID/AliHMPIDHit.h
HMPID/Hmenu.C

index 8504030a81036c994110ae9a399c4e874e1d5f22..9fbb3cde96f78f4e6fd572d7cddb3541368939b6 100644 (file)
@@ -24,25 +24,37 @@ void AliHMPIDCluster::CoG()
 // Calculates naive cluster position as a center of gravity of its digits.
 // Arguments: none 
 //   Returns: none
 // Calculates naive cluster position as a center of gravity of its digits.
 // Arguments: none 
 //   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
+  Int_t minPadX=999,minPadY=999,maxPadX=-1,maxPadY=-1;      //for box finding  
+  if(fDigs==0) return;                                      //no digits in this cluster
+  fX=fY=fQRaw=0;                                            //init summable parameters
+  Int_t maxQpad=-1,maxQ=-1;                                 //to calculate the pad with the highest charge
   AliHMPIDDigit *pDig;
   AliHMPIDDigit *pDig;
-  for(Int_t iDig=0;iDig<fDigs->GetEntriesFast();iDig++){//digits loop
-    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
+  for(Int_t iDig=0;iDig<fDigs->GetEntriesFast();iDig++){    //digits loop
+    pDig=(AliHMPIDDigit*)fDigs->At(iDig);                   //get pointer to next digit
+
+    if(pDig->PadPcX() > maxPadX) maxPadX = pDig->PadPcX();  // find the minimum box that contain the cluster  MaxX                            
+    if(pDig->PadPcY() > maxPadY) maxPadY = pDig->PadPcY();  //                                                MaxY
+    if(pDig->PadPcX() < minPadX) minPadX = pDig->PadPcX();  //                                                MinX   
+    if(pDig->PadPcY() < minPadY) minPadY = pDig->PadPcY();  //                                                MinY   
+    
+    Float_t q=pDig->Q();                                    //get QDC 
+    fX += pDig->LorsX()*q;fY +=pDig->LorsY()*q;             //add digit center weighted by QDC
+    fQRaw+=q;                                               //increment total charge 
+    if(q>maxQ) {maxQpad = pDig->Pad();maxQ=(Int_t)q;}       // to find pad with highest charge
   }//digits loop
   }//digits loop
-  if ( fQ != 0 )   fX/=fQ;fY/=fQ;                       //final center of gravity
+  
+  fBox=(maxPadX-minPadX+1)*100+maxPadY-minPadY+1;           // dimension of the box: format Xdim*100+Ydim
+  
+  if ( fQRaw != 0 )   fX/=fQRaw;fY/=fQRaw;                  //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 
+  CorrSin();                                               //correct it by sinoid   
+  
+  fQ  = fQRaw;                                             // Before starting fit procedure, Q and QRaw must be equal
+  fCh=pDig->Ch();                                          //initialize chamber number
+  fMaxQpad = maxQpad; fMaxQ=maxQ;                          //store max charge pad to the field
+  fChi2=0;                                                 // no Chi2 to find
+  fNlocMax=0;                                              // no maxima to find
   fSt=kCoG;
 }//CoG()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   fSt=kCoG;
 }//CoG()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -76,13 +88,13 @@ void AliHMPIDCluster::FitFunc(Int_t &iNpars, Double_t *, Double_t &chi2, Double_
   Int_t iNshape = iNpars/3;
     
   chi2 = 0;
   Int_t iNshape = iNpars/3;
     
   chi2 = 0;
-  for(Int_t i=0;i<pClu->Size();i++){                                       //loop on all pads of the cluster
-    Double_t dQpadMath = 0;                                                //pad charge collector  
-    for(Int_t j=0;j<iNshape;j++){                                          //Mathiesons loop as all of them may contribute to this pad
-      dQpadMath+=par[3*j+2]*pClu->Dig(i)->Mathieson(par[3*j],par[3*j+1]);  // par[3*j+2] is charge par[3*j] is x par[3*j+1] is y of current Mathieson
+  for(Int_t i=0;i<pClu->Size();i++){                                                   //loop on all pads of the cluster
+    Double_t dQpadMath = 0;                                                            //pad charge collector  
+    for(Int_t j=0;j<iNshape;j++){                                                      //Mathiesons loop as all of them may contribute to this pad
+      dQpadMath+=par[3*j+2]*pClu->Dig(i)->Mathieson(par[3*j],par[3*j+1]);              // par[3*j+2] is charge par[3*j] is x par[3*j+1] is y of current Mathieson
     }
     }
-    chi2 +=TMath::Power((pClu->Dig(i)->Q()-dQpadMath),2);                  //
-  }                                                                             //loop on all pads of the cluster     
+    if(dQpadMath>0)chi2 +=TMath::Power((pClu->Dig(i)->Q()-dQpadMath),2)/dQpadMath;     //
+  }                                                                                    //loop on all pads of the cluster     
 }//FitFunction()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDCluster::Print(Option_t* opt)const
 }//FitFunction()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDCluster::Print(Option_t* opt)const
@@ -94,7 +106,6 @@ void AliHMPIDCluster::Print(Option_t* opt)const
     case        kUnf  : status="unfolded (fit)"   ;break;
     case        kCoG  : status="coged         "   ;break;
     case        kLo1  : status="locmax 1 (fit)"   ;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        kMax  : status="exceeded (cog)"   ;break;
     case        kNot  : status="not done (cog)"   ;break;
     case        kEmp  : status="empty         "   ;break;
@@ -104,8 +115,10 @@ void AliHMPIDCluster::Print(Option_t* opt)const
     
     default:            status="??????"          ;break;   
   }
     
     default:            status="??????"          ;break;   
   }
-  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);
+  Double_t ratio=0;
+  if(Q()>0&&QRaw()>0) ratio = Q()/QRaw()*100;
+  Printf("%sCLU:(%7.3f,%7.3f) Qfitted=%8.3f QRaw=%8.3f(%3.0f%%) ch=%i  Npads=%2i  DimBox %i  NlocMax %i Chi2 %f   %s",
+         opt,X(),Y(),Q(),QRaw(),ratio,Ch(),Size(),fBox,fNlocMax,fChi2,status);
   if(fDigs) fDigs->Print();    
 }//Print()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   if(fDigs) fDigs->Print();    
 }//Print()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -120,38 +133,24 @@ Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
 //           isTryUnfold - flag to switch on/off unfolding   
 //  Returns: number of local maxima of original cluster
   CoG();
 //           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 
+  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;
   } 
     (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  
 //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;                                     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
+  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;                                     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. Also find the box contaning the cluster   
   fNlocMax=0;
 //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
-    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;
-    
+
+ 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 
     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 
@@ -160,78 +159,51 @@ Int_t AliHMPIDCluster::Solve(TClonesArray *pCluLst,Bool_t isTryUnfold)
       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(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&&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);
-      
+    if(iHowManyMoreCnt==0&&fNlocMax<kMaxLocMax){                                         //this pad has Q more then any neighbour so it's local maximum
+      pMinuit->mnparm(3*fNlocMax  ,Form("x%i",fNlocMax),pDig1->LorsX(),0.1,0,0,iErrFlg); // X,Y,Q initial values of the loc max pad w
+      pMinuit->mnparm(3*fNlocMax+1,Form("y%i",fNlocMax),pDig1->LorsY(),0.1,0,0,iErrFlg); //
+      pMinuit->mnparm(3*fNlocMax+2,Form("q%i",fNlocMax),pDig1->Q(),0.1,0,100000,iErrFlg);// constrained to be positive
       fNlocMax++;
     }//if this pad is local maximum
   }//first digits loop
   
       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
 //Phase 2. Fit loc max number of Mathiesons or add this current cluster to the list
-// 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);
+// case 1 -> no loc max found
+ if ( fNlocMax == 0) {                                                                   // case of no local maxima found: pads with same charge...
+   pMinuit->mnparm(3*fNlocMax  ,Form("x%i",fNlocMax),fX,0.1,0,0,iErrFlg);                // Init values taken from CoG() -> fX,fY,fQRaw
+   pMinuit->mnparm(3*fNlocMax+1,Form("y%i",fNlocMax),fY,0.1,0,0,iErrFlg);                //
+   pMinuit->mnparm(3*fNlocMax+2,Form("q%i",fNlocMax),fQRaw,0.1,0,100000,iErrFlg);        //
    fNlocMax = 1;
    fSt=kNoLoc;
  }
    fNlocMax = 1;
    fSt=kNoLoc;
  }
- if ( fNlocMax >= kMaxLocMax)
+
+// case 2 -> loc max found. Check # of loc maxima 
+ if ( fNlocMax >= kMaxLocMax)                                                            // if # of local maxima exceeds 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
+     fSt = kMax;   new ((*pCluLst)[iCluCnt++]) AliHMPIDCluster(*this);                   //...add this raw cluster  
+   }                                                                                     //or...
+ 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("SIMPLEX" ,arglist,2,iErrFlg);                                        //start fitting with Simplex
+   pMinuit->mnexcm("MIGRAD" ,arglist,2,iErrFlg);                                         //fitting improved by Migrad
    
    
-   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
-       }
-     }
- }
-    
+   Double_t dummy; TString sName;                                                        //vars to get results from Minuit
+   for(Int_t i=0;i<fNlocMax;i++){                                                        //store the local maxima parameters
+      pMinuit->mnpout(3*i   ,sName,  fX, fErrX , dummy, dummy, iErrFlg);                 // X 
+      pMinuit->mnpout(3*i+1 ,sName,  fY, fErrY , dummy, dummy, iErrFlg);                 // Y
+      pMinuit->mnpout(3*i+2 ,sName,  fQ, fErrQ , dummy, dummy, iErrFlg);                 // Q
+      pMinuit->mnstat(fChi2,dummy,dummy,iErrFlg,iErrFlg,iErrFlg);                        // Chi2 of the fit
 
 
-  
+      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 fNlocMax;
+ delete pMinuit;
+ return fNlocMax;
 }//Solve()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 }//Solve()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index bfbfcb8497619634cf4ff6d5af647282b7dd452a..63f14d2880b9daee9cdbc309c3046c3889741a2a 100644 (file)
@@ -10,19 +10,19 @@ class TClonesArray;        //Solve()
 class AliHMPIDCluster :public TObject
 {
 public:
 class AliHMPIDCluster :public TObject
 {
 public:
-  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
+  enum EClusterStatus {kFrm,kCoG,kLo1,kUnf,kMax,kNot,kEdg,kSi1,kNoLoc,kEmp=-1};      //status flags    
+  AliHMPIDCluster():TObject( ),fCh(-1),fSi(-1),fSt(kEmp),fBox(-1),fNlocMax(-1),fMaxQpad(-1),fMaxQ(-1),fQRaw(0),fQ(0),fErrQ(-1),fX(0),fErrX(-1),fY(0),fErrY(-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           (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),fQRaw(c.fQRaw),
+                                                                  fQ (c.fQ ),fErrQ(c.fErrQ),
+                                                                  fX (c.fX ),fErrX(c.fErrX),
+                                                                  fY (c.fY ),fErrY(c.fErrY),fChi2(c.fChi2),fDigs(0)                  {}//copy ctor
   AliHMPIDCluster &operator=(const AliHMPIDCluster &c) {if(this == &c)return *this;TObject::operator=(c);          
   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;}
+                                                        fSi=c.fSi;  fSt=c.fSt; fCh=c.fCh; fBox=c.fBox;fNlocMax=c.fNlocMax;fMaxQpad=c.fMaxQpad; fMaxQ=c.fMaxQ;fQRaw=c.fQRaw;
+                                                        fQ=c.fQ; fErrQ=c.fErrQ
+                                                        fX=c.fX; fErrX=c.fErrX;
+                                                        fY=c.fY; fErrY=c.fErrY; fChi2=c.fChi2;fDigs=c.fDigs ? new TObjArray(*c.fDigs):0; return *this;}
     
   virtual ~AliHMPIDCluster(                                     )                                                                        {if(fDigs) delete fDigs; fDigs=0;}
 //framework part                   
     
   virtual ~AliHMPIDCluster(                                     )                                                                        {if(fDigs) delete fDigs; fDigs=0;}
 //framework part                   
@@ -40,36 +40,32 @@ public:
          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                                  
          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       QRaw     (                                         )const{return fQRaw;                                  } //raw cluster charge in QDC channels 
+         Double_t       Q        (                                         )const{return fQ;                                     } //given cluster charge in QDC channels 
+         Double_t       Qe       (                                         )const{return fErrQ;                                  } //Error in cluster charge in QDC channels 
          Double_t       X        (                                         )const{return fX;                                     } //cluster x position in LRS
          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       Xe       (                                         )const{return fErrX;                                  } //cluster charge in QDC channels 
          Double_t       Y        (                                         )const{return fY;                                     } //cluster y position in LRS 
          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       Ye       (                                         )const{return fErrY;                                    } //cluster charge in QDC channels 
          Double_t       Chi2     (                                         )const{return fChi2;                                  } 
 protected:
   Int_t         fCh;          //chamber number
          Double_t       Chi2     (                                         )const{return fChi2;                                  } 
 protected:
   Int_t         fCh;          //chamber number
-  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         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
   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  
+  Double_t      fMaxQ;        //that max charge value
+  Double_t      fQRaw;        //QDC value of the raw cluster
+  Double_t      fQ;           //QDC value of the actual cluster
+  Double_t      fErrQ;        //error on Q
+  Double_t      fX;           //local x postion, [cm]
+  Double_t      fErrX;        //error on x postion, [cm]
+  Double_t      fY;           //local y postion, [cm]
+  Double_t      fErrY;        //error on y postion, [cm]
+  Double_t      fChi2;        //some estimator of the fit quality
   TObjArray    *fDigs;        //! list of digits forming this cluster
   TObjArray    *fDigs;        //! list of digits forming this cluster
-  ClassDef(AliHMPIDCluster,5) //HMPID cluster class       
+  ClassDef(AliHMPIDCluster,6) //HMPID cluster class
 };//class AliHMPIDCluster
 
 typedef AliHMPIDCluster AliRICHCluster; // for backward compatibility
 };//class AliHMPIDCluster
 
 typedef AliHMPIDCluster AliRICHCluster; // for backward compatibility
@@ -90,7 +86,11 @@ void AliHMPIDCluster::Reset()
 {
   if(fDigs) delete fDigs; 
   fDigs=0; 
 {
   if(fDigs) delete fDigs; 
   fDigs=0; 
-  fQ=fCh=fSi=-1;fX=fY=-1;fSt=kEmp;      
+  fSt=kEmp;
+  fQRaw=fQ=0;
+  fX=fY=0;
+  fCh=fSi=fBox=fNlocMax=fMaxQpad=-1;
+  fMaxQ=fErrQ=fErrX=fErrY=fChi2=-1; //empty ctor
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 Bool_t AliHMPIDCluster::IsInPc()
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 Bool_t AliHMPIDCluster::IsInPc()
index 3b14f0f7aa0cd5f2771d9c3797231bd23f37e90d..3000c92d1c79237ceb5ad9518008b60d2bc628d8 100644 (file)
@@ -15,7 +15,7 @@
 
 #include "AliHMPIDDigit.h" //class header
 #include <TClonesArray.h>  //Hit2Sdi() 
 
 #include "AliHMPIDDigit.h" //class header
 #include <TClonesArray.h>  //Hit2Sdi() 
-#include <TMarker.h>       //Draw() 
+#include <TBox.h>       //Draw() 
 ClassImp(AliHMPIDDigit)
 
 const Float_t AliHMPIDDigit::fMinPcX[]={0,SizePcX() + SizeDead(),0,SizePcX() + SizeDead(),       0,SizePcX() + SizeDead()};
 ClassImp(AliHMPIDDigit)
 
 const Float_t AliHMPIDDigit::fMinPcX[]={0,SizePcX() + SizeDead(),0,SizePcX() + SizeDead(),       0,SizePcX() + SizeDead()};
@@ -75,7 +75,10 @@ HMPID raw word is 32 bits with the structure:
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Draw(Option_t*)
 {
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Draw(Option_t*)
 {
-  TMarker *pMark=new TMarker(LorsX(),LorsY(),25); pMark->SetMarkerColor(kGreen); pMark->Draw();
+//  TMarker *pMark=new TMarker(LorsX(),LorsY(),25); pMark->SetMarkerColor(kGreen);pMark->Draw();
+  TBox *pad = new TBox(LorsX()-0.5*SizePadX(),LorsY()-0.5*SizePadY(),LorsX()+0.5*SizePadX(),LorsY()+0.5*SizePadY());
+  pad->SetFillStyle(0);pad->SetLineColor(kGreen);
+  pad->Draw();
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Hit2Sdi(AliHMPIDHit *pHit,TClonesArray *pSdiLst)
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDDigit::Hit2Sdi(AliHMPIDHit *pHit,TClonesArray *pSdiLst)
index df44e6e322a88a116274f0acdbeb22d55f108f63..52563a375678e2bf81a73a9a73ac1ea488091131 100644 (file)
@@ -18,7 +18,7 @@ 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    
   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(                                      ):AliDigit( ),fPad(Abs(-1,-1,-1,-1)),fQ(0)   {}                         //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    
            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    
@@ -39,7 +39,7 @@ public:
          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  
          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  IsOverTh    (Float_t q                      )     {return q >= 4;                                  }                  //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 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]
@@ -100,7 +100,7 @@ Bool_t AliHMPIDDigit::Set(AliHMPIDHit *pHit,Int_t pad)
 //            pad - for which pad to create
 //   Returns: none    
   
 //            pad - for which pad to create
 //   Returns: none    
   
-  fPad=Abs(-1,-1,-1,-1); fQ=-1; //reset
+  fPad=Abs(-1,-1,-1,-1); fQ=0; //reset
   Int_t pc,px,py;
   Float_t x=pHit->LorsX(),y=pHit->LorsY();
   
   Int_t pc,px,py;
   Float_t x=pHit->LorsX(),y=pHit->LorsY();
   
index 645249b3066dcec667d407bddb7a2b2ece09d14b..e84d8da21f9ab7b922f8bc491146cfb953e12095 100644 (file)
@@ -11,9 +11,9 @@
 class AliHMPIDHit : public AliHit //   TObject-AliHit-AliHMPIDHit
 {
 public:
 class AliHMPIDHit : public AliHit //   TObject-AliHit-AliHMPIDHit
 {
 public:
-  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 
+  AliHMPIDHit(                                                                         ):AliHit(    ),fCh(-1),fPid(-1),fQ(0),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(0),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(0),fLx(x),fLy(y){QdcTot(e);fTrack=tn;}//manual ctor 
   virtual ~AliHMPIDHit()                                                                                                                                 {}
 //framework part
          void    Print(Option_t *opt="")const;                                                 //from TObject to print current status
   virtual ~AliHMPIDHit()                                                                                                                                 {}
 //framework part
          void    Print(Option_t *opt="")const;                                                 //from TObject to print current status
index 4f47e58b115a136af863013159312fdb9df758f5..4a2dcf95f36779805d50baa76d693d9062fe6b41 100644 (file)
@@ -578,12 +578,11 @@ void hed()
 {//event display from files
   static TCanvas *pC=0;
   static Int_t iEvt=0;
 {//event display from files
   static TCanvas *pC=0;
   static Int_t iEvt=0;
-  static Int_t iEvtTot=0;
+  static Int_t iEvtTot=999;
   static TFile *pEsdFl=0;
   static TTree *pEsdTr=0;
   static AliESD *pEsd=0;
   static TFile *pEsdFl=0;
   static TTree *pEsdTr=0;
   static AliESD *pEsd=0;
-  
-  if(!pC){
+  if(!pC&&iEvt<iEvtTot){
     if(hl==0) {Printf("hed: no HMPID loader");return;}
     Printf("Opening session");
     pEsdFl=TFile::Open("AliESDs.root");     if(!pEsdFl || !pEsdFl->IsOpen()) return;//open AliESDs.root
     if(hl==0) {Printf("hed: no HMPID loader");return;}
     Printf("Opening session");
     pEsdFl=TFile::Open("AliESDs.root");     if(!pEsdFl || !pEsdFl->IsOpen()) return;//open AliESDs.root
@@ -597,13 +596,13 @@ void hed()
  
   if(iEvt<iEvtTot){
     pEsdTr->GetEntry(iEvt); al->GetEvent(iEvt); hl->TreeD()->GetEntry(0); hl->TreeR()->GetEntry(0);
  
   if(iEvt<iEvtTot){
     pEsdTr->GetEntry(iEvt); al->GetEvent(iEvt); hl->TreeD()->GetEntry(0); hl->TreeR()->GetEntry(0);
-    TLatex txt;   pC->cd(3);  txt.DrawLatex(0.2,0.2,Form("Event %i Total %i",iEvt,iEvtTot));
+    TLatex txt; pC->cd(3);  gPad->Clear(); 
+    txt.SetTextSize(0.1);txt.DrawLatex(0.2,0.2,Form("Event %i (total %i)",iEvt,iEvtTot));
     DrawEvt(pC,h->DigLst(),h->CluLst(),pEsd);
     iEvt++;
   }else{
     DrawEvt(pC,h->DigLst(),h->CluLst(),pEsd);
     iEvt++;
   }else{
-    Printf("Last event");
-    pC->Clear();
-    delete pC;pC=0x0;
+    Printf("--- No more events available...Bye.");
+    pC->Close();
   }
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
   }
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
@@ -639,10 +638,10 @@ void sed()
   AliHMPIDReconstructor::Dig2Clu(&ld,&lc);
 //        AliHMPIDTracker::Recon(&esd,&cl);
   
   AliHMPIDReconstructor::Dig2Clu(&ld,&lc);
 //        AliHMPIDTracker::Recon(&esd,&cl);
   
-  DrawEvt(pC1,&lh,&ld,&lc,&esd);  
+  DrawEvt(pC1,&ld,&lc,&esd);  
 }//SimEvt()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 }//SimEvt()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-void DrawEvt(TCanvas *pC,TClonesArray *pHitLst,TObjArray *pDigLst,TObjArray *pCluLst,AliESD *pEsd)
+void DrawEvt(TCanvas *pC,TObjArray *pDigLst,TObjArray *pCluLst,AliESD *pEsd)
 {//draws all the objects of current event
 
   AliHMPIDRecon rec;  
 {//draws all the objects of current event
 
   AliHMPIDRecon rec;  
@@ -681,10 +680,6 @@ void DrawEvt(TCanvas *pC,TClonesArray *pHitLst,TObjArray *pDigLst,TObjArray *pCl
     ((TClonesArray*)pCluLst->At(iCh))->Draw();  //draw clusters
                             pTxC[iCh]->Draw();  //draw intersections
                             pRin[iCh]->Draw();  //draw rings
     ((TClonesArray*)pCluLst->At(iCh))->Draw();  //draw clusters
                             pTxC[iCh]->Draw();  //draw intersections
                             pRin[iCh]->Draw();  //draw rings
-    for(Int_t iHit=0;iHit<pHitLst->GetEntries();iHit++) {
-      AliHMPIDHit *pHit=(AliHMPIDHit*)pHitLst->At(iHit);
-      if(pHit->Ch()==iCh)pHit->Draw();
-    }
     gPad->SetEditable(kFALSE);
   }//chambers loop
 //  TLatex txt; txt.SetTextSize(0.02);
     gPad->SetEditable(kFALSE);
   }//chambers loop
 //  TLatex txt; txt.SetTextSize(0.02);