Int_t markers[] = {20,20,21,22,23,28,29};
Int_t colors[] = {1,2,3,4,6,8,102};
+void loadlibs()
+{
+ gSystem->Load("libTree");
+ gSystem->Load("libVMC");
+
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libPWG0base");
+}
+
void SetRanges(TAxis* axis)
{
if (strcmp(axis->GetTitle(), "#eta") == 0)
for (Int_t i=0; i<folderCount; ++i)
{
- Track2Particle1DCreatePlots(fileNames[i], folderNames[i], upperPtLimit);
+ Correction1DCreatePlots(fileNames[i], folderNames[i], upperPtLimit);
- TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_x_div_meas_%s_nTrackToNPart_x", folderNames[i], folderNames[i])));
- TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_y_div_meas_%s_nTrackToNPart_y", folderNames[i], folderNames[i])));
- TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", folderNames[i], folderNames[i])));
+ TH1* corrX = dynamic_cast<TH1*> (gROOT->FindObject("generated_x_div_measured_x"));
+ TH1* corrY = dynamic_cast<TH1*> (gROOT->FindObject("generated_y_div_measured_y"));
+ TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject("generated_z_div_measured_z"));
Prepare1DPlot(corrX);
Prepare1DPlot(corrY);
return ptDists;
}
-void DrawpiKpAndCombinedZOnly(Float_t upperPtLimit=0.99)
-{
- gROOT->ProcessLine(".L drawPlots.C");
- gSystem->Load("libPWG0base");
-
- const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "correction_map.root" };
- const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "dndeta_correction" };
- const char* legendNames[] = { "#pi", "K", "p", "standard" };
- Int_t folderCount = 3;
-
- /*const char* fileNames[] = { "systematics.root", "systematics.root", "systematics.root", "systematics.root", "correction_map.root" };
- const char* folderNames[] = { "correction_0", "correction_1", "correction_2", "correction_3", "dndeta_correction" };
- const char* legendNames[] = { "#pi", "K", "p", "others", "standard" };
- Int_t folderCount = 5;*/
-
- TString canvasName;
- canvasName.Form("Track2Particle1DComposition");
- TCanvas* canvas = new TCanvas(canvasName, canvasName, 700, 500);
- canvas->SetGridx();
- canvas->SetGridy();
- canvas->SetBottomMargin(0.12);
- //InitPad();
-
- TLegend* legend = new TLegend(0.8, 0.7, 0.95, 0.95);
- legend->SetFillColor(0);
-
- Int_t mycolors[] = {1, 2, 4};
-
- for (Int_t i=0; i<folderCount; ++i)
- {
- Track2Particle1DCreatePlots(fileNames[i], folderNames[i], upperPtLimit);
-
- TH1* corrZ = dynamic_cast<TH1*> (gROOT->FindObject(Form("gene_%s_nTrackToNPart_z_div_meas_%s_nTrackToNPart_z", folderNames[i], folderNames[i])));
-
- Prepare1DPlot(corrZ);
-
- corrZ->SetTitle("");
- corrZ->GetXaxis()->SetRangeUser(0, upperPtLimit);
- corrZ->GetYaxis()->SetRangeUser(0.51, 6);
- corrZ->SetMarkerColor(mycolors[i]);
- corrZ->SetLineColor(mycolors[i]);
- corrZ->SetMarkerStyle(markers[i+1]);
- corrZ->GetYaxis()->SetTitle("correction factor");
-
- corrZ->DrawCopy(((i>0) ? "SAMEP" : "P"));
-
- legend->AddEntry(corrZ, legendNames[i]);
- }
-
- legend->Draw();
-
- canvas->SaveAs("ptcutoff_species.eps");
-}
-
void DrawCompareToReal()
{
gROOT->ProcessLine(".L drawPlots.C");
break;
// one species enhanced / reduced
- case 2: // + 50% pions
- case 3: // - 50% pions
- case 4: // + 50% kaons
- case 5: // - 50% kaons
- case 6: // + 50% protons
- case 7: // - 50% protons
- Int_t correctionIndex = (index - 2) / 2;
- Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5;
-
- fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Scale(scaleFactor);
- fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Scale(scaleFactor);
-
+ case 2: // + 50% kaons
+ case 3: // - 50% kaons
+ case 4: // + 50% protons
+ case 5: // - 50% protons
+ case 6: // + 50% kaons + 50% protons
+ case 7: // - 50% kaons - 50% protons
+ case 8: // + 50% kaons - 50% protons
+ case 9: // - 50% kaons + 50% protons
TString* str = new TString;
- str->Form("%s%s", (correctionIndex == 0) ? "Pi" : ((correctionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced");
+ if (index < 6)
+ {
+ Int_t correctionIndex = index / 2;
+ Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5;
+
+ fdNdEtaCorrection[correctionIndex]->GetTrack2ParticleCorrection()->GetTrackCorrection()->Scale(scaleFactor);
+ str->Form("%s%s", (correctionIndex == 0) ? "Pi" : ((correctionIndex == 1) ? "K" : (correctionIndex == 2) ? "p" : "others"), (index % 2 == 0) ? "Boosted" : "Reduced");
+ }
+ else
+ {
+ Double_t scaleFactor = (index % 2 == 0) ? 1.5 : 0.5;
+ fdNdEtaCorrection[1]->GetTrack2ParticleCorrection()->GetTrackCorrection()->Scale(scaleFactor);
+ str->Form("%s%s", "K", (scaleFactor > 1) ? "Boosted" : "Reduced");
+
+ if (index >= 8)
+ scaleFactor = (index % 2 == 0) ? 0.5 : 1.5;
+ fdNdEtaCorrection[2]->GetTrack2ParticleCorrection()->GetTrackCorrection()->Scale(scaleFactor);
+ *str += Form("%s%s", "p", (scaleFactor > 1) ? "Boosted" : "Reduced");
+ }
+
return str->Data();
break;
for (Int_t i=0; i<3; ++i)
{
- ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDists[i]);
- ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDists[i]);
+ ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram(), ptDists[i]);
+ ScalePtDependent(fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetGeneratedHistogram(), ptDists[i]);
}
TString* str = new TString;
str->Form("%s%s", (functionIndex == 0) ? "Pi" : ((functionIndex == 1) ? "K" : "p"), (index % 2 == 0) ? "Boosted" : "Reduced");
case 999:
TF1* ptDependence = new TF1("simple", "x", 0, 100);
- ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetGeneratedHistogram(), ptDependence);
- ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetMeasuredHistogram(), ptDependence);
+ ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetGeneratedHistogram(), ptDependence);
+ ScalePtDependent(fdNdEtaCorrection[0]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram(), ptDependence);
break;
}
void Composition()
{
- gSystem->Load("libPWG0base");
+ loadlibs();
gSystem->Unlink("new_compositions.root");
+ gSystem->Unlink("new_compositions_analysis.root");
+
+ const char* names[] = { "pi", "K", "p", "other" };
+ TH1* hRatios[20];
- Int_t nCompositions = 8;
- for (Int_t comp = 0; comp < nCompositions; ++comp)
+ Int_t nCompositions = 10;
+ Int_t counter = 0;
+ for (Int_t comp = 1; comp < nCompositions; ++comp)
{
AlidNdEtaCorrection* fdNdEtaCorrection[4];
- TFile::Open("systematics.root");
+ TFile::Open("correction_mapparticle-species.root");
for (Int_t i=0; i<4; ++i)
{
TString name;
- name.Form("correction_%d", i);
+ name.Form("dndeta_correction_%s", names[i]);
fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name);
- fdNdEtaCorrection[i]->LoadHistograms("systematics.root", name);
+ fdNdEtaCorrection[i]->LoadHistograms();
}
const char* newName = ChangeComposition(fdNdEtaCorrection, comp);
for (Int_t i=0; i<4; ++i)
{
- geneCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetGeneratedHistogram()->Integral();
- measCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetMeasuredHistogram()->Integral();
+ geneCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetGeneratedHistogram()->Integral();
+ measCount[i] = fdNdEtaCorrection[i]->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetMeasuredHistogram()->Integral();
geneCount[4] += geneCount[i];
measCount[4] += measCount[i];
TList* collection = new TList;
- for (Int_t i=0; i<4; ++i)
+ // skip "other" particle correction here
+ // with them has to be dealt differently, maybe just increasing the neutral particles...
+ for (Int_t i=1; i<3; ++i)
collection->Add(fdNdEtaCorrection[i]);
- AlidNdEtaCorrection* newComposition = new AlidNdEtaCorrection(newName, newName);
- newComposition->Merge(collection);
- newComposition->Finish();
+ fdNdEtaCorrection[0]->Merge(collection);
+ fdNdEtaCorrection[0]->Finish();
delete collection;
+ // save everything
TFile* file = TFile::Open("new_compositions.root", "UPDATE");
- newComposition->SaveHistograms();
+ fdNdEtaCorrection[0]->SetName(newName);
+ fdNdEtaCorrection[0]->SaveHistograms();
//file->Write();
file->Close();
- }
-
- gROOT->ProcessLine(".L drawPlots.C");
-
- const char* fileNames[] = {"new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root" };
- const char* folderNames[] = { "Pythia", "PythiaRatios", "PiBoosted", "PiReduced", "KBoosted", "KReduced", "pBoosted", "pReduced" };
-
- Track2Particle1DComposition(fileNames, nCompositions, folderNames);
-}
-
-Double_t ConvSigma1To2D(Double_t sigma)
-{
- return TMath::Sqrt( - TMath::Log( 1 - TMath::Erf(sigma / TMath::Sqrt(2)) ) * 2);
-}
-
-Double_t ConvDistance1DTo2D(Double_t distance)
-{
- return TMath::ErfInverse(1 - TMath::Exp(-distance * distance / 2)) * TMath::Sqrt(2);
-}
-
-Double_t Sigma2VertexCount(TH2F* tracks, Double_t nSigma)
-{
- Double_t count = 0;
-
- //nSigma = ConvSigma1To2D(nSigma);
-
- for (Int_t x=1; x<=tracks->GetNbinsX(); ++x)
- for (Int_t y=1; y<=tracks->GetNbinsY(); ++y)
- {
- Double_t impactX = tracks->GetXaxis()->GetBinCenter(x);
- Double_t impactY = tracks->GetYaxis()->GetBinCenter(y);
-
- Float_t d = TMath::Sqrt(impactX*impactX + impactY*impactY);
-
- d = ConvDistance1DTo2D(d);
-
- if (d < nSigma)
- count += tracks->GetBinContent(x, y);
- }
-
- return count;
-}
-
-TH2F* Sigma2VertexGaussianTracksHist()
-{
- TH2F* tracks = new TH2F("Sigma2Vertex_tracks", "Sigma2Vertex_tracks", 200, -5, 5, 200, -5, 5);
-
- TF2* gaussian2D = new TF2("gaussian2D", "xgausn(0) * ygausn(3)", -5, 5, -5, 5);
- gaussian2D->SetParameters(1, 0, 1, 1, 0, 1);
-
- for (Int_t x=1; x<=tracks->GetNbinsX(); ++x)
- for (Int_t y=1; y<=tracks->GetNbinsY(); ++y)
- tracks->SetBinContent(x, y, gaussian2D->Eval(tracks->GetXaxis()->GetBinCenter(x), tracks->GetYaxis()->GetBinCenter(y)));
-
- //normalize
- tracks->Scale(1.0 / tracks->Integral());
-
- return tracks;
-}
-
-TH1F* Sigma2VertexGaussian()
-{
- TH2F* tracks = Sigma2VertexGaussianTracksHist();
-
- TCanvas* canvas = new TCanvas("Sigma2VertexGaussian", "Sigma2VertexGaussian", 1000, 1000);
- canvas->Divide(2, 2);
-
- canvas->cd(1);
- tracks->Draw("COLZ");
-
- TH1F* ratio = new TH1F("Sigma2Vertex_ratio", "Sigma2Vertex_ratio;n sigma;included", 50, 0.05, 5.05);
- for (Double_t nSigma = 0.1; nSigma < 5.05; nSigma += 0.1)
- ratio->Fill(nSigma, Sigma2VertexCount(tracks, nSigma));
- ratio->SetMarkerStyle(21);
-
- canvas->cd(2);
- ratio->DrawCopy("P");
-
- TH1F* ratio2 = new TH1F("Sigma2Vertex_ratio2", "Sigma2Vertex_ratio2;nSigma;% included 3 sigma / % included n sigma", 50, 0.05, 5.05);
- Double_t sigma3 = Sigma2VertexCount(tracks, 3);
- for (Double_t nSigma = 0.1; nSigma < 5.05; nSigma += 0.1)
- ratio2->Fill(nSigma, sigma3 / ratio->GetBinContent(ratio->FindBin(nSigma)));
- ratio2->SetMarkerStyle(21);
+
+ // correct dNdeta distribution with modified correction map
+ TFile::Open("analysis_esd_raw.root");
- canvas->cd(3);
- ratio2->DrawCopy("P");
+ dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD");
+ fdNdEtaAnalysis->LoadHistograms();
- canvas->SaveAs("Sigma2Vertex.eps");
+ fdNdEtaAnalysis->Finish(fdNdEtaCorrection[0], 0.2, 3, newName);
+
+ hRatios[counter] = (TH1F*) fdNdEtaAnalysis->GetdNdEtaHistogram()->Clone(newName);
+ hRatios[counter]->SetTitle(newName);
+ hRatios[counter]->SetYTitle("dN_{ch}/d#eta ratio #frac{default composition}{modified composition}");
- return ratio2;
-}
+ if (counter > 0)
+ hRatios[counter]->Divide(hRatios[0],hRatios[counter],1,1);
-TH1F* Sigma2VertexSimulation()
-{
- TFile* file = TFile::Open("systematics.root");
+ file = TFile::Open("new_compositions_analysis.root", "UPDATE");
+ hRatios[counter]->Write();
+ file->Close();
+
+ delete fdNdEtaAnalysis;
- TH1F* sigmavertex = dynamic_cast<TH1F*> (file->Get("fSigmaVertex"));
- if (!sigmavertex)
- {
- printf("Could not read histogram\n");
- return;
+ counter++;
}
- TH1F* ratio = new TH1F("sigmavertexsimulation_ratio", "sigmavertexsimulation_ratio;N#sigma;% included in 3 #sigma / % included in N#sigma", sigmavertex->GetNbinsX(), sigmavertex->GetXaxis()->GetXmin(), sigmavertex->GetXaxis()->GetXmax());
-
- for (Int_t i=1; i<=sigmavertex->GetNbinsX(); ++i)
- ratio->SetBinContent(i, sigmavertex->GetBinContent(sigmavertex->GetXaxis()->FindBin(3)) / sigmavertex->GetBinContent(i));
-
- TCanvas* canvas = new TCanvas("Sigma2VertexSimulation", "Sigma2VertexSimulation", 1000, 500);
- canvas->Divide(2, 1);
-
- canvas->cd(1);
- sigmavertex->SetMarkerStyle(21);
- sigmavertex->Draw("P");
+ /*
+ gROOT->ProcessLine(".L drawPlots.C");
- canvas->cd(2);
- ratio->SetMarkerStyle(21);
- ratio->DrawCopy("P");
+ const char* fileNames[] = {"new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root", "new_compositions.root" };
+ const char* folderNames[] = { "PythiaRatios", "PiBoosted", "PiReduced", "KBoosted", "KReduced", "pBoosted", "pReduced" };
- return ratio;
+ Track2Particle1DComposition(fileNames, nCompositions, folderNames);
+ */
}
-void Sigma2VertexCompare()
-{
- TH1F* ratio1 = Sigma2VertexGaussian();
- TH1F* ratio2 = Sigma2VertexSimulation();
-
- ratio1->SetStats(kFALSE);
- ratio2->SetStats(kFALSE);
-
- ratio1->SetMarkerStyle(0);
- ratio2->SetMarkerStyle(0);
-
- ratio1->SetLineWidth(2);
- ratio2->SetLineWidth(2);
-
- TLegend* legend = new TLegend(0.7, 0.8, 0.95, 0.95);
- legend->SetFillColor(0);
- legend->AddEntry(ratio1, "Gaussian");
- legend->AddEntry(ratio2, "Simulation");
-
- ratio2->SetTitle("");
- ratio2->GetYaxis()->SetTitleOffset(1.5);
- ratio2->GetXaxis()->SetRangeUser(2, 4);
-
- TCanvas* canvas = new TCanvas("Sigma2VertexCompare", "Sigma2VertexCompare", 500, 500);
- InitPad();
-
- ratio2->SetMarkerStyle(21);
- ratio1->SetMarkerStyle(22);
-
- ratio2->GetYaxis()->SetRangeUser(0.8, 1.2);
- ratio2->SetLineColor(kRed);
- ratio2->SetMarkerColor(kRed);
- ratio2->Draw("PL");
- ratio1->Draw("SAMEPL");
-
- legend->Draw();
-
- canvas->SaveAs("Sigma2VertexCompare.eps");
-}
void drawSystematics()
{
canvas2->SaveAs("particlecomposition_result.eps");
}
-mergeCorrectionsWithDifferentCrosssections(Char_t* standardCorrectionFileName="correction_map.root",
- Char_t* systematicCorrectionFileName="systematics.root",
- Char_t* outputFileName="systematics_vtxtrigger_compositions.root") {
+void mergeCorrectionsWithDifferentCrosssections(Int_t correctionTarget = 3, Char_t* correctionFileName="correction_mapprocess-types.root", const char* analysisFileName = "analysis_esd_raw.root", const Char_t* outputFileName="systematics_vtxtrigger_compositions.root") {
//
// Function used to merge standard corrections with vertex
// reconstruction corrections obtained by a certain mix of ND, DD
// and SD events.
- //
+ //
+ // the dn/deta spectrum is corrected and the ratios
+ // (standard to changed x-section) of the different dN/deta
+ // distributions are saved to a file.
+ //
+ // correctionTarget is of type AlidNdEtaCorrection::CorrectionType
+ // kINEL = 3
+ // kNSD = 4
- gSystem->Load("libPWG0base");
+ loadlibs();
const Char_t* typeName[] = { "vertexreco", "trigger", "vtxtrigger" };
- const Char_t* changes[] = {"pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless"};
- Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.0, 1.0, 1.5, 0.5};
- Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5};
+
+
+ const Char_t* changes[] = { "pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless", "sdlessddmore", "sdmoreddless", "ddmore25","ddless25","sdmore25","sdless25", "dmore25", "dless25", "sdlessddmore25", "sdmoreddless25"};
+ Float_t scalesND[] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0 };
+ //Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5, 1.25, 0.75, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75};
+ //Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.5, 0.5, 1.5, 0.5, 0.5, 1.5, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75, 0.75, 1.25};
+ Float_t scalesDD[] = {1.0, 1.4, 0.6, 1.0, 1.0, 1.4, 0.6, 1.4, 0.6, 1.25, 0.75, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75};
+ Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.4, 0.6, 1.4, 0.6, 0.4, 1.6, 1.0, 1.0, 1.25, 0.75, 1.25, 0.75, 0.75, 1.25};
+ Int_t nChanges = 17;
+
+ /*
+ const Char_t* changes[] = { "pythia", "qgsm", "phojet"};
+ Float_t scalesND[] = {1.0, 1.10, 1.11};
+ Float_t scalesSD[] = {1.0, 0.69, 0.86};
+ Float_t scalesDD[] = {1.0, 0.98, 0.61};
+ Int_t nChanges = 3;
+ */
+
// cross section from Pythia
+ // 14 TeV!
Float_t sigmaND = 55.2;
Float_t sigmaDD = 9.78;
Float_t sigmaSD = 14.30;
- AlidNdEtaCorrection* corrections[21];
+ // standard correction
+ TFile::Open(correctionFileName);
+ AlidNdEtaCorrection* correctionStandard = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction");
+ correctionStandard->LoadHistograms();
+
+ // dont take vertexreco from this one
+ correctionStandard->GetVertexRecoCorrection()->Reset();
+ // dont take triggerbias from this one
+ correctionStandard->GetTriggerBiasCorrectionINEL()->Reset();
+ correctionStandard->GetTriggerBiasCorrectionNSD()->Reset();
+ correctionStandard->GetTriggerBiasCorrectionND()->Reset();
+
+ AlidNdEtaCorrection* corrections[100];
+ TH1F* hRatios[100];
+
+ Int_t counter = 0;
for (Int_t j=0; j<3; j++) { // j = 0 (change vtx), j = 1 (change trg), j = 2 (change both)
- AlidNdEtaCorrection* correctionStandard = new AlidNdEtaCorrection("dndeta_correction","dndeta_correction");
- correctionStandard->LoadHistograms(standardCorrectionFileName);
- // dont take vertexreco from this one
- correctionStandard->GetVertexRecoCorrection()->Reset();
- // dont take triggerbias from this one
- correctionStandard->GetTriggerBiasCorrectionINEL()->Reset();
+ for (Int_t i=0; i<nChanges; i++) {
+ TFile::Open(correctionFileName);
- for (Int_t i=0; i<7; i++) {
TString name;
name.Form("dndeta_correction_syst_%s_%s", typeName[j], changes[i]);
AlidNdEtaCorrection* current = new AlidNdEtaCorrection(name, name);
-
- name.Form("vertexRecoND");
- AlidNdEtaCorrection* dNdEtaCorrectionVtxND = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionVtxND->LoadHistograms(systematicCorrectionFileName, name);
- name.Form("vertexRecoDD");
- AlidNdEtaCorrection* dNdEtaCorrectionVtxDD = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionVtxDD->LoadHistograms(systematicCorrectionFileName, name);
- name.Form("vertexRecoSD");
- AlidNdEtaCorrection* dNdEtaCorrectionVtxSD = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionVtxSD->LoadHistograms(systematicCorrectionFileName, name);
-
- name.Form("triggerBiasND");
- AlidNdEtaCorrection* dNdEtaCorrectionTriggerND = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionTriggerND->LoadHistograms(systematicCorrectionFileName, name);
- name.Form("triggerBiasDD");
- AlidNdEtaCorrection* dNdEtaCorrectionTriggerDD = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionTriggerDD->LoadHistograms(systematicCorrectionFileName, name);
- name.Form("triggerBiasSD");
- AlidNdEtaCorrection* dNdEtaCorrectionTriggerSD = new AlidNdEtaCorrection(name,name);
- dNdEtaCorrectionTriggerSD->LoadHistograms(systematicCorrectionFileName, name);
+ current->LoadHistograms("dndeta_correction");
+ current->Reset();
+
+ name.Form("dndeta_correction_ND");
+ AlidNdEtaCorrection* dNdEtaCorrectionND = new AlidNdEtaCorrection(name,name);
+ dNdEtaCorrectionND->LoadHistograms();
+ name.Form("dndeta_correction_DD");
+ AlidNdEtaCorrection* dNdEtaCorrectionDD = new AlidNdEtaCorrection(name,name);
+ dNdEtaCorrectionDD->LoadHistograms();
+ name.Form("dndeta_correction_SD");
+ AlidNdEtaCorrection* dNdEtaCorrectionSD = new AlidNdEtaCorrection(name,name);
+ dNdEtaCorrectionSD->LoadHistograms();
// calculating relative
Float_t nd = 100 * sigmaND/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
Float_t dd = 100 * (scalesDD[i]*sigmaDD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
Float_t sd = 100 * (scalesSD[i]*sigmaSD)/(sigmaND + (scalesDD[i]*sigmaDD) + (scalesDD[i]*sigmaSD));
- printf(Form("%s : ND=%.1f\%, DD=%.1f\%, SD=%.1f\% \n",changes[i],nd,dd,sd));
+ printf(Form("%s : ND=%.2f\%, DD=%.2f\%, SD=%.2f\% \n",changes[i],nd,dd,sd));
current->SetTitle(Form("ND=%.2f\%,DD=%.2f\%,SD=%.2f\%",nd,dd,sd));
current->SetTitle(name);
// scale
if (j == 0 || j == 2)
{
- dNdEtaCorrectionVtxDD->GetVertexRecoCorrection()->GetMeasuredHistogram()->Scale(scalesDD[i]);
- dNdEtaCorrectionVtxSD->GetVertexRecoCorrection()->GetMeasuredHistogram()->Scale(scalesSD[i]);
- dNdEtaCorrectionVtxDD->GetVertexRecoCorrection()->GetGeneratedHistogram()->Scale(scalesDD[i]);
- dNdEtaCorrectionVtxSD->GetVertexRecoCorrection()->GetGeneratedHistogram()->Scale(scalesSD[i]);
+ dNdEtaCorrectionND->GetVertexRecoCorrection()->Scale(scalesND[i]);
+ dNdEtaCorrectionDD->GetVertexRecoCorrection()->Scale(scalesDD[i]);
+ dNdEtaCorrectionSD->GetVertexRecoCorrection()->Scale(scalesSD[i]);
}
if (j == 1 || j == 2)
{
- dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionINEL()->GetMeasuredHistogram()->Scale(scalesDD[i]);
- dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionINEL()->GetMeasuredHistogram()->Scale(scalesSD[i]);
- dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionINEL()->GetGeneratedHistogram()->Scale(scalesDD[i]);
- dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionINEL()->GetGeneratedHistogram()->Scale(scalesSD[i]);
+ dNdEtaCorrectionND->GetTriggerBiasCorrectionINEL()->Scale(scalesND[i]);
+ dNdEtaCorrectionDD->GetTriggerBiasCorrectionINEL()->Scale(scalesDD[i]);
+ dNdEtaCorrectionSD->GetTriggerBiasCorrectionINEL()->Scale(scalesSD[i]);
+
+ dNdEtaCorrectionND->GetTriggerBiasCorrectionNSD()->Scale(scalesND[i]);
+ dNdEtaCorrectionDD->GetTriggerBiasCorrectionNSD()->Scale(scalesDD[i]);
+ dNdEtaCorrectionSD->GetTriggerBiasCorrectionNSD()->Scale(scalesSD[i]);
+
+ dNdEtaCorrectionND->GetTriggerBiasCorrectionND()->Scale(scalesND[i]);
+ dNdEtaCorrectionDD->GetTriggerBiasCorrectionND()->Scale(scalesDD[i]);
+ dNdEtaCorrectionSD->GetTriggerBiasCorrectionND()->Scale(scalesSD[i]);
}
- //clear track, trigger in Vtx correction
- dNdEtaCorrectionVtxND->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionVtxND->GetTriggerBiasCorrectionINEL()->Reset();
- dNdEtaCorrectionVtxSD->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionVtxSD->GetTriggerBiasCorrectionINEL()->Reset();
- dNdEtaCorrectionVtxDD->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionVtxDD->GetTriggerBiasCorrectionINEL()->Reset();
-
- //clear track, vertexreco in trigger correction
- dNdEtaCorrectionTriggerND->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionTriggerND->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionTriggerND->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionTriggerND->GetVertexRecoCorrection()->Reset();
- dNdEtaCorrectionTriggerSD->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionTriggerSD->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionTriggerSD->GetVertexRecoCorrection()->Reset();
- dNdEtaCorrectionTriggerDD->GetTrack2ParticleCorrection()->Reset();
- dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionNSD()->Reset();
- dNdEtaCorrectionTriggerDD->GetTriggerBiasCorrectionND()->Reset();
- dNdEtaCorrectionTriggerDD->GetVertexRecoCorrection()->Reset();
+ //clear track in correction
+ dNdEtaCorrectionND->GetTrack2ParticleCorrection()->Reset();
+ dNdEtaCorrectionDD->GetTrack2ParticleCorrection()->Reset();
+ dNdEtaCorrectionSD->GetTrack2ParticleCorrection()->Reset();
TList collection;
collection.Add(correctionStandard);
- collection.Add(dNdEtaCorrectionVtxND);
- collection.Add(dNdEtaCorrectionVtxDD);
- collection.Add(dNdEtaCorrectionVtxSD);
- collection.Add(dNdEtaCorrectionTriggerND);
- collection.Add(dNdEtaCorrectionTriggerDD);
- collection.Add(dNdEtaCorrectionTriggerSD);
+ collection.Add(dNdEtaCorrectionND);
+ collection.Add(dNdEtaCorrectionDD);
+ collection.Add(dNdEtaCorrectionSD);
current->Merge(&collection);
current->Finish();
- corrections[i+j*7] = current;
- }
- }
+ corrections[counter] = current;
- TFile* fout = new TFile(outputFileName,"RECREATE");
+ // now correct dNdeta distribution with modified correction map
+ TFile::Open(analysisFileName);
- for (Int_t i=0; i<21; i++)
- corrections[i]->SaveHistograms();
+ dNdEtaAnalysis* fdNdEtaAnalysis = new dNdEtaAnalysis("fdNdEtaAnalysisESD", "fdNdEtaAnalysisESD");
+ fdNdEtaAnalysis->LoadHistograms();
- fout->Write();
- fout->Close();
-}
-
-
-DrawVertexRecoSyst() {
- // Draws the ratio of the dN/dEta obtained with changed SD and DD
- // cross-sections vertex reco corrections to the dN/dEta obtained
- // from the standard pythia cross-sections
- //
- // The files with the vertex reco corrections for different
- // processes (and with the other standard corrections) are expected
- // to have the names "analysis_esd_X.root", where the Xs are defined
- // in the array changes below.
-
- Char_t* changes[] = {"pythia","ddmore","ddless","sdmore","sdless", "dmore", "dless"};
- Char_t* descr[] = {"",
- "#sigma_{DD}' = 1.5#times#sigma_{DD}",
- "#sigma_{DD}' = 0.5#times#sigma_{DD}",
- "#sigma_{SD}' = 1.5#times#sigma_{SD}",
- "#sigma_{SD}' = 0.5#times#sigma_{SD}",
- "#sigma_{D}' = 1.5#times#sigma_{D}",
- "#sigma_{D}' = 0.5#times#sigma_{D}"};
-
- Float_t scalesDD[] = {1.0, 1.5, 0.5, 1.5, 0.5, 1.5, 0.5};
- Float_t scalesSD[] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.5, 0.5};
-
- Int_t markers[] = {20,20,21,22,23,28,29};
- Int_t colors[] = {1,2,3,4,6,8,102};
+ fdNdEtaAnalysis->Finish(current, 0.3, correctionTarget, Form("%d %d", j, i));
- // cross section from Pythia
- Float_t sigmaND = 55.2;
- Float_t sigmaDD = 9.78;
- Float_t sigmaSD = 14.30;
+ name = "ratio";
+ if (j==0) name.Append("_vetexReco_");
+ if (j==1) name.Append("_triggerBias_");
+ if (j==2) name.Append("_vertexReco_triggerBias_");
+ name.Append(changes[i]);
- TH1F* dNdEta[7];
- TH1F* ratios[7];
+ hRatios[counter] = (TH1F*) fdNdEtaAnalysis->GetdNdEtaHistogram()->Clone(name);
- TFile* fin;
+ name.Form("DD #times %0.2f, SD #times %0.2f",scalesDD[i],scalesSD[i]);
+ hRatios[counter]->SetTitle(name.Data());
+ hRatios[counter]->SetYTitle("dN_{ch}/d#eta ratio #frac{default cross-section}{modified cross-sections}");
- for (Int_t i=0; i<7; i++) {
- // calculating relative
- fin = TFile::Open(Form("analysis_esd_%s.root",changes[i]));
+ if (counter > 0)
+ hRatios[counter]->Divide(hRatios[0],hRatios[counter],1,1);
- dNdEta[i] = (TH1F*)(fin->Get("dndeta/dndeta_dNdEta_corrected_2"))->Clone();
+ delete fdNdEtaAnalysis;
- for (Int_t b=0; b<dNdEta[i]->GetNbinsX(); b++) {
- if (TMath::Abs(dNdEta[i]->GetBinCenter(b))>0.9) {
- dNdEta[i]->SetBinContent(b,0);
- dNdEta[i]->SetBinError(b,0);
- }
+ counter++;
}
-
- dNdEta[i]->Rebin(4);
-
- dNdEta[i]->SetLineWidth(2);
- dNdEta[i]->SetLineColor(colors[i]);
- dNdEta[i]->SetMarkerStyle(markers[i]);
- dNdEta[i]->SetMarkerSize(0.9);
- dNdEta[i]->SetMarkerColor(colors[i]);
-
- ratios[i] = (TH1F*)dNdEta[i]->Clone("ratio");
- ratios[i]->Divide(ratios[i],dNdEta[0],1,1,"B");
-
- ratios[i]->SetName(changes[i]);
- ratios[i]->SetTitle(changes[i]);
}
-
- //##########################################################
-
- gStyle->SetOptStat(0);
- gStyle->SetOptTitle(0);
- gStyle->SetOptFit(0);
-
- gStyle->SetTextSize(0.2);
- gStyle->SetTitleSize(0.05,"xyz");
- gStyle->SetLabelOffset(0.01, "xyz");
-
-
- gStyle->SetTitleOffset(1.2, "y");
- gStyle->SetTitleOffset(1.2, "x");
- gStyle->SetEndErrorSize(0.0);
-
- //##############################################
-
- //making canvas and pads
- TCanvas *c = new TCanvas(Form("vertex_reco_syst_%s", plot), Form("vertex_reco_syst_%s", plot),600,500);
- TPad* p1 = new TPad("pad1","", 0, 0.0, 1.0, 1.0, 0, 0, 0);
-
- p1->SetBottomMargin(0.15);
- p1->SetTopMargin(0.03);
- p1->SetLeftMargin(0.15);
- p1->SetRightMargin(0.03);
-
- p1->SetGridx();
- p1->SetGridy();
-
- p1->Draw();
- p1->cd();
-
-
- TH2F* null = new TH2F("","",100,-1.5,1.5,100,0.9601,1.099);
- null->SetXTitle("#eta");
- null->GetXaxis()->CenterTitle(kTRUE);
- null->SetYTitle("dN/d#eta / dN/d#eta_{pythia}");
- null->GetYaxis()->CenterTitle(kTRUE);
- null->Draw();
+ TFile* fout = new TFile(outputFileName,"RECREATE");
- for (Int_t i=1; i<7; i++)
- ratios[i]->Draw("same");
+ // to make everything consistent
+ hRatios[0]->Divide(hRatios[0],hRatios[0],1,1);
- TLegend* legend = new TLegend(0.6, 0.6, 0.95, 0.95);
- legend->SetFillColor(0);
-
- TLatex* text[7];
- for (Int_t i=1; i<7; i++) {
- legend->AddEntry(dNdEta[i], descr[i]);
+ for (Int_t i=0; i<nChanges * 3; i++)
+ {
+ corrections[i]->SaveHistograms();
+ hRatios[i]->Write();
}
- legend->Draw();
- //text(0.2,0.88,"Effect of changing",0.045,1,kTRUE);
- //text(0.2,0.83,"relative cross-sections",0.045,1,kTRUE);
- //text(0.2,0.78,"(vertex reconstruction corr.)",0.043,13,kTRUE);
-
- c->SaveAs(Form("%s.gif", c->GetName()));
- c->SaveAs(Form("%s.eps", c->GetName()));
+ fout->Write();
+ fout->Close();
}
-
-
DrawTriggerEfficiency(Char_t* fileName) {
gStyle->SetOptStat(0);
generatedPt[0]->Draw("same");
}
-void changePtSpectrum(const char* fileName = "analysis_mc.root")
+void changePtSpectrum(const char* fileName = "analysis_mc.root", Float_t ptCutOff = 0.2, const char* fileName2 = 0)
{
+ Float_t factor = 0.5;
+
TFile* file = TFile::Open(fileName);
- TH1F* hist = dynamic_cast<TH1F*> (file->Get("dndeta_check_pt"));
+ TH1F* hist = dynamic_cast<TH1F*> (file->Get("dndeta_check_pt")->Clone());
+
+ TH1* hist2 = 0;
+ if (fileName2)
+ {
+ file2 = TFile::Open(fileName2);
+ hist2 = dynamic_cast<TH1*> (file2->Get("dndeta_check_pt")->Clone());
+ hist2->Scale(hist->GetBinContent(hist->FindBin(ptCutOff)) / hist2->GetBinContent(hist2->FindBin(ptCutOff)));
+ }
+
+ //hist->Scale(1.0 / hist->Integral());
//hist->Rebin(3);
//hist->Scale(1.0/3);
TH1F* scale1 = dynamic_cast<TH1F*> (hist->Clone("scale1"));
TH1F* scale2 = dynamic_cast<TH1F*> (hist->Clone("scale2"));
- Float_t ptCutOff = 0.3;
-
for (Int_t i=1; i <= hist->GetNbinsX(); ++i)
{
if (hist->GetBinCenter(i) > ptCutOff)
else
{
// 90 % at pt = 0, 0% at pt = ptcutoff
- scale1->SetBinContent(i, 1 - (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3);
+ scale1->SetBinContent(i, 1 - (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * factor);
// 110% at pt = 0, ...
- scale2->SetBinContent(i, 1 + (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * 0.3);
+ scale2->SetBinContent(i, 1 + (ptCutOff - hist->GetBinCenter(i)) / ptCutOff * factor);
}
scale1->SetBinError(i, 0);
scale2->SetBinError(i, 0);
}
+ /*
new TCanvas;
-
scale1->Draw();
scale2->SetLineColor(kRed);
scale2->Draw("SAME");
+ */
clone1->Multiply(scale1);
clone2->Multiply(scale2);
clone2->SetMarkerStyle(markers[0]);*/
hist->SetTitle(";p_{T} in GeV/c;dN_{ch}/dp_{T} in c/GeV");
- hist->GetXaxis()->SetRangeUser(0, 0.7);
+ hist->GetXaxis()->SetRangeUser(0, 0.5);
hist->GetYaxis()->SetRangeUser(0.01, clone2->GetMaximum() * 1.1);
hist->GetYaxis()->SetTitleOffset(1);
- TCanvas* canvas = new TCanvas;
+ TCanvas* canvas = new TCanvas("c", "c", 600, 600);
+ gPad->SetGridx();
+ gPad->SetGridy();
+ gPad->SetTopMargin(0.05);
+ gPad->SetRightMargin(0.05);
gPad->SetBottomMargin(0.12);
hist->Draw("H");
clone1->SetLineColor(kRed);
clone2->SetLineColor(kBlue);
clone2->Draw("HSAME");
hist->Draw("HSAME");
+
+ if (hist2)
+ {
+ Prepare1DPlot(hist2);
+ hist2->SetLineStyle(2);
+ hist2->Draw("HSAME");
+ }
Float_t fraction = hist->Integral(hist->GetXaxis()->FindBin(ptCutOff), hist->GetNbinsX()) / hist->Integral(1, hist->GetNbinsX());
Float_t fraction1 = clone1->Integral(clone1->GetXaxis()->FindBin(ptCutOff), clone1->GetNbinsX()) / clone1->Integral(1, clone1->GetNbinsX());
printf("%f %f %f\n", fraction, fraction1, fraction2);
printf("Rel. %f %f\n", fraction1 / fraction, fraction2 / fraction);
- canvas->SaveAs("changePtSpectrum.gif");
+ //canvas->SaveAs("changePtSpectrum.gif");
canvas->SaveAs("changePtSpectrum.eps");
}
legend->Draw();
}
+
+void ChangePtInCorrection(const char* fileName = "correction_map.root", const char* dirName = "dndeta_correction")
+{
+ loadlibs();
+ if (!TFile::Open(fileName))
+ return;
+
+ AlidNdEtaCorrection* dNdEtaCorrection = new AlidNdEtaCorrection(dirName, dirName);
+ if (!dNdEtaCorrection->LoadHistograms())
+ return;
+
+ dNdEtaCorrection->Finish();
+
+ AliCorrection* corr = dNdEtaCorrection->GetCorrection(AlidNdEtaCorrection::kTrack2Particle);
+
+ Printf(">>>> Before");
+ corr->PrintInfo(0);
+
+ Float_t factor = 0.5;
+ Float_t ptCutOff = 0.2;
+
+ TH3* gene = corr->GetTrackCorrection()->GetGeneratedHistogram();
+ TH3* meas = corr->GetTrackCorrection()->GetMeasuredHistogram();
+
+ for (Int_t z = 1; z <= gene->GetZaxis()->FindBin(ptCutOff - 0.01); z++)
+ {
+ Float_t localFactor = 1 - (ptCutOff - gene->GetZaxis()->GetBinCenter(z)) / ptCutOff * factor;
+ Printf("%f %f", gene->GetZaxis()->GetBinCenter(z), localFactor);
+ for (Int_t x = 1; x <= gene->GetNbinsX(); x++)
+ for (Int_t y = 1; y <= gene->GetNbinsY(); y++)
+ {
+ gene->SetBinContent(x, y, z, gene->GetBinContent(x, y, z) * localFactor);
+ meas->SetBinContent(x, y, z, meas->GetBinContent(x, y, z) * localFactor);
+ }
+ }
+
+ dNdEtaCorrection->Finish();
+
+ Printf(">>>> After");
+ corr->PrintInfo(0);
+}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff