[u/mrichter/AliRoot.git] / ANALYSIS / AliHEPDataParser.cxx

//-------------------------------------------------------------------------
//           Implementation of Class AliHEPDataParser
//
//  This class is used to save the content of hisograms/graphs in the
//  HEP data format and viceversa. The HEP data format is not strictly
//  defined and there are many variants, the class only support a few
//  of them. More will be added as needed.  The input can be a set of
//  2 TH1, TGraphAsymmErrors or TGraphErrors (one for the stat and one
//  for the syst error). If the second one is a null pointer, only the
//  stat error is printed. The class can also import hepdata ascii
//  file (very preliminary)
// 
//  Author: Michele Floris, CERN
//-------------------------------------------------------------------------


#include "AliHEPDataParser.h"
#include "AliLog.h"
#include "TGraphAsymmErrors.h"
#include "TGraph.h"
#include "TGraphErrors.h"
#include "TH1.h"
#include "TObjArray.h"
#include "TObjString.h"
#include <fstream>
#include <iostream>
using namespace std;

ClassImp(AliHEPDataParser)

AliHEPDataParser::AliHEPDataParser() : TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("")
{
  // default ctor

}

AliHEPDataParser::AliHEPDataParser(TH1 * hStat, TH1 * hSyst): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("y")
{
  //ctor
  fHistStat = hStat;
  fHistSyst = hSyst;
  fHEPDataFileLines = new TObjArray;

}
AliHEPDataParser::AliHEPDataParser(TGraph * grStat, TGraph * grSyst): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("")
{
  // ctor
  fGraphStat = grStat;
  fGraphSyst = grSyst;
  fHEPDataFileLines = new TObjArray;
}

AliHEPDataParser::AliHEPDataParser(const char * hepfileName): TObject(), fHistStat(0),  fHistSyst(0),  fGraphStat(0),  fGraphSyst(0),  fHEPDataFileLines(0), fValueName("y")
{
  //ctor
  // Put results in graphs
  fGraphSyst = new TGraphAsymmErrors();
  fGraphStat = new TGraphAsymmErrors();
  ifstream infile;
  infile.open(hepfileName);
  TString line;
  Int_t ipoints = 0;
  while (line.ReadLine(infile)) {
    TObjArray * tokens = line.Tokenize(" ");
    if(tokens->GetEntries() < 1) {
      delete tokens;
      AliError("not enough columns");
      return;      
    }
    // TODO: Assumes the format
    // binmin binmax y +-stat +-syst. Try to make it smarter...
    TObjString * binMin = (TObjString*) tokens->At(0);
    TObjString * binMax = (TObjString*) tokens->At(1);
    TObjString * value  = (TObjString*) tokens->At(2);
    TObjString * stat   = (TObjString*) tokens->At(3);
    TObjString * syst   = (TObjString*) tokens->At(4);
    stat->String().ReplaceAll("+-","");
    if(syst) syst->String().ReplaceAll("+-","");
    if (!binMin->String().Atof()) {delete tokens; continue;} // skip headers
    Float_t binCenter = (binMax->String().Atof() + binMin->String().Atof())/2;
    Float_t binWidth  =  (binMax->String().Atof() - binMin->String().Atof())/2;
    cout << line.Data() << endl;//<< " " << binMin->String().Atof() <<" " <<  binCenter << " " << binWidth << endl;


    fGraphStat->SetPoint(ipoints, binCenter, value->String().Atof());
    fGraphSyst->SetPoint(ipoints, binCenter, value->String().Atof());
    ((TGraphAsymmErrors*)fGraphStat)->SetPointError(ipoints, 
						    binWidth,
						    binWidth,
						    stat->String().Atof(),
						    stat->String().Atof());
    if(syst) ((TGraphAsymmErrors*)fGraphSyst)->SetPointError(ipoints, 
							     binWidth,
							     binWidth,
							     syst->String().Atof(),
							     syst->String().Atof());
    ipoints++;
    delete tokens;
  }
  infile.close();
    

}

AliHEPDataParser::~AliHEPDataParser(){
  // dtor
  if(fHistStat) delete fHistStat;
  if(fHistSyst) delete fHistSyst;
  if(fGraphStat) delete fGraphStat;
  if(fGraphSyst) delete fGraphSyst;
  if(fHEPDataFileLines) delete fHEPDataFileLines;
}
  
