Fixes for reading zero-suppressed data. These should be propagated to
[u/mrichter/AliRoot.git] / FMD / scripts / DrawXsection.C
1 //____________________________________________________________________
2 //
3 // Script to draw a X-section, LOSS, or range made with MakeXsection
4 //
5 /** @ingroup xsec_script
6     @param scale 
7     @param filename 
8     @param var 
9     @param medName 
10     @param thick 
11     @param pdgName 
12 */
13 void
14 DrawXsection(Bool_t scale=kFALSE, 
15              const char* filename="xsec.root", 
16              const char* var="LOSS", 
17              const char* medName="FMD_Si$", 
18              Double_t thick=.03,
19              const char* pdgName="pi+")
20 {
21   TFile*   file = TFile::Open(filename, "READ");
22   TTree*   tree = static_cast<TTree*>(file->Get(Form("%s_%s",medName,
23                                                      pdgName)));
24   TLeaf* tb   = tree->GetLeaf("T");
25   TLeaf* vb   = tree->GetLeaf(var);
26   if (!vb) {
27     std::cerr << "Leaf " << var << " not found" << std::endl;
28     return;
29   }
30   Float_t tkine, value;
31   tb->SetAddress(&tkine);
32   vb->SetAddress(&value);
33   Int_t n = tree->GetEntries();
34
35   Float_t xscale = 1;
36   Float_t yscale = 1;
37   if (scale) {
38     TDatabasePDG* pdgDb = TDatabasePDG::Instance();
39     TParticlePDG* pdgP  = pdgDb->GetParticle(pdgName);
40     if (!pdgP) {
41       std::cerr << "Couldn't find particle " << pdgName << std::endl;
42       return;
43     }
44     Double_t m = pdgP->Mass();
45     Double_t q = pdgP->Charge() / 3;
46     if (m == 0 || q == 0) {
47       std::cerr  << "Mass is 0" << std::endl;
48       return;
49     }
50     xscale = 1 / m;
51     yscale = 1 / (q * q);
52   }
53   
54   TGraphErrors* graph = new TGraphErrors(n);
55   for (Int_t i = 0; i < n; i++) {
56     tree->GetEntry(i);
57     Double_t x = tkine*xscale;
58     Double_t y = value*yscale;
59     graph->SetPoint(i, x, y); 
60     // 5 sigma
61     graph->SetPointError(i, 0, 5 * .1 * y);
62   }
63   TCanvas* c = new TCanvas("c","c");
64   c->SetLogx();
65   c->SetLogy();
66   graph->SetLineWidth(2);
67   graph->SetFillStyle(3001);
68   graph->SetFillColor(6);
69   graph->Draw("L");
70   graph->DrawClone("AL3");
71   c->Modified();
72   c->Update();
73   c->cd();
74   c->SaveAs("xsec.C");
75   
76 }
77
78 //____________________________________________________________________
79 //
80 // EOF
81 //