]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG2/FLOW/AliFlowCommon/AliFlowAnalysisWithLeeYangZeros.cxx
fill the refmult
[u/mrichter/AliRoot.git] / PWG2 / FLOW / AliFlowCommon / AliFlowAnalysisWithLeeYangZeros.cxx
index 14ecb243339c8f417b61cc71f3f460a13d8c5a6c..ff29c7f61fd88a83b156dd5fdc2bf3058dd6f719 100644 (file)
 //Author: Naomi van der Kolk (kolk@nikhef.nl)
 /////////////////////////////////////////////////////////////////////////////////////////
 
-#include "Riostream.h"                 //needed as include
-#include "TObject.h"                   //needed as include
-#include "AliFlowCommonConstants.h"    //needed as include
-#include "AliFlowLYZConstants.h"       //needed as include
-#include "AliFlowAnalysisWithLeeYangZeros.h"
-#include "AliFlowLYZHist1.h"           //needed as include
-#include "AliFlowLYZHist2.h"           //needed as include
-#include "AliFlowCommonHist.h"         //needed as include
-#include "AliFlowCommonHistResults.h"  //needed as include
-#include "AliFlowEventSimple.h"        //needed as include
-#include "AliFlowTrackSimple.h"        //needed as include
-
-class AliFlowVector;
-
-#include "TMath.h" //needed as include
-#include "TFile.h" //needed as include
+#include "Riostream.h"
+#include "TObject.h" 
+#include "TMath.h" 
+#include "TFile.h"
 #include "TList.h"
-
-class TComplex;
-class TProfile;
-class TH1F;
-class TH1D;
-
+#include "TVector2.h"
+#include "TH1F.h"
+#include "TComplex.h"
+#include "TProfile.h"
+#include "TDirectoryFile.h"
+
+#include "AliFlowCommonConstants.h"
+#include "AliFlowLYZConstants.h"
+#include "AliFlowAnalysisWithLeeYangZeros.h"
+#include "AliFlowLYZHist1.h"
+#include "AliFlowLYZHist2.h"
+#include "AliFlowCommonHist.h"
+#include "AliFlowCommonHistResults.h"
+#include "AliFlowEventSimple.h"
+#include "AliFlowTrackSimple.h"
+#include "AliFlowVector.h"
 
 ClassImp(AliFlowAnalysisWithLeeYangZeros)
 
@@ -118,12 +116,14 @@ void AliFlowAnalysisWithLeeYangZeros::WriteHistograms(TString* outputFileName)
     //output->WriteObject(fHistList, "cobjLYZ1","SingleKey");
     if (fUseSum) { fHistList->SetName("cobjLYZ1SUM");}
     else {fHistList->SetName("cobjLYZ1PROD");}
+    fHistList->SetOwner(kTRUE);
     fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
   }
   else {
     //output->WriteObject(fHistList, "cobjLYZ2","SingleKey");
     if (fUseSum) { fHistList->SetName("cobjLYZ2SUM"); }
     else { fHistList->SetName("cobjLYZ2PROD"); }
+    fHistList->SetOwner(kTRUE);
     fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
   }
   delete output;
@@ -140,12 +140,14 @@ void AliFlowAnalysisWithLeeYangZeros::WriteHistograms(TString outputFileName)
     //output->WriteObject(fHistList, "cobjLYZ1","SingleKey");
     if (fUseSum) { fHistList->SetName("cobjLYZ1SUM");}
     else {fHistList->SetName("cobjLYZ1PROD");}
+    fHistList->SetOwner(kTRUE);
     fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
   }
   else {
     //output->WriteObject(fHistList, "cobjLYZ2","SingleKey");
     if (fUseSum) { fHistList->SetName("cobjLYZ2SUM"); }
     else { fHistList->SetName("cobjLYZ2PROD"); }
+    fHistList->SetOwner(kTRUE);
     fHistList->Write(fHistList->GetName(), TObject::kSingleKey);
   }
   delete output;
@@ -153,20 +155,47 @@ void AliFlowAnalysisWithLeeYangZeros::WriteHistograms(TString outputFileName)
 
 //-----------------------------------------------------------------------
 
