Import fmit upstream version 0.97.6
[fmit.git] / libs / CppAddons / CAMath.h
1 // Copyright 2003 "Gilles Degottex"
2
3 // This file is part of "CppAddons"
4
5 // "CppAddons" is free software; you can redistribute it and/or modify
6 // it under the terms of the GNU Lesser General Public License as published by
7 // the Free Software Foundation; either version 2.1 of the License, or
8 // (at your option) any later version.
9 //
10 // "CppAddons" is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU Lesser General Public License for more details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with this program; if not, write to the Free Software
17 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
18
19
20 #ifndef _Math_h_
21 #define _Math_h_
22
23 #include <math.h>
24 #include <cmath>
25 #include <complex>
26 using namespace std;
27
28 #undef min
29 #undef max
30
31 namespace Math
32 {
33         template<typename TypeData> inline TypeData sgn(TypeData a)                     {return (a<0)?-1:1;}
34
35         static const double Pi = 2*acos(0);
36         static const float fPi = 2*acos(0);
37
38         static const double E = exp(1);
39         static const float fE = exp(1);
40
41         // résoud une equation du 2ème degré
42         class SolOfEq2
43         {
44           public:
45                 enum ENError{NE_OK=0, NE_DISCRIMINENT_NEG, NE_A_AND_B_EQ_ZERO, NE_RACINE_NEG, NE_X1_AND_X2_NEG, NE_X1_AND_X2_POS};
46
47           private:
48                 ENError m_err;
49
50                 double  x1;
51                 double  x2;
52
53           public:
54                 double getX1(){return x1;}
55                 double getX2(){return x2;}
56                 double getPosSol();
57
58                 SolOfEq2(double a, double b, double c);
59         };
60
61         // calcul l'intérale de f sur [a;b] avec un pas de h
62         // méhode de Simpson
63         template<class Function>
64         double Simpson(double a, double b, Function f, double h)
65         {
66                 double I4=f(a+h/2.0), I2=0;
67                 for(double x4=a+(h/2.0)+h, x2=a+h; x4<b; x4+=h, x2+=h)
68                 {
69                         I4+=f(x4);
70                         I2+=f(x2);
71                 }
72                 return (h/6.0)*(f(a)+4*I4+2*I2+f(b));
73         }
74
75         template<typename Type> std::complex<Type> make_complex(Type value[]){return std::complex<Type>(value[0], value[1]);}
76
77         inline double modulo(double d1, double d2)
78         {
79                     return d1-int(d1/d2)*d2;
80         }
81         inline double mod2(const double c[2])
82         {
83                 return c[0]*c[0]+c[1]*c[1];
84         }
85         inline double mod(const double c[2])
86         {
87                 return sqrt(mod2(c));
88         }
89         inline double mod2(const std::complex<double>& c)
90         {
91                 return c.real()*c.real()+c.imag()*c.imag();
92         }
93         inline double mod(const std::complex<double>& c)
94         {
95                 return sqrt(mod2(c));
96         }
97
98         inline double mod_equal(double& d1, double d2)
99         {
100                     return d1 -= int(d1/d2)*d2;
101         }
102
103         //! gauss fonction
104         /*!
105                 * \param x \f$\in ]-\infty,\infty[\f$
106         */
107         inline double gauss(double x)
108         {
109                 return exp(-Math::Pi*x*x);
110         }
111
112         inline double sinc(double t)
113         {
114                 if(t==0.0)      return 1.0;
115
116                 return sin(Math::Pi*t)/(Math::Pi*t);
117         }
118
119 /*      doesn't need for Linux for sure, seems to be probelatic under win32 (macro ambiquity)
120         template<typename TypeData1, typename TypeData2>
121         TypeData1 max(const TypeData1& a, const TypeData2& b)
122         {
123         return (a>b)?a:b;
124 }
125
126         template<typename TypeData1, typename TypeData2>
127         TypeData1 min(const TypeData1& a, const TypeData2& b)
128         {
129         return (a<b)?a:b;
130 }*/
131
132 //      template<typename TypeData> TypeData abs(const TypeData& a) {return (a<0)?-a:a;}        // include <cmath> instead
133 //      template<typename TypeData> TypeData abs(TypeData a)            {return (a<0)?-a:a;}    // include <cmath> instead
134
135 }
136
137 #endif
138