Logo ROOT   6.08/07
Reference Guide
TF1NormSum.cxx
Go to the documentation of this file.
1 // @(#)root/hist:$Id$
2 // Authors: L. Moneta, A. Flandi 08/2014
3 //
4 //
5 /**********************************************************************
6  * *
7  * Copyright (c) 2015 ROOT Team, CERN/PH-SFT *
8  * *
9  * *
10  **********************************************************************/
11 //
12 // TF1NormSum.cxx
13 //
14 //
15 //
16 #include "TROOT.h"
17 #include "TClass.h"
18 #include "TMath.h"
19 #include "TF1NormSum.h"
20 #include "Math/WrappedFunction.h"
21 #include "Math/WrappedTF1.h"
22 
23 
24 /** \class TF1NormSum
25  \ingroup Hist
26 Class adding two functions: c1*f1+c2*f2
27 */
28 
29 //ClassImp(TF1NormSum)
30 
31 
32 // function to find and rename duplicate parameters with the same name
33 
34 template<class Iterator>
35 void FixDuplicateNames(Iterator begin, Iterator end) {
36 
37 
38  // make a map of values
39 
40  std::multimap<TString, int > parMap;
41  for (Iterator it = begin; it != end; ++it) {
42  parMap.insert( std::make_pair( *it, std::distance(begin,it) ) );
43  }
44  for ( auto & elem : parMap) {
45  TString name = elem.first;
46  int n = parMap.count( name);
47  if (n > 1 ) {
48  std::pair <std::multimap<TString,int>::iterator, std::multimap<TString,int>::iterator> ret;
49  ret = parMap.equal_range(name);
50  int i = 0;
51  for (std::multimap<TString,int>::iterator it=ret.first; it!=ret.second; ++it) {
52  *(begin+it->second) = TString::Format("%s%d",name.Data(),++i);
53  }
54  }
55  }
56 
57  // for (Iterator it = begin; it != end; ++it)
58  // std::cout << *it << " ";
59  // std::cout << std::endl;
60 
61 }
62 
63 
64 
65 void TF1NormSum::InitializeDataMembers(const std::vector <std::shared_ptr < TF1 >> &functions, const std::vector <Double_t> &coeffs, Double_t scale)
66 {
67 
68  fScale = scale;
69  fFunctions = functions;
70  fCoeffs = coeffs;
71  fNOfFunctions = functions.size();
72  fCstIndexes = std::vector < Int_t > (fNOfFunctions);
73  fParNames = std::vector<TString> (fNOfFunctions);
74  fParNames.reserve(3*fNOfFunctions); // enlarge capacity for function parameters
75 
76  for (unsigned int n=0; n < fNOfFunctions; n++)
77  {
78  int npar = fFunctions[n] -> GetNpar();
79  fCstIndexes[n] = fFunctions[n] -> GetParNumber("Constant");//return -1 if there is no constant parameter
80  //std::cout << " cst index of function " << n << " : " << fCstIndexes[n] << std::endl;
81  //std::cout << "nofparam of function " << n <<" : " << fNOfParams[n] << std::endl;
82  fParNames[n] = TString::Format("Coeff%d",n);
83  //printf("examing function %s \n",fFunctions[n]->GetName() );
84  if (fCstIndexes[n]!= -1) //if there exists a constant parameter
85  {
86  fFunctions[n] -> FixParameter(fCstIndexes[n], 1.); //fixes the parameters called "Constant" to 1
87  int k = 0; //index for the temp arry, k wil go form 0 until fNofNonCstParameter
88  for (int i=0; i<npar; i++) //go through all the parameter to
89  {
90  if (i==fCstIndexes[n]) continue; //go to next step if this is the constant parameter
91  fParNames.push_back( fFunctions[n] -> GetParName(i) );
92  k++;
93  }
94  }
95  else {
96  for (int i=0; i < npar; i++) //go through all the parameter to
97  {
98  fParNames.push_back( fFunctions[n] -> GetParName(i) );
99  }
100  }
101  //normalize the functions if it is not already done (do at the end so constant parameter is not zero)
102  if (!fFunctions[n] -> IsEvalNormalized()) fFunctions[n] -> SetNormalized(true);
103  }
104 
105  FixDuplicateNames(fParNames.begin()+fNOfFunctions, fParNames.end());
106 
107 }
109 {
110  fNOfFunctions = 0;
111  fScale = 1.;
112  fFunctions = std::vector< std::shared_ptr < TF1 >>(0) ; // Vector of size fNOfFunctions containing TF1 functions
113  fCoeffs = std::vector < Double_t >(0) ; // Vector of size fNOfFunctions containing coefficients in front of each function
114  fCstIndexes = std::vector < Int_t > (0);
115 }
116 
117 ////////////////////////////////////////////////////////////////////////////////
118 
119 TF1NormSum::TF1NormSum(const std::vector <TF1*> &functions, const std::vector <Double_t> &coeffs, Double_t scale)
120 {
121  std::vector <std::shared_ptr < TF1 > >f;
122  for (unsigned int i = 0; i<functions.size(); i++)
123  {
124  f[i] = std::shared_ptr < TF1 >((TF1*)functions[i]->Clone());
125  }
126 
127  InitializeDataMembers(f,coeffs,scale);
128 }
129 
130 ////////////////////////////////////////////////////////////////////////////////
131 /// TF1NormSum constructor taking 2 functions, and 2 coefficients (if not equal to 1)
132 
133 TF1NormSum::TF1NormSum(TF1* function1, TF1* function2, Double_t coeff1, Double_t coeff2, Double_t scale)
134 {
135  std::vector < std::shared_ptr < TF1 > > functions(2);
136  std::vector < Double_t > coeffs(2);
137  TF1 * fnew1 = 0;
138  TF1 * fnew2 = 0;
139  // need to use Copy because clone does not work for functor-based functions
140  if (function1) {
141  fnew1 = (TF1*) function1->IsA()->New();
142  function1->Copy(*fnew1);
143  }
144  if (function2) {
145  fnew2 = (TF1*) function2->IsA()->New();
146  function2->Copy(*fnew2);
147  }
148  if (fnew1 == nullptr || fnew2 == nullptr)
149  Fatal("TF1NormSum","Invalid input functions - Abort");
150 
151  std::shared_ptr < TF1 > f1( fnew1);
152  std::shared_ptr < TF1 > f2( fnew2);
153 
154  functions = {f1, f2};
155  coeffs = {coeff1, coeff2};
156 
157  InitializeDataMembers(functions, coeffs,scale);
158 }
159 
160 ////////////////////////////////////////////////////////////////////////////////
161 /// TF1NormSum constructor taking 3 functions, and 3 coefficients (if not equal to 1)
162 
163 TF1NormSum::TF1NormSum(TF1* function1, TF1* function2, TF1* function3, Double_t coeff1, Double_t coeff2, Double_t coeff3, Double_t scale)
164 {
165  std::vector < std::shared_ptr < TF1 > > functions(3);
166  std::vector < Double_t > coeffs(3);
167  TF1 * fnew1 = 0;
168  TF1 * fnew2 = 0;
169  TF1 * fnew3 = 0;
170  if (function1) {
171  fnew1 = (TF1*) function1->IsA()->New();
172  function1->Copy(*fnew1);
173  }
174  if (function2) {
175  fnew2 = (TF1*) function2->IsA()->New();
176  function2->Copy(*fnew2);
177  }
178  if (function3) {
179  fnew3 = (TF1*) function3->IsA()->New();
180  function3->Copy(*fnew2);
181  }
182  if (!fnew1 || !fnew2 || !fnew3 )
183  Fatal("TF1NormSum","Invalid input functions - Abort");
184 
185  std::shared_ptr < TF1 > f1( fnew1);
186  std::shared_ptr < TF1 > f2( fnew2);
187  std::shared_ptr < TF1 > f3( fnew3);
188 
189 
190  functions = {f1, f2, f3};
191  coeffs = {coeff1, coeff2, coeff3};
192 
193  InitializeDataMembers(functions, coeffs,scale);
194 }
195 
196 ////////////////////////////////////////////////////////////////////////////////
197 /// TF1NormSum constructor taking any addition of formulas with coefficient or not
198 ///
199 /// - example 1 : 2.*expo + gauss + 0.5* gauss
200 /// - example 2 : expo + 0.3*f1 if f1 is defined in the list of fucntions
201 
203 {
205 
206  TObjArray *arrayall = formula.Tokenize("*+");
207  TObjArray *arraytimes = formula.Tokenize("*") ;
208  Int_t noffunctions = (formula.Tokenize("+")) -> GetEntries();
209  Int_t nofobj = arrayall -> GetEntries();
210  Int_t nofcoeffs = nofobj - noffunctions;
211 
212  std::vector < std::shared_ptr < TF1 > > functions(noffunctions);
213  std::vector < Double_t > coeffs(noffunctions);
214  std::vector < TString > funcstringall(nofobj);
215  std::vector < Int_t > indexsizetimes(nofcoeffs+1);
216  std::vector < Bool_t > isacoeff(nofobj);//1 is it is a coeff, 0 if it is a functions
217 
218  for (int i=0; i<nofobj; i++)
219  {
220  funcstringall[i] = ((TObjString*)((*arrayall)[i])) -> GetString();
221  funcstringall[i].ReplaceAll(" ","");
222  }
223  //algorithm to determine which object is a coefficient and which is a function
224  //uses the fact that the last item of funcstringtimes[i].Tokenize("+") is always a coeff.
225  Int_t j = 0;
226  Int_t k = 1;
227  for (int i=0; i<nofcoeffs+1; i++)
228  {
229  indexsizetimes[i] = ( ( ( (TObjString*)(*arraytimes)[i] ) -> GetString() ).Tokenize("+") ) -> GetEntries();
230  while (k < indexsizetimes[i])
231  {
232  isacoeff[k+j-1] = 0;
233  k++;
234  }
235  j = j+indexsizetimes[i];
236  if (j==nofobj) isacoeff[j-1] = 0; //the last one is never a coeff
237  else isacoeff[j-1] = 1;
238  k = 1;
239  }
240  k = 0;
241  for (int i=0; i<noffunctions; i++)
242  {
243  if (isacoeff[k]==0)
244  {
245  coeffs[i] = 1.;
246  TF1* f = (TF1*)(gROOT -> GetListOfFunctions() -> FindObject(funcstringall[k]));
247  if (!f) Error("TF1NormSum", "Function %s does not exist", funcstringall[k].Data());
248  functions[i] = std::shared_ptr < TF1 > ((TF1*)f->Clone(TString::Format("function_%s_%d",funcstringall[k].Data(), i)));
249  functions[i]->SetRange(xmin,xmax);
250  k++;
251  }
252  else
253  {
254  coeffs[i] = funcstringall[k].Atof();
255  TF1* f = (TF1*)(gROOT -> GetListOfFunctions() -> FindObject(funcstringall[k+1]));
256  if (!f) Error("TF1NormSum", "Function %s does not exist", funcstringall[k+1].Data());
257  functions[i] = std::shared_ptr < TF1 >((TF1*)f->Clone(TString::Format("function_%s_%d",funcstringall[k+1].Data(), i) ));
258  functions[i]->SetRange(xmin,xmax);
259  k=k+2;
260  }
261  }
262  InitializeDataMembers(functions, coeffs,1.);
263 
264  /*for (auto f : functions)
265  {
266  f->Print();
267  }
268  for (auto c : coeffs)
269  {
270  std::cout << "coeff " << c << std::endl;
271  }*/
272 }
273 
274 ////////////////////////////////////////////////////////////////////////////////
275 /// Overload the parenthesis to add the functions
276 
277 double TF1NormSum::operator()(double* x, double* p)
278 {
279  if (p!=0) TF1NormSum::SetParameters(p); // first refresh the parameters
280 
281  Double_t sum = 0.;
282  for (unsigned int n=0; n<fNOfFunctions; n++)
283  {
284  sum += fCoeffs[n]*(fFunctions[n] -> EvalPar(x,0));
285  }
286  // normalize by a scale parameter (typically the bin width)
287  return fScale * sum;
288 }
289 
290 ////////////////////////////////////////////////////////////////////////////////
291 /// Return array of parameters
292 
293 std::vector<double> TF1NormSum::GetParameters() const {
294  std::vector<double> params(GetNpar() );
295  int offset = 0;
296  int nOfNonCstParams = 0;
297  for (unsigned int n=0; n<fNOfFunctions; n++)
298  {
299  params[n] = fCoeffs[n]; // copy the coefficients
300  offset += nOfNonCstParams; // offset to go along the list of parameters
301  int k = 0;
302  for (int j = 0; j < fFunctions[n]->GetNpar(); ++j) {
303  if (j != fCstIndexes[n]) {
304  params[k+fNOfFunctions+offset] = fFunctions[n]->GetParameter(j);
305  k++;
306  }
307  }
308  nOfNonCstParams = k;
309  }
310  return params;
311 }
312 ////////////////////////////////////////////////////////////////////////////////
313 /// Initialize array of all parameters.
314 ///
315 /// double *params must contains first an array of the coefficients, then an array of the parameters.
316 
317 void TF1NormSum::SetParameters(const double* params)//params should have the size [fNOfFunctions][fNOfNonCstParams]
318 {
319  for (unsigned int n=0; n<fNOfFunctions; n++) //initialization of the coefficients
320  {
321  fCoeffs[n] = params[n];
322  }
323  Int_t offset = 0;
324  int k = 0; // k indicates the nnumber of non-constant parameter per function
325  for (unsigned int n=0; n<fNOfFunctions; n++)
326  {
327  bool equalParams = true;
328  Double_t * funcParams = fFunctions[n]->GetParameters();
329  int npar = fFunctions[n]->GetNpar();
330  offset += k; // offset to go along the list of parameters
331  k = 0; // reset k value for next function
332  for (int i = 0; i < npar; ++i) {
333  // constant parameters can be only one
334  if (i != fCstIndexes[n])
335  {
336  // check if they are equal
337  equalParams &= (funcParams[i] == params[k+fNOfFunctions+offset] );
338  funcParams[i] = params[k+fNOfFunctions+offset];
339  k++;
340  }
341  }
342  // update function integral if not equal
343  if (!equalParams) fFunctions[n]->Update();
344 
345  }
346 }
347 
348 ////////////////////////////////////////////////////////////////////////////////
349 /// Initialize array of all parameters.
350 ///
351 /// Overload the TF1::SetParameters() method.
352 /// A maximum of 10 parameters must be used, with first the coefficients, then the parameters
353 
355  Double_t p5, Double_t p6, Double_t p7, Double_t p8, Double_t p9, Double_t p10)
356 {
357  const double params[] = {p0,p1,p2,p3,p4,p5,p6,p7,p8,p9,p10};
359 
360 }
361 //return the number of (non constant) paramters including the coefficients: for 2 functions: c1,c2,p0,p1,p2,p3...
363 {
364  Int_t nofparams = 0;
365  for (unsigned int n=0; n<fNOfFunctions; ++n)
366  {
367  nofparams += fFunctions[n]->GetNpar();
368  if (fCstIndexes[n] >= 0) nofparams -= 1;
369  }
370  return nofparams + fNOfFunctions; //fNOfFunctions for the coefficientws
371 }
static long int sum(long int i)
Definition: Factory.cxx:1786
An array of TObjects.
Definition: TObjArray.h:39
Int_t GetNpar() const
Definition: TF1NormSum.cxx:362
float xmin
Definition: THbookFile.cxx:93
std::vector< TString > fParNames
Definition: TF1NormSum.h:42
static double p3(double t, double a, double b, double c, double d)
virtual void Copy(TObject &f1) const
Copy this F1 to a new F1.
Definition: TF1.cxx:770
void Fatal(const char *location, const char *msgfmt,...)
Collectable string class.
Definition: TObjString.h:32
const char * GetParName(Int_t ipar) const
Definition: TF1NormSum.h:70
#define gROOT
Definition: TROOT.h:364
Basic string class.
Definition: TString.h:137
int Int_t
Definition: RtypesCore.h:41
unsigned int fNOfFunctions
Definition: TF1NormSum.h:32
Double_t x[n]
Definition: legend1.C:17
std::vector< Int_t > fCstIndexes
Definition: TF1NormSum.h:41
static TString Format(const char *fmt,...)
Static method which formats a string using a printf style format descriptor and return a TString...
Definition: TString.cxx:2335
static double p2(double t, double a, double b, double c)
std::vector< std::vector< double > > Data
void InitializeDataMembers(const std::vector< std::shared_ptr< TF1 >> &functions, const std::vector< Double_t > &coeffs, Double_t scale)
Definition: TF1NormSum.cxx:65
void Error(const char *location, const char *msgfmt,...)
void FixDuplicateNames(Iterator begin, Iterator end)
Definition: TF1NormSum.cxx:35
static void InitStandardFunctions()
Create the basic function objects.
Definition: TF1.cxx:2249
static double p1(double t, double a, double b)
float xmax
Definition: THbookFile.cxx:93
std::vector< Double_t > fCoeffs
Definition: TF1NormSum.h:36
TObjArray * Tokenize(const TString &delim) const
This function is used to isolate sequential tokens in a TString.
Definition: TString.cxx:2241
double f(double x)
double Double_t
Definition: RtypesCore.h:55
Double_t fScale
Number of functions to add.
Definition: TF1NormSum.h:34
std::vector< double > GetParameters() const
Return array of parameters.
Definition: TF1NormSum.cxx:293
void SetParameters(const double *params)
Initialize array of all parameters.
Definition: TF1NormSum.cxx:317
std::vector< std::shared_ptr< TF1 > > fFunctions
Definition: TF1NormSum.h:35
virtual TObject * Clone(const char *newname="") const
Make a clone of an object using the Streamer facility.
Definition: TNamed.cxx:65
double f2(const double *x)
1-Dim function class
Definition: TF1.h:149
TF1 * f1
Definition: legend1.C:11
double operator()(double *x, double *p)
Overload the parenthesis to add the functions.
Definition: TF1NormSum.cxx:277
const Int_t n
Definition: legend1.C:16
char name[80]
Definition: TGX11.cxx:109
const char * Data() const
Definition: TString.h:349