Logo ROOT   6.10/09
Reference Guide
SparseDataComparer.cxx
Go to the documentation of this file.
1 #include "TH2.h"
2 #include "TF2.h"
3 #include "TCanvas.h"
4 #include "TApplication.h"
5 
6 #include "TMath.h"
7 #include "Fit/SparseData.h"
8 #include "HFitInterface.h"
9 #include "Fit/Fitter.h"
10 #include "Math/WrappedMultiTF1.h"
11 
12 #include <iostream>
13 #include <iterator>
14 #include <algorithm>
15 #include <functional>
16 
17 using namespace std;
18 
19 
20 double minRange[3] = { -5., -5., -5.};
21 double maxRange[3] = { 5., 5., 5.};
22 int nbins[3] = {10 , 10 , 100 };
23 
24 bool showGraphics = false;
25 
26 ostream& operator << (ostream& out, ROOT::Fit::BinData& bd)
27 {
28  const unsigned int ndim( bd.NDim() );
29  const unsigned int npoints( bd.NPoints() );
30  for ( unsigned int i = 0; i < npoints; ++i )
31  {
32  double value = 0, error = 0;
33  const double *x = bd.GetPoint(i, value, error);
34  for ( unsigned int j = 0; j < ndim; ++j )
35  {
36  out << " x[" << j << "]: " << x[j];
37  }
38  out << " value: " << value;
39  out << " error: " << error;
40  out << endl;
41  }
42  return out;
43 }
44 
45 int findBin(ROOT::Fit::BinData& bd, const double *x)
46 {
47  const unsigned int ndim = bd.NDim();
48  const unsigned int npoints = bd.NPoints();
49 
50  for ( unsigned int i = 0; i < npoints; ++i )
51  {
52  double value1 = 0, error1 = 0;
53  const double *x1 = bd.GetPoint(i, value1, error1);
54  bool thisIsIt = true;
55  for ( unsigned int j = 0; j < ndim; ++j )
56  {
57  thisIsIt &= fabs(x1[j] - x[j]) < 1E-15;
58  }
59  if ( thisIsIt ) return i;
60  }
61 
62  cout << "ERROR FINDING BIN!" << endl;
63  return -1;
64 }
65 
67 {
68  const unsigned int ndim = bd1.NDim();
69  const unsigned int npoints = bd1.NPoints();
70 // unsigned int npoints2 = 0;
71 
72  bool equals = true;
73 
74  cout << "Equals" << endl;
75 
76  for ( unsigned int i = 0; i < npoints && equals; ++i )
77  {
78  double value1 = 0, error1 = 0 ;
79  const double *x1 = bd1.GetPoint(i, value1, error1);
80 
81  int bin = findBin(bd2, x1);
82 
83  double value2 = 0, error2 = 0;
84  const double *x2 = bd2.GetPoint(bin, value2, error2);
85 
86  equals &= ( value1 == value2 );
87  cout << " v: " << equals;
88  equals &= ( error1 == error2 );
89  cout << " e: " << equals;
90  for ( unsigned int j = 0; j < ndim; ++j )
91  {
92  equals &= fabs(x1[j] - x2[j]) < 1E-15;
93  cout << " x[" << j << "]: " << equals;
94  }
95 
96  cout << " bd1: ";
97  std::copy(x1, &x1[ndim], ostream_iterator<double>(cout, " "));
98  cout << " value:" << value1;
99  cout << " error:" << error1;
100 
101  cout << " bd2: ";
102  std::copy(x2, &x2[ndim], ostream_iterator<double>(cout, " "));
103  cout << " value:" << value2;
104  cout << " error:" << error2;
105 
106  cout << " equals: " << equals;
107 
108  cout << endl;
109  }
110 
111  return equals;
112 }
113 
114 
116 {
117  TH1D* h = new TH1D("h1", "h1-title", nbins[0], minRange[0], maxRange[0]);
118  h->FillRandom("gaus", 40);
119 
120 
122  ROOT::Fit::FillData(bd, h);
123 
124  cout << bd << endl;
125 
126  double min[] = { minRange[0] };
127  double max[] = { maxRange[0] };
128  ROOT::Fit::SparseData sd(1, min,max);
129  ROOT::Fit::FillData(sd, h);
130  ROOT::Fit::BinData bd2;
131  sd.GetBinData(bd2);
132 
133  cout << bd2 << endl;
134 
135  cout << " equals : ";
136  bool ok = (bd == bd2);
137 
138  cout << "One Dimension test ............\t";
139  if (ok)
140  cout << "OK\n";
141  else
142  cout << "FAILED\n";
143 
144 
145  if (showGraphics) {
146  h->Draw();
147  gPad->Update();
148  }
149 }
150 
151 double gaus2D(double *x, double *p)
152 {
153  return p[0]*TMath::Gaus(x[0],p[1],p[2]) * TMath::Gaus(x[1],p[3],p[4]);
154 }
155 
156 
158 {
159  TH2D* h = new TH2D("h2", "h2-title",
160  nbins[0], minRange[0], maxRange[0],
161  nbins[1], minRange[1], maxRange[1]);
162 
163  TF2* f2 = new TF2("gaus2D", gaus2D,
164  minRange[0],maxRange[0], minRange[1], maxRange[1], 5);
165  double initialPars[] = {300,0.,2.,0.,3.};
166  f2->SetParameters(initialPars);
167 
168  h->FillRandom("gaus2D",20);
169 
171  ROOT::Fit::FillData(bd, h);
172 
173  cout << bd << endl;
174 
175  double min[] = { minRange[0], minRange[1] };
176  double max[] = { maxRange[0], maxRange[1] };
177  ROOT::Fit::SparseData sd(2, min, max);
178 
179  ROOT::Fit::FillData(sd, h);
180  ROOT::Fit::BinData bd2;
181  sd.GetBinData(bd2);
182 
183  cout << bd2 << endl;
184 
185  cout << " equals : ";
186  bool ok = (bd == bd2);
187 
188  if (showGraphics) {
189  new TCanvas();
190  h->Draw("lego2");
191  gPad->Update();
192  }
193 
194  cout << "Two Dimension test............\t";
195  if (ok)
196  cout << "OK\n";
197  else
198  cout << "FAILED\n";
199 }
200 
201 int main(int argc, char** argv)
202 {
203 
204  // Parse command line arguments
205  for (Int_t i=1 ; i<argc ; i++) {
206  std::string arg = argv[i] ;
207  if (arg == "-g") {
208  showGraphics = true;
209  }
210  if (arg == "-v") {
211  showGraphics = true;
212  //verbose = true;
213  }
214  if (arg == "-h") {
215  cerr << "Usage: " << argv[0] << " [-g] [-v]\n";
216  cerr << " where:\n";
217  cerr << " -g : graphics mode\n";
218  cerr << " -v : verbose mode";
219  cerr << endl;
220  return -1;
221  }
222  }
223 
224  TApplication* theApp = 0;
225 
226  if (showGraphics)
227  theApp = new TApplication("App",&argc,argv);
228 
229  cout << "\nONE DIMENSION" << endl;
230  OneDimension();
231  cout << "\nTWO DIMENSIONS" << endl;
232  TwoDimensions();
233 
234  if (showGraphics) {
235  theApp->Run();
236  delete theApp;
237  theApp = 0;
238  }
239 
240  return 0;
241 }
242 
virtual void SetParameters(const Double_t *params)
Definition: TF1.h:588
virtual void FillRandom(const char *fname, Int_t ntimes=5000)
Fill histogram following distribution in function fname.
Definition: TH2.cxx:597
TH1 * h
Definition: legend2.C:5
int Int_t
Definition: RtypesCore.h:41
int nbins[3]
int equals(Double_t n1, Double_t n2, double ERRORLIMIT=1.E-10)
Definition: UnitTesting.cxx:24
void OneDimension()
STL namespace.
const double * GetPoint(unsigned int ipoint, double &value) const
retrieve at the same time a pointer to the coordinate data and the fit value More efficient than call...
Definition: BinData.h:347
static const double x2[5]
Double_t x[n]
Definition: legend1.C:17
virtual void Run(Bool_t retrn=kFALSE)
Main application eventloop. Calls system dependent eventloop via gSystem.
double gaus2D(double *x, double *p)
int main(int argc, char **argv)
VecExpr< UnaryOp< Fabs< T >, VecExpr< A, T, D >, T >, T, D > fabs(const VecExpr< A, T, D > &rhs)
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
ostream & operator<<(ostream &out, ROOT::Fit::BinData &bd)
virtual void FillRandom(const char *fname, Int_t ntimes=5000)
Fill histogram following distribution in function fname.
Definition: TH1.cxx:3294
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2851
int findBin(ROOT::Fit::BinData &bd, const double *x)
bool showGraphics
double maxRange[3]
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition: BinData.h:53
A 2-Dim function with parameters.
Definition: TF2.h:29
tomato 1-D histogram with a double per channel (see TH1 documentation)}
Definition: TH1.h:594
constexpr Double_t E()
Definition: TMath.h:74
Double_t Gaus(Double_t x, Double_t mean=0, Double_t sigma=1, Bool_t norm=kFALSE)
Calculate a gaussian function with mean and sigma.
Definition: TMath.cxx:452
The Canvas class.
Definition: TCanvas.h:31
static const double x1[5]
void GetBinData(BinData &) const
Definition: SparseData.cxx:298
bool operator==(ROOT::Fit::BinData &bd1, ROOT::Fit::BinData &bd2)
THist< 2, double, THistStatContent, THistStatUncertainty > TH2D
Definition: THist.hxx:316
double f2(const double *x)
THist< 1, double, THistStatContent, THistStatUncertainty > TH1D
Definition: THist.hxx:310
#define gPad
Definition: TVirtualPad.h:284
This class creates the ROOT Application Environment that interfaces to the windowing system eventloop...
Definition: TApplication.h:39
unsigned int NPoints() const
return number of fit points
Definition: FitData.h:294
unsigned int NDim() const
return coordinate data dimension
Definition: FitData.h:310
void TwoDimensions()
double minRange[3]
tomato 2-D histogram with a double per channel (see TH1 documentation)}
Definition: TH2.h:290