adding phi, etaphi plots

adding function to determine acceptance map
adding process type statistics

diff --git a/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx b/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx
index 1dcc6dbac433f9ad761f9e63e61e49b7dddb6d4a..9966f1482570eacce1cbfcd55de85ff8c3682e03 100644 (file)
--- a/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx
+++ b/PWG0/dNdEta/AlidNdEtaCorrectionTask.cxx
@@ -205,10 +205,12 @@ void AlidNdEtaCorrectionTask::CreateOutputObjects()
     fOutput->Add(fDeltaPhi[i]);
   }
 
-  fEventStats = new TH2F("fEventStats", "fEventStats;event type;status;count", 2, -0.5, 1.5, 4, -0.5, 3.5);
+  fEventStats = new TH2F("fEventStats", "fEventStats;event type;status;count", 104, -3.5, 100.5, 4, -0.5, 3.5);
   fOutput->Add(fEventStats);
-  fEventStats->GetXaxis()->SetBinLabel(1, "INEL");
-  fEventStats->GetXaxis()->SetBinLabel(2, "NSD");
+  fEventStats->GetXaxis()->SetBinLabel(1, "INEL"); // x = -3
+  fEventStats->GetXaxis()->SetBinLabel(2, "NSD");  // x = -2
+  for (Int_t i=0; i<100; i++)
+    fEventStats->GetXaxis()->SetBinLabel(4+i, Form("%d", i)); 
 
   fEventStats->GetYaxis()->SetBinLabel(1, "nothing");
   fEventStats->GetYaxis()->SetBinLabel(2, "trg");
@@ -271,7 +273,7 @@ void AlidNdEtaCorrectionTask::Exec(Option_t*)
   AliDebug(AliLog::kDebug+1, Form("Found process type %d", processType));
 
   if (processType == AliPWG0Helper::kInvalidProcess)
-    AliDebug(AliLog::kError, Form("Unknown process type %d.", processType));
+    AliDebug(AliLog::kError, "Unknown process type.");
 
   // get the MC vertex
   AliGenEventHeader* genHeader = header->GenEventHeader();
@@ -298,10 +300,10 @@ void AlidNdEtaCorrectionTask::Exec(Option_t*)
   // fill process type
   Int_t biny = (Int_t) eventTriggered + 2 * (Int_t) eventVertex;
   // INEL
-  fEventStats->Fill(0.0, biny);
+  fEventStats->Fill(-3, biny);
   // NSD
-  if (processType != AliPWG0Helper::kSD )
-    fEventStats->Fill(1.0, biny);
+  if (processType != AliPWG0Helper::kSD)
+    fEventStats->Fill(-2, biny);
   
   // create list of (label, eta, pt) tuples
   Int_t inputCount = 0;
@@ -446,10 +448,13 @@ void AlidNdEtaCorrectionTask::Exec(Option_t*)
         fdNdEtaAnalysisMC->FillTrack(vtxMC[2], eta, pt);
 
     // TODO this value might be needed lower for the SPD study (only used for control histograms anyway)
-    if (TMath::Abs(eta) < 1 && pt > 0.2)
+    if (TMath::Abs(eta) < 1.5) // && pt > 0.2)
       nAccepted++;
   }
 
+  if (nAccepted == 0)
+    fEventStats->Fill(AliPWG0Helper::GetLastProcessType(), biny);
+
   fMultAll->Fill(nAccepted);
   if (eventTriggered) {
     fMultTr->Fill(nAccepted);

diff --git a/PWG0/dNdEta/AlidNdEtaTask.cxx b/PWG0/dNdEta/AlidNdEtaTask.cxx
index 73204dd21845d8d57558e18aeb9fe8c8139f1864..46e1a02a309ad5732b59b44a35740c51d486c70b 100644 (file)
--- a/PWG0/dNdEta/AlidNdEtaTask.cxx
+++ b/PWG0/dNdEta/AlidNdEtaTask.cxx
@@ -60,6 +60,8 @@ AlidNdEtaTask::AlidNdEtaTask(const char* opt) :
   fdNdEtaAnalysisTracks(0),
   fVertex(0),
   fPartPt(0),
+  fPhi(0),
+  fEtaPhi(0),
   fDeltaPhi(0)
 {
   //
@@ -168,6 +170,15 @@ void AlidNdEtaTask::CreateOutputObjects()
   fVertexResolution = new TH1F("dndeta_vertex_resolution_z", "dndeta_vertex_resolution_z", 1000, 0, 10);
   fOutput->Add(fVertexResolution);
 
+  fPhi = new TH1F("fPhi", "fPhi;#phi in rad.;count", 720, 0, 2 * TMath::Pi());
+  fOutput->Add(fPhi);
+
+  fEtaPhi = new TH2F("fEtaPhi", "fEtaPhi;#eta;#phi in rad.;count", 80, -4, 4, 80, 0, 2 * TMath::Pi());
+  fOutput->Add(fEtaPhi);
+
+  if (fAnalysisMode == AliPWG0Helper::kSPD)
+    fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 1000, -3.14, 3.14);
+
   if (fReadMC)
   {
     fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta", fAnalysisMode);
@@ -192,9 +203,6 @@ void AlidNdEtaTask::CreateOutputObjects()
     fPartPt->Sumw2();
     fOutput->Add(fPartPt);
   }
-
-  if (fAnalysisMode == AliPWG0Helper::kSPD)
-    fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 1000, -3.14, 3.14);
 }
 
 void AlidNdEtaTask::Exec(Option_t*)
