1 /* 2 Example implementation of the PADsynth basic algorithm 3 By: Nasca O. Paul, Tg. Mures, Romania 4 5 Ported to pure C by Paul Batchelor 6 7 This implementation and the algorithm are released under Public Domain 8 Feel free to use it into your projects or your products ;-) 9 10 This implementation is tested under GCC/Linux, but it's 11 very easy to port to other compiler/OS. 12 */ 13 14 #include <stdlib.h> 15 #include <math.h> 16 #include "soundpipe.h" sp_osc_create(sp_osc ** osc)17 18 #ifndef M_PI 19 #define M_PI 3.14159265358979323846 20 #endif 21 22 int sp_gen_padsynth(sp_data *sp, sp_ftbl *ps, sp_ftbl *amps, 23 SPFLOAT f, SPFLOAT bw) 24 { 25 26 int i, nh; 27 int N = (int) ps->size; 28 int number_harmonics = (int) amps->size; 29 SPFLOAT *A = amps->tbl; 30 SPFLOAT *smp = ps->tbl; 31 32 SPFLOAT *freq_amp = malloc((N / 2) * sizeof(SPFLOAT)); 33 SPFLOAT *freq_phase = malloc((N / 2) * sizeof(SPFLOAT)); 34 35 for (i=0;i<N/2;i++) freq_amp[i]=0.0; 36 37 for (nh=1;nh<number_harmonics;nh++) { 38 SPFLOAT bw_Hz; 39 SPFLOAT bwi; 40 SPFLOAT fi; 41 bw_Hz = (pow(2.0, bw/1200.0) - 1.0) * f * nh; 42 bwi = bw_Hz/(2.0*ps->size); 43 fi = f*nh/ps->size; 44 for (i = 0; i < N/2 ; i++) { 45 SPFLOAT hprofile; 46 hprofile = sp_padsynth_profile((i / (SPFLOAT) N) - fi, bwi); 47 freq_amp[i] += hprofile*A[nh]; 48 } 49 } 50 51 for (i=0;i<N/2;i++) { 52 freq_phase[i]= (sp_rand(sp) / (SP_RANDMAX + 1.0)) * 2.0 * M_PI; 53 }; 54 55 sp_padsynth_ifft(N,freq_amp,freq_phase,smp); 56 sp_padsynth_normalize(N,smp); 57 58 free(freq_amp); 59 free(freq_phase); 60 return SP_OK; 61 } 62 63 /* This is the profile of one harmonic 64 In this case is a Gaussian distribution (e^(-x^2)) 65 The amplitude is divided by the bandwidth to ensure that the harmonic 66 keeps the same amplitude regardless of the bandwidth */ 67 68 SPFLOAT sp_padsynth_profile(SPFLOAT fi, SPFLOAT bwi) 69 { 70 SPFLOAT x =fi/bwi; 71 x *= x; 72 73 /* 74 * this avoids computing the e^(-x^2) where 75 * it's results are very close to zero 76 */ 77 if (x>14.71280603) return 0.0; 78 79 return exp(-x)/bwi; 80 } 81 82 int sp_padsynth_ifft(int N, SPFLOAT *freq_amp, 83 SPFLOAT *freq_phase, SPFLOAT *smp) 84 { 85 int i; 86 FFTwrapper *fft; 87 FFTwrapper_create(&fft, N); 88 FFTFREQS fftfreqs; 89 newFFTFREQS(&fftfreqs,N/2); 90 91 for (i=0; i<N/2; i++){ 92 fftfreqs.c[i]=freq_amp[i]*cos(freq_phase[i]); 93 fftfreqs.s[i]=freq_amp[i]*sin(freq_phase[i]); 94 }; 95 freqs2smps(fft, &fftfreqs,smp); 96 deleteFFTFREQS(&fftfreqs); 97 FFTwrapper_destroy(&fft); 98 return SP_OK; 99 } 100 101 /* 102 Simple normalization function. It normalizes the sound to 1/sqrt(2) 103 */ 104 105 int sp_padsynth_normalize(int N, SPFLOAT *smp) 106 { 107 int i; 108 SPFLOAT max=0.0; 109 for (i=0;i<N;i++) if (fabs(smp[i])>max) max=fabs(smp[i]); 110 if (max<1e-5) max=1e-5; 111 for (i=0;i<N;i++) smp[i]/=max*1.4142; 112 return SP_OK; 113 } 114