+void AliFlowAnalysisWithLeeYangZeros::WriteHistograms(TDirectoryFile *outputFileName)
+{
+ //store the final results in output .root file
+ if (GetFirstRun()) {
+   if (fUseSum) { fHistList->SetName("cobjLYZ1SUM");}
+   else {fHistList->SetName("cobjLYZ1PROD");}
+   fHistList->SetOwner(kTRUE);
+   outputFileName->Add(fHistList);
+   outputFileName->Write(outputFileName->GetName(), TObject::kSingleKey); 
+ }
+ else {
+  if (fUseSum) { fHistList->SetName("cobjLYZ2SUM"); }
+  else { fHistList->SetName("cobjLYZ2PROD"); }
+  fHistList->SetOwner(kTRUE);
+  outputFileName->Add(fHistList);
+  outputFileName->Write(outputFileName->GetName(), TObject::kSingleKey); 
+ }
+}
+
+//-----------------------------------------------------------------------
+
 Bool_t AliFlowAnalysisWithLeeYangZeros::Init() 
 {
   //init method 
   if (fDebug) cout<<"****AliFlowAnalysisWithLeeYangZeros::Init()****"<<endl;
 
+  //save old value and prevent histograms from being added to directory
+  //to avoid name clashes in case multiple analaysis objects are used
+  //in an analysis
+  Bool_t oldHistAddStatus = TH1::AddDirectoryStatus();
+  TH1::AddDirectory(kFALSE);
   // Book histograms
-  Int_t iNtheta = AliFlowLYZConstants::kTheta;
-  Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
-  Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
-
-  Double_t  dPtMin = AliFlowCommonConstants::GetPtMin();            
-  Double_t  dPtMax = AliFlowCommonConstants::GetPtMax();
-  Double_t  dEtaMin = AliFlowCommonConstants::GetEtaMin();          
-  Double_t  dEtaMax = AliFlowCommonConstants::GetEtaMax();
+  Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
+  Int_t iNbinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
+  Int_t iNbinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
+
+  Double_t  dPtMin = AliFlowCommonConstants::GetMaster()->GetPtMin();       
+  Double_t  dPtMax = AliFlowCommonConstants::GetMaster()->GetPtMax();
+  Double_t  dEtaMin = AliFlowCommonConstants::GetMaster()->GetEtaMin();             
+  Double_t  dEtaMax = AliFlowCommonConstants::GetMaster()->GetEtaMax();
     
   //for control histograms
   if (fFirstRun){ 
@@ -311,7 +340,10 @@ Bool_t AliFlowAnalysisWithLeeYangZeros::Init()
   if (fDebug) cout<<"****Histograms initialised****"<<endl;
     
   fEventNumber = 0; //set event counter to zero
+  
+  //resore old stuff
+  TH1::AddDirectory(oldHistAddStatus);
+
   return kTRUE; 
 }
  
@@ -349,7 +381,7 @@ Bool_t AliFlowAnalysisWithLeeYangZeros::Make(AliFlowEventSimple* anEvent)
 void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHistos) {
  // get the pointers to all output histograms before calling Finish()
  
-  const Int_t iNtheta = AliFlowLYZConstants::kTheta;
+  const Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
   
   if (outputListHistos) {
 
@@ -365,9 +397,18 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
     TProfile* pHistProVPtRP   = NULL;
     TProfile* pHistProVPtPOI  = NULL;
     TH1F* pHistQsumforChi = NULL;
-    AliFlowLYZHist1 *pLYZHist1[iNtheta] = {NULL};      //array of pointers to AliFlowLYZHist1
-    AliFlowLYZHist2 *pLYZHist2RP[iNtheta] = {NULL};    //array of pointers to AliFlowLYZHist2
-    AliFlowLYZHist2 *pLYZHist2POI[iNtheta] = {NULL};   //array of pointers to AliFlowLYZHist2
+    //AliFlowLYZHist1 *pLYZHist1[iNtheta] = {NULL};      //array of pointers to AliFlowLYZHist1
+    //AliFlowLYZHist2 *pLYZHist2RP[iNtheta] = {NULL};    //array of pointers to AliFlowLYZHist2
+    //AliFlowLYZHist2 *pLYZHist2POI[iNtheta] = {NULL};   //array of pointers to AliFlowLYZHist2
+    AliFlowLYZHist1 **pLYZHist1 = new AliFlowLYZHist1*[iNtheta];      //array of pointers to AliFlowLYZHist1
+    AliFlowLYZHist2 **pLYZHist2RP = new AliFlowLYZHist2*[iNtheta];    //array of pointers to AliFlowLYZHist2
+    AliFlowLYZHist2 **pLYZHist2POI = new AliFlowLYZHist2*[iNtheta];   //array of pointers to AliFlowLYZHist2
+    for (Int_t i=0; i<iNtheta; i++)
+    {
+      pLYZHist1[i] = NULL;
+      pLYZHist2RP[i] = NULL;
+      pLYZHist2POI[i] = NULL;
+    }
 
     if (GetFirstRun()) { //first run
       //Get the common histograms from the output list
@@ -571,6 +612,9 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
       }
     } //secondrun
     //outputListHistos->Print(); 
+    delete [] pLYZHist1;
+    delete [] pLYZHist2RP;
+    delete [] pLYZHist2POI;
   } //listhistos
   else { 
     cout << "histogram list pointer is empty in method AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms() " << endl;}
