]>
Commit | Line | Data |
---|---|---|
a1d7e11f | 1 | /** |
56199f2b | 2 | * A task to do a comparison between tracklets and clusers in the SPD |
a1d7e11f | 3 | * |
56199f2b | 4 | * Since the class SPDComparisonTask derives from a compiled class |
a1d7e11f | 5 | * (AliAnalysisTaskSE) we need to compile that code. The script will, |
6 | * when executed in the AliROOT prompt load it self again and byte | |
7 | * compile it with the preprocessor flag BUILD defined. \ | |
8 | * | |
9 | * The preprocessor define BUILD is not defined at first, when the | |
10 | * script is loaded using | |
11 | * | |
12 | * @verbatim | |
56199f2b | 13 | * Root> .x SPDComparison.C |
a1d7e11f | 14 | * @endverbatim |
15 | * | |
56199f2b | 16 | * which means that CINT will only see the function SPDComparison. |
a1d7e11f | 17 | * In that function, we define the BUILD preprocessor symbol |
18 | * | |
19 | * @code | |
20 | * gSystem->AddIncludePath("-DBUILD=1 ..."); | |
21 | * @endcode | |
22 | * | |
23 | * and then ACLic compile ourselves | |
24 | * | |
25 | * @code | |
26 | * gROOT->LoadMacro("./SPDComparison.C++"); | |
27 | * @endcode | |
28 | * | |
29 | * But since BUILD is now defined, it means that ACLic will only see | |
30 | * the class and not the function SPDComparison. | |
31 | * | |
32 | * This trick hinges on that when you initially load the script and | |
33 | * when it is done inside the script it is done using two distinct | |
34 | * paths - otherwise ROOT will try to unload the script first, and | |
35 | * that fails. | |
36 | * | |
37 | */ | |
38 | ||
39 | #ifdef BUILD | |
40 | #include <AliAnalysisTaskSE.h> | |
41 | #include <AliESDEvent.h> | |
42 | #include <AliLog.h> | |
43 | #include <AliMultiplicity.h> | |
44 | #include "AliFMDEventInspector.h" | |
45 | #include "AliAODForwardMult.h" | |
46 | #include <TH1D.h> | |
47 | #include <TH2D.h> | |
48 | #include <TGraphErrors.h> | |
49 | #include <TList.h> | |
50 | #include <TObjArray.h> | |
51 | #include <TMath.h> | |
52 | ||
53 | ||
54 | /** | |
55 | * @class SPDComparisonTask | |
56 | * | |
57 | * Task to do forward/backward correlations using FMD and SPD | |
58 | * (cluster) data. | |
59 | * | |
60 | * The class contains the sub-structure Bin that represents an | |
61 | * @f$\eta@f$ range. One can add (possibly overlapping) @f$\eta@f$ | |
62 | * ranges by calling the member function AddBin | |
63 | * | |
bd6f5206 | 64 | * @ingroup pwglf_forward_scripts |
a1d7e11f | 65 | */ |
66 | class SPDComparisonTask : public AliAnalysisTaskSE | |
67 | { | |
68 | public: | |
69 | /** | |
70 | * A single vertex bin - contains summed histogram of clusters and | |
71 | * tracklets | |
72 | * | |
73 | */ | |
74 | struct VtxBin : public TNamed | |
75 | { | |
76 | /** | |
77 | * Static function to get a bin name | |
78 | * | |
79 | * @param low Low limit | |
80 | * @param high High limit | |
81 | * | |
82 | * @return Pointer to static array | |
83 | */ | |
84 | static const char* BinName(Double_t low, Double_t high) | |
85 | { | |
86 | static TString n; | |
87 | n = (Form("vtx%+05.1f_%+05.1f", low, high)); | |
88 | n.ReplaceAll("+", "p"); | |
89 | n.ReplaceAll("-", "m"); | |
90 | n.ReplaceAll(".", "d"); | |
91 | return n.Data(); | |
92 | } | |
93 | /** | |
94 | * Constructor | |
95 | * | |
96 | * @param low Low limit | |
97 | * @param high High limit | |
98 | */ | |
99 | VtxBin(Double_t low, Double_t high) | |
100 | : TNamed(BinName(low,high), ""), | |
101 | fClusters(0), | |
102 | fTracklets(0) | |
103 | { | |
104 | } | |
105 | /** | |
106 | * Constructor | |
107 | */ | |
108 | VtxBin() : TNamed(), fClusters(0), fTracklets(0) {} | |
109 | /** | |
110 | * Copy constructor | |
111 | * | |
112 | * @param o Object to copy from | |
113 | */ | |
114 | VtxBin(const VtxBin& o) : | |
115 | TNamed(o), fClusters(o.fClusters), fTracklets(o.fTracklets) | |
116 | {} | |
117 | /** | |
118 | * Assignment operator | |
119 | * | |
120 | * @param o Object to assign from | |
121 | * | |
122 | * @return Reference to this. | |
123 | */ | |
124 | VtxBin& operator=(const VtxBin& o) | |
125 | { | |
126 | TNamed::operator=(o); | |
127 | fClusters = o.fClusters; | |
128 | fTracklets = o.fTracklets; | |
129 | return *this; | |
130 | } | |
131 | /** | |
132 | * Define outputs | |
133 | * | |
134 | * @param l Output list | |
135 | * @param eta Eta axis to use | |
136 | */ | |
137 | void DefineOutput(TList* l, const TAxis& eta) | |
138 | { | |
139 | TList* ll = new TList; | |
140 | ll->SetName(GetName()); | |
141 | ll->SetOwner(); | |
142 | l->Add(ll); | |
143 | ||
144 | fClusters = new TH1D("clusters", "Clusters", | |
145 | eta.GetNbins(), | |
146 | eta.GetXmin(), | |
147 | eta.GetXmax()); | |
148 | fClusters->SetXTitle("#eta"); | |
149 | fClusters->SetYTitle("dN/d#eta"); | |
150 | fClusters->SetDirectory(0); | |
151 | fClusters->Sumw2(); | |
152 | ll->Add(fClusters); | |
153 | ||
154 | fTracklets = static_cast<TH1D*>(fClusters->Clone("tracklets")); | |
155 | fTracklets->SetTitle("Tracklets"); | |
156 | fTracklets->SetDirectory(0); | |
157 | ll->Add(fTracklets); | |
158 | } | |
159 | /** | |
160 | * Process the input | |
161 | * | |
162 | * @param spdmult Input | |
163 | */ | |
164 | void Process(const AliMultiplicity* spdmult) | |
165 | { | |
166 | //Filling clusters in layer 1 used for tracklets... | |
167 | for(Int_t j = 0; j< spdmult->GetNumberOfTracklets();j++) { | |
168 | Double_t eta = spdmult->GetEta(j); | |
169 | fClusters->Fill(eta); | |
170 | fTracklets->Fill(eta); | |
171 | } | |
172 | ||
173 | //...and then the unused ones in layer 1 | |
174 | for(Int_t j = 0; j< spdmult->GetNumberOfSingleClusters();j++) { | |
175 | Double_t eta = -TMath::Log(TMath::Tan(spdmult->GetThetaSingle(j)/2.)); | |
176 | fClusters->Fill(eta); | |
177 | } | |
178 | } | |
179 | /** | |
180 | * End of job processing | |
181 | * | |
182 | * @param l Output list | |
183 | * @param nEvents Number of events to normalise to | |
184 | */ | |
185 | void Finish(TList* l, Int_t nEvents) | |
186 | { | |
187 | TList* ll = static_cast<TList*>(l->FindObject(GetName())); | |
188 | fClusters = static_cast<TH1D*>(ll->FindObject("clusters")); | |
189 | fTracklets = static_cast<TH1D*>(ll->FindObject("tracklets")); | |
190 | fClusters->Scale(1. / nEvents, "width"); | |
191 | fTracklets->Scale(1. / nEvents, "width"); | |
192 | TH1* ratio = static_cast<TH1D*>(fClusters->Clone("ratio")); | |
193 | ratio->SetDirectory(0); | |
194 | ratio->SetTitle("Clusters/Tracklets"); | |
195 | ratio->Divide(fTracklets); | |
196 | ll->Add(ratio); | |
197 | } | |
198 | TH1D* fClusters; // Cluster histogram | |
199 | TH1D* fTracklets; // Tracklet histogram | |
200 | }; | |
201 | ||
202 | //__________________________________________________________________ | |
203 | /** | |
204 | * Default constructor - do not use | |
205 | */ | |
206 | SPDComparisonTask() | |
207 | : AliAnalysisTaskSE(), | |
208 | fInspector(), | |
209 | fBins(0), | |
210 | fEvents(0), | |
211 | fList(0), | |
212 | fFirstEvent(true), | |
213 | fVertexAxis(20, -20, 20) | |
214 | {} | |
215 | /** | |
216 | * Constructor with a single argument. | |
217 | */ | |
218 | SPDComparisonTask(const char*) | |
219 | : AliAnalysisTaskSE("spdComparision"), | |
220 | fInspector("eventInspector"), | |
221 | fBins(0), | |
222 | fEvents(0), | |
223 | fList(0), | |
224 | fFirstEvent(true), | |
225 | fVertexAxis(20, -20, 20) | |
226 | { | |
227 | // Declare our output container | |
228 | DefineOutput(1,TList::Class()); | |
229 | } | |
230 | /** | |
231 | * Copy constructor | |
232 | * | |
233 | * @param o Object to copy from | |
234 | */ | |
235 | SPDComparisonTask(const SPDComparisonTask& o) | |
236 | : AliAnalysisTaskSE(o), | |
237 | fInspector(o.fInspector), | |
238 | fBins(o.fBins), | |
239 | fEvents(o.fEvents), | |
240 | fList(o.fList), | |
241 | fFirstEvent(o.fFirstEvent), | |
242 | fVertexAxis(20, -20, 20) | |
243 | { | |
244 | SetVertexAxis(o.fVertexAxis); | |
245 | } | |
246 | /** | |
247 | * Assigment operator | |
248 | * | |
249 | * @param o Object to assign from | |
250 | * | |
251 | * @return Reference to this | |
252 | */ | |
253 | SPDComparisonTask& operator=(const SPDComparisonTask& o) | |
254 | { | |
255 | AliAnalysisTaskSE::operator=(o); | |
256 | fInspector = o.fInspector; | |
257 | fEvents = o.fEvents; | |
258 | fList = o.fList; | |
259 | fFirstEvent = o.fFirstEvent; | |
260 | SetVertexAxis(o.fVertexAxis); | |
261 | return *this; | |
262 | } | |
263 | /** | |
264 | * Destructor | |
265 | */ | |
266 | ~SPDComparisonTask() {} | |
267 | /** | |
268 | * Set the vertex axis to use | |
269 | * | |
270 | * @param a Axis to set from | |
271 | */ | |
272 | void SetVertexAxis(const TAxis& a) | |
273 | { | |
274 | SetVertexAxis(a.GetNbins(), | |
275 | a.GetXmin(), | |
276 | a.GetXmax()); | |
277 | } | |
278 | /** | |
279 | * Set the vertex axis to use | |
280 | * | |
281 | * @param n Number of bins | |
282 | * @param xmin Least @f$v_z@f$ | |
283 | * @param xmax Most @f$v_z@f$ | |
284 | */ | |
285 | void SetVertexAxis(Int_t n, Double_t xmin, Double_t xmax) | |
286 | { | |
287 | fVertexAxis.Set(n, xmin, xmax); | |
288 | } | |
289 | /** | |
290 | * Create output objects | |
291 | * | |
292 | */ | |
293 | void UserCreateOutputObjects() | |
294 | { | |
295 | fList = new TList; | |
296 | fList->SetName(GetName()); | |
297 | fList->SetOwner(); | |
298 | ||
299 | ||
300 | fEvents = new TH1D("events", "Events", | |
301 | fVertexAxis.GetNbins(), | |
302 | fVertexAxis.GetXmin(), | |
303 | fVertexAxis.GetXmax()); | |
304 | fEvents->SetDirectory(0); | |
305 | fEvents->SetXTitle("v_{z} [cm]"); | |
306 | fList->Add(fEvents); | |
307 | ||
308 | TAxis& vtxAxis = fVertexAxis; | |
309 | fBins = new TObjArray(vtxAxis.GetNbins(),1); | |
310 | ||
311 | TAxis etaAxis(120, -3, 3); | |
312 | for (Int_t i = 1; i <= vtxAxis.GetNbins(); i++) { | |
313 | VtxBin* bin = new VtxBin(vtxAxis.GetBinLowEdge(i), | |
314 | vtxAxis.GetBinUpEdge(i)); | |
315 | bin->DefineOutput(fList, etaAxis); | |
316 | fBins->AddAt(bin,i); | |
317 | } | |
318 | ||
319 | ||
320 | fInspector.DefineOutput(fList); | |
321 | fInspector.Init(*(fEvents->GetXaxis())); | |
322 | ||
323 | PostData(1, fList); | |
324 | } | |
325 | /** | |
326 | * Process one event | |
327 | * | |
328 | * @param option Not used | |
329 | */ | |
330 | void UserExec(Option_t* option="") | |
331 | { | |
332 | // Get the input data - ESD event | |
333 | AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent()); | |
334 | if (!esd) { | |
335 | AliWarning("No ESD event found for input event"); | |
336 | return; | |
337 | } | |
338 | ||
339 | //--- Read run information ----------------------------------------- | |
340 | if (fFirstEvent && esd->GetESDRun()) { | |
341 | fInspector.ReadRunDetails(esd); | |
342 | ||
343 | AliInfo(Form("Initializing with parameters from the ESD:\n" | |
344 | " AliESDEvent::GetBeamEnergy() ->%f\n" | |
345 | " AliESDEvent::GetBeamType() ->%s\n" | |
346 | " AliESDEvent::GetCurrentL3() ->%f\n" | |
347 | " AliESDEvent::GetMagneticField()->%f\n" | |
348 | " AliESDEvent::GetRunNumber() ->%d\n", | |
349 | esd->GetBeamEnergy(), | |
350 | esd->GetBeamType(), | |
351 | esd->GetCurrentL3(), | |
352 | esd->GetMagneticField(), | |
353 | esd->GetRunNumber())); | |
354 | ||
355 | Print(); | |
356 | fFirstEvent = false; | |
357 | } | |
358 | ||
359 | // Some variables | |
360 | UInt_t triggers; // Trigger bits | |
361 | Bool_t lowFlux; // Low flux flag | |
362 | UShort_t iVz; // Vertex bin from ESD | |
363 | Double_t vZ; // Z coordinate from ESD | |
364 | Double_t cent; // Centrality | |
365 | UShort_t nClusters;// Number of SPD clusters | |
366 | // Process the data | |
367 | UInt_t retESD = fInspector.Process(esd, triggers, lowFlux, iVz, vZ, | |
368 | cent, nClusters); | |
369 | Bool_t isInel = triggers & AliAODForwardMult::kB; | |
370 | Bool_t hasVtx = retESD == AliFMDEventInspector::kOk; | |
371 | ||
372 | if (!isInel || !hasVtx) return; | |
373 | fEvents->Fill(vZ); | |
374 | ||
375 | const AliMultiplicity* spdmult = esd->GetMultiplicity(); | |
376 | ||
377 | VtxBin* bin = static_cast<VtxBin*>(fBins->At(iVz)); | |
378 | if (!bin) { | |
379 | AliError(Form("No bin @ %d (%fcm)", iVz, vZ)); | |
380 | return; | |
381 | } | |
382 | bin->Process(spdmult); | |
383 | ||
384 | // Post data to output container | |
385 | PostData(1,fList); | |
386 | } | |
387 | /** | |
388 | * Called at the end of the processing | |
389 | * | |
390 | * @param option | |
391 | */ | |
392 | void Terminate(Option_t* option="") | |
393 | { | |
394 | fList = dynamic_cast<TList*>(GetOutputData(1)); | |
395 | if (!fList) { | |
396 | AliError("No output list defined"); | |
397 | return; | |
398 | } | |
399 | ||
400 | fEvents = static_cast<TH1D*>(fList->FindObject("events")); | |
401 | ||
402 | Int_t nEvents = fEvents->GetEntries(); | |
403 | AliInfo(Form("Got a total of %d events", nEvents)); | |
404 | ||
405 | TH1* clusters = 0; | |
406 | TH1* tracklets = 0; | |
407 | TIter next(fBins); | |
408 | VtxBin* bin = 0; | |
409 | Int_t i = 1; | |
410 | while ((bin = static_cast<VtxBin*>(next()))) { | |
411 | bin->Finish(fList, fEvents->GetBinContent(i++)); | |
412 | if (!clusters) clusters = static_cast<TH1D*>(bin->fClusters->Clone()); | |
413 | else clusters->Add(bin->fClusters); | |
414 | if (!tracklets) tracklets = static_cast<TH1D*>(bin->fTracklets->Clone()); | |
415 | else tracklets->Add(bin->fTracklets); | |
416 | } | |
417 | clusters->SetDirectory(0); | |
418 | tracklets->SetDirectory(0); | |
419 | clusters->Scale(1. / i); | |
420 | tracklets->Scale(1. / i); | |
421 | ||
422 | ||
423 | TH1D* ratio = static_cast<TH1D*>(clusters->Clone("ratio")); | |
424 | ratio->SetDirectory(0); | |
425 | ratio->SetTitle("Clusters/Tracklets"); | |
426 | ratio->Divide(tracklets); | |
427 | ||
428 | fList->Add(clusters); | |
429 | fList->Add(tracklets); | |
430 | fList->Add(ratio); | |
431 | } | |
432 | protected: | |
433 | AliFMDEventInspector fInspector; // Inspector | |
434 | TObjArray* fBins; // Vertex bins | |
435 | TH1D* fEvents; // Events | |
436 | TList* fList; // Output list | |
437 | Bool_t fFirstEvent;// First event flag | |
438 | TAxis fVertexAxis; | |
439 | ||
440 | ClassDef(SPDComparisonTask,1); // BF analysis | |
441 | }; | |
442 | #else | |
443 | ||
444 | //==================================================================== | |
445 | /** | |
446 | * Run the analysis | |
447 | * | |
c6115ede | 448 | * @param esddir Input directory |
a1d7e11f | 449 | * @param nEvents Number of events, negative means all |
450 | */ | |
451 | void | |
452 | SPDComparison(const char* esddir, Int_t nEvents=-1) | |
453 | { | |
454 | // --- Libraries to load ------------------------------------------- | |
bd6f5206 | 455 | gROOT->Macro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/LoadLibs.C"); |
a1d7e11f | 456 | |
457 | // --- Our data chain ---------------------------------------------- | |
bd6f5206 | 458 | gROOT->LoadMacro("$ALICE_ROOT/PWGLF/FORWARD/analysis2/scripts/MakeChain.C"); |
a1d7e11f | 459 | TChain* chain = MakeChain("ESD", esddir, true); |
460 | // If 0 or less events is select, choose all | |
461 | if (nEvents <= 0) nEvents = chain->GetEntries(); | |
462 | ||
463 | // --- Manager ----------------------------------------------------- | |
464 | AliAnalysisManager* mgr = new AliAnalysisManager("FB", "FB train"); | |
465 | AliAnalysisManager::SetCommonFileName("spd_comps.root"); | |
466 | ||
467 | // --- AOD input handler ------------------------------------------- | |
468 | AliESDInputHandler *inputHandler = new AliESDInputHandler(); | |
469 | mgr->SetInputEventHandler(inputHandler); | |
470 | ||
471 | // --- compile our code -------------------------------------------- | |
bd6f5206 | 472 | gSystem->AddIncludePath("-I${ALICE_ROOT}/PWGLF/FORWARD/analysis2 " |
a1d7e11f | 473 | "-I${ALICE_ROOT}/ANALYSIS " |
474 | "-I${ALICE_ROOT}/include -DBUILD=1"); | |
475 | gROOT->LoadMacro("./SPDComparison.C++g"); | |
476 | ||
477 | // --- Make our object --------------------------------------------- | |
478 | SPDComparisonTask* task = new SPDComparisonTask("SPD_COMP"); | |
479 | task->SetVertexAxis(10, -10, 10); | |
480 | mgr->AddTask(task); | |
481 | ||
482 | // --- create containers for input/output -------------------------- | |
483 | AliAnalysisDataContainer *sums = | |
484 | mgr->CreateContainer("spdComp", TList::Class(), | |
485 | AliAnalysisManager::kOutputContainer, | |
486 | AliAnalysisManager::GetCommonFileName()); | |
487 | // --- connect input/output ---------------------------------------- | |
488 | mgr->ConnectInput(task, 0, mgr->GetCommonInputContainer()); | |
489 | mgr->ConnectOutput(task, 1, sums); | |
490 | ||
491 | // --- Run the analysis -------------------------------------------- | |
492 | TStopwatch t; | |
493 | if (!mgr->InitAnalysis()) { | |
494 | Error("SPDComparison", "Failed to initialize analysis train!"); | |
495 | return; | |
496 | } | |
497 | // Some informative output | |
498 | mgr->PrintStatus(); | |
499 | // mgr->SetDebugLevel(3); | |
500 | if (mgr->GetDebugLevel() < 1) mgr->SetUseProgressBar(kTRUE,100); | |
501 | ||
502 | // Run the train | |
503 | t.Start(); | |
504 | Printf("=== RUNNING ANALYSIS on %9d events ==================", nEvents); | |
505 | mgr->StartAnalysis("local", chain, nEvents); | |
506 | t.Stop(); | |
507 | t.Print(); | |
508 | ||
509 | DrawSPDComparison(); | |
510 | } | |
511 | ||
512 | //==================================================================== | |
513 | /** | |
514 | * Draw results | |
515 | * | |
516 | * @param filename | |
517 | */ | |
518 | void | |
519 | DrawSPDComparison(const char* filename="spd_comps.root") | |
520 | { | |
521 | // --- Open the file ----------------------------------------------- | |
522 | TFile* file = TFile::Open(filename, "READ"); | |
523 | if (!file) { | |
524 | Error("DrawSPDComparison", "Failed to open file %s", filename); | |
525 | return; | |
526 | } | |
527 | ||
528 | // --- Get our list ------------------------------------------------ | |
529 | TList* spd = static_cast<TList*>(file->Get("spdComp")); | |
530 | if (!spd) { | |
531 | Error("DrawSPDComparison", "Failed to get list SPD_COMP from file %s", | |
532 | filename); | |
533 | return; | |
534 | } | |
535 | ||
536 | // --- Loop over list and get correlation plots -------------------- | |
537 | TH1* clusters = static_cast<TH1*>(spd->FindObject("clusters")); | |
538 | TH1* tracklets = static_cast<TH1*>(spd->FindObject("tracklets")); | |
539 | TH1* ratio = static_cast<TH1*>(spd->FindObject("ratio")); | |
540 | TH1* events = static_cast<TH1*>(spd->FindObject("events")); | |
541 | ||
542 | // --- Set style parameters ---------------------------------------- | |
543 | gStyle->SetPalette(1); | |
544 | gStyle->SetTitleFillColor(0); | |
545 | gStyle->SetTitleStyle(0); | |
546 | gStyle->SetTitleBorderSize(0); | |
547 | gStyle->SetOptStat(0); | |
548 | gStyle->SetTitleX(0.69); | |
549 | gStyle->SetTitleY(0.99); | |
550 | gStyle->SetTitleW(0.30); | |
551 | gStyle->SetTitleH(0.10); | |
552 | ||
553 | // --- Make canvas for correlation plots and plot them ------------- | |
554 | TCanvas* c1 = new TCanvas("c1", "c1", 900, 600); | |
555 | c1->SetTopMargin(0.05); | |
556 | c1->SetRightMargin(0.05); | |
557 | c1->SetFillColor(0); | |
558 | c1->SetFillStyle(0); | |
559 | c1->SetBorderSize(0); | |
560 | c1->SetBorderMode(0); | |
561 | ||
562 | c1->cd(); | |
563 | TPad* p1 = new TPad("p1","p1", 0.0, 0.3, 1.0, 1.0, 0, 0, 0); | |
564 | p1->SetFillColor(0); | |
565 | p1->SetFillStyle(0); | |
566 | p1->SetBorderSize(0); | |
567 | p1->SetBorderMode(0); | |
568 | p1->SetTopMargin(0.05); | |
569 | p1->SetBottomMargin(0.001); | |
570 | p1->SetRightMargin(0.05); | |
571 | p1->Draw(); | |
572 | p1->cd(); | |
573 | clusters->SetMarkerStyle(20); | |
574 | clusters->SetMarkerColor(kRed+1); | |
575 | clusters->Draw(); | |
576 | ||
577 | tracklets->SetMarkerStyle(20); | |
578 | tracklets->SetMarkerColor(kBlue+1); | |
579 | tracklets->Draw("same"); | |
580 | ||
581 | TString v(Form("%+5.1f<|v_{z}|<%+5.1f", | |
582 | events->GetXaxis()->GetXmin(), | |
583 | events->GetXaxis()->GetXmax())); | |
584 | TLatex* ltx = new TLatex(.2,.80, v.Data()); | |
585 | ltx->Draw(); | |
586 | ||
587 | TLegend* l = p1->BuildLegend(.6,.6,.94,.94); | |
588 | l->SetFillColor(0); | |
589 | l->SetFillStyle(0); | |
590 | l->SetBorderSize(0); | |
591 | ||
592 | c1->cd(); | |
593 | TPad* p2 = new TPad("p2","p2", 0.0, 0.0, 1.0, 0.3, 0, 0, 0); | |
594 | p2->SetFillColor(0); | |
595 | p2->SetFillStyle(0); | |
596 | p2->SetBorderSize(0); | |
597 | p2->SetBorderMode(0); | |
598 | p2->SetTopMargin(0.001); | |
599 | p2->SetBottomMargin(0.2); | |
600 | p2->SetRightMargin(0.05); | |
601 | p2->Draw(); | |
602 | p2->cd(); | |
603 | ratio->SetMarkerStyle(20); | |
604 | ratio->SetMarkerColor(kGray+1); | |
605 | ratio->Draw(); | |
606 | ratio->GetYaxis()->SetRangeUser(0,2); | |
607 | ||
608 | c1->SaveAs("SPDComparison.png"); | |
609 | c1->SaveAs("SPDComparison.root"); | |
610 | } | |
611 | ||
612 | ||
613 | ||
614 | ||
615 | ||
616 | #endif | |
617 | // | |
618 | // EOF | |
619 | // | |
620 |