ROOT   Reference Guide
ChebyshevPol.h
Go to the documentation of this file.
1// @(#)root/mathcore:$Id$
2// Author: L. Moneta, 11/2012
3
4/*************************************************************************
5 * Copyright (C) 1995-2012, Rene Brun and Fons Rademakers. *
6 * All rights reserved. *
7 * *
8 * For the licensing terms see $ROOTSYS/LICENSE. * 9 * For the list of contributors see$ROOTSYS/README/CREDITS. *
10 *************************************************************************/
11
12//////////////////////////////////////////////////////////////////////////
13// //
14// Header file declaring functions for the evaluation of the Chebyshev //
15// polynomials and the ChebyshevPol class which can be used for //
16// creating a TF1. //
17// //
18//////////////////////////////////////////////////////////////////////////
19
20#ifndef ROOT_Math_ChebyshevPol
21#define ROOT_Math_ChebyshevPol
22
23#include <sys/types.h>
24#include <cstring>
25
26namespace ROOT {
27
28 namespace Math {
29
30 /// template recursive functions for defining evaluation of Chebyshev polynomials
31 /// T_n(x) and the series S(x) = Sum_i c_i* T_i(x)
32 namespace Chebyshev {
33
34 template<int N> double T(double x) {
35 return (2.0 * x * T<N-1>(x)) - T<N-2>(x);
36 }
37
38 template<> double T<0> (double );
39 template<> double T<1> (double x);
40 template<> double T<2> (double x);
41 template<> double T<3> (double x);
42
43 template<int N> double Eval(double x, const double * c) {
44 return c[N]*T<N>(x) + Eval<N-1>(x,c);
45 }
46
47 template<> double Eval<0> (double , const double *c);
48 template<> double Eval<1> (double x, const double *c);
49 template<> double Eval<2> (double x, const double *c);
50 template<> double Eval<3> (double x, const double *c);
51
52 } // end namespace Chebyshev
53
54
55 // implementation of Chebyshev polynomials using all coefficients
56 // needed for creating TF1 functions
57 inline double Chebyshev0(double , double c0) {
58 return c0;
59 }
60 inline double Chebyshev1(double x, double c0, double c1) {
61 return c0 + c1*x;
62 }
63 inline double Chebyshev2(double x, double c0, double c1, double c2) {
64 return c0 + c1*x + c2*(2.0*x*x - 1.0);
65 }
66 inline double Chebyshev3(double x, double c0, double c1, double c2, double c3) {
67 return c3*Chebyshev::T<3>(x) + Chebyshev2(x,c0,c1,c2);
68 }
69 inline double Chebyshev4(double x, double c0, double c1, double c2, double c3, double c4) {
70 return c4*Chebyshev::T<4>(x) + Chebyshev3(x,c0,c1,c2,c3);
71 }
72 inline double Chebyshev5(double x, double c0, double c1, double c2, double c3, double c4, double c5) {
73 return c5*Chebyshev::T<5>(x) + Chebyshev4(x,c0,c1,c2,c3,c4);
74 }
75 inline double Chebyshev6(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6) {
76 return c6*Chebyshev::T<6>(x) + Chebyshev5(x,c0,c1,c2,c3,c4,c5);
77 }
78 inline double Chebyshev7(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7) {
79 return c7*Chebyshev::T<7>(x) + Chebyshev6(x,c0,c1,c2,c3,c4,c5,c6);
80 }
81 inline double Chebyshev8(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8) {
82 return c8*Chebyshev::T<8>(x) + Chebyshev7(x,c0,c1,c2,c3,c4,c5,c6,c7);
83 }
84 inline double Chebyshev9(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8, double c9) {
85 return c9*Chebyshev::T<9>(x) + Chebyshev8(x,c0,c1,c2,c3,c4,c5,c6,c7,c8);
86 }
87 inline double Chebyshev10(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8, double c9, double c10) {
88 return c10*Chebyshev::T<10>(x) + Chebyshev9(x,c0,c1,c2,c3,c4,c5,c6,c7,c8,c9);
89 }
90
91
92 // implementation of Chebyshev polynomial with run time parameter
93 inline double ChebyshevN(unsigned int n, double x, const double * c) {
94
95 if (n == 0) return Chebyshev0(x,c[0]);
96 if (n == 1) return Chebyshev1(x,c[0],c[1]);
97 if (n == 2) return Chebyshev2(x,c[0],c[1],c[2]);
98 if (n == 3) return Chebyshev3(x,c[0],c[1],c[2],c[3]);
99 if (n == 4) return Chebyshev4(x,c[0],c[1],c[2],c[3],c[4]);
100 if (n == 5) return Chebyshev5(x,c[0],c[1],c[2],c[3],c[4],c[5]);
101
102 /* do not use recursive formula
103 (2.0 * x * Tn(n - 1, x)) - Tn(n - 2, x) ;
104 which is too slow for large n
105 */
106
107 size_t i;
108 double d1 = 0.0;
109 double d2 = 0.0;
110
111 // if not in range [-1,1]
112 //double y = (2.0 * x - a - b) / (b - a);
113 //double y = x;
114 double y2 = 2.0 * x;
115
116 for (i = n; i >= 1; i--)
117 {
118 double temp = d1;
119 d1 = y2 * d1 - d2 + c[i];
120 d2 = temp;
121 }
122
123 return x * d1 - d2 + c[0];
124 }
125
126
127 // implementation of Chebyshev Polynomial class
128 // which can be used for building TF1 classes
130 public:
131 ChebyshevPol(unsigned int n) : fOrder(n) {}
132
133 double operator() (const double *x, const double * coeff) {
134 return ChebyshevN(fOrder, x[0], coeff);
135 }
136 private:
137 unsigned int fOrder;
138 };
139
140
141
142 } // end namespace Math
143
144} // end namespace ROOT
145
146
147
148#endif // ROOT_Math_Chebyshev
#define c(i)
Definition: RSha256.hxx:101
#define N
ChebyshevPol(unsigned int n)
Definition: ChebyshevPol.h:131
double operator()(const double *x, const double *coeff)
Definition: ChebyshevPol.h:133
return c1
Definition: legend1.C:41
Double_t x[n]
Definition: legend1.C:17
const Int_t n
Definition: legend1.C:16
return c2
Definition: legend2.C:14
return c3
Definition: legend3.C:15
Namespace for new Math classes and functions.
double Eval< 3 >(double x, const double *c)
double Eval(double x, const double *c)
Definition: ChebyshevPol.h:43
double Eval< 2 >(double x, const double *c)
double T< 1 >(double x)
double T< 2 >(double x)
double T< 3 >(double x)
double Eval< 0 >(double, const double *c)
double T< 0 >(double)
double Eval< 1 >(double x, const double *c)
double T(double x)
Definition: ChebyshevPol.h:34
double Chebyshev10(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8, double c9, double c10)
Definition: ChebyshevPol.h:87
double ChebyshevN(unsigned int n, double x, const double *c)
Definition: ChebyshevPol.h:93
double Chebyshev0(double, double c0)
Definition: ChebyshevPol.h:57
double Chebyshev2(double x, double c0, double c1, double c2)
Definition: ChebyshevPol.h:63
double Chebyshev3(double x, double c0, double c1, double c2, double c3)
Definition: ChebyshevPol.h:66
double Chebyshev4(double x, double c0, double c1, double c2, double c3, double c4)
Definition: ChebyshevPol.h:69
double Chebyshev1(double x, double c0, double c1)
Definition: ChebyshevPol.h:60
double Chebyshev6(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6)
Definition: ChebyshevPol.h:75
double Chebyshev7(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7)
Definition: ChebyshevPol.h:78
double Chebyshev5(double x, double c0, double c1, double c2, double c3, double c4, double c5)
Definition: ChebyshevPol.h:72
double Chebyshev8(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8)
Definition: ChebyshevPol.h:81
double Chebyshev9(double x, double c0, double c1, double c2, double c3, double c4, double c5, double c6, double c7, double c8, double c9)
Definition: ChebyshevPol.h:84
tbb::task_arena is an alias of tbb::interface7::task_arena, which doesn't allow to forward declare tb...