end-of-line normalization
[u/mrichter/AliRoot.git] / ANALYSIS / AliHEPDataParser.cxx
1 //-------------------------------------------------------------------------
2 //           Implementation of Class AliHEPDataParser
3 //
4 //  This class is used to save the content of hisograms/graphs in the
5 //  HEP data format and viceversa. The HEP data format is not strictly
6 //  defined and there are many variants, the class only support a few
7 //  of them. More will be added as needed.  The input can be a set of
8 //  2 TH1, TGraphAsymmErrors or TGraphErrors (one for the stat and one
9 //  for the syst error). If the second one is a null pointer, only the
10 //  stat error is printed. The class can also import hepdata ascii
11 //  file (very preliminary)
12 // 
13 //  Author: Michele Floris, CERN
14 //-------------------------------------------------------------------------
15
16
17 #include "AliHEPDataParser.h"
18 #include "AliLog.h"
19 #include "TGraphAsymmErrors.h"
20 #include "TGraph.h"
21 #include "TGraphErrors.h"
22 #include "TH1.h"
23 #include "TObjArray.h"
24 #include "TObjString.h"
25 #include <fstream>
26 #include <iostream>
27 using namespace std;
28
29 ClassImp(AliHEPDataParser)
30
31 AliHEPDataParser::AliHEPDataParser() : TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("")
32 {
33   // default ctor
34
35 }
36
37 AliHEPDataParser::AliHEPDataParser(TH1 * hStat, TH1 * hSyst): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("y")
38 {
39   //ctor
40   fHistStat = hStat;
41   fHistSyst = hSyst;
42   fHEPDataFileLines = new TObjArray;
43
44 }
45 AliHEPDataParser::AliHEPDataParser(TGraph * grStat, TGraph * grSyst): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("")
46 {
47   // ctor
48   fGraphStat = grStat;
49   fGraphSyst = grSyst;
50   fHEPDataFileLines = new TObjArray;
51 }
52
53 AliHEPDataParser::AliHEPDataParser(const char * hepfileName): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("y")
54 {
55   //ctor
56   // Put results in graphs
57   fGraphSyst = new TGraphAsymmErrors();
58   fGraphStat = new TGraphAsymmErrors();
59   ifstream infile;
60   infile.open(hepfileName);
61   TString line;
62   Int_t ipoints = 0;
63   while (line.ReadLine(infile)) {
64     TObjArray * tokens = line.Tokenize(" ");
65     if(tokens->GetEntries() < 1) {
66       delete tokens;
67       AliError("not enough columns");
68       return;      
69     }
70     // TODO: Assumes the format
71     // binmin binmax y +-stat +-syst. Try to make it smarter...
72     TObjString * binMin = (TObjString*) tokens->At(0);
73     TObjString * binMax = (TObjString*) tokens->At(1);
74     TObjString * value  = (TObjString*) tokens->At(2);
75     TObjString * stat   = (TObjString*) tokens->At(3);
76     TObjString * syst   = (TObjString*) tokens->At(4);
77     stat->String().ReplaceAll("+-","");
78     if(syst) syst->String().ReplaceAll("+-","");
79     if (!binMin->String().Atof()) {delete tokens; continue;} // skip headers
80     Float_t binCenter = (binMax->String().Atof() + binMin->String().Atof())/2;
81     Float_t binWidth  =  (binMax->String().Atof() - binMin->String().Atof())/2;
82     cout << line.Data() << endl;//<< " " << binMin->String().Atof() <<" " <<  binCenter << " " << binWidth << endl;
83
84
85     fGraphStat->SetPoint(ipoints, binCenter, value->String().Atof());
86     fGraphSyst->SetPoint(ipoints, binCenter, value->String().Atof());
87     ((TGraphAsymmErrors*)fGraphStat)->SetPointError(ipoints, 
88                                                     binWidth,
89                                                     binWidth,
90                                                     stat->String().Atof(),
91                                                     stat->String().Atof());
92     if(syst) ((TGraphAsymmErrors*)fGraphSyst)->SetPointError(ipoints, 
93                                                              binWidth,
94                                                              binWidth,
95                                                              syst->String().Atof(),
96                                                              syst->String().Atof());
97     ipoints++;
98     delete tokens;
99   }
100   infile.close();
101     
102
103 }
104
105 AliHEPDataParser::~AliHEPDataParser(){
106   // dtor
107   if(fHistStat) delete fHistStat;
108   if(fHistSyst) delete fHistSyst;
109   if(fGraphStat) delete fGraphStat;
110   if(fGraphSyst) delete fGraphSyst;
111   if(fHEPDataFileLines) delete fHEPDataFileLines;
112 }
113   
114 void AliHEPDataParser::SaveHEPDataFile(const char * hepfileName, Bool_t trueUseGraphFalesUseHisto) {
115   // Fills fHEPDataFileLines and saves its content to a file
116   if(!fHEPDataFileLines)   fHEPDataFileLines = new TObjArray;
117   if(trueUseGraphFalesUseHisto) {
118     if(!fGraphStat) {
119       AliError("Graph not set");
120       return;
121     }
122     Bool_t asym = kFALSE; // check if this has asymmetric errors
123     if       (!strcmp(fGraphStat->ClassName(), "TGraphErrors")) asym = kFALSE;
124     else if  (!strcmp(fGraphStat->ClassName(), "TGraphAsymmErrors")) asym = kTRUE;
125     else     {AliError("Unsupported graph type"); return;}
126     Int_t npoint = fGraphStat->GetN();
127     if(asym) AliInfo("Assymmetric errors");
128     for(Int_t ipoint = 0; ipoint < npoint; ipoint++){            
129       if(ipoint == 0) {
130         if(fGraphSyst) {
131           if (asym)
132             fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +stat -stat +syst -syst", fValueName.Data()))); 
133           else 
134             fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +-stat +-syst", fValueName.Data()))); 
135         }
136         else {
137           if(asym)
138             fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +stat -stat", fValueName.Data()))); 
139           else 
140             fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +-stat", fValueName.Data()))); 
141         }
142       }
143       // Skip empty bins
144       if(!fGraphStat->GetY()[ipoint]) continue;
145       TObjString * line = new TObjString;      
146       if(fGraphSyst) {
147         if (asym)
148           line->String().Form("%f %f +%f -%f +%f -%f", 
149                               fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
150                               ((TGraphAsymmErrors*)fGraphStat)->GetEYhigh()[ipoint], 
151                               ((TGraphAsymmErrors*)fGraphStat)->GetEYlow()[ipoint], 
152                               ((TGraphAsymmErrors*)fGraphSyst)->GetEYhigh()[ipoint], 
153                               ((TGraphAsymmErrors*)fGraphSyst)->GetEYlow()[ipoint]);
154         else 
155           line->String().Form("%f %f +-%f +-%f", 
156                               fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
157                               ((TGraphErrors*)fGraphStat)->GetEY()[ipoint], 
158                               ((TGraphErrors*)fGraphSyst)->GetEY()[ipoint]);
159
160         fHEPDataFileLines->Add(line);
161       } else {
162         if (asym)
163           line->String().Form("%f %f +%f -%f", 
164                               fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
165                               ((TGraphAsymmErrors*)fGraphStat)->GetEYhigh()[ipoint], ((TGraphAsymmErrors*)fGraphStat)->GetEYlow()[ipoint]);
166         else { 
167           line->String().Form("%f %f +-%f", 
168                               fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
169                               ((TGraphErrors*)fGraphStat)->GetEY()[ipoint]);
170         }
171
172         fHEPDataFileLines->Add(line);
173       }
174     }    
175   }
176   else {
177     if(!fHistStat) {
178       AliError("Hist not set");
179       return;
180     }
181     Int_t nbin = fHistStat->GetNbinsX();
182     
183     for(Int_t ibin = 1; ibin <= nbin; ibin++){
184       if(ibin == 1) {
185         if(fHistSyst)
186           fHEPDataFileLines->Add(new TObjString(Form("BinLow BinHigh %s +-stat +-syst", fValueName.Data()))); 
187         else 
188           fHEPDataFileLines->Add(new TObjString(Form("BinLow BinHigh %s +-stat", fValueName.Data())));          
189       }
190       // Skip empty bins
191       if(!fHistStat->GetBinContent(ibin)) continue;
192       TObjString * line = new TObjString;      
193       if(fHistSyst) {
194         line->String().Form("%f %f %f +-%f +-%f", 
195                             fHistStat->GetBinLowEdge(ibin), 
196                             fHistStat->GetBinLowEdge(ibin)+fHistStat->GetBinWidth(ibin), 
197                             fHistStat->GetBinContent(ibin), fHistStat->GetBinError(ibin), fHistSyst->GetBinError(ibin));
198         fHEPDataFileLines->Add(line);
199       } else {
200         line->String().Form("%f %f %f +-%f", 
201                             fHistStat->GetBinLowEdge(ibin), 
202                             fHistStat->GetBinLowEdge(ibin)+fHistStat->GetBinWidth(ibin), 
203                             fHistStat->GetBinContent(ibin), fHistStat->GetBinError(ibin));
204         fHEPDataFileLines->Add(line);
205       }
206       //      delete line;      
207     }           
208   }
209
210   TIterator * lineIter = fHEPDataFileLines->MakeIterator();
211   TObjString * obj = 0;
212   ofstream outfile;
213   outfile.open (hepfileName);
214   cout << "Saving HEP File " << hepfileName << endl;
215   
216   while ((obj = (TObjString*) lineIter->Next())) {    
217     cout << obj->String().Data() << endl;    
218     outfile << obj->String().Data() << endl;
219   }
220   outfile.close();
221 }
222
223