New class AliHMPIDPid added. Now HMPID could run PID even on ESD within the same...
authordibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 13 Dec 2007 16:06:29 +0000 (16:06 +0000)
committerdibari <dibari@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 13 Dec 2007 16:06:29 +0000 (16:06 +0000)
HMPID/AliHMPIDParam.cxx
HMPID/AliHMPIDParam.h
HMPID/AliHMPIDPid.cxx [new file with mode: 0644]
HMPID/AliHMPIDPid.h [new file with mode: 0644]
HMPID/AliHMPIDRecon.cxx
HMPID/AliHMPIDRecon.h
HMPID/AliHMPIDReconstructor.cxx
HMPID/HMPIDrecLinkDef.h
HMPID/Hdisp.C
HMPID/libHMPIDrec.pkg

index 78ebef3..ac5475a 100644 (file)
@@ -195,3 +195,119 @@ Int_t AliHMPIDParam::StackCount(Int_t pid,Int_t evt)
   return iCnt;
 }
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Double_t AliHMPIDParam::Sigma2(Double_t trkTheta,Double_t trkPhi,Double_t ckovTh, Double_t ckovPh)
+{
+// Analithical calculation of total error (as a sum of localization, geometrical and chromatic errors) on Cerenkov angle for a given Cerenkov photon 
+// created by a given MIP. Fromulae according to CERN-EP-2000-058 
+// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
+//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
+//            MIP beta
+//   Returns: absolute error on Cerenkov angle, [radians]    
+  
+  TVector3 v(-999,-999,-999);
+  Double_t trkBeta = 1./(TMath::Cos(ckovTh)*GetRefIdx());
+  
+  if(trkBeta > 1) trkBeta = 1;                 //protection against bad measured thetaCer  
+  if(trkBeta < 0) trkBeta = 0.0001;            //
+
+  v.SetX(SigLoc (trkTheta,trkPhi,ckovTh,ckovPh,trkBeta));
+  v.SetY(SigGeom(trkTheta,trkPhi,ckovTh,ckovPh,trkBeta));
+  v.SetZ(SigCrom(trkTheta,trkPhi,ckovTh,ckovPh,trkBeta));
+
+  return v.Mag2();
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Double_t AliHMPIDParam::SigLoc(Double_t trkTheta,Double_t trkPhi,Double_t thetaC, Double_t phiC,Double_t betaM)
+{
+// Analitical calculation of localization error (due to finite segmentation of PC) on Cerenkov angle for a given Cerenkov photon 
+// created by a given MIP. Fromulae according to CERN-EP-2000-058 
+// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
+//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
+//            MIP beta
+//   Returns: absolute error on Cerenkov angle, [radians]    
+  
+  Double_t phiDelta = phiC - trkPhi;
+
+  Double_t sint     = TMath::Sin(trkTheta);
+  Double_t cost     = TMath::Cos(trkTheta);
+  Double_t sinf     = TMath::Sin(trkPhi);
+  Double_t cosf     = TMath::Cos(trkPhi);
+  Double_t sinfd    = TMath::Sin(phiDelta);
+  Double_t cosfd    = TMath::Cos(phiDelta);
+  Double_t tantheta = TMath::Tan(thetaC);
+  
+  Double_t alpha =cost-tantheta*cosfd*sint;                                                 // formula (11)
+  Double_t k = 1.-GetRefIdx()*GetRefIdx()+alpha*alpha/(betaM*betaM);        // formula (after 8 in the text)
+  if (k<0) return 1e10;
+  Double_t mu =sint*sinf+tantheta*(cost*cosfd*sinf+sinfd*cosf);                             // formula (10)
+  Double_t e  =sint*cosf+tantheta*(cost*cosfd*cosf-sinfd*sinf);                             // formula (9)
+
+  Double_t kk = betaM*TMath::Sqrt(k)/(GapThick()*alpha);                            // formula (6) and (7)
+  Double_t dtdxc = kk*(k*(cosfd*cosf-cost*sinfd*sinf)-(alpha*mu/(betaM*betaM))*sint*sinfd); // formula (6)           
+  Double_t dtdyc = kk*(k*(cosfd*sinf+cost*sinfd*cosf)+(alpha* e/(betaM*betaM))*sint*sinfd); // formula (7)            pag.4
+
+  Double_t errX = 0.2,errY=0.25;                                                            //end of page 7
+  return  TMath::Sqrt(errX*errX*dtdxc*dtdxc + errY*errY*dtdyc*dtdyc);
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Double_t AliHMPIDParam::SigCrom(Double_t trkTheta,Double_t trkPhi,Double_t thetaC, Double_t phiC,Double_t betaM)
+{
+// Analitical calculation of chromatic error (due to lack of knowledge of Cerenkov photon energy) on Cerenkov angle for a given Cerenkov photon 
+// created by a given MIP. Fromulae according to CERN-EP-2000-058 
+// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
+//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
+//            MIP beta
+//   Returns: absolute error on Cerenkov angle, [radians]    
+  
+  Double_t phiDelta = phiC - trkPhi;
+
+  Double_t sint     = TMath::Sin(trkTheta);
+  Double_t cost     = TMath::Cos(trkTheta);
+  Double_t cosfd    = TMath::Cos(phiDelta);
+  Double_t tantheta = TMath::Tan(thetaC);
+  
+  Double_t alpha =cost-tantheta*cosfd*sint;                                                 // formula (11)
+  Double_t dtdn = cost*GetRefIdx()*betaM*betaM/(alpha*tantheta);                    // formula (12)
+            
+//  Double_t f = 0.00928*(7.75-5.635)/TMath::Sqrt(12.);
+  Double_t f = 0.0172*(7.75-5.635)/TMath::Sqrt(24.);
+
+  return f*dtdn;
+}//SigCrom()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+Double_t AliHMPIDParam::SigGeom(Double_t trkTheta,Double_t trkPhi,Double_t thetaC, Double_t phiC,Double_t betaM)
+{
+// Analitical calculation of geometric error (due to lack of knowledge of creation point in radiator) on Cerenkov angle for a given Cerenkov photon 
+// created by a given MIP. Formulae according to CERN-EP-2000-058 
+// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
+//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
+//            MIP beta
+//   Returns: absolute error on Cerenkov angle, [radians]    
+
+  Double_t phiDelta = phiC - trkPhi;
+
+  Double_t sint     = TMath::Sin(trkTheta);
+  Double_t cost     = TMath::Cos(trkTheta);
+  Double_t sinf     = TMath::Sin(trkPhi);
+  Double_t cosfd    = TMath::Cos(phiDelta);
+  Double_t costheta = TMath::Cos(thetaC);
+  Double_t tantheta = TMath::Tan(thetaC);
+  
+  Double_t alpha =cost-tantheta*cosfd*sint;                                                  // formula (11)
+  
+  Double_t k = 1.-GetRefIdx()*GetRefIdx()+alpha*alpha/(betaM*betaM);         // formula (after 8 in the text)
+  if (k<0) return 1e10;
+
+  Double_t eTr = 0.5*RadThick()*betaM*TMath::Sqrt(k)/(GapThick()*alpha);     // formula (14)
+  Double_t lambda = 1.-sint*sint*sinf*sinf;                                                  // formula (15)
+
+  Double_t c1 = 1./(1.+ eTr*k/(alpha*alpha*costheta*costheta));                              // formula (13.a)
+  Double_t c2 = betaM*TMath::Power(k,1.5)*tantheta*lambda/(GapThick()*alpha*alpha);  // formula (13.b)
+  Double_t c3 = (1.+eTr*k*betaM*betaM)/((1+eTr)*alpha*alpha);                                // formula (13.c)
+  Double_t c4 = TMath::Sqrt(k)*tantheta*(1-lambda)/(GapThick()*betaM);               // formula (13.d)
+  Double_t dtdT = c1 * (c2+c3*c4);
+  Double_t trErr = RadThick()/(TMath::Sqrt(12.)*cost);
+
+  return trErr*dtdT;
+}//SigGeom()
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
index bb54fe8..d66b278 100644 (file)
@@ -90,6 +90,11 @@ public:
 
   void     SetRefIdx   (Double_t refRadIdx                                  ) {fRadNmean = refRadIdx;}             //set refractive index of freon
     
+  //For PID
+  Double_t SigLoc      (Double_t trkTheta,Double_t trkPhi,Double_t ckovTh,Double_t ckovPh,Double_t beta);//error due to cathode segmetation
+  Double_t SigGeom     (Double_t trkTheta,Double_t trkPhi,Double_t ckovTh,Double_t ckovPh,Double_t beta);//error due to unknown photon origin
+  Double_t SigCrom     (Double_t trkTheta,Double_t trkPhi,Double_t ckovTh,Double_t ckovPh,Double_t beta);//error due to unknonw photon energy
+  Double_t Sigma2      (Double_t trkTheta,Double_t trkPhi,Double_t ckovTh,Double_t ckovPh              );//photon candidate sigma^2
   
   enum EPlaneId {kPc,kRad,kAnod};            //3 planes in chamber 
   enum ETrackingFlags {kMipDistCut=-9,kMipQdcCut=-5,kNoPhotAccept=-11};     //flags for Reconstruction
diff --git a/HMPID/AliHMPIDPid.cxx b/HMPID/AliHMPIDPid.cxx
new file mode 100644 (file)
index 0000000..0715c35
--- /dev/null
@@ -0,0 +1,74 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+//////////////////////////////////////////////////////////////////////////
+//                                                                      //
+// AliHMPIDPid                                                          //
+//                                                                      //
+// HMPID class to perfom particle identification                        //
+//                                                                      //
+//////////////////////////////////////////////////////////////////////////
+
+#include "AliHMPIDPid.h"       //class header
+#include "AliHMPIDParam.h"     //class header
+#include <AliESDEvent.h>       //FindPid()
+
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+AliHMPIDPid::AliHMPIDPid():TTask("HMPIDrec","HMPIDPid")
+{
+//..
+//init of data members
+//..
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
+void AliHMPIDPid::FindPid(AliESDEvent *pESD)
+{
+// Calculates probability to be a electron-muon-pion-kaon-proton
+// from the given Cerenkov angle and momentum assuming no initial particle composition
+// (i.e. apriory probability to be the particle of the given sort is the same for all sorts)
+
+  AliPID ppp; //needed
+  Double_t pid[AliPID::kSPECIES],h[AliPID::kSPECIES];
+   
+  for(Int_t iTrk=0;iTrk<pESD->GetNumberOfTracks();iTrk++){//ESD tracks loop
+    AliESDtrack *pTrk = pESD->GetTrack(iTrk);// get next reconstructed track
+    if(pTrk->GetHMPIDsignal()<=0){//HMPID does not find anything reasonable for this track, assign 0.2 for all species
+      for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++) pid[iPart]=1.0/AliPID::kSPECIES;
+      pTrk->SetHMPIDpid(pid);
+      continue;
+    } 
+    Double_t pmod = pTrk->GetP();
+    Double_t hTot=0;
+    for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++){
+      Double_t mass = AliPID::ParticleMass(iPart);
+      Double_t cosThetaTh = TMath::Sqrt(mass*mass+pmod*pmod)/(AliHMPIDParam::Instance()->MeanIdxRad()*pmod);
+      if(cosThetaTh<1) //calculate the height of theoretical theta ckov on the gaus of experimental one
+        h[iPart] =TMath::Gaus(TMath::ACos(cosThetaTh),pTrk->GetHMPIDsignal(),TMath::Sqrt(pTrk->GetHMPIDchi2()),kTRUE);      
+      else             //beta < 1/ref. idx. => no light at all  
+        h[iPart] =0 ;       
+      hTot    +=h[iPart]; //total height of all theoretical heights for normalization
+    }//species loop
+     
+    Double_t hMin=TMath::Gaus(pTrk->GetHMPIDsignal()-4*TMath::Sqrt(pTrk->GetHMPIDchi2()),pTrk->GetHMPIDsignal(),TMath::Sqrt(pTrk->GetHMPIDchi2()),kTRUE);//5 sigma protection
+    
+    for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++)//species loop to assign probabilities
+      if(hTot>hMin)  
+        pid[iPart]=h[iPart]/hTot;
+      else                               //all theoretical values are far away from experemental one
+        pid[iPart]=1.0/AliPID::kSPECIES; 
+    pTrk->SetHMPIDpid(pid);
+  }//ESD tracks loop
+}
+//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
diff --git a/HMPID/AliHMPIDPid.h b/HMPID/AliHMPIDPid.h
new file mode 100644 (file)
index 0000000..b6b43b8
--- /dev/null
@@ -0,0 +1,39 @@
+#ifndef AliHMPIDPid_h
+#define AliHMPIDPid_h
+
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+//////////////////////////////////////////////////////////////////////////
+//                                                                      //
+// AliHMPIDPid                                                          //
+//                                                                      //
+// HMPID class to perfom pattern recognition based on Hough transfrom   //
+//                                                                      //
+//////////////////////////////////////////////////////////////////////////
+
+
+#include <TTask.h>        //base class
+
+class AliESDEvent;
+
+class AliHMPIDPid : public TTask 
+{
+public : 
+             AliHMPIDPid();    //ctor
+    virtual ~AliHMPIDPid() {;} //dtor
+    
+    void     FindPid(AliESDEvent *pESD);             //Find PID for tracks
+
+//
+protected:
+  
+private:
+  AliHMPIDPid(const AliHMPIDPid& r);                //dummy copy constructor
+  AliHMPIDPid &operator=(const AliHMPIDPid& r);     //dummy assignment operator
+//
+  ClassDef(AliHMPIDPid,0)
+};
+
+#endif // #ifdef AliHMPIDPid_cxx
+
index cf23978..0f3257c 100644 (file)
@@ -324,7 +324,7 @@ Double_t AliHMPIDRecon::FindRingCkov(Int_t)
       weightThetaCerenkov += fPhotCkov[i]*fPhotWei[i];
       wei += fPhotWei[i];                                                    //collect weight as sum of all candidate weghts   
       
-      sigma2 += 1./Sigma2(fPhotCkov[i],fPhotPhi[i]);
+      sigma2 += 1./fParam->Sigma2(fTrkDir.Theta(),fTrkDir.Phi(),fPhotCkov[i],fPhotPhi[i]);
     }
   }//candidates loop
   
@@ -456,118 +456,3 @@ Double_t AliHMPIDRecon::HoughResponse()
   return (Double_t)(locMax*fDTheta+0.5*fDTheta); //final most probable track theta ckov   
 }//HoughResponse()
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Double_t AliHMPIDRecon::Sigma2(Double_t ckovTh, Double_t ckovPh)const
-{
-// Analithical calculation of total error (as a sum of localization, geometrical and chromatic errors) on Cerenkov angle for a given Cerenkov photon 
-// created by a given MIP. Fromulae according to CERN-EP-2000-058 
-// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
-//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
-//            MIP beta
-//   Returns: absolute error on Cerenkov angle, [radians]    
-  
-  TVector3 v(-999,-999,-999);
-  Double_t trkBeta = 1./(TMath::Cos(ckovTh)*fParam->GetRefIdx());
-  
-  if(trkBeta > 1) trkBeta = 1;                 //protection against bad measured thetaCer  
-  if(trkBeta < 0) trkBeta = 0.0001;            //
-
-  v.SetX(SigLoc (ckovTh,ckovPh,trkBeta));
-  v.SetY(SigGeom(ckovTh,ckovPh,trkBeta));
-  v.SetZ(SigCrom(ckovTh,ckovPh,trkBeta));
-
-  return v.Mag2();
-}
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Double_t AliHMPIDRecon::SigLoc(Double_t thetaC, Double_t phiC,Double_t betaM)const
-{
-// Analithical calculation of localization error (due to finite segmentation of PC) on Cerenkov angle for a given Cerenkov photon 
-// created by a given MIP. Fromulae according to CERN-EP-2000-058 
-// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
-//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
-//            MIP beta
-//   Returns: absolute error on Cerenkov angle, [radians]    
-  
-  Double_t phiDelta = phiC - fTrkDir.Phi();
-
-  Double_t sint     = TMath::Sin(fTrkDir.Theta());
-  Double_t cost     = TMath::Cos(fTrkDir.Theta());
-  Double_t sinf     = TMath::Sin(fTrkDir.Phi());
-  Double_t cosf     = TMath::Cos(fTrkDir.Phi());
-  Double_t sinfd    = TMath::Sin(phiDelta);
-  Double_t cosfd    = TMath::Cos(phiDelta);
-  Double_t tantheta = TMath::Tan(thetaC);
-  
-  Double_t alpha =cost-tantheta*cosfd*sint;                                                 // formula (11)
-  Double_t k = 1.-fParam->GetRefIdx()*fParam->GetRefIdx()+alpha*alpha/(betaM*betaM);        // formula (after 8 in the text)
-  if (k<0) return 1e10;
-  Double_t mu =sint*sinf+tantheta*(cost*cosfd*sinf+sinfd*cosf);                             // formula (10)
-  Double_t e  =sint*cosf+tantheta*(cost*cosfd*cosf-sinfd*sinf);                             // formula (9)
-
-  Double_t kk = betaM*TMath::Sqrt(k)/(fParam->GapThick()*alpha);                            // formula (6) and (7)
-  Double_t dtdxc = kk*(k*(cosfd*cosf-cost*sinfd*sinf)-(alpha*mu/(betaM*betaM))*sint*sinfd); // formula (6)           
-  Double_t dtdyc = kk*(k*(cosfd*sinf+cost*sinfd*cosf)+(alpha* e/(betaM*betaM))*sint*sinfd); // formula (7)            pag.4
-
-  Double_t errX = 0.2,errY=0.25;                                                            //end of page 7
-  return  TMath::Sqrt(errX*errX*dtdxc*dtdxc + errY*errY*dtdyc*dtdyc);
-}
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Double_t AliHMPIDRecon::SigCrom(Double_t thetaC, Double_t phiC,Double_t betaM)const
-{
-// Analithical calculation of chromatic error (due to lack of knowledge of Cerenkov photon energy) on Cerenkov angle for a given Cerenkov photon 
-// created by a given MIP. Fromulae according to CERN-EP-2000-058 
-// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
-//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
-//            MIP beta
-//   Returns: absolute error on Cerenkov angle, [radians]    
-  
-  Double_t phiDelta = phiC - fTrkDir.Phi();
-
-  Double_t sint     = TMath::Sin(fTrkDir.Theta());
-  Double_t cost     = TMath::Cos(fTrkDir.Theta());
-  Double_t cosfd    = TMath::Cos(phiDelta);
-  Double_t tantheta = TMath::Tan(thetaC);
-  
-  Double_t alpha =cost-tantheta*cosfd*sint;                                                 // formula (11)
-  Double_t dtdn = cost*fParam->GetRefIdx()*betaM*betaM/(alpha*tantheta);                    // formula (12)
-            
-//  Double_t f = 0.00928*(7.75-5.635)/TMath::Sqrt(12.);
-  Double_t f = 0.0172*(7.75-5.635)/TMath::Sqrt(24.);
-
-  return f*dtdn;
-}//SigCrom()
-//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
-Double_t AliHMPIDRecon::SigGeom(Double_t thetaC, Double_t phiC,Double_t betaM)const
-{
-// Analithical calculation of geometric error (due to lack of knowledge of creation point in radiator) on Cerenkov angle for a given Cerenkov photon 
-// created by a given MIP. Formulae according to CERN-EP-2000-058 
-// Arguments: Cerenkov and azimuthal angles for Cerenkov photon, [radians]
-//            dip and azimuthal angles for MIP taken at the entrance to radiator, [radians]        
-//            MIP beta
-//   Returns: absolute error on Cerenkov angle, [radians]    
-
-  Double_t phiDelta = phiC - fTrkDir.Phi();
-
-  Double_t sint     = TMath::Sin(fTrkDir.Theta());
-  Double_t cost     = TMath::Cos(fTrkDir.Theta());
-  Double_t sinf     = TMath::Sin(fTrkDir.Phi());
-  Double_t cosfd    = TMath::Cos(phiDelta);
-  Double_t costheta = TMath::Cos(thetaC);
-  Double_t tantheta = TMath::Tan(thetaC);
-  
-  Double_t alpha =cost-tantheta*cosfd*sint;                                                  // formula (11)
-  
-  Double_t k = 1.-fParam->GetRefIdx()*fParam->GetRefIdx()+alpha*alpha/(betaM*betaM);         // formula (after 8 in the text)
-  if (k<0) return 1e10;
-
-  Double_t eTr = 0.5*fParam->RadThick()*betaM*TMath::Sqrt(k)/(fParam->GapThick()*alpha);     // formula (14)
-  Double_t lambda = 1.-sint*sint*sinf*sinf;                                                  // formula (15)
-
-  Double_t c1 = 1./(1.+ eTr*k/(alpha*alpha*costheta*costheta));                              // formula (13.a)
-  Double_t c2 = betaM*TMath::Power(k,1.5)*tantheta*lambda/(fParam->GapThick()*alpha*alpha);  // formula (13.b)
-  Double_t c3 = (1.+eTr*k*betaM*betaM)/((1+eTr)*alpha*alpha);                                // formula (13.c)
-  Double_t c4 = TMath::Sqrt(k)*tantheta*(1-lambda)/(fParam->GapThick()*betaM);               // formula (13.d)
-  Double_t dtdT = c1 * (c2+c3*c4);
-  Double_t trErr = fParam->RadThick()/(TMath::Sqrt(12.)*cost);
-
-  return trErr*dtdT;
-}//SigGeom()
index ba0a12c..185d634 100644 (file)
@@ -53,10 +53,6 @@ public :
                                 {fPc.Set(xPc,yPc);}                                                //set track impact to PC 
   void     SetMip       (Double_t xmip,Double_t ymip                                        )
                                 {fMipPos.Set(xmip,ymip);}                                          //set track impact to PC
-  Double_t SigLoc       (Double_t ckovTh,Double_t ckovPh,Double_t beta                      )const;//error due to cathode segmetation
-  Double_t SigGeom      (Double_t ckovTh,Double_t ckovPh,Double_t beta                      )const;//error due to unknown photon origin
-  Double_t SigCrom      (Double_t ckovTh,Double_t ckovPh,Double_t beta                      )const;//error due to unknonw photon energy
-  Double_t Sigma2       (Double_t ckovTh,Double_t ckovPh                                    )const;//photon candidate sigma^2
   enum ETrackingFlags {kMipDistCut=-9,kMipQdcCut=-5,kNoPhotAccept=-11};
 //
 protected:
index 2265284..e040af2 100644 (file)
@@ -20,6 +20,7 @@
 #include "AliHMPIDReconstructor.h" //class header
 #include "AliHMPID.h"              //Reconstruct() 
 #include "AliHMPIDCluster.h"       //Dig2Clu()
+#include "AliHMPIDPid.h"           //FillEsd() 
 #include "AliHMPIDParam.h"         //FillEsd() 
 #include <AliCDBEntry.h>           //ctor
 #include <AliCDBManager.h>         //ctor
@@ -193,39 +194,12 @@ void AliHMPIDReconstructor::ConvertDigits(AliRawReader *pRR,TTree *pDigTree)cons
 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 void AliHMPIDReconstructor::FillESD(TTree */*digitsTree*/, TTree */*clustersTree*/, AliESDEvent *pESD) const
 {
-// Calculates probability to be a electron-muon-pion-kaon-proton
-// from the given Cerenkov angle and momentum assuming no initial particle composition
-// (i.e. apriory probability to be the particle of the given sort is the same for all sorts)
+// Fill ESD with all the infos from HMPID
+// Probability vector from AliHMPIDPid
+//...
+  
+  AliHMPIDPid *pPid = new AliHMPIDPid();
+  pPid->FindPid(pESD);
+  delete pPid;
 
-  AliPID ppp; //needed
-  Double_t pid[AliPID::kSPECIES],h[AliPID::kSPECIES];
-   
-  for(Int_t iTrk=0;iTrk<pESD->GetNumberOfTracks();iTrk++){//ESD tracks loop
-    AliESDtrack *pTrk = pESD->GetTrack(iTrk);// get next reconstructed track
-    if(pTrk->GetHMPIDsignal()<=0){//HMPID does not find anything reasonable for this track, assign 0.2 for all species
-      for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++) pid[iPart]=1.0/AliPID::kSPECIES;
-      pTrk->SetHMPIDpid(pid);
-      continue;
-    } 
-    Double_t pmod = pTrk->GetP();
-    Double_t hTot=0;
-    for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++){
-      Double_t mass = AliPID::ParticleMass(iPart);
-      Double_t cosThetaTh = TMath::Sqrt(mass*mass+pmod*pmod)/(AliHMPIDParam::Instance()->MeanIdxRad()*pmod);
-      if(cosThetaTh<1) //calculate the height of theoretical theta ckov on the gaus of experimental one
-        h[iPart] =TMath::Gaus(TMath::ACos(cosThetaTh),pTrk->GetHMPIDsignal(),TMath::Sqrt(pTrk->GetHMPIDchi2()),kTRUE);      
-      else             //beta < 1/ref. idx. => no light at all  
-        h[iPart] =0 ;       
-      hTot    +=h[iPart]; //total height of all theoretical heights for normalization
-    }//species loop
-     
-    Double_t hMin=TMath::Gaus(pTrk->GetHMPIDsignal()-4*TMath::Sqrt(pTrk->GetHMPIDchi2()),pTrk->GetHMPIDsignal(),TMath::Sqrt(pTrk->GetHMPIDchi2()),kTRUE);//5 sigma protection
-    
-    for(Int_t iPart=0;iPart<AliPID::kSPECIES;iPart++)//species loop to assign probabilities
-      if(hTot>hMin)  
-        pid[iPart]=h[iPart]/hTot;
-      else                               //all theoretical values are far away from experemental one
-        pid[iPart]=1.0/AliPID::kSPECIES; 
-    pTrk->SetHMPIDpid(pid);
-  }//ESD tracks loop
 }//FillESD()