@@ -319,6 +327,12 @@ void AlidNdEtaTask::Exec(Option_t*)
       if (TMath::Abs(deltaPhi) > 1)
         printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), deltaPhi);
 
+      Float_t phi = mult->GetPhi(i);
+      if (phi < 0)
+        phi += TMath::Pi() * 2;
+      fPhi->Fill(phi);
+      fEtaPhi->Fill(mult->GetEta(i), phi);
+
       fDeltaPhi->Fill(deltaPhi);
 
       etaArr[inputCount] = mult->GetEta(i);
@@ -355,6 +369,12 @@ void AlidNdEtaTask::Exec(Option_t*)
         continue;
       }
 
+      Float_t phi = esdTrack->Phi();
+      if (phi < 0)
+        phi += TMath::Pi() * 2;
+      fPhi->Fill(phi);
+      fEtaPhi->Fill(esdTrack->Eta(), phi);
+
       etaArr[inputCount] = esdTrack->Eta();
       labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel());
       ptArr[inputCount] = esdTrack->Pt();
@@ -473,7 +493,7 @@ void AlidNdEtaTask::Exec(Option_t*)
       Float_t pt = particle->Pt();
 
        // make a rough eta cut (so that nAcceptedParticles is not too far off the true measured value (NB: this histograms are only gathered for comparison))
-      if (TMath::Abs(eta) < 1.5 && pt > 0.3)
+      if (TMath::Abs(eta) < 1.5) // && pt > 0.3)
         nAcceptedParticles++;
 
       fdNdEtaAnalysis->FillTrack(vtxMC[2], eta, pt);
@@ -494,7 +514,7 @@ void AlidNdEtaTask::Exec(Option_t*)
     }
 
     fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles);
-    if (processType != AliPWG0Helper::kSD )
+    if (processType != AliPWG0Helper::kSD)
       fdNdEtaAnalysisNSD->FillEvent(vtxMC[2], nAcceptedParticles);
 
     if (eventTriggered)
@@ -552,10 +572,14 @@ void AlidNdEtaTask::Terminate(Option_t *)
   fdNdEtaAnalysisESD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("fdNdEtaAnalysisESD"));
   fMult = dynamic_cast<TH1F*> (fOutput->FindObject("fMult"));
   fMultVtx = dynamic_cast<TH1F*> (fOutput->FindObject("fMultVtx"));
-  fEvents = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_vertex"));
-  fVertexResolution = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_vertex_resolution_z"));
   for (Int_t i=0; i<3; ++i)
     fPartEta[i] = dynamic_cast<TH1F*> (fOutput->FindObject(Form("dndeta_check_%d", i)));
