1 //____________________________________________________________________
12 * @ingroup pwglf_forward_scripts_corr
15 MakeSecMap(TList* list, Double_t low, Double_t high,
16 AliFMDCorrSecondaryMap* m)
18 // --- Get the list ------------------------------------------------
19 TString lName(AliForwardMCCorrectionsTask::VtxBin::BinName(low, high));
20 TList* vl = static_cast<TList*>(list->FindObject(lName));
22 Error("MakeSecMap", "List %s not found in %s",
23 lName.Data(), list->GetName());
27 // --- Get primary distribution ------------------------------------
28 const char* primaryName = "primary";
29 TH2D* primary = static_cast<TH2D*>(vl->FindObject(primaryName));
31 Error("MakeSecMap", "Couldn't not find histogram %s in %s",
32 primaryName, lName.Data());
35 TH2D* primaryI = static_cast<TH2D*>(primary->Clone("primaryI"));
36 TH2D* primaryO = static_cast<TH2D*>(primary->Clone("primaryO"));
37 primaryI->SetDirectory(0);
38 primaryO->SetDirectory(0);
42 // --- Calculate vertex --------------------------------------------
43 Double_t vz = (high+low) / 2;
45 // --- Loop over rings ---------------------------------------------
46 for (UShort_t d = 1; d <= 3; d++) {
47 UShort_t nr = (d == 1 ? 1 : 2);
48 for (UShort_t q = 0; q < nr; q++) {
49 Char_t r = (q == 0 ? 'I' : 'O');
51 const char* ringName = Form("FMD%d%c_cache", d, r);
52 TH2D* ring = static_cast<TH2D*>(vl->FindObject(ringName));
54 Error("MakeSecMap", "Didn't find histogram %s in %s",
55 ringName, vl->GetName());
60 TH2D* tmp = static_cast<TH2D*>(ring->Clone("tmp"));
61 TString tmpName = tmp->GetName();
62 tmpName.ReplaceAll("_cache", "");
63 tmp->SetName(tmpName);
65 tmp->Divide(q == 0 ? primaryI : primaryO);
68 m->SetCorrection(d, r, vz, tmp);
74 //____________________________________________________________________
83 * @ingroup pwglf_forward_scripts_corr
86 MakeCorrSecMap(const char* filename,
91 // --- Load code ---------------------------------------------------
92 // gROOT->Macro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/LoadLibs.C");
95 // --- Get the file ------------------------------------------------
96 TFile* file = TFile::Open(filename, "READ");
98 Error("MakeCorrSecMap", "Couldn't open file %s", filename);
102 // --- Get the parent list -----------------------------------------
103 const char* forwardName = "ForwardSums";
104 TList* forward = static_cast<TList*>(file->Get(forwardName));
106 Error("MakeCorrSecMap", "Couldn't get list %s from %s",
107 forwardName, filename);
111 // --- Get the vertex axis -----------------------------------------
112 const char* vtxName = "vtxAxis";
113 TH1* vtxHist = static_cast<TH1*>(forward->FindObject(vtxName));
115 Error("MakeCorrSecMap", "Couldn't get histogram %s from %s",
116 vtxName, forwardName);
119 const TAxis& vtxAxis = *(vtxHist->GetXaxis());
121 // --- Get the vertex axis -----------------------------------------
122 const char* etaName = "etaAxis";
123 TH1* etaHist = static_cast<TH1*>(forward->FindObject(etaName));
125 Error("MakeCorrSecMap", "Couldn't get histogram %s from %s",
126 etaName, forwardName);
129 const TAxis& etaAxis = *(etaHist->GetXaxis());
131 // --- Fill correction object --------------------------------------
132 AliFMDCorrSecondaryMap* corr = new AliFMDCorrSecondaryMap;
133 corr->SetVertexAxis(vtxAxis);
134 corr->SetEtaAxis(etaAxis);
136 for (Int_t i = 1; i <= vtxAxis.GetNbins(); i++) {
137 Double_t low = vtxAxis.GetBinLowEdge(i);
138 Double_t high = vtxAxis.GetBinUpEdge(i);
140 MakeSecMap(forward, low, high, corr);
143 // --- Get output filename and open --------------------------------
144 UShort_t isys = AliForwardUtil::ParseCollisionSystem(sys);
145 AliForwardCorrectionManager& mgr = AliForwardCorrectionManager::Instance();
146 TString fname(mgr.GetFileName(AliForwardCorrectionManager::kSecondaryMap,
147 isys, cms, field, false));
148 TFile* output = TFile::Open(fname.Data(), "RECREATE");
150 Warning("Run", "Failed to open output file %s", fname.Data());
154 // --- Write to output ---------------------------------------------
155 corr->Write(mgr.GetObjectName(AliForwardCorrectionManager::kSecondaryMap));
158 Info("Run", "File %s created. It should be copied to %s and stored in SVN",
159 fname.Data(),mgr.GetFileDir(AliForwardCorrectionManager::kSecondaryMap));
162 //____________________________________________________________________