index cd6ae28..c7e9a8a 100644 (file)
@@ -6,6 +6,7 @@
 #pragma link C++ class  AliHMPIDReconstructor+;
 #pragma link C++ class  AliHMPIDTracker+;
 #pragma link C++ class  AliHMPIDRecon+;
+#pragma link C++ class  AliHMPIDPid+;
 #pragma link C++ class  AliHMPIDRecoParam+;
 #pragma link C++ class  AliHMPIDReconHTA+;
 #endif
index d661f4b..169ac7b 100644 (file)
@@ -303,10 +303,11 @@ void RenderEsd(AliESDEvent *pEsd)
     Float_t thRa,phRa,xRa,yRa; pTrk->GetHMPIDtrk(xRa,yRa,thRa,phRa);            //get info on current track
     ch/=1000000;                            
     Float_t xPc=0,yPc=0; AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc);              //find again intersection of track with PC--> it is not stored in ESD!
-    fRenTxC[ch]->SetNextPoint(xPc,yPc);                                         //add this intersection point
+    Int_t npTrk = fRenTxC[ch]->SetNextPoint(xPc,yPc);                           //add this intersection point
     Float_t ckov=pTrk->GetHMPIDsignal();                                        //get ckov angle stored for this track  
     if(ckov>0){
       rec.SetTrack(xRa,yRa,thRa,phRa);
+      fRenRin[ch]->SetUniqueID(npTrk);
       for(Int_t j=0;j<500;j++){
         TVector2 pos; pos=rec.TracePhot(ckov,j*6.28/500.);
        if(!AliHMPIDParam::IsInDead(pos.X(),pos.Y())) fRenRin[ch]->SetNextPoint(pos.X(),pos.Y());
@@ -440,7 +441,7 @@ void DisplayInfo(Int_t evt,Int_t px, Int_t py)
     text0.Append(Form("Pad(%i,%i)-LORS(%6.2f,%6.2f)-MARS(%7.2f,%7.2f,%7.2f)",padX,padY,x,y,xyz.X(),xyz.Y(),xyz.Z()));
     text1.Append(Form("Module %i Sector %i",ch,pc));
     text2="";
-    text3.Append(Form("Pads = %i - Clusters = %i - Multiplicity %5.2f%%",fTotPads[ch],fTotClus[ch],100.*fTotPads[ch]/(144.*160.)));
+    text3.Append(Form("Pads = %i - Clusters = %i - Occupancy %5.2f%%",fTotPads[ch],fTotClus[ch],100.*fTotPads[ch]/(144.*160.)));
   }
   
   TObject *obj = fCanvas->GetSelected();
@@ -543,16 +544,20 @@ void DisplayInfo(Int_t evt,Int_t px, Int_t py)
     text0="";text0.Append(Form("CLUSTER: x %6.2f y %6.2f",pClu->X(),pClu->Y()));
     text2="";text2.Append(Form("charge = %i ADC",(Int_t)pClu->Q()));
     }