void AliHEPDataParser::SaveHEPDataFile(const char * hepfileName, Bool_t trueUseGraphFalesUseHisto) {
  // Fills fHEPDataFileLines and saves its content to a file
  if(!fHEPDataFileLines)   fHEPDataFileLines = new TObjArray;
  if(trueUseGraphFalesUseHisto) {
    if(!fGraphStat) {
      AliError("Graph not set");
      return;
    }
    Bool_t asym = kFALSE; // check if this has asymmetric errors
    if       (!strcmp(fGraphStat->ClassName(), "TGraphErrors")) asym = kFALSE;
    else if  (!strcmp(fGraphStat->ClassName(), "TGraphAsymmErrors")) asym = kTRUE;
    else     {AliError("Unsupported graph type"); return;}
    Int_t npoint = fGraphStat->GetN();
    if(asym) AliInfo("Assymmetric errors");
    for(Int_t ipoint = 0; ipoint < npoint; ipoint++){            
      if(ipoint == 0) {
	if(fGraphSyst) {
	  if (asym)
	    fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +stat -stat +syst -syst", fValueName.Data()))); 
	  else 
	    fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +-stat +-syst", fValueName.Data()))); 
	}
	else {
	  if(asym)
	    fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +stat -stat", fValueName.Data()))); 
	  else 
	    fHEPDataFileLines->Add(new TObjString(Form("BinCenter %s +-stat", fValueName.Data()))); 
	}
      }
      // Skip empty bins
      if(!fGraphStat->GetY()[ipoint]) continue;
      TObjString * line = new TObjString;      
      if(fGraphSyst) {
	if (asym)
	  line->String().Form("%f %f +%f -%f +%f -%f", 
			      fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
			      ((TGraphAsymmErrors*)fGraphStat)->GetEYhigh()[ipoint], 
			      ((TGraphAsymmErrors*)fGraphStat)->GetEYlow()[ipoint], 
			      ((TGraphAsymmErrors*)fGraphSyst)->GetEYhigh()[ipoint], 
			      ((TGraphAsymmErrors*)fGraphSyst)->GetEYlow()[ipoint]);
	else 
	  line->String().Form("%f %f +-%f +-%f", 
			      fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
			      ((TGraphErrors*)fGraphStat)->GetEY()[ipoint], 
			      ((TGraphErrors*)fGraphSyst)->GetEY()[ipoint]);

	fHEPDataFileLines->Add(line);
      } else {
	if (asym)
	  line->String().Form("%f %f +%f -%f", 
			      fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
			      ((TGraphAsymmErrors*)fGraphStat)->GetEYhigh()[ipoint], ((TGraphAsymmErrors*)fGraphStat)->GetEYlow()[ipoint]);
	else { 
	  line->String().Form("%f %f +-%f", 
			      fGraphStat->GetX()[ipoint], fGraphStat->GetY()[ipoint],
			      ((TGraphErrors*)fGraphStat)->GetEY()[ipoint]);
	}

	fHEPDataFileLines->Add(line);
      }
    }    
  }
  else {
    if(!fHistStat) {
      AliError("Hist not set");
      return;
    }
    Int_t nbin = fHistStat->GetNbinsX();
    
    for(Int_t ibin = 1; ibin <= nbin; ibin++){
      if(ibin == 1) {
	if(fHistSyst)
	  fHEPDataFileLines->Add(new TObjString(Form("BinLow BinHigh %s +-stat +-syst", fValueName.Data()))); 
	else 
	  fHEPDataFileLines->Add(new TObjString(Form("BinLow BinHigh %s +-stat", fValueName.Data())));       	
      }
      // Skip empty bins
      if(!fHistStat->GetBinContent(ibin)) continue;
      TObjString * line = new TObjString;      
      if(fHistSyst) {
	line->String().Form("%f %f %f +-%f +-%f", 
			    fHistStat->GetBinLowEdge(ibin), 
			    fHistStat->GetBinLowEdge(ibin)+fHistStat->GetBinWidth(ibin), 
			    fHistStat->GetBinContent(ibin), fHistStat->GetBinError(ibin), fHistSyst->GetBinError(ibin));
	fHEPDataFileLines->Add(line);
      } else {
	line->String().Form("%f %f %f +-%f", 
			    fHistStat->GetBinLowEdge(ibin), 
			    fHistStat->GetBinLowEdge(ibin)+fHistStat->GetBinWidth(ibin), 
			    fHistStat->GetBinContent(ibin), fHistStat->GetBinError(ibin));
	fHEPDataFileLines->Add(line);
      }
      //      delete line;      
    }           
  }

  TIterator * lineIter = fHEPDataFileLines->MakeIterator();
  TObjString * obj = 0;
  ofstream outfile;
  outfile.open (hepfileName);
  cout << "Saving HEP File " << hepfileName << endl;
  
  while ((obj = (TObjString*) lineIter->Next())) {    
    cout << obj->String().Data() << endl;    
    outfile << obj->String().Data() << endl;
  }
  outfile.close();
}
Commit	Line	Data
a5ae512a	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("+-","");
09d5920f	79	if (!binMin->String().Atof()) {delete tokens; continue;} // skip headers
a5ae512a	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++;
09d5920f	98	delete tokens;
a5ae512a	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