@@ -584,7 +628,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
   
   //define variables for both runs
   Double_t  dJ01 = 2.405; 
-  Int_t iNtheta = AliFlowLYZConstants::kTheta;
+  Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
   //set the event number
   SetEventNumber((int)fCommonHists->GetHistMultOrig()->GetEntries());
   //cout<<"number of events processed is "<<fEventNumber<<endl; 
@@ -609,15 +653,17 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
        dR0 = fHist1[theta]->GetR0();
                                   
        //calculate integrated flow
-       if (dR0!=0.) { dVtheta = dJ01/dR0; }
+       if (!TMath::AreEqualAbs(dR0, 0., 1e-100)) { dVtheta = dJ01/dR0; }
        else { cout<<"r0 is not found! Leaving LYZ analysis."<<endl; return kFALSE; }
 
        //for estimating systematic error resulting from d0
-       Double_t dBinsize = (AliFlowLYZConstants::fgMax)/(AliFlowLYZConstants::kNbins);
+       Double_t dBinsize =0.;
+       if (fUseSum){ dBinsize = (AliFlowLYZConstants::GetMaster()->GetMaxSUM())/(AliFlowLYZConstants::GetMaster()->GetNbins());}
+       else { dBinsize = (AliFlowLYZConstants::GetMaster()->GetMaxPROD())/(AliFlowLYZConstants::GetMaster()->GetNbins());}
        Double_t dVplus = -1.;
        Double_t dVmin  = -1.;
-       if (dR0+dBinsize!=0.) {dVplus = dJ01/(dR0+dBinsize);}
-       if (dR0-dBinsize!=0.) {dVmin = dJ01/(dR0-dBinsize);}
+       if (!TMath::AreEqualAbs(dR0+dBinsize, 0., 1e-100)) {dVplus = dJ01/(dR0+dBinsize);}
+       if (!TMath::AreEqualAbs(dR0-dBinsize, 0., 1e-100)) {dVmin = dJ01/(dR0-dBinsize);}
        //convert V to v 
        Double_t dvplus = -1.;
        Double_t dvmin= -1.;
