src/common/RTMath.h

/***************************************************************************
 *                                                                         *
 *   LinuxSampler - modular, streaming capable sampler                     *
 *                                                                         *
 *   Copyright (C) 2003, 2004 by Benno Senoner and Christian Schoenebeck   *
 *                                                                         *
 *   This program is free software; you can redistribute it and/or modify  *
 *   it under the terms of the GNU General Public License as published by  *
 *   the Free Software Foundation; either version 2 of the License, or     *
 *   (at your option) any later version.                                   *
 *                                                                         *
 *   This program is distributed in the hope that it will be useful,       *
 *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
 *   GNU General Public License for more details.                          *
 *                                                                         *
 *   You should have received a copy of the GNU General Public License     *
 *   along with this program; if not, write to the Free Software           *
 *   Foundation, Inc., 59 Temple Place, Suite 330, Boston,                 *
 *   MA  02111-1307  USA                                                   *
 ***************************************************************************/

#ifndef __RT_MATH_H__
#define __RT_MATH_H__

#include <math.h>
#include <stdint.h>
#include "global.h"

/// Needed for calculating frequency ratio used to pitch a sample
#define TWELVEHUNDREDTH_ROOT_OF_TWO     1.000577789506555

enum implementation_t {
    CPP
    #if ARCH_X86
    ,ASM_X86_MMX_SSE
    #endif // ARCH_X86
};

class RTMathBase {
    public:
        /**
         * Highly accurate time stamp.
         */
        typedef uint32_t time_stamp_t;

        /**
         * We read the processor's cycle count register as a reference
         * for the real time. These are of course only abstract values
         * with arbitrary time entity, but that's not a problem as long
         * as we calculate relatively.
         */
        static time_stamp_t CreateTimeStamp();

        /**
         * Calculates the frequency ratio for a pitch value given in cents
         * (assuming equal tempered scale of course, divided into 12
         * semitones per octave and 100 cents per semitone).
         *
         * Note: MAX_PITCH (defined in global.h) has to be defined to an
         * appropriate value, otherwise the behavior of this function is
         * undefined, but most probably if MAX_PITCH is too small, the
         * application will crash due to segmentation fault here.
         *
         * @param cents - pitch value in cents (+1200 cents means +1 octave)
         * @returns  frequency ratio (e.g. +2.0 for +1 octave)
         */
        inline static double CentsToFreqRatio(double Cents) {
            int   index_int   = (int) (Cents);      // integer index
            float index_fract = Cents - index_int;  // fractional part of index
            return pCentsToFreqTable[index_int] + index_fract * (pCentsToFreqTable[index_int+1] - pCentsToFreqTable[index_int]);
        }

    private:
        static float  CentsToFreqTable[MAX_PITCH * 1200 * 2 + 1];
        static float* pCentsToFreqTable;

        static float* InitCentsToFreqTable();
};

/** Real Time Math
 *
 * Math functions for real time operation.
 */
template<implementation_t IMPL = CPP>
class __RTMath : public RTMathBase {
    public:
        // conversion using truncate
        inline static int Int(const float a) {
            switch (IMPL) {
                case CPP: {
                    return (int) a;
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    int ret;
                    asm (
                        "cvttss2si %1, %0  # convert to int\n\t"
                        : "=r" (ret)
                        : "m" (a)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        //for doubles and everything else except floats
        template<class T_a> inline static int Int(const T_a a) {
            return (int) a;
        }

        inline static float Float(const int a) {
            switch (IMPL) {
                case CPP: {
                    return (float) a;
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "cvtsi2ss %1, %%xmm0  # convert to float\n\t"
                        "movss    %%xmm0,%0   # output\n\t"
                        : "=m" (ret)
                        : "r" (a)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

#if 0
        //for everything except ints
        template<class T_a> inline static float Float(T_a a) {
            return (float) a;
        }
#endif

        inline static float Sum(const float& a, const float& b) {
            switch (IMPL) {
                case CPP: {
                    return (a + b);
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "addss    %2, %%xmm0  # a + b\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Sum(const T_a a, const T_b b) {
            return (a + b);
        }

        inline static float Sub(const float& a, const float& b) {
            switch (IMPL) {
                case CPP: {
                    return (a - b);
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "subss    %2, %%xmm0  # a - b\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Sub(const T_a a, const T_b b) {
            return (a - b);
        }

        inline static float Mul(const float a, const float b) {
            switch (IMPL) {
                case CPP: {
                    return (a * b);
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "mulss    %2, %%xmm0  # a * b\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Mul(const T_a a, const T_b b) {
            return (a * b);
        }

        inline static float Div(const float a, const float b) {
            switch (IMPL) {
                case CPP: {
                    return (a / b);
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "divss    %2, %%xmm0  # a / b\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Div(const T_a a, const T_b b) {
            return (a / b);
        }

        inline static float Min(const float a, const float b) {
            switch (IMPL) {
                case CPP: {
                    return (b < a) ? b : a;
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "minss    %2, %%xmm0  # Minimum(a, b)\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Min(const T_a a, const T_b b) {
            return (b < a) ? b : a;
        }

        inline static float Max(const float a, const float b) {
            switch (IMPL) {
                case CPP: {
                    return (b > a) ? b : a;
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "maxss    %2, %%xmm0  # Maximum(a, b)\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }

        template<class T_a, class T_b> inline static T_a Max(const T_a a, const T_b b) {
            return (b > a) ? b : a;
        }

        inline static float Fmodf(const float &a, const float &b) {
            switch (IMPL) {
                case CPP: {
                    return fmodf(a, b);
                }
                #if ARCH_X86
                case ASM_X86_MMX_SSE: {
                    float ret;
                    asm (
                        "movss    %1, %%xmm0  # load a\n\t"
                        "movss    %2, %%xmm1  # load b\n\t"
                        "movss    %%xmm0,%%xmm2\n\t"
                        "divss    %%xmm1, %%xmm2  # xmm2 = a / b\n\t"
                        "cvttss2si %%xmm2, %%ecx  #convert to int\n\t"
                        "cvtsi2ss %%ecx, %%xmm2  #convert back to float\n\t"
                        "mulss    %%xmm1, %%xmm2  # xmm2 = b * int(a/b)\n\t"
                        "subss    %%xmm2, %%xmm0  #sub a\n\t"
                        "movss    %%xmm0, %0  # output\n\t"
                        : "=m" (ret)
                        : "m" (a), "m" (b)
                        : "%ecx"
                    );
                    return ret;
                }
                #endif // ARCH_X86
            }
        }
};

/// convenience typedef for using the default implementation (which is CPP)
typedef __RTMath<> RTMath;

#endif // __RT_MATH_H__
1	schoenebeck	53	/***************************************************************************
2			* *
3			* LinuxSampler - modular, streaming capable sampler *
4			* *
5	schoenebeck	56	* Copyright (C) 2003, 2004 by Benno Senoner and Christian Schoenebeck *
6	schoenebeck	53	* *
7			* This program is free software; you can redistribute it and/or modify *
8			* it under the terms of the GNU General Public License as published by *
9			* the Free Software Foundation; either version 2 of the License, or *
10			* (at your option) any later version. *
11			* *
12			* This program is distributed in the hope that it will be useful, *
13			* but WITHOUT ANY WARRANTY; without even the implied warranty of *
14			* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15			* GNU General Public License for more details. *
16			* *
17			* You should have received a copy of the GNU General Public License *
18			* along with this program; if not, write to the Free Software *
19			* Foundation, Inc., 59 Temple Place, Suite 330, Boston, *
20			* MA 02111-1307 USA *
21			***************************************************************************/
22
23			#ifndef __RT_MATH_H__
24			#define __RT_MATH_H__
25
26			#include <math.h>
27	schoenebeck	328	#include <stdint.h>
28	schoenebeck	53	#include "global.h"
29
30			/// Needed for calculating frequency ratio used to pitch a sample
31			#define TWELVEHUNDREDTH_ROOT_OF_TWO 1.000577789506555
32
33	schoenebeck	319	enum implementation_t {
34	schoenebeck	328	CPP
35			#if ARCH_X86
36			,ASM_X86_MMX_SSE
37			#endif // ARCH_X86
38	schoenebeck	319	};
39
40			class RTMathBase {
41	schoenebeck	53	public:
42			/**
43	schoenebeck	328	* Highly accurate time stamp.
44			*/
45			typedef uint32_t time_stamp_t;
46
47			/**
48			* We read the processor's cycle count register as a reference
49			* for the real time. These are of course only abstract values
50			* with arbitrary time entity, but that's not a problem as long
51			* as we calculate relatively.
52			*/
53			static time_stamp_t CreateTimeStamp();
54
55			/**
56	schoenebeck	53	* Calculates the frequency ratio for a pitch value given in cents
57			* (assuming equal tempered scale of course, divided into 12
58			* semitones per octave and 100 cents per semitone).
59			*
60			* Note: MAX_PITCH (defined in global.h) has to be defined to an
61			* appropriate value, otherwise the behavior of this function is
62			* undefined, but most probably if MAX_PITCH is too small, the
63			* application will crash due to segmentation fault here.
64			*
65			* @param cents - pitch value in cents (+1200 cents means +1 octave)
66			* @returns frequency ratio (e.g. +2.0 for +1 octave)
67			*/
68			inline static double CentsToFreqRatio(double Cents) {
69	schoenebeck	319	int index_int = (int) (Cents); // integer index
70	schoenebeck	53	float index_fract = Cents - index_int; // fractional part of index
71			return pCentsToFreqTable[index_int] + index_fract * (pCentsToFreqTable[index_int+1] - pCentsToFreqTable[index_int]);
72			}
73
74	schoenebeck	319	private:
75			static float CentsToFreqTable[MAX_PITCH * 1200 * 2 + 1];
76			static float* pCentsToFreqTable;
77
78			static float* InitCentsToFreqTable();
79			};
80
81			/** Real Time Math
82			*
83			* Math functions for real time operation.
84			*/
85			template<implementation_t IMPL = CPP>
86			class __RTMath : public RTMathBase {
87			public:
88			// conversion using truncate
89			inline static int Int(const float a) {
90			switch (IMPL) {
91			case CPP: {
92			return (int) a;
93			}
94	schoenebeck	328	#if ARCH_X86
95	schoenebeck	319	case ASM_X86_MMX_SSE: {
96			int ret;
97			asm (
98			"cvttss2si %1, %0 # convert to int\n\t"
99			: "=r" (ret)
100			: "m" (a)
101			);
102			return ret;
103			}
104	schoenebeck	328	#endif // ARCH_X86
105	schoenebeck	319	}
106			}
107
108			//for doubles and everything else except floats
109			template<class T_a> inline static int Int(const T_a a) {
110			return (int) a;
111			}
112
113			inline static float Float(const int a) {
114			switch (IMPL) {
115			case CPP: {
116			return (float) a;
117			}
118	schoenebeck	328	#if ARCH_X86
119	schoenebeck	319	case ASM_X86_MMX_SSE: {
120			float ret;
121			asm (
122			"cvtsi2ss %1, %%xmm0 # convert to float\n\t"
123			"movss %%xmm0,%0 # output\n\t"
124			: "=m" (ret)
125			: "r" (a)
126			);
127			return ret;
128			}
129	schoenebeck	328	#endif // ARCH_X86
130	schoenebeck	319	}
131			}
132
133			#if 0
134			//for everything except ints
135			template<class T_a> inline static float Float(T_a a) {
136			return (float) a;
137			}
138			#endif
139
140			inline static float Sum(const float& a, const float& b) {
141			switch (IMPL) {
142			case CPP: {
143			return (a + b);
144			}
145	schoenebeck	328	#if ARCH_X86
146	schoenebeck	319	case ASM_X86_MMX_SSE: {
147			float ret;
148			asm (
149			"movss %1, %%xmm0 # load a\n\t"
150			"addss %2, %%xmm0 # a + b\n\t"
151			"movss %%xmm0, %0 # output\n\t"
152			: "=m" (ret)
153			: "m" (a), "m" (b)
154			);
155			return ret;
156			}
157	schoenebeck	328	#endif // ARCH_X86
158	schoenebeck	319	}
159			}
160
161			template<class T_a, class T_b> inline static T_a Sum(const T_a a, const T_b b) {
162			return (a + b);
163			}
164
165			inline static float Sub(const float& a, const float& b) {
166			switch (IMPL) {
167			case CPP: {
168			return (a - b);
169			}
170	schoenebeck	328	#if ARCH_X86
171	schoenebeck	319	case ASM_X86_MMX_SSE: {
172			float ret;
173			asm (
174			"movss %1, %%xmm0 # load a\n\t"
175			"subss %2, %%xmm0 # a - b\n\t"
176			"movss %%xmm0, %0 # output\n\t"
177			: "=m" (ret)
178			: "m" (a), "m" (b)
179			);
180			return ret;
181			}
182	schoenebeck	328	#endif // ARCH_X86
183	schoenebeck	319	}
184			}
185
186			template<class T_a, class T_b> inline static T_a Sub(const T_a a, const T_b b) {
187			return (a - b);
188			}
189
190			inline static float Mul(const float a, const float b) {
191			switch (IMPL) {
192			case CPP: {
193			return (a * b);
194			}
195	schoenebeck	328	#if ARCH_X86
196	schoenebeck	319	case ASM_X86_MMX_SSE: {
197			float ret;
198			asm (
199			"movss %1, %%xmm0 # load a\n\t"
200			"mulss %2, %%xmm0 # a * b\n\t"
201			"movss %%xmm0, %0 # output\n\t"
202			: "=m" (ret)
203			: "m" (a), "m" (b)
204			);
205			return ret;
206			}
207	schoenebeck	328	#endif // ARCH_X86
208	schoenebeck	319	}
209			}
210
211			template<class T_a, class T_b> inline static T_a Mul(const T_a a, const T_b b) {
212			return (a * b);
213			}
214
215			inline static float Div(const float a, const float b) {
216			switch (IMPL) {
217			case CPP: {
218			return (a / b);
219			}
220	schoenebeck	328	#if ARCH_X86
221	schoenebeck	319	case ASM_X86_MMX_SSE: {
222			float ret;
223			asm (
224			"movss %1, %%xmm0 # load a\n\t"
225			"divss %2, %%xmm0 # a / b\n\t"
226			"movss %%xmm0, %0 # output\n\t"
227			: "=m" (ret)
228			: "m" (a), "m" (b)
229			);
230			return ret;
231			}
232	schoenebeck	328	#endif // ARCH_X86
233	schoenebeck	319	}
234			}
235
236			template<class T_a, class T_b> inline static T_a Div(const T_a a, const T_b b) {
237			return (a / b);
238			}
239
240			inline static float Min(const float a, const float b) {
241			switch (IMPL) {
242			case CPP: {
243			return (b < a) ? b : a;
244			}
245	schoenebeck	328	#if ARCH_X86
246	schoenebeck	319	case ASM_X86_MMX_SSE: {
247			float ret;
248			asm (
249			"movss %1, %%xmm0 # load a\n\t"
250			"minss %2, %%xmm0 # Minimum(a, b)\n\t"
251			"movss %%xmm0, %0 # output\n\t"
252			: "=m" (ret)
253			: "m" (a), "m" (b)
254			);
255			return ret;
256			}
257	schoenebeck	328	#endif // ARCH_X86
258	schoenebeck	319	}
259			}
260
261			template<class T_a, class T_b> inline static T_a Min(const T_a a, const T_b b) {
262	schoenebeck	53	return (b < a) ? b : a;
263			}
264
265	schoenebeck	319	inline static float Max(const float a, const float b) {
266			switch (IMPL) {
267			case CPP: {
268			return (b > a) ? b : a;
269			}
270	schoenebeck	328	#if ARCH_X86
271	schoenebeck	319	case ASM_X86_MMX_SSE: {
272			float ret;
273			asm (
274			"movss %1, %%xmm0 # load a\n\t"
275			"maxss %2, %%xmm0 # Maximum(a, b)\n\t"
276			"movss %%xmm0, %0 # output\n\t"
277			: "=m" (ret)
278			: "m" (a), "m" (b)
279			);
280			return ret;
281			}
282	schoenebeck	328	#endif // ARCH_X86
283	schoenebeck	319	}
284			}
285
286			template<class T_a, class T_b> inline static T_a Max(const T_a a, const T_b b) {
287	schoenebeck	53	return (b > a) ? b : a;
288			}
289
290	schoenebeck	319	inline static float Fmodf(const float &a, const float &b) {
291			switch (IMPL) {
292			case CPP: {
293			return fmodf(a, b);
294			}
295	schoenebeck	328	#if ARCH_X86
296	schoenebeck	319	case ASM_X86_MMX_SSE: {
297			float ret;
298			asm (
299			"movss %1, %%xmm0 # load a\n\t"
300			"movss %2, %%xmm1 # load b\n\t"
301			"movss %%xmm0,%%xmm2\n\t"
302			"divss %%xmm1, %%xmm2 # xmm2 = a / b\n\t"
303			"cvttss2si %%xmm2, %%ecx #convert to int\n\t"
304			"cvtsi2ss %%ecx, %%xmm2 #convert back to float\n\t"
305			"mulss %%xmm1, %%xmm2 # xmm2 = b * int(a/b)\n\t"
306			"subss %%xmm2, %%xmm0 #sub a\n\t"
307			"movss %%xmm0, %0 # output\n\t"
308			: "=m" (ret)
309			: "m" (a), "m" (b)
310			: "%ecx"
311			);
312			return ret;
313			}
314	schoenebeck	328	#endif // ARCH_X86
315	schoenebeck	319	}
316			}
317	schoenebeck	53	};
318
319	schoenebeck	319	/// convenience typedef for using the default implementation (which is CPP)
320			typedef __RTMath<> RTMath;
321
322	schoenebeck	53	#endif // __RT_MATH_H__