3 #ifndef SPECTMORPH_MATH_HH 4 #define SPECTMORPH_MATH_HH 10 #include <bse/bseieee754.hh> 12 #include <xmmintrin.h> 45 template<
class Iterator,
int MODE>
47 internal_fast_vector_sin (
const VectorSinParams& params, Iterator begin, Iterator end)
49 g_return_if_fail (params.
mix_freq > 0 && params.
freq > 0 && params.
phase > -99 && params.
mag > 0);
51 const double phase_inc = params.
freq / params.
mix_freq * 2 * M_PI;
52 const double inc_re = cos (phase_inc);
53 const double inc_im = sin (phase_inc);
59 sincos (params.
phase, &state_im, &state_re);
60 state_re *= params.
mag;
61 state_im *= params.
mag;
63 for (Iterator x = begin; x != end; x++)
69 if ((n++ & 255) == 255)
71 sincos (phase_inc * n + params.
phase, &state_im, &state_re);
72 state_re *= params.
mag;
73 state_im *= params.
mag;
81 const double re = state_re * inc_re - state_im * inc_im;
82 const double im = state_re * inc_im + state_im * inc_re;
89 template<
class Iterator,
int MODE>
91 internal_fast_vector_sincos (
const VectorSinParams& params, Iterator sin_begin, Iterator sin_end, Iterator cos_begin)
93 g_return_if_fail (params.
mix_freq > 0 && params.
freq > 0 && params.
phase > -99 && params.
mag > 0);
95 const double phase_inc = params.
freq / params.
mix_freq * 2 * M_PI;
96 const double inc_re = cos (phase_inc);
97 const double inc_im = sin (phase_inc);
103 sincos (params.
phase, &state_im, &state_re);
104 state_re *= params.
mag;
105 state_im *= params.
mag;
107 for (Iterator x = sin_begin, y = cos_begin; x != sin_end; x++, y++)
119 if ((n++ & 255) == 255)
121 sincos (phase_inc * n + params.
phase, &state_im, &state_re);
122 state_re *= params.
mag;
123 state_im *= params.
mag;
131 const double re = state_re * inc_re - state_im * inc_im;
132 const double im = state_re * inc_im + state_im * inc_re;
139 template<
class Iterator>
141 fast_vector_sin (
const VectorSinParams& params, Iterator sin_begin, Iterator sin_end)
145 internal_fast_vector_sin<Iterator, VectorSinParams::ADD> (params, sin_begin, sin_end);
149 internal_fast_vector_sin<Iterator, VectorSinParams::REPLACE> (params, sin_begin, sin_end);
153 g_assert_not_reached();
157 template<
class Iterator>
159 fast_vector_sincos (
const VectorSinParams& params, Iterator sin_begin, Iterator sin_end, Iterator cos_begin)
163 internal_fast_vector_sincos<Iterator, VectorSinParams::ADD> (params, sin_begin, sin_end, cos_begin);
167 internal_fast_vector_sincos<Iterator, VectorSinParams::REPLACE> (params, sin_begin, sin_end, cos_begin);
171 g_assert_not_reached();
187 template<
bool NEED_COS,
int MODE>
189 internal_fast_vector_sincosf (
const VectorSinParams& params,
float *sin_begin,
float *sin_end,
float *cos_begin)
192 g_return_if_fail (params.
mix_freq > 0 && params.
freq > 0 && params.
phase > -99 && params.
mag > 0);
194 const int TABLE_SIZE = 32;
196 const double phase_inc = params.
freq / params.
mix_freq * 2 * M_PI;
197 const double inc_re16 = cos (phase_inc * TABLE_SIZE * 4);
198 const double inc_im16 = sin (phase_inc * TABLE_SIZE * 4);
204 sincos (params.
phase, &state_im, &state_re);
205 state_re *= params.
mag;
206 state_im *= params.
mag;
208 F4Vector incf_re[TABLE_SIZE];
209 F4Vector incf_im[TABLE_SIZE];
213 table_params.
phase = 0;
214 table_params.
mag = 1;
216 fast_vector_sincos (table_params, incf_im[0].f, incf_im[0].f + (TABLE_SIZE * 4), incf_re[0].f);
219 int todo = sin_end - sin_begin;
220 while (todo >= 4 * TABLE_SIZE)
224 sf_re.f[0] = state_re;
225 sf_re.f[1] = state_re;
226 sf_re.f[2] = state_re;
227 sf_re.f[3] = state_re;
228 sf_im.f[0] = state_im;
229 sf_im.f[1] = state_im;
230 sf_im.f[2] = state_im;
231 sf_im.f[3] = state_im;
239 F4Vector *new_im =
reinterpret_cast<F4Vector *
> (sin_begin + n);
240 F4Vector *new_re =
reinterpret_cast<F4Vector *
> (cos_begin + n);
241 for (
int k = 0; k < TABLE_SIZE; k++)
247 new_re[k].v = _mm_add_ps (new_re[k].v, _mm_sub_ps (_mm_mul_ps (sf_re.v, incf_re[k].v),
248 _mm_mul_ps (sf_im.v, incf_im[k].v)));
250 new_im[k].v = _mm_add_ps (new_im[k].v, _mm_add_ps (_mm_mul_ps (sf_re.v, incf_im[k].v),
251 _mm_mul_ps (sf_im.v, incf_re[k].v)));
257 new_re[k].v = _mm_sub_ps (_mm_mul_ps (sf_re.v, incf_re[k].v),
258 _mm_mul_ps (sf_im.v, incf_im[k].v));
260 new_im[k].v = _mm_add_ps (_mm_mul_ps (sf_re.v, incf_im[k].v),
261 _mm_mul_ps (sf_im.v, incf_re[k].v));
271 const double re = state_re * inc_re16 - state_im * inc_im16;
272 const double im = state_re * inc_im16 + state_im * inc_re16;
276 todo -= 4 * TABLE_SIZE;
281 rest_params.
phase += n * phase_inc;
283 fast_vector_sincos (rest_params, sin_begin + n, sin_end, cos_begin + n);
285 fast_vector_sin (rest_params, sin_begin + n, sin_end);
288 fast_vector_sincos (params, sin_begin, sin_end, cos_begin);
290 fast_vector_sin (params, sin_begin, sin_end);
295 fast_vector_sincosf (
const VectorSinParams& params,
float *sin_begin,
float *sin_end,
float *cos_begin)
299 internal_fast_vector_sincosf<true, VectorSinParams::ADD> (params, sin_begin, sin_end, cos_begin);
303 internal_fast_vector_sincosf<true, VectorSinParams::REPLACE> (params, sin_begin, sin_end, cos_begin);
307 g_assert_not_reached();
312 fast_vector_sinf (
const VectorSinParams& params,
float *sin_begin,
float *sin_end)
316 internal_fast_vector_sincosf<false, VectorSinParams::ADD> (params, sin_begin, sin_end, NULL);
320 internal_fast_vector_sincosf<false, VectorSinParams::REPLACE> (params, sin_begin, sin_end, NULL);
324 g_assert_not_reached();
329 zero_float_block (
size_t n_values,
float *values)
331 memset (values, 0, n_values *
sizeof (
float));
337 extern float *int_sincos_table;
339 return int_sincos_table[i];
345 extern float *int_sincos_table;
348 return int_sincos_table[i];
354 extern float *int_sincos_table;
356 int_sincos_table = (
float *) malloc (
sizeof (
float) * 256);
357 for (
int i = 0; i < 256; i++)
358 int_sincos_table[i] = sin (
double (i / 256.0) * 2 * M_PI);
362 window_blackman_harris_92 (
double x)
367 const double a0 = 0.35875, a1 = 0.48829, a2 = 0.14128, a3 = 0.01168;
369 return a0 + a1 * cos (M_PI * x) + a2 * cos (2.0 * M_PI * x) + a3 * cos (3.0 * M_PI * x);
372 #if defined (__i386__) && defined (__GNUC__) 374 sm_round_positive (
double d)
380 sm_round_positive (
float f)
386 sm_round_positive (
double d)
388 return int (d + 0.5);
392 sm_round_positive (
float f)
394 return int (f + 0.5);
400 static float idb2f_high[256];
401 static float idb2f_low[256];
403 static float ifreq2f_high[256];
404 static float ifreq2f_low[256];
407 #define SM_IDB_CONST_M96 uint16_t ((512 - 96) * 64) 409 int sm_factor2delta_idb (
double factor);
410 double sm_idb2factor_slow (uint16_t idb);
414 uint16_t sm_freq2ifreq (
double freq);
415 double sm_ifreq2freq_slow (uint16_t ifreq);
418 sm_idb2factor (uint16_t idb)
420 return MathTables::idb2f_high[idb >> 8] * MathTables::idb2f_low[idb & 0xff];
424 sm_ifreq2freq (uint16_t ifreq)
426 return MathTables::ifreq2f_high[ifreq >> 8] * MathTables::ifreq2f_low[ifreq & 0xff];
430 sm_factor2idb (
double factor)
438 const double db = 20 * log10 (std::max (factor, 1e-25));
440 return sm_round_positive (db * 64 + 512 * 64);
Definition: smmath.hh:398
double freq
the frequency of the sin (and cos) wave to be created
Definition: smmath.hh:25
add computed values to the values that are in the output array
Definition: smmath.hh:31
double phase
the start phase of the wave
Definition: smmath.hh:26
double mag
the magnitude (amplitude) of the wave
Definition: smmath.hh:27
enum SpectMorph::VectorSinParams::@1 mode
whether to overwrite or add output
parameter structure for the various optimized vector sine functions
Definition: smmath.hh:21
Definition: smaudio.hh:15
double mix_freq
the mix freq (sampling rate) of the sin (and cos) wave to be created
Definition: smmath.hh:24
replace values in the output array with computed values
Definition: smmath.hh:32