#include <TNtuple.h>
#include <TObjString.h>
#include <TF1.h>
+#include <TGraph.h>
#include <AliLog.h>
#include <AliESDVertex.h>
#include <AliMCEventHandler.h>
#include <AliMCEvent.h>
#include <AliESDInputHandler.h>
+#include <AliESDHeader.h>
#include "AliESDtrackCuts.h"
#include "AliPWG0Helper.h"
fOption(opt),
fAnalysisMode(AliPWG0Helper::kTPC),
fTrigger(AliPWG0Helper::kMB1),
+ fFillPhi(kFALSE),
+ fDeltaPhiCut(-1),
fReadMC(kFALSE),
fUseMCVertex(kFALSE),
fOnlyPrimaries(kFALSE),
fEvents(0),
fVertexResolution(0),
fdNdEtaAnalysis(0),
+ fdNdEtaAnalysisND(0),
fdNdEtaAnalysisNSD(0),
fdNdEtaAnalysisTr(0),
fdNdEtaAnalysisTrVtx(0),
fPartPt(0),
fVertex(0),
fPhi(0),
+ fRawPt(0),
fEtaPhi(0),
- fDeltaPhi(0)
+ fDeltaPhi(0),
+ fFiredChips(0),
+ fTriggerVsTime(0),
+ fStats(0)
{
//
// Constructor. Initialization of pointers
// Define input and output slots here
DefineInput(0, TChain::Class());
DefineOutput(0, TList::Class());
+
+ fZPhi[0] = 0;
+ fZPhi[1] = 0;
AliLog::SetClassDebugLevel("AlidNdEtaTask", AliLog::kWarning);
}
{
static Int_t count = 0;
count++;
- Printf("Processing %d. file", count);
+ Printf("Processing %d. file: %s", count, ((TTree*) GetInputData(0))->GetCurrentFile()->GetName());
return kTRUE;
}
fEvents = new TH1F("dndeta_check_vertex", "dndeta_check_vertex", 800, -40, 40);
fOutput->Add(fEvents);
- fVertexResolution = new TH1F("dndeta_vertex_resolution_z", "dndeta_vertex_resolution_z", 1000, 0, 10);
+ Float_t resMax = (fAnalysisMode == AliPWG0Helper::kSPD) ? 0.2 : 2;
+ fVertexResolution = new TH1F("dndeta_vertex_resolution_z", "dndeta_vertex_resolution_z", 1000, 0, resMax);
fOutput->Add(fVertexResolution);
fPhi = new TH1F("fPhi", "fPhi;#phi in rad.;count", 720, 0, 2 * TMath::Pi());
fEtaPhi = new TH2F("fEtaPhi", "fEtaPhi;#eta;#phi in rad.;count", 80, -4, 4, 18*5, 0, 2 * TMath::Pi());
fOutput->Add(fEtaPhi);
+ fTriggerVsTime = new TGraph; //TH1F("fTriggerVsTime", "fTriggerVsTime;trigger time;count", 100, 0, 100);
+ fTriggerVsTime->SetName("fTriggerVsTime");
+ fTriggerVsTime->GetXaxis()->SetTitle("trigger time");
+ fTriggerVsTime->GetYaxis()->SetTitle("count");
+ fOutput->Add(fTriggerVsTime);
+
+ fStats = new TH1F("fStats", "fStats", 2, 0.5, 2.5);
+ fStats->GetXaxis()->SetBinLabel(1, "vertexer 3d");
+ fStats->GetXaxis()->SetBinLabel(2, "vertexer z");
+ fOutput->Add(fStats);
+
if (fAnalysisMode == AliPWG0Helper::kSPD)
- fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 18*50, -3.14, 3.14);
+ {
+ fDeltaPhi = new TH1F("fDeltaPhi", "fDeltaPhi;#Delta #phi;Entries", 500, -0.2, 0.2);
+ fOutput->Add(fDeltaPhi);
+ fFiredChips = new TH2F("fFiredChips", "fFiredChips;Chips L1 + L2;tracklets", 1201, -0.5, 1201, 50, -0.5, 49.5);
+ fOutput->Add(fFiredChips);
+ for (Int_t i=0; i<2; i++)
+ {
+ fZPhi[i] = new TH2F(Form("fZPhi_%d", i), Form("fZPhi Layer %d;z (cm);#phi (rad.)", i), 200, -20, 20, 180, 0, TMath::Pi() * 2);
+ fOutput->Add(fZPhi[i]);
+ }
+ }
+
+ if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS)
+ {
+ fRawPt = new TH1F("fRawPt", "raw pt;p_{T};Count", 2000, 0, 100);
+ fOutput->Add(fRawPt);
+ }
fVertex = new TH3F("vertex_check", "vertex_check", 100, -1, 1, 100, -1, 1, 100, -30, 30);
fOutput->Add(fVertex);
fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta", fAnalysisMode);
fOutput->Add(fdNdEtaAnalysis);
+ fdNdEtaAnalysisND = new dNdEtaAnalysis("dndetaND", "dndetaND", fAnalysisMode);
+ fOutput->Add(fdNdEtaAnalysisND);
+
fdNdEtaAnalysisNSD = new dNdEtaAnalysis("dndetaNSD", "dndetaNSD", fAnalysisMode);
fOutput->Add(fdNdEtaAnalysisNSD);
// ESD analysis
if (fESD)
{
+ // check event type (should be PHYSICS = 7)
+ AliESDHeader* esdHeader = fESD->GetHeader();
+ if (!esdHeader)
+ {
+ Printf("ERROR: esdHeader could not be retrieved");
+ return;
+ }
+
+ /*
+ UInt_t eventType = esdHeader->GetEventType();
+ if (eventType != 7)
+ {
+ Printf("Skipping event because it is of type %d", eventType);
+ return;
+ }
+ */
+
// trigger definition
- eventTriggered = AliPWG0Helper::IsEventTriggered(fESD->GetTriggerMask(), fTrigger);
+ eventTriggered = AliPWG0Helper::IsEventTriggered(fESD, fTrigger);
// get the ESD vertex
vtxESD = AliPWG0Helper::GetVertex(fESD, fAnalysisMode);
if (AliPWG0Helper::TestVertex(vtxESD, fAnalysisMode))
{
vtxESD->GetXYZ(vtx);
+
+ // vertex stats
+ if (strcmp(vtxESD->GetTitle(), "vertexer: 3D") == 0)
+ {
+ fStats->Fill(1);
+ }
+ else if (strcmp(vtxESD->GetTitle(), "vertexer: Z") == 0)
+ fStats->Fill(2);
}
else
vtxESD = 0;
Int_t inputCount = 0;
Int_t* labelArr = 0;
Float_t* etaArr = 0;
- Float_t* ptArr = 0;
+ Float_t* thirdDimArr = 0;
if (fAnalysisMode == AliPWG0Helper::kSPD)
{
// get tracklets
labelArr = new Int_t[mult->GetNumberOfTracklets()];
etaArr = new Float_t[mult->GetNumberOfTracklets()];
- ptArr = new Float_t[mult->GetNumberOfTracklets()];
+ thirdDimArr = new Float_t[mult->GetNumberOfTracklets()];
- // get multiplicity from ITS tracklets
+ // get multiplicity from SPD tracklets
for (Int_t i=0; i<mult->GetNumberOfTracklets(); ++i)
{
//printf("%d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), mult->GetDeltaPhi(i));
if (TMath::Abs(deltaPhi) > 1)
printf("WARNING: Very high Delta Phi: %d %f %f %f\n", i, mult->GetTheta(i), mult->GetPhi(i), deltaPhi);
+ Int_t label = mult->GetLabel(i, 0);
+ Float_t eta = mult->GetEta(i);
+
+ // control histograms
Float_t phi = mult->GetPhi(i);
if (phi < 0)
phi += TMath::Pi() * 2;
fPhi->Fill(phi);
- fEtaPhi->Fill(mult->GetEta(i), phi);
+ fEtaPhi->Fill(eta, phi);
+
+ if (deltaPhi < 0.01)
+ {
+ // layer 0
+ Float_t z = vtx[2] + 3.9 / TMath::Tan(2 * TMath::ATan(TMath::Exp(- eta)));
+ fZPhi[0]->Fill(z, phi);
+ // layer 1
+ z = vtx[2] + 7.6 / TMath::Tan(2 * TMath::ATan(TMath::Exp(- eta)));
+ fZPhi[1]->Fill(z, phi);
+ }
fDeltaPhi->Fill(deltaPhi);
- Int_t label = mult->GetLabel(i, 0);
- Float_t eta = mult->GetEta(i);
+ if (fDeltaPhiCut > 0 && TMath::Abs(deltaPhi) > fDeltaPhiCut)
+ continue;
if (fUseMCKine)
{
{
TParticle* particle = stack->Particle(label);
eta = particle->Eta();
+ phi = particle->Phi();
}
else
Printf("WARNING: fUseMCKine set without fOnlyPrimaries and no label found");
}
-
+
etaArr[inputCount] = eta;
labelArr[inputCount] = label;
- ptArr[inputCount] = 0; // no pt for tracklets
+ thirdDimArr[inputCount] = phi;
++inputCount;
}
- // fill multiplicity in pt bin
- for (Int_t i=0; i<inputCount; ++i)
- ptArr[i] = inputCount;
+ if (!fFillPhi)
+ {
+ // fill multiplicity in third axis
+ for (Int_t i=0; i<inputCount; ++i)
+ thirdDimArr[i] = inputCount;
+ }
- Printf("Accepted %d tracklets", inputCount);
+ Int_t firedChips = mult->GetNumberOfFiredChips(0) + mult->GetNumberOfFiredChips(1);
+ fFiredChips->Fill(firedChips, inputCount);
+ Printf("Accepted %d tracklets (%d fired chips)", inputCount, firedChips);
}
else if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS)
{
return;
}
- // get multiplicity from ESD tracks
- TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD, (fAnalysisMode == AliPWG0Helper::kTPC));
- Int_t nGoodTracks = list->GetEntries();
-
- labelArr = new Int_t[nGoodTracks];
- etaArr = new Float_t[nGoodTracks];
- ptArr = new Float_t[nGoodTracks];
-
- // loop over esd tracks
- for (Int_t i=0; i<nGoodTracks; i++)
+ if (vtxESD)
{
- AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i));
- if (!esdTrack)
+ // get multiplicity from ESD tracks
+ TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD, (fAnalysisMode == AliPWG0Helper::kTPC));
+ Int_t nGoodTracks = list->GetEntries();
+ Printf("Accepted %d tracks", nGoodTracks);
+
+ labelArr = new Int_t[nGoodTracks];
+ etaArr = new Float_t[nGoodTracks];
+ thirdDimArr = new Float_t[nGoodTracks];
+
+ // loop over esd tracks
+ for (Int_t i=0; i<nGoodTracks; i++)
{
- AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i));
- continue;
- }
-
- Float_t phi = esdTrack->Phi();
- if (phi < 0)
- phi += TMath::Pi() * 2;
- fPhi->Fill(phi);
- fEtaPhi->Fill(esdTrack->Eta(), phi);
-
- Float_t eta = esdTrack->Eta();
- Int_t label = TMath::Abs(esdTrack->GetLabel());
- Float_t pT = esdTrack->Pt();
-
- if (fOnlyPrimaries && label == 0)
- continue;
-
- if (fUseMCKine)
- {
- if (label > 0)
+ AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(i));
+ if (!esdTrack)
{
- TParticle* particle = stack->Particle(label);
- eta = particle->Eta();
- pT = particle->Pt();
+ AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", i));
+ continue;
}
- else
- Printf("WARNING: fUseMCKine set without fOnlyPrimaries and no label found");
+
+ Float_t phi = esdTrack->Phi();
+ if (phi < 0)
+ phi += TMath::Pi() * 2;
+
+ Float_t eta = esdTrack->Eta();
+ Int_t label = TMath::Abs(esdTrack->GetLabel());
+ Float_t pT = esdTrack->Pt();
+
+ fPhi->Fill(phi);
+ fEtaPhi->Fill(eta, phi);
+ if (eventTriggered && vtxESD)
+ fRawPt->Fill(pT);
+
+ if (fOnlyPrimaries && label == 0)
+ continue;
+
+ if (fUseMCKine)
+ {
+ if (label > 0)
+ {
+ TParticle* particle = stack->Particle(label);
+ eta = particle->Eta();
+ pT = particle->Pt();
+ }
+ else
+ Printf("WARNING: fUseMCKine set without fOnlyPrimaries and no label found");
+ }
+
+ etaArr[inputCount] = eta;
+ labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel());
+ thirdDimArr[inputCount] = pT;
+ ++inputCount;
}
-
- etaArr[inputCount] = eta;
- labelArr[inputCount] = TMath::Abs(esdTrack->GetLabel());
- ptArr[inputCount] = pT;
- ++inputCount;
+
+ // TODO restrict inputCount used as measure for the multiplicity to |eta| < 1
+
+ delete list;
}
-
- Printf("Accepted %d tracks", nGoodTracks);
-
- delete list;
}
else
return;
fMult->Fill(inputCount);
fdNdEtaAnalysisESD->FillTriggeredEvent(inputCount);
+ fTriggerVsTime->SetPoint(fTriggerVsTime->GetN(), fESD->GetTimeStamp(), 1);
+
if (vtxESD)
{
// control hist
for (Int_t i=0; i<inputCount; ++i)
{
Float_t eta = etaArr[i];
- Float_t pt = ptArr[i];
+ Float_t thirdDim = thirdDimArr[i];
- fdNdEtaAnalysisESD->FillTrack(vtx[2], eta, pt);
+ fdNdEtaAnalysisESD->FillTrack(vtx[2], eta, thirdDim);
- if (TMath::Abs(vtx[2]) < 20)
+ if (TMath::Abs(vtx[2]) < 10)
{
fPartEta[0]->Fill(eta);
continue;
}
- fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], particle->Eta(), (fAnalysisMode == AliPWG0Helper::kSPD) ? inputCount : particle->Pt());
+ Float_t thirdDim = -1;
+ if (fAnalysisMode == AliPWG0Helper::kSPD)
+ {
+ if (fFillPhi)
+ {
+ thirdDim = particle->Phi();
+ }
+ else
+ thirdDim = inputCount;
+ }
+ else
+ thirdDim = particle->Pt();
+
+ fdNdEtaAnalysisTracks->FillTrack(vtxMC[2], particle->Eta(), thirdDim);
} // end of track loop
// for event count per vertex
}
}
- delete[] etaArr;
- delete[] labelArr;
- delete[] ptArr;
+ if (etaArr)
+ delete[] etaArr;
+ if (labelArr)
+ delete[] labelArr;
+ if (thirdDimArr)
+ delete[] thirdDimArr;
}
if (fReadMC) // Processing of MC information (optional)
continue;
Float_t eta = particle->Eta();
- Float_t thirdDim = (fAnalysisMode == AliPWG0Helper::kSPD) ? nAcceptedParticles : particle->Pt();
+ Float_t thirdDim = -1;
+
+ if (fAnalysisMode == AliPWG0Helper::kSPD)
+ {
+ if (fFillPhi)
+ {
+ thirdDim = particle->Phi();
+ }
+ else
+ thirdDim = nAcceptedParticles;
+ }
+ else
+ thirdDim = particle->Pt();
fdNdEtaAnalysis->FillTrack(vtxMC[2], eta, thirdDim);
- if (processType != AliPWG0Helper::kSD )
+ if (processType != AliPWG0Helper::kSD)
fdNdEtaAnalysisNSD->FillTrack(vtxMC[2], eta, thirdDim);
+ if (processType == AliPWG0Helper::kND)
+ fdNdEtaAnalysisND->FillTrack(vtxMC[2], eta, thirdDim);
+
if (eventTriggered)
{
fdNdEtaAnalysisTr->FillTrack(vtxMC[2], eta, thirdDim);
fdNdEtaAnalysisTrVtx->FillTrack(vtxMC[2], eta, thirdDim);
}
- if (TMath::Abs(eta) < 0.8)
+ if (TMath::Abs(eta) < 1.0)
fPartPt->Fill(particle->Pt());
}
fdNdEtaAnalysis->FillEvent(vtxMC[2], nAcceptedParticles);
if (processType != AliPWG0Helper::kSD)
fdNdEtaAnalysisNSD->FillEvent(vtxMC[2], nAcceptedParticles);
+ if (processType == AliPWG0Helper::kND)
+ fdNdEtaAnalysisND->FillEvent(vtxMC[2], nAcceptedParticles);
if (eventTriggered)
{
fVertexResolution = dynamic_cast<TH1F*> (fOutput->FindObject("dndeta_vertex_resolution_z"));
fPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fPhi"));
+ fRawPt = dynamic_cast<TH1F*> (fOutput->FindObject("fRawPt"));
fEtaPhi = dynamic_cast<TH2F*> (fOutput->FindObject("fEtaPhi"));
+ for (Int_t i=0; i<2; ++i)
+ fZPhi[i] = dynamic_cast<TH2F*> (fOutput->FindObject(Form("fZPhi_%d", i)));
fDeltaPhi = dynamic_cast<TH1F*> (fOutput->FindObject("fDeltaPhi"));
+ fFiredChips = dynamic_cast<TH2F*> (fOutput->FindObject("fFiredChips"));
+ fTriggerVsTime = dynamic_cast<TGraph*> (fOutput->FindObject("fTriggerVsTime"));
+ fStats = dynamic_cast<TH1F*> (fOutput->FindObject("fStats"));
fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fOutput->FindObject("fEsdTrackCuts"));
}
fMultVtx->Draw("SAME");
}
+ if (fFiredChips)
+ {
+ new TCanvas;
+ fFiredChips->Draw("COLZ");
+ }
+
if (fPartEta[0])
{
Int_t events1 = (Int_t) fEvents->Integral(fEvents->GetXaxis()->FindBin(-19.9), fEvents->GetXaxis()->FindBin(-0.001));
if (fPhi)
fPhi->Write();
+ if (fRawPt)
+ fRawPt->Write();
+
if (fEtaPhi)
fEtaPhi->Write();
+ for (Int_t i=0; i<2; ++i)
+ if (fZPhi[i])
+ fZPhi[i]->Write();
+
+ if (fFiredChips)
+ fFiredChips->Write();
+
+ if (fTriggerVsTime)
+ fTriggerVsTime->Write();
+
+ if (fStats)
+ fStats->Write();
+
fVertex = dynamic_cast<TH3F*> (fOutput->FindObject("vertex_check"));
if (fVertex)
fVertex->Write();
if (fOutput)
{
fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
+ fdNdEtaAnalysisND = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaND"));
fdNdEtaAnalysisNSD = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaNSD"));
fdNdEtaAnalysisTr = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTr"));
fdNdEtaAnalysisTrVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndetaTrVtx"));
}
fdNdEtaAnalysis->Finish(0, -1, AlidNdEtaCorrection::kNone);
+ fdNdEtaAnalysisND->Finish(0, -1, AlidNdEtaCorrection::kNone);
fdNdEtaAnalysisNSD->Finish(0, -1, AlidNdEtaCorrection::kNone);
fdNdEtaAnalysisTr->Finish(0, -1, AlidNdEtaCorrection::kNone);
fdNdEtaAnalysisTrVtx->Finish(0, -1, AlidNdEtaCorrection::kNone);
TFile* fout = new TFile("analysis_mc.root","RECREATE");
fdNdEtaAnalysis->SaveHistograms();
+ fdNdEtaAnalysisND->SaveHistograms();
fdNdEtaAnalysisNSD->SaveHistograms();
fdNdEtaAnalysisTr->SaveHistograms();
fdNdEtaAnalysisTrVtx->SaveHistograms();
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff