-Added base Calo QA class.

[u/mrichter/AliRoot.git] / PHOS / AliPHOSAnalyze.cxx

diff --git a/PHOS/AliPHOSAnalyze.cxx b/PHOS/AliPHOSAnalyze.cxx
index 2f8da5d4a73540f9e11c0a8d7a4afa92369ede9d..69c72ca7c8edf8400a2bc7a08d0b73228d5359c7 100644 (file)
--- a/PHOS/AliPHOSAnalyze.cxx
+++ b/PHOS/AliPHOSAnalyze.cxx
@@ -68,6 +68,7 @@
 #include "TH2.h"
 #include "TH2.h"
 #include "TParticle.h"
+#include "TDatabasePDG.h"
 #include "TClonesArray.h"
 #include "TMath.h"
 #include "TROOT.h"
@@ -92,19 +93,23 @@
 ClassImp(AliPHOSAnalyze)
 
 //____________________________________________________________________________
-  AliPHOSAnalyze::AliPHOSAnalyze()
+AliPHOSAnalyze::AliPHOSAnalyze():
+  fCorrection(1.2),  //Value calculated for default parameters of reconstruction
+  fEvt(0),
+  ffileName(),
+  fRunLoader(0)
 {
   // default ctor (useless)
-  fCorrection = 1.2 ;  //Value calculated for default parameters of reconstruction  
-  fRunLoader = 0x0;
 }
 
 //____________________________________________________________________________
-AliPHOSAnalyze::AliPHOSAnalyze(Text_t * fileName)
+AliPHOSAnalyze::AliPHOSAnalyze(Text_t * fileName):
+  fCorrection(1.05),  //Value calculated for default parameters of reconstruction   
+  fEvt(0),
+  ffileName(fileName),
+  fRunLoader(0)
 {
   // ctor: analyze events from root file "name"
-  ffileName = fileName;
-  fCorrection = 1.05 ;  //Value calculated for default parameters of reconstruction   
   fRunLoader = AliRunLoader::Open(fileName,"AliPHOSAnalyze");
   if (fRunLoader == 0x0)
    {
@@ -113,8 +118,12 @@ AliPHOSAnalyze::AliPHOSAnalyze(Text_t * fileName)
 }
 
 //____________________________________________________________________________
-AliPHOSAnalyze::AliPHOSAnalyze(const AliPHOSAnalyze & ana)
-  : TObject(ana)
+AliPHOSAnalyze::AliPHOSAnalyze(const AliPHOSAnalyze & ana): 
+  TObject(ana),
+  fCorrection(0.),
+  fEvt(0),
+  ffileName(),
+  fRunLoader(0)
 {
   // copy ctor
   ( (AliPHOSAnalyze &)ana ).Copy(*this) ;
@@ -151,12 +160,12 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
     AliError(Form("There is no event %d only %d events available", Nevent, fRunLoader->GetNumberOfEvents() )) ;
     return ;
   }
-  const AliPHOSGeometry * phosgeom = gime->PHOSGeometry() ; 
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
   fRunLoader->GetEvent(Nevent);
 
   Int_t nx = phosgeom->GetNPhi() ;
   Int_t nz = phosgeom->GetNZ() ;
-  Float_t * cri= phosgeom->GetEMCAGeometry()->GetCrystalHalfSize() ;
+  const Float_t * cri= phosgeom->GetEMCAGeometry()->GetCrystalHalfSize() ;
   Float_t x = nx*cri[0] ;
   Float_t z = nz*cri[2] ;
   Int_t nxCPV = (Int_t) (nx*phosgeom->GetPadSizePhi()/(2.*cri[0])) ;
@@ -196,6 +205,10 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
     recPhot->Delete() ;
   recPhot = new TH2F("recPhot","RecParticles with primary Photon",nx,-x,x,nz,-z,z);
   
+  //Get Vertex
+  Double_t vtx[3]={0.,0.,0.} ;  
+//DP: extract vertex either from Generator or from data
+
   
   //Plot Primary Particles
   
@@ -220,7 +233,7 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
       if( primaryType == 22 ) {
         Int_t moduleNumber ;
         Double_t primX, primZ ;
-        phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+        phosgeom->ImpactOnEmc(vtx,primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
         if(moduleNumber==Nmod) 
           phot->Fill(primZ,primX,primary->Energy()) ;
       }
@@ -246,16 +259,15 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
        sdigit = (AliPHOSDigit *) sdigits->At(iSDigit) ;
        Int_t relid[4];
        phosgeom->AbsToRelNumbering(sdigit->GetId(), relid) ;
-       Float_t x,z ;
-       phosgeom->RelPosInModule(relid,x,z);
-       AliPHOSSDigitizer* sd = dynamic_cast<AliPHOSSDigitizer*>(gime->SDigitizer());
-       Float_t e = sd->Calibrate(sdigit->GetAmp()) ;
+       Float_t xd,zd ;
+       phosgeom->RelPosInModule(relid,xd,zd);
+       Float_t e = sdigit->GetEnergy() ;
        nsdig[relid[0]-1]++ ;
        if(relid[0]==Nmod){
          if(relid[1]==0)  //EMC
-           emcSdigits->Fill(x,z,e) ;
+           emcSdigits->Fill(xd,zd,e) ;
          if( relid[1]!=0 )
-           cpvSdigits->Fill(x,z,e) ;
+           cpvSdigits->Fill(xd,zd,e) ;
        }
       }
   }