-  else if (symbol==kTrack) {
+  else if (symbol==kTrack || symbol==kRing) {
+    if(symbol==kRing) index = b->GetUniqueID();
     AliESDtrack *pTrk=fEsd->GetTrack(index);
     Float_t thRa,phRa,xRa,yRa; pTrk->GetHMPIDtrk(xRa,yRa,thRa,phRa);
     Float_t xPc,yPc; AliHMPIDTracker::IntTrkCha(pTrk,xPc,yPc);
-    text0="";text0.Append(Form("TRACK: x %6.2f y %6.2f at PC plane",xPc,yPc));
+    text0="";text0.Append(Form("TRACK n.%d: x %6.2f y %6.2f at PC plane",index,xPc,yPc));
     text2="";text2.Append(Form("p = %7.2f GeV/c",pTrk->GetP()));
     Float_t ckov=pTrk->GetHMPIDsignal();                             
+    Double_t prob[5];
+    pTrk->GetHMPIDpid(prob);
     if(ckov>0){
       Float_t x,y;Int_t q,nacc;   pTrk->GetHMPIDmip(x,y,q,nacc);
-      text3="";text3.Append(Form("Theta Cherenkov %5.3f with %i photons",pTrk->GetHMPIDsignal(),nacc));
+      text3="";text3.Append(Form("Theta Cherenkov %5.3f with %i photons |  e %5.1f%% | u %5.1f%% | K %5.1f%% | pi %5.1f%% | p %5.1f%%",
+          pTrk->GetHMPIDsignal(),nacc,prob[0]*100,prob[1]*100,prob[2]*100,prob[3]*100,prob[4]*100));
     }      
   }
 //Update toolbar status  
index 9ec328a..da6cbc6 100644 (file)
@@ -1,4 +1,6 @@
-SRCS:=  AliHMPIDReconstructor.cxx AliHMPIDTracker.cxx AliHMPIDRecon.cxx AliHMPIDRecoParam.cxx AliHMPIDReconHTA.cxx
+SRCS:=  AliHMPIDReconstructor.cxx AliHMPIDTracker.cxx \
+        AliHMPIDRecon.cxx AliHMPIDRecoParam.cxx       \
+        AliHMPIDReconHTA.cxx AliHMPIDPid.cxx
 
 HDRS:= $(SRCS:.cxx=.h)
 DHDR:= HMPIDrecLinkDef.h