@@ -628,7 +674,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
          dvmin = dVmin; }
        else {
          //for PRODUCT: v=V/M
-         if (dMultRP!=0.){
+         if (!TMath::AreEqualAbs(dMultRP, 0., 1e-100)){
            dv = dVtheta/dMultRP;
            dvplus = dVplus/dMultRP;
            dvmin = dVmin/dMultRP; }}
@@ -645,14 +691,17 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
     //get average value of fVtheta = fV
     dV /=iNtheta;
     if (!fUseSum) { if (dMultRP!=0.){dV /=dMultRP;}} //scale with multiplicity for PRODUCT
-
+    
     //sigma2 and chi 
     Double_t  dSigma2 = 0;
     Double_t  dChi= 0;
     if (fEventNumber!=0) {
       *fQsum /= fEventNumber;
+      //cout<<"fQsum is "<<fQsum->X()<<" "<<fQsum->Y()<<endl; 
       fQ2sum /= fEventNumber;
+      //cout<<"fQ2sum is "<<fQ2sum<<endl; 
       dSigma2 = fQ2sum - TMath::Power(fQsum->X(),2.) - TMath::Power(fQsum->Y(),2.) - TMath::Power(dV,2.);  //BP eq. 62
+      //cout<<"dSigma2 is "<<dSigma2<<endl; 
       if (dSigma2>0) dChi = dV/TMath::Sqrt(dSigma2);
       else dChi = -1.;
       fCommonHistsRes->FillChiRP(dChi);
@@ -714,8 +763,8 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
     Int_t m = 1;
     TComplex i = TComplex::I();
     Double_t dBesselRatio[3] = {1., 1.202, 2.69};
-    Int_t iNbinsPt = AliFlowCommonConstants::GetNbinsPt();
-    Int_t iNbinsEta = AliFlowCommonConstants::GetNbinsEta();
+    Int_t iNbinsPt = AliFlowCommonConstants::GetMaster()->GetNbinsPt();
+    Int_t iNbinsEta = AliFlowCommonConstants::GetMaster()->GetNbinsEta();
 
     Double_t dEtaRP, dPtRP, dReRatioRP, dVetaRP, dVPtRP, dEtaPOI, dPtPOI, dReRatioPOI, dVetaPOI, dVPtPOI;
          
@@ -730,7 +779,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
       }        else {
        dR0 = fHistProR0theta->GetBinContent(theta+1); //histogram starts at bin 1
        if (fDebug) cerr<<"dR0 = "<<dR0<<endl;
-       if (dR0!=0) dVtheta = dJ01/dR0;
+       if (!TMath::AreEqualAbs(dR0, 0., 1e-100)) dVtheta = dJ01/dR0;
        dV += dVtheta; 
        // BP Eq. 9  -> Vn^theta = j01/r0^theta
        if (!fHistProReDenom && !fHistProImDenom) {
@@ -758,7 +807,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
              if (fDebug) cerr<<"WARNING: modulus of cNumerRP is zero in Finish()"<<endl;
              dReRatioRP = 0;
            }
-           else if (cDenom.Rho()==0) {
+           else if (TMath::AreEqualAbs(cDenom.Rho(), 0, 1e-100)) {
              if (fDebug) cerr<<"WARNING: modulus of cDenom is zero"<<endl;
              dReRatioRP = 0;
            }
@@ -838,7 +887,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
     //calculate the average of fVtheta = fV
     dV /= iNtheta;
     if (!fUseSum) { if (dMultRP!=0.) { dV /=dMultRP; }} //scale by the multiplicity for PRODUCT
-    if (dV==0.) { cout<<"dV = 0! Leaving LYZ analysis."<<endl; return kFALSE; }
+    if (TMath::AreEqualAbs(dV, 0., 1e-100)) { cout<<"dV = 0! Leaving LYZ analysis."<<endl; return kFALSE; }
 
     //sigma2 and chi (for statistical error calculations)
     Double_t  dSigma2 = 0;
@@ -898,7 +947,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
       Double_t dErrdifcomb = 0.;  //set error to zero
       Double_t dErr2difcomb = 0.; //set error to zero
       //calculate error
-      if (dNprime!=0.) { 
+      if (!TMath::AreEqualAbs(dNprime, 0., 1e-100)) { 
        for (Int_t theta=0;theta<iNtheta;theta++) {
          Double_t dTheta = ((double)theta/iNtheta)*TMath::Pi(); 
          Double_t dApluscomb = TMath::Exp((dJ01*dJ01)/(2*dChi*dChi)*
@@ -912,7 +961,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
        } //loop over theta
       } 
       
-      if (dErr2difcomb!=0.) {
+      if (!TMath::AreEqualAbs(dErr2difcomb, 0., 1e-100)) {
        dErr2difcomb /= iNtheta;
        dErrdifcomb = TMath::Sqrt(dErr2difcomb);
       }
@@ -997,7 +1046,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
       } else { cout<<"fHistPtRP is NULL"<<endl; }
     } //loop over bins b
 
-    if (dSum != 0.) {
+    if (!TMath::AreEqualAbs(dSum, 0., 1e-100)) {
       dVRP /= dSum; //the pt distribution should be normalised
       dErrV /= (dSum*dSum);
       dErrV = TMath::Sqrt(dErrV);
@@ -1085,8 +1134,8 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
    
   //define variables
   TComplex cExpo, cGtheta, cGthetaNew, cZ;
-  Int_t iNtheta = AliFlowLYZConstants::kTheta;
-  Int_t iNbins = AliFlowLYZConstants::kNbins;
+  Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
+  Int_t iNbins = AliFlowLYZConstants::GetMaster()->GetNbins();
    
 
   //calculate flow
@@ -1097,7 +1146,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
   //weight by the multiplicity
   Double_t dQX = 0;
   Double_t dQY = 0;
-  if (vQ.GetMult() != 0) {
+  if (!TMath::AreEqualAbs(vQ.GetMult(), 0., 1e-100)) {
     dQX = vQ.X()/vQ.GetMult(); 
     dQY = vQ.Y()/vQ.GetMult();
   }
@@ -1166,7 +1215,7 @@ void AliFlowAnalysisWithLeeYangZeros::GetOutputHistograms(TList *outputListHisto
   Double_t dCosTermPOI = 0.;
   Double_t m = 1.;
   Double_t dOrder = 2.;
-  Int_t iNtheta = AliFlowLYZConstants::kTheta;
+  Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
   
    
   //get the Q vector 
@@ -1328,7 +1377,7 @@ TComplex AliFlowAnalysisWithLeeYangZeros::GetDiffFlow(AliFlowEventSimple* const
 
   Int_t iNumberOfTracks = anEvent->NumberOfTracks();
   
-  Int_t iNtheta = AliFlowLYZConstants::kTheta;
+  Int_t iNtheta = AliFlowLYZConstants::GetMaster()->GetNtheta();
   Double_t dTheta = ((double)theta/iNtheta)*TMath::Pi()/dOrder;
   
   //for the denominator (use all RP selected particles)