+  fEvents = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_vertex"));
+  fVertexResolution = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_vertex_resolution_z"));
+
+  fPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fPhi"));
+  fEtaPhi = dynamic_cast<TH2F*> (fOutput->FindObject("fEtaPhi"));
+  fDeltaPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fDeltaPhi"));
 
   if (!fdNdEtaAnalysisESD)
   {
@@ -629,6 +653,12 @@ void AlidNdEtaTask::Terminate(Option_t *)
   if (fDeltaPhi)
     fDeltaPhi->Write();
 
+  if (fPhi)
+    fPhi->Write();
+
+  if (fEtaPhi)
+    fEtaPhi->Write();
+
   fout->Write();
   fout->Close();
 
@@ -642,6 +672,7 @@ void AlidNdEtaTask::Terminate(Option_t *)
     fdNdEtaAnalysisTrVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTrVtx"));
     fdNdEtaAnalysisTracks = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTracks"));
     fPartPt = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_check_pt"));
+    fVertex = dynamic_cast<TH3F*> (fOutput->FindObject("vertex_check"));
 
     if (!fdNdEtaAnalysis || !fdNdEtaAnalysisTr || !fdNdEtaAnalysisTrVtx || !fPartPt)
     {
@@ -666,7 +697,12 @@ void AlidNdEtaTask::Terminate(Option_t *)
     fdNdEtaAnalysisTr->SaveHistograms();
     fdNdEtaAnalysisTrVtx->SaveHistograms();
     fdNdEtaAnalysisTracks->SaveHistograms();
-    fPartPt->Write();
+
+    if (fPartPt)
+      fPartPt->Write();
+
+    if (fVertex)
+      fVertex->Write();
 
     fout->Write();
     fout->Close();

diff --git a/PWG0/dNdEta/AlidNdEtaTask.h b/PWG0/dNdEta/AlidNdEtaTask.h
index 8c40dd006a843bb6e4fd65b9494fd546f1c08653..96fd24af684bf25f14eb6851929a5a66196f04c0 100644 (file)
--- a/PWG0/dNdEta/AlidNdEtaTask.h
+++ b/PWG0/dNdEta/AlidNdEtaTask.h
@@ -10,6 +10,7 @@
 class AliESDtrackCuts;
 class dNdEtaAnalysis;
 class TH1F;
+class TH2F;
 class TH3F;
 class AliESDEvent;
 
@@ -62,9 +63,14 @@ class AlidNdEtaTask : public AliAnalysisTask {
     dNdEtaAnalysis* fdNdEtaAnalysisTr;      //! contains the dndeta from the triggered events
     dNdEtaAnalysis* fdNdEtaAnalysisTrVtx;   //! contains the dndeta from the triggered events with vertex
     dNdEtaAnalysis* fdNdEtaAnalysisTracks;  //! contains the dndeta from the triggered events with vertex counted from the mc particles associated to the tracks (comparing this to the raw values from the esd shows the effect of the detector resolution)
-    // the following are control histograms to check the dNdEtaAnalysis class
+
+    // control histograms (MC)
     TH3F* fVertex;                //! vertex of counted particles
     TH1F* fPartPt;                //! counted particles as function of pt
+
+    // control histograms (ESD)
+    TH1F* fPhi;                   //! raw phi distribution
+    TH2F* fEtaPhi;                //! raw eta - phi distribution
     TH1F* fDeltaPhi;              //! histogram of delta_phi values for tracklets (only for SPD analysis)
 
  private:

diff --git a/PWG0/dNdEta/dNdEtaAnalysis.cxx b/PWG0/dNdEta/dNdEtaAnalysis.cxx
index 6614bc2cbf4f904ebc06e86d720391ec8edc1598..e0956e3e3c2a7e129410409f85341a447ceb9477 100644 (file)
--- a/PWG0/dNdEta/dNdEtaAnalysis.cxx
+++ b/PWG0/dNdEta/dNdEtaAnalysis.cxx
@@ -390,29 +390,42 @@ void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, Alid
       Int_t vertexBinBegin = vertexHist->GetXaxis()->FindBin(vertexRangeBegin[vertexPos]);
       Int_t vertexBinEnd   = vertexHist->GetXaxis()->FindBin(vertexRangeEnd[vertexPos]);
 
+      const Int_t *binBegin = 0;
       // adjust acceptance range for SPD
       if (fAnalysisMode == AliPWG0Helper::kSPD)
       {
         //const Int_t binBegin[30] = { 18, 16, 15, 14, 13, 13, 12, 11, 9, 7, 6, 5, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 1 };
-        const Int_t binBegin[30] = { -1, -1, -1, -1, 16, 14, 12, 10, 9, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 2, 2, 2, 2, 2, -1, -1, -1, -1}; // limit in correction map is 5
+        const Int_t binBeginSPD[30] = { -1, -1, -1, -1, 16, 14, 12, 10, 9, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 2, 2, 2, 2, 2, -1, -1, -1, -1}; // limit in correction map is 5
         //const Int_t binBegin[30] = { -1, -1, -1, 17, 15, 14, 12, 10, 8, 7, 6, 5, 4, 4, 4, 4, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, -1, -1, -1};  // limit in correction map is 10
+        binBegin = binBeginSPD;
 
-        Int_t vtxBegin = binBegin[iEta - 1];
-        Int_t vtxEnd   = 18 + 1 - binBegin[30 - iEta];
-
-        // eta range not accessible
-        if (vtxBegin == -1)
-          continue;
-
-        //Float_t eta = vtxVsEta->GetYaxis()->GetBinCenter(iEta);
-        //Printf("%d %d | %d %d", vtxBegin, vertexHist->GetXaxis()->FindBin(GetVtxMin(eta)), vtxEnd, vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta)));
-        //vtxBegin = vertexHist->GetXaxis()->FindBin(GetVtxMin(eta));
-        //vtxEnd = vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta));
-        //printf("Eta bin: %d (%f) Vertex range: %d; before: %d %d ", iEta, eta, vertexPos, vertexBinBegin, vertexBinEnd);
-        vertexBinBegin = TMath::Max(vertexBinBegin, vtxBegin);
-        vertexBinEnd =   TMath::Min(vertexBinEnd, vtxEnd);
-        //Printf("after:  %d %d", vertexBinBegin, vertexBinEnd);
       }
+      else if (fAnalysisMode == AliPWG0Helper::kTPC)
+      {
+//        const Int_t binBegin[30] = { -1, -1, -1, -1, 16, 14, 12, 10, 9, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 2, 2, 2, 2, 2, -1, -1, -1, -1}; // limit in correction map is 5
+        const Int_t binBeginTPC[30] = {-1, -1, -1, -1, -1, -1, -1, -1, 9, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, -1, -1, -1, -1, -1, -1, -1, -1};
+        binBegin = binBeginTPC;
+      }
+      else if (fAnalysisMode == AliPWG0Helper::kTPCITS)
+      {
+        // TODO create map
+      }
+
+      Int_t vtxBegin = binBegin[iEta - 1];
+      Int_t vtxEnd   = 18 + 1 - binBegin[30 - iEta];
+      
+      // eta range not accessible
+      if (vtxBegin == -1)
+        continue;
+      
+      //Float_t eta = vtxVsEta->GetYaxis()->GetBinCenter(iEta);
+      //Printf("%d %d | %d %d", vtxBegin, vertexHist->GetXaxis()->FindBin(GetVtxMin(eta)), vtxEnd, vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta)));
+      //vtxBegin = vertexHist->GetXaxis()->FindBin(GetVtxMin(eta));
+      //vtxEnd = vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta));
+      //printf("Eta bin: %d (%f) Vertex range: %d; before: %d %d ", iEta, eta, vertexPos, vertexBinBegin, vertexBinEnd);
+      vertexBinBegin = TMath::Max(vertexBinBegin, vtxBegin);
+      vertexBinEnd =   TMath::Min(vertexBinEnd, vtxEnd);
+      //Printf("after:  %d %d", vertexBinBegin, vertexBinEnd);
 
       // no data for this bin
       if (vertexBinBegin > vertexBinEnd)

diff --git a/PWG0/dNdEta/drawPlots.C b/PWG0/dNdEta/drawPlots.C
index 2b944a279075ca2f2489a38af8ea803050265c62..73a28c68199e3265503dcf5d4e221b6456f1246b 100644 (file)
--- a/PWG0/dNdEta/drawPlots.C
+++ b/PWG0/dNdEta/drawPlots.C
@@ -1729,7 +1729,7 @@ void DrawTrackletOrigin()
   }
 }
 
-void DetermineAcceptance(const char* fileName = "correction_map.root", const char* dirName = "dndeta_correction")
+void DetermineAcceptance(const char* fileName = "correction_map.root", const char* dirName = "dndeta_correction", Double_t ptmin=0.2)
 {
   loadlibs();
 
@@ -1739,12 +1739,36 @@ void DetermineAcceptance(const char* fileName = "correction_map.root", const cha
   if (!dNdEtaCorrection->LoadHistograms())
     return;
 
-  hist = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetCorrectionHistogram();
+  //  hist = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetCorrectionHistogram();
 
-  proj = hist->Project3D("yx");
+  gener = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetGeneratedHistogram();
+  measu = dNdEtaCorrection->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram();
+
+  gener->GetZaxis()->SetRange(gener->GetZaxis()->FindBin(ptmin), gener->GetNbinsZ()+1);
+  TH2D *gener_xy = gener->Project3D("yx");
+
+  measu->GetZaxis()->SetRange(measu->GetZaxis()->FindBin(ptmin), measu->GetNbinsZ()+1);
+  TH2D *measu_xy = measu->Project3D("yx");
+
+  cout << measu->GetZaxis()->FindBin(ptmin) << " " << measu->GetNbinsZ()+1 << endl;
+
+  TCanvas *canp = new TCanvas("canp","canp",600,1000);
+  canp->Divide(1,2,0.0001,0.0001);
+  canp->cd(1);
+  gener_xy->Draw("COLZ");
+  canp->cd(2);
+  measu_xy->Draw("COLZ");
+
+
+  TCanvas *canpr = new TCanvas("canpr","canpr",700,500);
+  canpr->cd();
+  TH2D *proj = new TH2D(*gener_xy);
+  proj->Divide(measu_xy);
+
+//   proj = hist->Project3D("yx");
   proj->Draw("COLZ");
 
-  const Float_t limit = 10;
+  const Float_t limit = 5;
 
   TString array = "{";
   TString arrayEnd = "}";