}
//____________________________________________________________________
-dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::AnalysisMode analysisMode) :
+dNdEtaAnalysis::dNdEtaAnalysis(const Char_t* name, const Char_t* title, AliPWG0Helper::AnalysisMode analysisMode) :
TNamed(name, title),
fData(0),
fMult(0),
// the measured result is not used up to a multiplicity of multCut (the bin at multCut is the first that is used!)
//
- fTag.Form("Correcting dN/deta spectrum >>> %s <<<. Correction type: %d, pt cut: %.2f.", tag, (Int_t) correctionType, ptCut);
+ fTag.Form("Correcting dN/deta spectrum (data: %d) >>> %s <<<. Correction type: %d, pt cut: %.2f.", (Int_t) fAnalysisMode, tag, (Int_t) correctionType, ptCut);
Printf("\n\n%s", fTag.Data());
// set corrections to 1
if (correction && correctionType != AlidNdEtaCorrection::kNone)
{
- TH3F* trackCorr = fData->GetTrackCorrection()->GetCorrectionHistogram();
- TH2F* eventCorr = fData->GetEventCorrection()->GetCorrectionHistogram();
+ TH3* trackCorr = fData->GetTrackCorrection()->GetCorrectionHistogram();
+ TH2* eventCorr = fData->GetEventCorrection()->GetCorrectionHistogram();
if (correctionType >= AlidNdEtaCorrection::kTrack2Particle)
trackCorr->Multiply(correction->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetCorrectionHistogram());
TH2* eTrig = correction->GetVertexRecoCorrection()->GetEventCorrection()->GetGeneratedHistogram();
TH2* eTrigVtx = correction->GetVertexRecoCorrection()->GetEventCorrection()->GetMeasuredHistogram();
+ //TH1* eTrigVtx_projx = eTrigVtx->ProjectionX("eTrigVtx_projx", 2, rawMeasured->GetNbinsY()+1);
TH1* eTrigVtx_projx = eTrigVtx->ProjectionX("eTrigVtx_projx", 2, rawMeasured->GetNbinsY()+1);
//new TCanvas; correctedEvents->DrawCopy("TEXT");
Printf("%d triggered events with 0 mult. -- %d triggered events with vertex", triggeredEventsWith0Mult, allEventsWithVertex);
+ TH1* kineBias = (TH1*) vertexDist->Clone("kineBias");
+ kineBias->Reset();
+
for (Int_t i = 1; i <= rawMeasured->GetNbinsX(); i++)
{
Double_t alpha = (Double_t) vertexDist->GetBinContent(i) / allEventsWithVertex;
Double_t fZ = eTrigVtx_projx->Integral(0, eTrigVtx_projx->GetNbinsX()+1) / eTrigVtx_projx->GetBinContent(i) *
eTrig->GetBinContent(i, 1) / eTrig->Integral(0, eTrig->GetNbinsX()+1, 1, 1);
+ kineBias->SetBinContent(i, fZ);
events *= fZ;
}
//new TCanvas; correctedEvents->DrawCopy("TEXT");
+ //new TCanvas; kineBias->DrawCopy();
}
fData->PrintInfo(ptCut);
- TH3F* dataHist = fData->GetTrackCorrection()->GetGeneratedHistogram();
+ TH3* dataHist = fData->GetTrackCorrection()->GetGeneratedHistogram();
// integrate multiplicity axis out (include under/overflow bins!!!)
- TH2F* tmp = fData->GetEventCorrection()->GetGeneratedHistogram();
+ TH2* tmp = fData->GetEventCorrection()->GetGeneratedHistogram();
TH1D* vertexHist = (TH1D*) tmp->ProjectionX("_px", 0, tmp->GetNbinsY() + 1, "e");
// create pt hist
+ if (fAnalysisMode == AliPWG0Helper::kTPC || fAnalysisMode == AliPWG0Helper::kTPCITS)
{
// reset all ranges
dataHist->GetXaxis()->SetRange(0, 0);
dataHist->GetYaxis()->SetRange(0, 0);
dataHist->GetZaxis()->SetRange(0, 0);
- // integrate over pt (with pt cut)
+ // integrate over pt (with pt cut) (TPC, TPCITS) or multiplicity (SPD)
Int_t ptLowBin = 1;
- if (ptCut > 0)
+ if (ptCut > 0 && fAnalysisMode != AliPWG0Helper::kSPD)
ptLowBin = dataHist->GetZaxis()->FindBin(ptCut);
+
+ //new TCanvas; dataHist->DrawCopy();
+ //dataHist->Sumw2();
dataHist->GetZaxis()->SetRange(ptLowBin, dataHist->GetZaxis()->GetNbins()+1);
- printf("pt range %d %d\n", ptLowBin, dataHist->GetZaxis()->GetNbins()+1);
+ printf("pt/multiplicity range %d %d\n", ptLowBin, dataHist->GetZaxis()->GetNbins()+1);
TH2D* vtxVsEta = dynamic_cast<TH2D*> (dataHist->Project3D("yx2e"));
+
+ //new TCanvas; vtxVsEta->Draw("COLZ");
dataHist->GetZaxis()->SetRange(0, 0);
vtxVsEta->GetXaxis()->SetTitle(dataHist->GetXaxis()->GetTitle());
if (vtxVsEta == 0)
{
- printf("ERROR: pt integration failed\n");
+ printf("ERROR: pt/multiplicity integration failed\n");
return;
}
+ //new TCanvas(tag, tag, 800, 600); vtxVsEta->DrawCopy("COLZ");
+
// clear result histograms
for (Int_t vertexPos=0; vertexPos<kVertexBinning; ++vertexPos)
{
Int_t vertexBinEnd = vertexHist->GetXaxis()->FindBin(vertexRangeEnd[vertexPos]);
const Int_t *binBegin = 0;
- // adjust acceptance range for SPD
+ const Int_t maxBins = 60;
+
+ // adjust acceptance range
+ // produce with drawPlots.C: DetermineAcceptance(...)
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 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 binBeginSPD[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 }; // by eye
+ //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;
+ //const Int_t binBeginSPD[30] = { -1, -1, -1, -1, 16, 15, 13, 11, 9, 8, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 3, -1, -1, -1, -1}; // limit 2
+ const Int_t binBeginSPD[maxBins] = {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, 15, 14, 13, 12, 11, 10, 9, 9, 8, 7, 7, 6, 6, 5, 5, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}; //limit 5
+
+ binBegin = binBeginSPD;
}
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};
+ //const Int_t binBeginTPC[30] = {-1, -1, -1, -1, -1, -1, -1, -1, -1, 4, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, -1, -1, -1, -1, -1, -1, -1, -1, -1}; // limit 5, pt cut off 0.2 mev/c
+ const Int_t binBeginTPC[maxBins] = {-1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, 9, 5, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1}; // limit 5
+
binBegin = binBeginTPC;
}
else if (fAnalysisMode == AliPWG0Helper::kTPCITS)
// TODO create map
}
- Int_t vtxBegin = binBegin[iEta - 1];
- Int_t vtxEnd = 18 + 1 - binBegin[30 - iEta];
+ Int_t vtxBegin = 1;
+ Int_t vtxEnd = maxBins;
+
+ if (binBegin)
+ {
+ vtxBegin = binBegin[iEta - 1];
+ vtxEnd = 18 + 1 - binBegin[maxBins - iEta];
+ }
+ else
+ Printf("WARNING: No acceptance applied!");
// 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);
+
+ //Float_t eta = vtxVsEta->GetYaxis()->GetBinCenter(iEta);
+ //printf("Eta bin: %d (%f) Vertex range: %d; before: %d %d (range) %d %d (acceptance)", iEta, eta, vertexPos, vertexBinBegin, vertexBinEnd, vtxBegin, vtxEnd);
vertexBinBegin = TMath::Max(vertexBinBegin, vtxBegin);
vertexBinEnd = TMath::Min(vertexBinEnd, vtxEnd);
- //Printf("after: %d %d", vertexBinBegin, vertexBinEnd);
+ //Printf(" after: %d %d", vertexBinBegin, vertexBinEnd);
// no data for this bin
if (vertexBinBegin > vertexBinEnd)
continue;
}
- Float_t totalEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd);
- if (totalEvents == 0)
- {
- printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
- continue;
- }
-
+ Float_t totalEvents = 0;
Float_t sum = 0;
Float_t sumError2 = 0;
- for (Int_t iVtx = vertexBinBegin; iVtx <= vertexBinEnd; iVtx++)
+ Float_t unusedTracks = 0;
+ Float_t unusedEvents = 0;
+ for (Int_t iVtx = 1; iVtx <= vtxVsEta->GetNbinsX(); iVtx++)
{
- if (vtxVsEta->GetBinContent(iVtx, iEta) != 0)
+ if (iVtx >= vertexBinBegin && iVtx <= vertexBinEnd)
{
- sum = sum + vtxVsEta->GetBinContent(iVtx, iEta);
+ if (vtxVsEta->GetBinContent(iVtx, iEta) != 0)
+ {
+ sum += vtxVsEta->GetBinContent(iVtx, iEta);
- if (sumError2 > 10e30)
- Printf("WARNING: sum of error2 is dangerously large - be prepared for crash... ");
+ if (sumError2 > 10e30)
+ Printf("WARNING: sum of error2 is dangerously large - be prepared for crash... ");
- sumError2 = sumError2 + TMath::Power(vtxVsEta->GetBinError(iVtx, iEta),2);
+ sumError2 = sumError2 + TMath::Power(vtxVsEta->GetBinError(iVtx, iEta),2);
+ }
+ totalEvents += vertexHist->GetBinContent(iVtx);
}
+ else
+ {
+ unusedTracks += vtxVsEta->GetBinContent(iVtx, iEta);
+ unusedEvents += vertexHist->GetBinContent(iVtx);
+ }
+ }
+
+ if (totalEvents == 0)
+ {
+ printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
+ continue;
}
Float_t ptCutOffCorrection = 1;
- if (correction && ptCut > 0)
- ptCutOffCorrection = correction->GetMeasuredFraction(correctionType, ptCut, vtxVsEta->GetYaxis()->GetBinCenter(iEta));
- if (ptCutOffCorrection <= 0)
+ // find pt cut off correction factor
+ if (fAnalysisMode != AliPWG0Helper::kSPD)
{
- printf("UNEXPECTED: ptCutOffCorrection is %f for hist %d %d %d\n", ptCutOffCorrection, vertexPos, vertexBinBegin, vertexBinEnd);
- continue;
+ if (correction && ptCut > 0)
+ ptCutOffCorrection = correction->GetMeasuredFraction(correctionType, ptCut, vtxVsEta->GetYaxis()->GetBinCenter(iEta), vertexBinBegin, vertexBinEnd);
+
+ if (ptCutOffCorrection <= 0)
+ {
+ printf("UNEXPECTED: ptCutOffCorrection is %f for hist %d %d %d\n", ptCutOffCorrection, vertexPos, vertexBinBegin, vertexBinEnd);
+ continue;
+ }
}
- //printf("Eta: %d Vertex Range: %d %d, Event Count %f, Track Sum: %f, Track Sum corrected: %f\n", iEta, vertexBinBegin, vertexBinEnd, totalEvents, sum, sum / ptCutOffCorrection);
+ //printf("Eta: %d (%f) Vertex Range: %d %d, Event Count %f, Track Sum: %f, Track Sum corrected: %f \n", iEta, vtxVsEta->GetYaxis()->GetBinCenter(iEta), vertexBinBegin, vertexBinEnd, totalEvents, sum, sum / ptCutOffCorrection);
Int_t bin = fdNdEta[vertexPos]->FindBin(vtxVsEta->GetYaxis()->GetBinCenter(iEta));
if (bin > 0 && bin <= fdNdEta[vertexPos]->GetNbinsX())
fdNdEtaPtCutOffCorrected[vertexPos]->SetBinContent(bin, dndeta);
fdNdEtaPtCutOffCorrected[vertexPos]->SetBinError(bin, error);
- //Printf("Bin %d has dN/deta = %f +- %f; %f +- %f", bin, dndeta, error, fdNdEta[vertexPos]->GetBinContent(bin), fdNdEta[vertexPos]->GetBinError(bin));
+ //Printf("Bin %d has dN/deta = %f +- %f; %.2f tracks %.2f events (outside acceptance: %.2f tracks, %.2f events)", bin, dndeta, error, sum, totalEvents, unusedTracks, unusedEvents);
}
}
}
fPtDist = dynamic_cast<TH1F*> (gDirectory->Get(fPtDist->GetName()));
- if (dynamic_cast<TNamed*> (gFile->Get("fTag")))
- fTag = (dynamic_cast<TNamed*> (gFile->Get("fTag")))->GetTitle();
+ if (dynamic_cast<TNamed*> (gDirectory->Get("fTag")))
+ fTag = (dynamic_cast<TNamed*> (gDirectory->Get("fTag")))->GetTitle();
- if (dynamic_cast<TParameter<Int_t>*> (gFile->Get("fAnalysisMode")))
- fAnalysisMode = (AliPWG0Helper::AnalysisMode) (dynamic_cast<TParameter<Int_t>*> (gFile->Get("fAnalysisMode")))->GetVal();
+ if (dynamic_cast<TParameter<Int_t>*> (gDirectory->Get("fAnalysisMode")))
+ fAnalysisMode = (AliPWG0Helper::AnalysisMode) (dynamic_cast<TParameter<Int_t>*> (gDirectory->Get("fAnalysisMode")))->GetVal();
gDirectory->cd("../");
}