Initial import.
[Faustine.git] / dsp_files / music.lib
1 /************************************************************************
2 ************************************************************************
3 FAUST library file
4 Copyright (C) 2003-2012 GRAME, Centre National de Creation Musicale
5 ---------------------------------------------------------------------
6 This program is free software; you can redistribute it and/or modify
7 it under the terms of the GNU Lesser General Public License as
8 published by the Free Software Foundation; either version 2.1 of the
9 License, or (at your option) any later version.
10
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 Lesser General Public License for more details.
15
16 You should have received a copy of the GNU Lesser General Public
17 License along with the GNU C Library; if not, write to the Free
18 Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
19 02111-1307 USA.
20
21 EXCEPTION TO THE LGPL LICENSE : As a special exception, you may create a
22 larger FAUST program which directly or indirectly imports this library
23 file and still distribute the compiled code generated by the FAUST
24 compiler, or a modified version of this compiled code, under your own
25 copyright and license. This EXCEPTION TO THE LGPL LICENSE explicitly
26 grants you the right to freely choose the license for the resulting
27 compiled code. In particular the resulting compiled code has no obligation
28 to be LGPL or GPL. For example you are free to choose a commercial or
29 closed source license or any other license if you decide so.
30
31 ************************************************************************
32 ************************************************************************/
33
34 declare name "Music Library";
35 declare author "GRAME";
36 declare copyright "GRAME";
37 declare version "1.0";
38 declare license "LGPL with exception";
39
40 import("math.lib");
41
42 //-----------------------------------------------
43 // DELAY LINE
44 //-----------------------------------------------
45 frac(n) = n-int(n);
46 index(n) = &(n-1) ~ +(1); // n = 2**i
47 //delay(n,d,x) = rwtable(n, 0.0, index(n), x, (index(n)-int(d)) & (n-1));
48 delay(n,d,x) = x@(int(d)&(n-1));
49 fdelay(n,d,x) = delay(n,int(d),x)*(1 - frac(d)) + delay(n,int(d)+1,x)*frac(d);
50
51
52 delay1s(d) = delay(65536,d);
53 delay2s(d) = delay(131072,d);
54 delay5s(d) = delay(262144,d);
55 delay10s(d) = delay(524288,d);
56 delay21s(d) = delay(1048576,d);
57 delay43s(d) = delay(2097152,d);
58
59 fdelay1s(d) = fdelay(65536,d);
60 fdelay2s(d) = fdelay(131072,d);
61 fdelay5s(d) = fdelay(262144,d);
62 fdelay10s(d) = fdelay(524288,d);
63 fdelay21s(d) = fdelay(1048576,d);
64 fdelay43s(d) = fdelay(2097152,d);
65
66 millisec = SR/1000.0;
67
68 time1s = hslider("time", 0, 0, 1000, 0.1)*millisec;
69 time2s = hslider("time", 0, 0, 2000, 0.1)*millisec;
70 time5s = hslider("time", 0, 0, 5000, 0.1)*millisec;
71 time10s = hslider("time", 0, 0, 10000, 0.1)*millisec;
72 time21s = hslider("time", 0, 0, 21000, 0.1)*millisec;
73 time43s = hslider("time", 0, 0, 43000, 0.1)*millisec;
74
75
76 echo1s = vgroup("echo 1000", +~(delay(65536, int(hslider("millisecond", 0, 0, 1000, 0.10)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
77 echo2s = vgroup("echo 2000", +~(delay(131072, int(hslider("millisecond", 0, 0, 2000, 0.25)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
78 echo5s = vgroup("echo 5000", +~(delay(262144, int(hslider("millisecond", 0, 0, 5000, 0.50)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
79 echo10s = vgroup("echo 10000", +~(delay(524288, int(hslider("millisecond", 0, 0, 10000, 1.00)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
80 echo21s = vgroup("echo 21000", +~(delay(1048576, int(hslider("millisecond", 0, 0, 21000, 1.00)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
81 echo43s = vgroup("echo 43000", +~(delay(2097152, int(hslider("millisecond", 0, 0, 43000, 1.00)*millisec)-1) * (hslider("feedback", 0, 0, 100, 0.1)/100.0)));
82
83
84 //--------------------------sdelay(N,it,dt)----------------------------
85 // s(mooth)delay : a mono delay that doesn't click and doesn't
86 // transpose when the delay time is changed. It takes 4 input signals
87 // and produces a delayed output signal
88 //
89 // USAGE : ... : sdelay(N,it,dt) : ...
90 //
91 // Where :
92 // <N> = maximal delay in samples (must be a constant power of 2, for example 65536)
93 // <it> = interpolation time (in samples) for example 1024
94 // <dt> = delay time (in samples)
95 // < > = input signal we want to delay
96 //--------------------------------------------------------------------------
97
98 sdelay(N, it, dt) = ctrl(it,dt),_ : ddi(N)
99
100 with {
101
102 //---------------------------ddi(N,i,d0,d1)-------------------------------
103 // DDI (Double Delay with Interpolation) : the input signal is sent to two
104 // delay lines. The outputs of these delay lines are crossfaded with
105 // an interpolation stage. By acting on this interpolation value one
106 // can move smoothly from one delay to another. When <i> is 0 we can
107 // freely change the delay time <d1> of line 1, when it is 1 we can freely change
108 // the delay time <d0> of line 0.
109 //
110 // <N> = maximal delay in samples (must be a power of 2, for example 65536)
111 // <i> = interpolation value between 0 and 1 used to crossfade the outputs of the
112 // two delay lines (0.0: first delay line, 1.0: second delay line)
113 // <d0> = delay time of delay line 0 in samples between 0 and <N>-1
114 // <d1> = delay time of delay line 1 in samples between 0 and <N>-1
115 // < > = the input signal we want to delay
116 //-------------------------------------------------------------------------
117 ddi(N, i, d0, d1) = _ <: delay(N,d0), delay(N,d1) : interpolate(i);
118
119
120 //----------------------------ctrl(it,dt)------------------------------------
121 // Control logic for a Double Delay with Interpolation according to two
122 //
123 // USAGE : ctrl(it,dt)
124 // where :
125 // <it> an interpolation time (in samples, for example 256)
126 // <dt> a delay time (in samples)
127 //
128 // ctrl produces 3 outputs : an interpolation value <i> and two delay
129 // times <d0> and <d1>. These signals are used to control a ddi (Double Delay with Interpolation).
130 // The principle is to detect changes in the input delay time dt, then to
131 // change the delay time of the delay line currently unused and then by a
132 // smooth crossfade to remove the first delay line and activate the second one.
133 //
134 // The control logic has an internal state controlled by 4 elements
135 // <v> : the interpolation variation (0, 1/it, -1/it)
136 // <i> : the interpolation value (between 0 and 1)
137 // <d0>: the delay time of line 0
138 // <d1>: the delay time of line 1
139 //
140 // Please note that the last stage (!,_,_,_) cut <v> because it is only
141 // used internally.
142 //-------------------------------------------------------------------------
143 ctrl(it, dt) = \(v,ip,d0,d1).( (nv, nip, nd0, nd1)
144 with {
145
146 // interpolation variation
147 nv = if (v!=0.0, // if variation we are interpolating
148 if( (ip>0.0) & (ip<1.0), v , 0), // should we continue or not ?
149 if ((ip==0.0) & (dt!=d0), 1.0/it, // if true xfade from dl0 to dl1
150 if ((ip==1.0) & (dt!=d1), -1.0/it, // if true xfade from dl1 to dl0
151 0))); // nothing to change
152 // interpolation value
153 nip = ip+nv : min(1.0) : max(0.0);
154
155 // update delay time of line 0 if needed
156 nd0 = if ((ip >= 1.0) & (d1!=dt), dt, d0);
157
158 // update delay time of line 0 if needed
159 nd1 = if ((ip <= 0.0) & (d0!=dt), dt, d1);
160
161 } ) ~ (_,_,_,_) : (!,_,_,_);
162 };
163
164
165
166
167 //-----------------------------------------------
168 // Tempo, beats and pulses
169 //-----------------------------------------------
170
171 tempo(t) = (60*SR)/t; // tempo(t) -> samples
172
173 period(p) = %(int(p))~+(1); // signal en dent de scie de periode p
174 pulse(t) = period(t)==0; // pulse (10000...) de periode p
175 pulsen(n,t) = period(t)<n; // pulse (1110000...) de taille n et de periode p
176 beat(t) = pulse(tempo(t)); // pulse au tempo t
177
178
179
180 //-----------------------------------------------
181 // conversions between db and linear values
182 //-----------------------------------------------
183
184 db2linear(x) = pow(10, x/20.0);
185 linear2db(x) = 20*log10(x);
186
187
188
189
190
191 //===============================================
192 // Random and Noise generators
193 //===============================================
194
195
196 //-----------------------------------------------
197 // noise : Noise generator
198 //-----------------------------------------------
199
200 random = +(12345) ~ *(1103515245);
201 RANDMAX = 2147483647.0;
202
203 noise = random / RANDMAX;
204
205
206 //-----------------------------------------------
207 // Generates multiple decorrelated random numbers
208 // in parallel. Expects n>0.
209 //-----------------------------------------------
210
211 multirandom(n) = randomize(n) ~_
212 with {
213 randomize (1) = +(12345) : *(1103515245);
214 randomize (n) = randomize(1) <: randomize(n-1),_;
215 };
216
217
218 //-----------------------------------------------
219 // Generates multiple decorrelated noises
220 // in parallel. Expects n>0.
221 //-----------------------------------------------
222
223 multinoise(n) = multirandom(n) : par(i,n,/(RANDMAX))
224 with {
225 RANDMAX = 2147483647.0;
226 };
227
228
229 //------------------------------------------------
230
231 noises(N,i) = multinoise(N) : selector(i,N);
232
233
234 //-----------------------------------------------
235 // osc(freq) : Sinusoidal Oscillator
236 //-----------------------------------------------
237
238 tablesize = 1 << 16;
239 samplingfreq = SR;
240
241 time = (+(1)~_ ) - 1; // 0,1,2,3,...
242 sinwaveform = float(time)*(2.0*PI)/float(tablesize) : sin;
243
244 decimal(x) = x - floor(x);
245 phase(freq) = freq/float(samplingfreq) : (+ : decimal) ~ _ : *(float(tablesize));
246 osc(freq) = rdtable(tablesize, sinwaveform, int(phase(freq)) );
247 osci(freq) = s1 + d * (s2 - s1)
248 with {
249 i = int(phase(freq));
250 d = decimal(phase(freq));
251 s1 = rdtable(tablesize+1,sinwaveform,i);
252 s2 = rdtable(tablesize+1,sinwaveform,i+1);};
253
254
255 //-----------------------------------------------
256 // ADSR envelop
257 //-----------------------------------------------
258
259 // a,d,s,r = attack (sec), decay (sec), sustain (percentage of t), release (sec)
260 // t = trigger signal ( >0 for attack, then release is when t back to 0)
261
262 adsr(a,d,s,r,t) = env ~ (_,_) : (!,_) // the 2 'state' signals are fed back
263 with {
264 env (p2,y) =
265 (t>0) & (p2|(y>=1)), // p2 = decay-sustain phase
266 (y + p1*u - (p2&(y>s))*v*y - p3*w*y) // y = envelop signal
267 *((p3==0)|(y>=eps)) // cut off tails to prevent denormals
268 with {
269 p1 = (p2==0) & (t>0) & (y<1); // p1 = attack phase
270 p3 = (t<=0) & (y>0); // p3 = release phase
271 // #samples in attack, decay, release, must be >0
272 na = SR*a+(a==0.0); nd = SR*d+(d==0.0); nr = SR*r+(r==0.0);
273 // correct zero sustain level
274 z = s+(s==0.0)*db2linear(-60);
275 // attack, decay and (-60dB) release rates
276 u = 1/na; v = 1-pow(z, 1/nd); w = 1-1/pow(z*db2linear(60), 1/nr);
277 // values below this threshold are considered zero in the release phase
278 eps = db2linear(-120);
279 };
280 };
281
282
283 //-----------------------------------------------
284 // Spatialisation
285 //-----------------------------------------------
286
287 panner(c) = _ <: *(1-c), *(c);
288
289 bus2 = _,_;
290 bus3 = _,_,_;
291 bus4 = _,_,_,_;
292 bus5 = _,_,_,_,_;
293 bus6 = _,_,_,_,_,_;
294 bus7 = _,_,_,_,_,_,_;
295 bus8 = _,_,_,_,_,_,_,_;
296
297 gain2(g) = *(g),*(g);
298 gain3(g) = *(g),*(g),*(g);
299 gain4(g) = *(g),*(g),*(g),*(g);
300 gain5(g) = *(g),*(g),*(g),*(g),*(g);
301 gain6(g) = *(g),*(g),*(g),*(g),*(g),*(g);
302 gain7(g) = *(g),*(g),*(g),*(g),*(g),*(g),*(g);
303 gain8(g) = *(g),*(g),*(g),*(g),*(g),*(g),*(g),*(g);
304
305
306 //------------------------------------------------------
307 //
308 // GMEM SPAT
309 // n-outputs spatializer
310 // implementation of L. Pottier
311 //
312 //------------------------------------------------------
313 //
314 // n = number of outputs
315 // r = rotation (between 0 et 1)
316 // d = distance of the source (between 0 et 1)
317 //
318 //------------------------------------------------------
319 spat(n,a,d) = _ <: par(i, n, *( scaler(i, n, a, d) : smooth(0.9999) ))
320 with {
321 scaler(i,n,a,d) = (d/2.0+0.5)
322 * sqrt( max(0.0, 1.0 - abs(fmod(a+0.5+float(n-i)/n, 1.0) - 0.5) * n * d) );
323 smooth(c) = *(1-c) : +~*(c);
324 };
325
326
327
328 //--------------- Second Order Generic Transfert Function -------------------------
329 // TF2(b0,b1,b2,a1,a2)
330 //
331 //---------------------------------------------------------------------------------
332
333 TF2(b0,b1,b2,a1,a2) = sub ~ conv2(a1,a2) : conv3(b0,b1,b2)
334 with {
335 conv3(k0,k1,k2,x) = k0*x + k1*x' + k2*x'';
336 conv2(k0,k1,x) = k0*x + k1*x';
337 sub(x,y) = y-x;
338 };
339
340
341 /*************************** Break Point Functions ***************************
342
343 bpf is an environment (a group of related definitions) tha can be used to
344 create break-point functions. It contains three functions :
345 - start(x,y) to start a break-point function
346 - end(x,y) to end a break-point function
347 - point(x,y) to add intermediate points to a break-point function
348
349 A minimal break-point function must contain at least a start and an end point :
350
351 f = bpf.start(x0,y0) : bpf.end(x1,y1);
352
353 A more involved break-point function can contains any number of intermediate
354 points
355
356 f = bpf.start(x0,y0) : bpf.point(x1,y1) : bpf.point(x2,y2) : bpf.end(x3,y3);
357
358 In any case the x_{i} must be in increasing order (for all i, x_{i} < x_{i+1})
359
360 For example the following definition :
361
362 f = bpf.start(x0,y0) : ... : bpf.point(xi,yi) : ... : bpf.end(xn,yn);
363
364 implements a break-point function f such that :
365
366 f(x) = y_{0} when x < x_{0}
367 f(x) = y_{n} when x > x_{n}
368 f(x) = y_{i} + (y_{i+1}-y_{i})*(x-x_{i})/(x_{i+1}-x_{i}) when x_{i} <= x and x < x_{i+1}
369
370 ******************************************************************************/
371
372 bpf = environment
373 {
374 // Start a break-point function
375 start(x0,y0) = \(x).(x0,y0,x,y0);
376
377 // Add a break-point
378 point(x1,y1) = \(x0,y0,x,y).(x1, y1, x , if (x < x0, y, if (x < x1, y0 + (x-x0)*(y1-y0)/(x1-x0), y1)));
379
380 // End a break-point function
381 end (x1,y1) = \(x0,y0,x,y).(if (x < x0, y, if (x < x1, y0 + (x-x0)*(y1-y0)/(x1-x0), y1)));
382
383 // definition of if
384 if (c,t,e) = select2(c,e,t);
385 };
386
387
388