@@ -274,15 +286,14 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
        digit = (AliPHOSDigit *) digits->At(iDigit) ;
        Int_t relid[4];
        phosgeom->AbsToRelNumbering(digit->GetId(), relid) ;
-       Float_t x,z ;
-       phosgeom->RelPosInModule(relid,x,z) ;
-       AliPHOSSDigitizer* sd = dynamic_cast<AliPHOSSDigitizer*>(gime->SDigitizer());
-       Float_t e = sd->Calibrate(digit->GetAmp()) ;
+       Float_t xd,zd ;
+       phosgeom->RelPosInModule(relid,xd,zd) ;
+       Float_t e = digit->GetEnergy() ;
        if(relid[0]==Nmod){
          if(relid[1]==0)  //EMC
-           emcDigits->Fill(x,z,e) ;
+           emcDigits->Fill(xd,zd,e) ;
          if( relid[1]!=0 )
-           cpvDigits->Fill(x,z,e) ;
+           cpvDigits->Fill(xd,zd,e) ;
        }
       }
   }
@@ -324,7 +335,7 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
        recParticle = (AliPHOSRecParticle *) rp->At(iRecParticle) ;
        Int_t moduleNumberRec ;
        Double_t recX, recZ ;
-       phosgeom->ImpactOnEmc(recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
+       phosgeom->ImpactOnEmc(vtx,recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
        if(moduleNumberRec == Nmod){
          
          Double_t minDistance = 5. ;
@@ -335,14 +346,14 @@ void AliPHOSAnalyze::DrawRecon(Int_t Nevent,Int_t Nmod){
          Int_t numberofprimaries ;
          Int_t * listofprimaries  = ((AliPHOSRecPoint*) emcrp->At(emcIndex))->GetPrimaries(numberofprimaries)  ;
          Int_t index ;
-         const TParticle * primary ;
+         const TParticle * primPart ;
          Double_t distance = minDistance ;
          
          for ( index = 0 ; index < numberofprimaries ; index++){
-           primary = fRunLoader->Stack()->Particle(listofprimaries[index]) ;
+           primPart = fRunLoader->Stack()->Particle(listofprimaries[index]) ;
            Int_t moduleNumber ;
            Double_t primX, primZ ;
-           phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+           phosgeom->ImpactOnEmc(vtx,primPart->Theta(), primPart->Phi(), moduleNumber, primX, primZ) ;
            if(moduleNumberRec == moduleNumber)
              distance = TMath::Sqrt((recX-primX)*(recX-primX)+(recZ-primZ)*(recZ-primZ) ) ;
            if(minDistance > distance)
@@ -558,7 +569,6 @@ void AliPHOSAnalyze::Ls(){
     nRecParticles[mevent] = iRecPhot-1 ;  
     
     //check, if it is time to calculate invariant mass?
-    Int_t maxevent = (Int_t)gAlice->TreeE()->GetEntries() ; 
     if((mevent == 0) && (event +1 == maxevent)){
       
       //   if((mevent == 0) && (event +1 == gime->MaxEvent())){
@@ -669,7 +679,7 @@ void AliPHOSAnalyze::Ls(){
    } 
 
 
-  const AliPHOSGeometry * phosgeom = gime->PHOSGeometry();
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
 
   Int_t ievent;
   Int_t maxevent = (Int_t)fRunLoader->TreeE()->GetEntries();
@@ -682,6 +692,8 @@ void AliPHOSAnalyze::Ls(){
     //read the current event
     fRunLoader->GetEvent(ievent) ;
 
+    Double_t vtx[3]={0.,0.,0.} ;  
+
     const AliPHOSRecParticle * recParticle ;
     Int_t iRecParticle ;
     TClonesArray * rp = gime->RecParticles() ;
@@ -706,7 +718,7 @@ void AliPHOSAnalyze::Ls(){
       //find the closest primary
       Int_t moduleNumberRec ;
       Double_t recX, recZ ;
-      phosgeom->ImpactOnEmc(recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
+      phosgeom->ImpactOnEmc(vtx,recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
       
       Double_t minDistance  = 100. ;
       Int_t closestPrimary = -1 ;
@@ -728,7 +740,7 @@ void AliPHOSAnalyze::Ls(){
 
        Int_t moduleNumber ;
        Double_t primX, primZ ;
-       phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+       phosgeom->ImpactOnEmc(vtx,primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
        if(moduleNumberRec == moduleNumber) {
          dX = recX - primX;
          dZ = recZ - primZ;
@@ -744,7 +756,7 @@ void AliPHOSAnalyze::Ls(){
 
       //if found primary, fill histograms
       if(closestPrimary >=0 ){
-       const TParticle * primary = fRunLoader->Stack()->Particle(closestPrimary) ;
+       primary = fRunLoader->Stack()->Particle(closestPrimary) ;
        if(primary->GetPdgCode() == 22){
          hAllEnergy->Fill(primary->Energy(), recParticle->Energy()) ;
          if(recParticle->GetType() == AliPHOSFastRecParticle::kNEUTRALEMFAST){
@@ -823,7 +835,7 @@ void AliPHOSAnalyze::PositionResolution()
   
   if (fRunLoader->TreeE() == 0x0) fRunLoader->LoadHeader();
 
-  const AliPHOSGeometry * phosgeom = gime->PHOSGeometry() ; 
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
 
   Int_t ievent;
   Int_t maxevent = (Int_t)fRunLoader->TreeE()->GetEntries() ; 
@@ -831,6 +843,10 @@ void AliPHOSAnalyze::PositionResolution()
     
     //read the current event
     fRunLoader->GetEvent(ievent) ;
+
+    //DP:Extract vertex position
+    Double_t vtx[3]={0.,0.,0.} ;  
+
     TClonesArray * rp = gime->RecParticles() ;
     if(!rp) {
       AliError(Form("Event %d,  Can't find RecParticles", ievent)) ;
@@ -856,7 +872,7 @@ void AliPHOSAnalyze::PositionResolution()
       //find the closest primary
       Int_t moduleNumberRec ;
       Double_t recX, recZ ;
-      phosgeom->ImpactOnEmc(recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
+      phosgeom->ImpactOnEmc(vtx,recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
       
       Double_t minDistance  = 100. ;
       Int_t closestPrimary = -1 ;
@@ -877,7 +893,7 @@ void AliPHOSAnalyze::PositionResolution()
        primary = fRunLoader->Stack()->Particle(listofprimaries[index]) ;
        Int_t moduleNumber ;
        Double_t primX, primZ ;
-       phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+       phosgeom->ImpactOnEmc(vtx,primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
        if(moduleNumberRec == moduleNumber) {
          dX = recX - primX;
          dZ = recZ - primZ;
@@ -893,7 +909,7 @@ void AliPHOSAnalyze::PositionResolution()
       
       //if found primary, fill histograms
       if(closestPrimary >=0 ){
-       const TParticle * primary = fRunLoader->Stack()->Particle(closestPrimary) ;
+       primary = fRunLoader->Stack()->Particle(closestPrimary) ;
        if(primary->GetPdgCode() == 22){
          hAllPosition->Fill(primary->Energy(), minDistance) ;
          hAllPositionX->Fill(primary->Energy(), dX) ;
@@ -1044,7 +1060,7 @@ void AliPHOSAnalyze::Contamination(){
    } 
   
   if (fRunLoader->TreeE() == 0x0) fRunLoader->LoadHeader();
-  const AliPHOSGeometry * phosgeom = gime->PHOSGeometry() ; 
+  AliPHOSGeometry * phosgeom =  AliPHOSGeometry::GetInstance() ;
   
   Int_t ievent;
   Int_t maxevent = (Int_t)fRunLoader->TreeE()->GetEntries() ; 
@@ -1052,6 +1068,9 @@ void AliPHOSAnalyze::Contamination(){
     
     fRunLoader->GetEvent(ievent) ;
     
+    //DP:Extract vertex position
+    Double_t vtx[3]={0.,0.,0.} ;
+
     TClonesArray * rp = gime->RecParticles() ;
     if(!rp) {
       AliError(Form("Event %d,  Can't find RecParticles", ievent)) ;
@@ -1079,7 +1098,7 @@ void AliPHOSAnalyze::Contamination(){
        //check, if photons folls onto PHOS
        Int_t moduleNumber ;
        Double_t primX, primZ ;
-       phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+       phosgeom->ImpactOnEmc(vtx,primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
        if(moduleNumber)
          hPrimary->Fill(primary->Energy()) ;
        
@@ -1098,7 +1117,7 @@ void AliPHOSAnalyze::Contamination(){
       //==========find the closest primary       
       Int_t moduleNumberRec ;
       Double_t recX, recZ ;
-      phosgeom->ImpactOnEmc(recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
+      phosgeom->ImpactOnEmc(vtx,recParticle->Theta(), recParticle->Phi(), moduleNumberRec, recX, recZ) ;
       
       Double_t minDistance  = 100. ;
       Int_t closestPrimary = -1 ;
@@ -1108,7 +1127,6 @@ void AliPHOSAnalyze::Contamination(){
       Int_t numberofprimaries ;
       Int_t * listofprimaries  = ((AliPHOSEmcRecPoint *) emcrp->At(emcIndex))->GetPrimaries(numberofprimaries)  ;
       Int_t index ;
-      const TParticle * primary ;
       Double_t distance = minDistance ;
       Double_t dX, dZ; 
       Double_t dXmin = 0.; 
@@ -1117,7 +1135,7 @@ void AliPHOSAnalyze::Contamination(){
        primary = fRunLoader->Stack()->Particle(listofprimaries[index]) ;
        Int_t moduleNumber ;
        Double_t primX, primZ ;
-       phosgeom->ImpactOnEmc(primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
+       phosgeom->ImpactOnEmc(vtx,primary->Theta(), primary->Phi(), moduleNumber, primX, primZ) ;
        if(moduleNumberRec == moduleNumber) {
          dX = recX - primX;
          dZ = recZ - primZ;
@@ -1134,7 +1152,7 @@ void AliPHOSAnalyze::Contamination(){
       //===========define the "type" of closest primary
       if(closestPrimary >=0 ){
        Int_t primaryCode = -1;
-       const TParticle * primary = fRunLoader->Stack()->Particle(closestPrimary) ;
+       primary = fRunLoader->Stack()->Particle(closestPrimary) ;
        Int_t primaryType = primary->GetPdgCode() ;
        if(primaryType == 22) // photon ?
          primaryCode = 0 ;
@@ -1243,7 +1261,7 @@ void AliPHOSAnalyze::Contamination(){
  
   
   //print Final Table
-  maxevent = (Int_t)gAlice->TreeE()->GetEntries() ; 
+  maxevent = (Int_t)AliRunLoader::Instance()->TreeE()->GetEntries() ; 
 
   TString message ; 
   message  = "Resolutions: Analyzed %d event(s)\n" ;