2 SuperCollider real time audio synthesis system
3 Copyright (c) 2002 James McCartney. All rights reserved.
4 http://www.audiosynth.com
6 This program is free software; you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation; either version 2 of the License, or
9 (at your option) any later version.
11 This program is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
16 You should have received a copy of the GNU General Public License
17 along with this program; if not, write to the Free Software
18 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
25 #include "SC_Constants.h"
30 #include <smmintrin.h>
33 #if _XOPEN_SOURCE >= 600 || _ISOC99_SOURCE /* c99 compliant compiler */
38 /* needed for msvc compiler at least */
39 template <typename float_type
>
40 inline float_type
trunc(float_type arg
)
42 if(arg
> float_type(0))
43 return std::floor(arg
);
45 return std::ceil(arg
);
49 ///////////////////////////////////////////////////////////////////////////////////////
51 inline bool sc_isnan(float x
)
53 #if defined(__cplusplus) && defined(__GNUC__) && _GLIBCXX_HAVE_ISNAN
56 return (!(x
>= 0.f
|| x
<= 0.f
));
60 ///////////////////////////////////////////////////////////////////////////////////////
62 // versions provided for float32 and float64
63 // did not supply template because do not want to instantiate for integers.
64 // all constants explicitly cast to prevent PowerPC frsp instruction generation.
66 ///////////////////////////////////////////////////////////////////////////////////////
68 // this is a function for preventing pathological math operations in ugens.
69 // can be used at the end of a block to fix any recirculating filter values.
70 inline float32
zapgremlins(float32 x
)
72 float32 absx
= std::abs(x
);
73 // very small numbers fail the first test, eliminating denormalized numbers
74 // (zero also fails the first test, but that is OK since it returns zero.)
75 // very large numbers fail the second test, eliminating infinities
76 // Not-a-Numbers fail both tests and are eliminated.
77 return (absx
> (float32
)1e-15 && absx
< (float32
)1e15
) ? x
: (float32
)0.;
80 inline float32
sc_log2(float32 x
)
83 return ::log2f(std::abs(x
));
85 return static_cast<float32
>(std::log(std::abs(x
)) * (float32
)rlog2
);
89 inline float32
sc_log10(float32 x
)
91 return std::log10(std::abs(x
));
94 inline float32
sc_midicps(float32 note
)
96 return (float32
)440. * std::pow((float32
)2., (note
- (float32
)69.) * (float32
)0.083333333333);
99 inline float32
sc_cpsmidi(float32 freq
)
101 return sc_log2(freq
* (float32
)0.0022727272727) * (float32
)12. + (float32
)69.;
104 inline float32
sc_midiratio(float32 midi
)
106 return std::pow((float32
)2. , midi
* (float32
)0.083333333333);
109 inline float32
sc_ratiomidi(float32 ratio
)
111 return (float32
)12. * sc_log2(ratio
);
114 inline float32
sc_octcps(float32 note
)
116 return (float32
)440. * std::pow((float32
)2., note
- (float32
)4.75);
119 inline float32
sc_cpsoct(float32 freq
)
121 return sc_log2(freq
* (float32
)0.0022727272727) + (float32
)4.75;
124 inline float32
sc_ampdb(float32 amp
)
126 return std::log10(amp
) * (float32
)20.;
129 inline float32
sc_dbamp(float32 db
)
131 return std::pow((float32
)10., db
* (float32
).05);
134 inline float32
sc_squared(float32 x
)
139 inline float32
sc_cubed(float32 x
)
144 inline float32
sc_sqrt(float32 x
)
146 return x
< (float32
)0. ? -sqrt(-x
) : sqrt(x
);
150 inline float32
sc_hanwindow(float32 x
)
152 if (x
< (float32
)0. || x
> (float32
)1.) return (float32
)0.;
153 return (float32
)0.5 - (float32
)0.5 * static_cast<float32
>(cos(x
* (float32
)twopi
));
156 inline float32
sc_welwindow(float32 x
)
158 if (x
< (float32
)0. || x
> (float32
)1.) return (float32
)0.;
159 return static_cast<float32
>(sin(x
* pi
));
162 inline float32
sc_triwindow(float32 x
)
164 if (x
< (float32
)0. || x
> (float32
)1.) return (float32
)0.;
165 if (x
< (float32
)0.5) return (float32
)2. * x
;
166 else return (float32
)-2. * x
+ (float32
)2.;
169 inline float32
sc_bitriwindow(float32 x
)
171 float32 ax
= (float32
)1. - std::abs(x
);
172 if (ax
<= (float32
)0.) return (float32
)0.;
176 inline float32
sc_rectwindow(float32 x
)
178 if (x
< (float32
)0. || x
> (float32
)1.) return (float32
)0.;
182 inline float32
sc_scurve(float32 x
)
184 if (x
<= (float32
)0.) return (float32
)0.;
185 if (x
>= (float32
)1.) return (float32
)1.;
186 return x
* x
* ((float32
)3. - (float32
)2. * x
);
189 inline float32
sc_scurve0(float32 x
)
191 // assumes that x is in range
192 return x
* x
* ((float32
)3. - (float32
)2. * x
);
195 inline float32
sc_ramp(float32 x
)
197 if (x
<= (float32
)0.) return (float32
)0.;
198 if (x
>= (float32
)1.) return (float32
)1.;
202 inline float32
sc_sign(float32 x
)
204 return x
< (float32
)0. ? (float32
)-1. : (x
> (float32
)0. ? (float32
)1.f
: (float32
)0.f
);
207 inline float32
sc_distort(float32 x
)
209 return x
/ ((float32
)1. + std::abs(x
));
212 inline float32
sc_distortneg(float32 x
)
215 return x
/((float32
)1. - x
);
220 inline float32
sc_softclip(float32 x
)
222 float32 absx
= std::abs(x
);
223 if (absx
<= (float32
)0.5) return x
;
224 else return (absx
- (float32
)0.25) / x
;
227 // Taylor expansion out to x**9/9! factored into multiply-adds
229 inline float32
taylorsin(float32 x
)
231 // valid range from -pi/2 to +3pi/2
232 x
= static_cast<float32
>((float32
)pi2
- std::abs(pi2
- x
));
234 return static_cast<float32
>(x
*(x2
*(x2
*(x2
*(x2
*(1.0/362880.0)
241 inline float32
sc_trunc(float32 x
)
254 inline float32
sc_ceil(float32 x
)
257 __m128 a
= _mm_set_ss(x
);
258 const int cntrl
= _MM_FROUND_TO_POS_INF
;
259 __m128 b
= _mm_round_ss(a
, a
, cntrl
);
260 return _mm_cvtss_f32(b
);
266 inline float32
sc_floor(float32 x
)
269 __m128 a
= _mm_set_ss(x
);
270 const int cntrl
= _MM_FROUND_TO_NEG_INF
;
271 __m128 b
= _mm_round_ss(a
, a
, cntrl
);
272 return _mm_cvtss_f32(b
);
274 return std::floor(x
);
278 inline float32
sc_frac(float32 x
)
280 return x
- sc_floor(x
);
284 inline float32
sc_lg3interp(float32 x1
, float32 a
, float32 b
, float32 c
, float32 d
)
286 // cubic lagrange interpolator
287 float32 x0
= x1
+ 1.f
;
288 float32 x2
= x1
- 1.f
;
289 float32 x3
= x1
- 2.f
;
291 float32 x03
= x0
* x3
* 0.5f
;
292 float32 x12
= x1
* x2
* 0.16666666666666667f
;
294 return x12
* (d
* x0
- a
* x3
) + x03
* (b
* x2
- c
* x1
);
297 inline float32
sc_CalcFeedback(float32 delaytime
, float32 decaytime
)
299 if (delaytime
== 0.f
|| decaytime
== 0.f
)
303 float32 absret
= static_cast<float32
>(exp(log001
* delaytime
/ sc_abs(decaytime
)));
304 float32 ret
= copysignf(absret
, decaytime
);
308 return static_cast<float32
>(exp(log001
* delaytime
/ decaytime
));
310 return static_cast<float32
>(-exp(log001
* delaytime
/ -decaytime
));
315 inline float32
sc_wrap1(float32 x
)
317 if (x
>= (float32
) 1.) return x
+ (float32
)-2.;
318 if (x
< (float32
)-1.) return x
+ (float32
) 2.;
322 inline float32
sc_fold1(float32 x
)
324 if (x
>= (float32
) 1.) return (float32
) 2. - x
;
325 if (x
< (float32
)-1.) return (float32
)-2. - x
;
330 ///////////////////////////////////////////////////////////////////////////////////////
332 inline float64
zapgremlins(float64 x
)
334 float64 absx
= std::abs(x
);
335 // very small numbers fail the first test, eliminating denormalized numbers
336 // (zero also fails the first test, but that is OK since it returns zero.)
337 // very large numbers fail the second test, eliminating infinities
338 // Not-a-Numbers fail both tests and are eliminated.
339 return (absx
> (float64
)1e-15 && absx
< (float64
)1e15
) ? x
: (float64
)0.;
342 inline float64
sc_log2(float64 x
)
345 return ::log2(std::abs(x
));
347 return std::log(std::abs(x
)) * rlog2
;
351 inline float64
sc_log10(float64 x
)
353 return std::log10(std::abs(x
));
356 inline float64
sc_midicps(float64 note
)
358 return (float64
)440. * std::pow((float64
)2., (note
- (float64
)69.) * (float64
)0.08333333333333333333333333);
361 inline float64
sc_cpsmidi(float64 freq
)
363 return sc_log2(freq
* (float64
)0.002272727272727272727272727) * (float64
)12. + (float64
)69.;
366 inline float64
sc_midiratio(float64 midi
)
368 return std::pow((float64
)2. , midi
* (float64
)0.083333333333);
371 inline float64
sc_ratiomidi(float64 ratio
)
373 return (float64
)12. * sc_log2(ratio
);
376 inline float64
sc_octcps(float64 note
)
378 return (float64
)440. * std::pow((float64
)2., note
- (float64
)4.75);
381 inline float64
sc_cpsoct(float64 freq
)
383 return sc_log2(freq
* (float64
)0.0022727272727) + (float64
)4.75;
386 inline float64
sc_ampdb(float64 amp
)
388 return std::log10(amp
) * (float64
)20.;
391 inline float64
sc_dbamp(float64 db
)
393 return std::pow((float64
)10., db
* (float64
).05);
396 inline float64
sc_squared(float64 x
)
401 inline float64
sc_cubed(float64 x
)
406 inline float64
sc_sqrt(float64 x
)
408 return x
< (float64
)0. ? -sqrt(-x
) : sqrt(x
);
411 inline float64
sc_hanwindow(float64 x
)
413 if (x
< (float64
)0. || x
> (float64
)1.) return (float64
)0.;
414 return (float64
)0.5 - (float64
)0.5 * cos(x
* twopi
);
417 inline float64
sc_welwindow(float64 x
)
419 if (x
< (float64
)0. || x
> (float64
)1.) return (float64
)0.;
423 inline float64
sc_triwindow(float64 x
)
425 if (x
< (float64
)0. || x
> (float64
)1.) return (float64
)0.;
426 if (x
< (float64
)0.5) return (float64
)2. * x
;
427 else return (float64
)-2. * x
+ (float64
)2.;
430 inline float64
sc_bitriwindow(float64 x
)
432 float64 ax
= std::abs(x
);
433 if (ax
> (float64
)1.) return (float64
)0.;
434 return (float64
)1. - ax
;
437 inline float64
sc_rectwindow(float64 x
)
439 if (x
< (float64
)0. || x
> (float64
)1.) return (float64
)0.;
443 inline float64
sc_scurve(float64 x
)
445 if (x
<= (float64
)0.) return (float64
)0.;
446 if (x
>= (float64
)1.) return (float64
)1.;
447 return x
* x
* ((float64
)3. - (float64
)2. * x
);
450 inline float64
sc_scurve0(float64 x
)
452 // assumes that x is in range
453 return x
* x
* ((float64
)3. - (float64
)2. * x
);
456 inline float64
sc_ramp(float64 x
)
458 if (x
<= (float64
)0.) return (float64
)0.;
459 if (x
>= (float64
)1.) return (float64
)1.;
463 inline float64
sc_sign(float64 x
)
465 return x
< (float64
)0. ? (float64
)-1. : (x
> (float64
)0. ? (float64
)1.f
: (float64
)0.f
);
468 inline float64
sc_distort(float64 x
)
470 return x
/ ((float64
)1. + std::abs(x
));
473 inline float64
sc_distortneg(float64 x
)
476 return x
/((float64
)1. - x
);
482 inline float64
sc_softclip(float64 x
)
484 float64 absx
= std::abs(x
);
485 if (absx
<= (float64
)0.5) return x
;
486 else return (absx
- (float64
)0.25) / x
;
489 // Taylor expansion out to x**9/9! factored into multiply-adds
491 inline float64
taylorsin(float64 x
)
493 x
= pi2
- std::abs(pi2
- x
);
495 return x
*(x2
*(x2
*(x2
*(x2
*(1.0/362880.0)
502 inline float64
sc_trunc(float64 x
)
510 inline float64
sc_ceil(float64 x
)
513 __m128d a
= _mm_set_sd(x
);
514 const int cntrl
= _MM_FROUND_TO_POS_INF
;
515 __m128d b
= _mm_round_sd(a
, a
, cntrl
);
516 return _mm_cvtsd_f64(b
);
522 inline float64
sc_floor(float64 x
)
525 __m128d a
= _mm_set_sd(x
);
526 const int cntrl
= _MM_FROUND_TO_NEG_INF
;
527 __m128d b
= _mm_round_sd(a
, a
, cntrl
);
528 return _mm_cvtsd_f64(b
);
530 return std::floor(x
);
534 inline float64
sc_frac(float64 x
)
536 return x
- sc_floor(x
);
539 inline float64
sc_wrap1(float64 x
)
541 if (x
>= (float64
) 1.) return x
+ (float64
)-2.;
542 if (x
< (float64
)-1.) return x
+ (float64
) 2.;
546 inline float64
sc_fold1(float64 x
)
548 if (x
>= (float64
) 1.) return (float64
) 2. - x
;
549 if (x
< (float64
)-1.) return (float64
)-2. - x
;
553 inline int32
sc_grayCode(int32 x
)
560 ///////////////////////////////////////////////////////////////////////////////////////