File: physmodels.lib

package info (click to toggle)
faust 2.81.10%2Bds-1
links: PTS, VCS
area: main
in suites: forky, sid
size: 431,496 kB
sloc: cpp: 283,941; ansic: 116,215; javascript: 18,529; sh: 14,356; vhdl: 14,052; java: 5,900; python: 5,091; objc: 3,852; makefile: 2,725; cs: 1,672; lisp: 1,146; ruby: 954; yacc: 586; xml: 471; lex: 247; awk: 111; tcl: 26
file content (4116 lines) | stat: -rw-r--r-- 176,404 bytes
//##################################### physmodels.lib ###################################
// Faust physical modeling library. Its official prefix is `pm`.
//
// This library provides an environment to facilitate physical modeling of musical
// instruments. It contains dozens of functions implementing low and high level
// elements going from a simple waveguide to fully operational models with
// built-in UI, etc.
//
// It is organized as follows:
//
// * [Global Variables](#global-variables): useful pre-defined variables for
// physical modeling (e.g., speed of sound, etc.).
// * [Conversion Tools](#conversion-tools-1): conversion functions specific
// to physical modeling (e.g., length to frequency, etc.).
// * [Bidirectional Utilities](#bidirectional-utilities): functions to create
// bidirectional block diagrams for physical modeling.
// * [Basic Elements](#basic-elements-1): waveguides, specific types of filters, etc.
// * [String Instruments](#string-instruments): various types of strings
// (e.g., steel, nylon, etc.), bridges, guitars, etc.
// * [Bowed String Instruments](#bowed-string-instruments): parts and models
// specific to bowed string instruments (e.g., bows, bridges, violins, etc.).
// * [Wind Instrument](#wind-instruments): parts and models specific to wind
// instruments (e.g., reeds, mouthpieces, flutes, clarinets, etc.).
// * [Exciters](#exciters): pluck generators, "blowers", etc.
// * [Modal Percussions](#modal-percussions): percussion instruments based on
// modal models.
// * [Vocal Synthesis](#vocal-synthesis): functions for various vocal synthesis
// techniques (e.g., fof, source/filter, etc.) and vocal synthesizers.
// * [Misc Functions](#misc-functions): any other functions that don't fit in
// the previous category (e.g., nonlinear filters, etc.).
//
// This library is part of the Faust Physical Modeling ToolKit.
// More information on how to use this library can be found on this [page](https://ccrma.stanford.edu/~rmichon/pmFaust) or this [video](https://faust.grame.fr/community/events/#introduction-to-the-faust-physical-modeling-toolkit-romain-michon). Tutorials on how to make
// physical models of musical instruments using Faust can be found
// [here](https://ccrma.stanford.edu/~rmichon/faustTutorials/#making-physical-models-of-musical-instruments-with-faust) as well.
//
// #### References
// * <https://github.com/grame-cncm/faustlibraries/blob/master/physmodels.lib>
//########################################################################################
// Authors: Romain Michon, Pierre-Amaury Grumiaux, and Yann Orlarey

import("stdfaust.lib");

declare name "Faust Physical Models Library";
declare version "1.2.0";

/*
TODO:
    - It'd be cool to have a version of the block diagram generator that automatically flips
    things based on the use of chains, etc.
    - When setting pole of filters by hand (e.g. smooth, should adjust pole in function of SR)
    - Probably need a single resonator function / see how to integrate that with "mode"
    - Need a non-linear function and see how this can be integrated with modal synthesis
    - See how bowed modal models could be integrated to this
    - Currently still missing keyboard instruments
    - Currently still missing vocal synth: easy to fix (create a formant filter function)
    - Real polyphonic instruments should be designated with some kind of prefix (e.g.,
    full)
*/

//=============================Global Variables===========================================
// Useful pre-defined variables for physical modeling.
//========================================================================================

//--------------`(pm.)speedOfSound`----------
// Speed of sound in meters per second (340m/s).
//--------------------------------------
speedOfSound = 340;

//--------------`(pm.)maxLength`----------
// The default maximum length (3) in meters of strings and tubes used in this
// library. This variable should be overriden to allow longer strings or tubes.
//--------------------------------------
maxLength = 3;

//================================Conversion Tools=======================================
// Useful conversion tools for physical modeling.
//========================================================================================

//--------------`(pm.)f2l`----------
// Frequency to length in meters.
//
// #### Usage
//
// ```
// f2l(freq) : distanceInMeters
// ```
//
// Where:
//
// * `freq`: the frequency
//-------------------------------
f2l(freq) = speedOfSound/freq;

//--------------`(pm.)l2f`----------
// Length in meters to frequency.
//
// #### Usage
//
// ```
// l2f(length) : freq
// ```
//
// Where:
//
// * `length`: length/distance in meters
//-------------------------------
l2f(length) = speedOfSound/length;

//--------------`(pm.)l2s`----------
// Length in meters to number of samples.
//
// #### Usage
//
// ```
// l2s(l) : numberOfSamples
// ```
//
// Where:
//
// * `l`: length in meters
//-------------------------------
l2s(l) = l*ma.SR/speedOfSound;

//=============================Bidirectional Utilities====================================
// Set of fundamental functions to create bi-directional block diagrams in Faust.
// These elements are used as the basis of this library to connect high level
// elements (e.g., mouthpieces, strings, bridge, instrument body, etc.). Each
// block has 3 inputs and 3 outputs. The first input/output carry left going
// waves, the second input/output carry right going waves, and the third
// input/output is used to carry any potential output signal to the end of the
// algorithm.
//========================================================================================

//--------------`(pm.)basicBlock`----------
// Empty bidirectional block to be used with [`chain`](#chain): 3 signals ins
// and 3 signals out.
//
// #### Usage
//
// ```
// chain(basicBlock : basicBlock : etc.)
// ```
//-------------------------------
basicBlock = _,_,_;

//-------`(pm.)chain`----------
// Creates a chain of bidirectional blocks.
// Blocks must have 3 inputs and outputs. The first input/output carry left
// going waves, the second input/output carry right going waves, and the third
// input/output is used to carry any potential output signal to the end of the
// algorithm. The implied one sample delay created by the `~` operator is
// generalized to the left and right going waves. Thus, `n` blocks in `chain()`
// will add an `n` samples delay to both left and right going waves.
//
// #### Usage
//
// ```
// leftGoingWaves,rightGoingWaves,mixedOutput : chain( A : B ) : leftGoingWaves,rightGoingWaves,mixedOutput
// with {
//    A = _,_,_;
//    B = _,_,_;
// };
// ```
//-----------------------------
chain(A:As) = ((ro.crossnn(1),_',_ : _,A : ro.crossnn(1),_,_ : _,chain(As) : ro.crossnn(1),_,_)) ~ _ : !,_,_,_;
chain(A) = A;

//-------`(pm.)inLeftWave`--------------
// Adds a signal to left going waves anywhere in a [`chain`](#chain) of blocks.
//
// #### Usage
//
// ```
// model(x) = chain(A : inLeftWave(x) : B)
// ```
//
// Where `A` and `B` are bidirectional blocks and `x` is the signal added to left
// going waves in that chain.
//--------------------------------
inLeftWave(x) = +(x),_,_;

//-------`(pm.)inRightWave`--------------
// Adds a signal to right going waves anywhere in a [`chain`](#chain) of blocks.
//
// #### Usage
//
// ```
// model(x) = chain(A : inRightWave(x) : B)
// ```
//
// Where `A` and `B` are bidirectional blocks and `x` is the signal added to right
// going waves in that chain.
//--------------------------------
inRightWave(x) = _,+(x),_;

//-------`(pm.)in`--------------
// Adds a signal to left and right going waves anywhere in a [`chain`](#chain)
// of blocks.
//
// #### Usage
//
// ```
// model(x) = chain(A : in(x) : B)
// ```
//
// Where `A` and `B` are bidirectional blocks and `x` is the signal added to
// left and right going waves in that chain.
//--------------------------------
in(x) = +(x),+(x),_;

//-------`(pm.)outLeftWave`--------------
// Sends the signal of left going waves to the output channel of the [`chain`](#chain).
//
// #### Usage
//
// ```
// chain(A : outLeftWave : B)
// ```
//
// Where `A` and `B` are bidirectional blocks.
//--------------------------------
outLeftWave(x,y,s) = x,y,x+s;

//-------`(pm.)outRightWave`--------------
// Sends the signal of right going waves to the output channel of the [`chain`](#chain).
//
// #### Usage
//
// ```
// chain(A : outRightWave : B)
// ```
//
// Where `A` and `B` are bidirectional blocks.
//--------------------------------
outRightWave(x,y,s) = x,y,y+s;

//-------`(pm.)out`--------------
// Sends the signal of right and left going waves to the output channel of the
// [`chain`](#chain).
//
// #### Usage
//
// ```
// chain(A : out : B)
// ```
//
// Where `A` and `B` are bidirectional blocks.
//--------------------------------
out(x,y,s) = x,y,x+y+s;

//-------`(pm.)terminations`--------------
// Creates terminations on both sides of a [`chain`](#chain) without closing
// the inputs and outputs of the bidirectional signals chain. As for
// [`chain`](#chain), this function adds a 1 sample delay to the bidirectional
// signal, both ways. Of course, this function can be nested within a
// [`chain`](#chain).
//
// #### Usage
//
// ```
// terminations(a,b,c)
// with {
//    a = *(-1); // left termination
//    b = chain(D : E : F); // bidirectional chain of blocks (D, E, F, etc.)
//    c = *(-1); // right termination
// };
// ```
//----------------------------------------
terminations(a,b,c) = (_,ro.crossnn(1),_,_ : +,+,_ : b) ~ (a,c : ro.crossnn(1));

//-------`(pm.)lTermination`----------
// Creates a termination on the left side of a [`chain`](#chain) without
// closing the inputs and outputs of the bidirectional signals chain. This
// function adds a 1 sample delay near the termination and can be nested
// within another [`chain`](#chain).
//
// #### Usage
//
// ```
// lTerminations(a,b)
// with {
//    a = *(-1); // left termination
//    b = chain(D : E : F); // bidirectional chain of blocks (D, E, F, etc.)
// };
// ```
//----------------------------------------
lTermination(a,b) = (ro.crossnn(1),_,_ : _,+,_ : b) ~ a;

//-------`(pm.)rTermination`----------
// Creates a termination on the right side of a [`chain`](#chain) without
// closing the inputs and outputs of the bidirectional signals chain. This
// function adds a 1 sample delay near the termination and can be nested
// within another [`chain`](#chain).
//
// #### Usage
//
// ```
// rTerminations(b,c)
// with {
//    b = chain(D : E : F); // bidirectional chain of blocks (D, E, F, etc.)
//    c = *(-1); // right termination
// };
// ```
//----------------------------------------
rTermination(b,c) = (_,_,_,_ : +,_,_ : b) ~ (!,c);

//-------`(pm.)closeIns`----------
// Closes the inputs of a bidirectional chain in all directions.
//
// #### Usage
//
// ```
// closeIns : chain(...) : _,_,_
// ```
//----------------------------------------
closeIns = 0,0,0;

//-------`(pm.)closeOuts`----------
// Closes the outputs of a bidirectional chain in all directions except for the
// main signal output (3d output).
//
// #### Usage
//
// ```
// _,_,_ : chain(...) : _
// ```
//----------------------------------------
closeOuts = !,!,_;

//-------`(pm.)endChain`----------
// Closes the inputs and outputs of a bidirectional chain in all directions
// except for the main signal output (3d output).
//
// #### Usage
//
// ```
// endChain(chain(...)) : _
// ```
//----------------------------------------
endChain(b) = closeIns : b : closeOuts;


//==================================Basic Elements========================================
// Basic elements for physical modeling (e.g., waveguides, specific filters,
// etc.).
//========================================================================================

//-------`(pm.)waveguideN`----------
// A series of waveguide functions based on various types of delays (see
// [`fdelay[n]`](#fdelayn)).
//
// #### List of functions
//
// * `waveguideUd`: unit delay waveguide
// * `waveguideFd`: fractional delay waveguide
// * `waveguideFd2`: second order fractional delay waveguide
// * `waveguideFd4`: fourth order fractional delay waveguide
//
// #### Usage
//
// ```
// chain(A : waveguideUd(nMax,n) : B)
// ```
//
// Where:
//
// * `nMax`: the maximum length of the delays in the waveguide
// * `n`: the length of the delay lines in samples.
//----------------------------------
waveguideUd(nMax,n) = par(i,2,de.delay(nMax,n)),_;
waveguideFd(nMax,n) = par(i,2,de.fdelay(nMax,n)),_;
waveguideFd2(nMax,n) = par(i,2,de.fdelay2(nMax,n)),_;
waveguideFd4(nMax,n) = par(i,2,de.fdelay4(nMax,n)),_;


//-------`(pm.)waveguide`----------
// Standard `pm.lib` waveguide (based on [`waveguideFd4`](#waveguiden)).
//
// #### Usage
//
// ```
// chain(A : waveguide(nMax,n) : B)
// ```
//
// Where:
//
// * `nMax`: the maximum length of the delays in the waveguide
// * `n`: the length of the delay lines in samples.
//----------------------------------
waveguide(nMax,n) = waveguideFd4(nMax,n);


//-------`(pm.)bridgeFilter`----------
// Generic two zeros bridge FIR filter (as implemented in the
// [STK](https://ccrma.stanford.edu/software/stk/)) that can be used to
// implement the reflectance violin, guitar, etc. bridges.
//
// #### Usage
//
// ```
// _ : bridge(brightness,absorption) : _
// ```
//
// Where:
//
// * `brightness`: controls the damping of high frequencies (0-1)
// * `absorption`: controls the absorption of the brige and thus the t60 of
// the string plugged to it (0-1) (1 = 20 seconds)
//----------------------------------
// TODO: perhaps, the coefs of this filter should be adapted in function of SR
bridgeFilter(brightness,absorption,x) = rho * (h0 * x' + h1*(x+x''))
with {
    freq = 320;
    t60 = (1-absorption)*20;
    h0 = (1.0 + brightness)/2;
    h1 = (1.0 - brightness)/4;
    rho = pow(0.001,1.0/(freq*t60));
};


//-------`(pm.)modeFilter`----------
// Resonant bandpass filter that can be used to implement a single resonance
// (mode).
//
// #### Usage
//
// ```
// _ : modeFilter(freq,t60,gain) : _
// ```
//
// Where:
//
// * `freq`: mode frequency
// * `t60`: mode resonance duration (in seconds)
// * `gain`: mode gain (0-1)
//----------------------------------
modeFilter(freq,t60,gain) = fi.tf2(b0,b1,b2,a1,a2)*gain
with {
    b0 = 1;
    b1 = 0;
    b2 = -1;
    w = 2*ma.PI*freq/ma.SR;
    r = pow(0.001,1/float(t60*ma.SR));
    a1 = -2*r*cos(w);
    a2 = r^2;
};


//================================String Instruments======================================
// Low and high level string instruments parts. Most of the elements in
// this section can be used in a bidirectional chain.
//========================================================================================

//-------`(pm.)stringSegment`----------
// A string segment without terminations (just a simple waveguide).
//
// #### Usage
//
// ```
// chain(A : stringSegment(maxLength,length) : B)
// ```
//
// Where:
//
// * `maxLength`: the maximum length of the string in meters (should be static)
// * `length`: the length of the string in meters
//----------------------------------
stringSegment(maxLength,length) = waveguide(nMax,n)
with {
    nMax = maxLength : l2s;
    n = length : l2s/2;
};


//-------`(pm.)openString`----------
// A bidirectional block implementing a basic "generic" string with a
// selectable excitation position. Lowpass filters are built-in and
// allow to simulate the effect of dispersion on the sound and thus
// to change the "stiffness" of the string.
//
// #### Usage
//
// ```
// chain(... : openString(length,stiffness,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `stiffness`: the stiffness of the string (0-1) (1 for max stiffness)
// * `pluckPosition`: excitation position (0-1) (1 is bottom)
// * `excitation`: the excitation signal
//----------------------------------
openString(length,stiffness,pluckPosition,excitation) = chain(stringSegment(maxStringLength,ntbd) : in(excitation) : dispersionFilters : stringSegment(maxStringLength,btbd))
with {
    dispersionFilters = par(i,2,si.smooth(stiffness)),_; // one pole filters
    maxStringLength = maxLength;
    ntbd = length*pluckPosition; // length of the upper portion of the string
    btbd = length*(1-pluckPosition); // length of the lower portion of the string
};


//-------`(pm.)nylonString`----------
// A bidirectional block implementing a basic nylon string with selectable
// excitation position. This element is based on [`openString`](#openstring)
// and has a fix stiffness corresponding to that of a nylon string.
//
// #### Usage
//
// ```
// chain(... : nylonString(length,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: excitation position (0-1) (1 is bottom)
// * `excitation`: the excitation signal
//----------------------------------
nylonString(length,pluckPosition,excitation) =
openString(length,stiffness,pluckPosition,excitation)
with {
    stiffness = 0.4; // empirically set but it sounds good ;)
};


//-------`(pm.)steelString`----------
// A bidirectional block implementing a basic steel string with selectable
// excitation position. This element is based on [`openString`](#openstring)
// and has a fix stiffness corresponding to that of a steel string.
//
// #### Usage
//
// ```
// chain(... : steelString(length,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: excitation position (0-1) (1 is bottom)
// * `excitation`: the excitation signal
//----------------------------------
steelString(length,pluckPosition,excitation) =
openString(length,stiffness,pluckPosition,excitation)
with {
    stiffness = 0.05; // empirically set but it sounds good ;)
    // in fact, we could almost get rid of the filters in that case,
    // but I think it's good to keep them for consistency
};


//-------`(pm.)openStringPick`----------
// A bidirectional block implementing a "generic" string with selectable
// excitation position. It also has a built-in pickup whose position is the
// same as the excitation position. Thus, moving the excitation position
// will also move the pickup.
//
// #### Usage
//
// ```
// chain(... : openStringPick(length,stiffness,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `stiffness`: the stiffness of the string (0-1) (1 for max stiffness)
// * `pluckPosition`: excitation position (0-1) (1 is bottom)
// * `excitation`: the excitation signal
//----------------------------------
openStringPick(length,stiffness,pluckPosition,excitation) = strChain
with {
    dispersionFilters = par(i,2,si.smooth(stiffness)),_;
    maxStringLength = maxLength;
    nti = length*pluckPosition; // length of the upper portion of the string
    itb = length*(1-pluckPosition); // length of the lower portion of the string
    strChain = chain(stringSegment(maxStringLength,nti) : in(excitation) : out :
        dispersionFilters : stringSegment(maxStringLength,itb));
};


//-------`(pm.)openStringPickUp`----------
// A bidirectional block implementing a "generic" string with selectable
// excitation position and stiffness. It also has a built-in pickup whose
// position can be independenly selected. The only constraint is that the
// pickup has to be placed after the excitation position.
//
// #### Usage
//
// ```
// chain(... : openStringPickUp(length,stiffness,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `stiffness`: the stiffness of the string (0-1) (1 for max stiffness)
// * `pluckPosition`: pluck position between the top of the string and the
// pickup (0-1) (1 for same as pickup position)
// * `pickupPosition`: position of the pickup on the string (0-1) (1 is bottom)
// * `excitation`: the excitation signal
//----------------------------------
openStringPickUp(length,stiffness,pluckPosition,pickupPosition,excitation) = strChain
with {
    dispersionFilters = par(i,2,si.smooth(stiffness)),_;
    maxStringLength = maxLength;
    nti = length*pluckPosition; // top to excitation length
    nto = nti*pickupPosition; // nuts to pickup length
    oti = nti*(1-pickupPosition); // pickup to excitation length
    itb = length*(1-pluckPosition); // pickup to bottom length
    strChain = chain(stringSegment(maxStringLength,nto) : out :
    stringSegment(maxStringLength,oti) : in(excitation) : dispersionFilters :
    stringSegment(maxStringLength,itb));
};


//-------`(pm.)openStringPickDown`----------
// A bidirectional block implementing a "generic" string with selectable
// excitation position and stiffness. It also has a built-in pickup whose
// position can be independenly selected. The only constraint is that the
// pickup has to be placed before the excitation position.
//
// #### Usage
//
// ```
// chain(... : openStringPickDown(length,stiffness,pluckPosition,excitation) : ...)
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `stiffness`: the stiffness of the string (0-1) (1 for max stiffness)
// * `pluckPosition`: pluck position on the string (0-1) (1 is bottom)
// * `pickupPosition`: position of the pickup between the top of the string
// and the excitation position (0-1) (1 is excitation position)
// * `excitation`: the excitation signal
//----------------------------------
openStringPickDown(length,stiffness,pluckPosition,pickupPosition,excitation) =
strChain
with {
    dispersionFilters = par(i,2,si.smooth(stiffness)),_;
    maxStringLength = maxLength;
    nto = length*pickupPosition; // top to pickup length
    nti = nto*pluckPosition; // top to excitation length
    ito = nto*(1-pluckPosition); // excitation to pickup length
    otb = length*(1-pickupPosition); // pickup to bottom length
    strChain = chain(stringSegment(maxStringLength,nti) : in(excitation) :
    stringSegment(maxStringLength,ito) : out : dispersionFilters :
    stringSegment(maxStringLength,otb));
};


// TODO: eventually, we'd want to implement a generic function here that
// automatically switches the position of elements in the algorithm
// depending on the position of the pick. Even though this is currently
// possible, it will pose optimization issues (we'd want the new mute
// feature of Faust to be generalized in order to do that)

//-------`(pm.)ksReflexionFilter`----------
// The "typical" one-zero Karplus-strong feedforward reflexion filter. This
// filter will be typically used in a termination (see below).
//
// #### Usage
//
// ```
// terminations(_,chain(...),ksReflexionFilter)
// ```
//----------------------------------
ksReflexionFilter = _ <: (_+_')/2;


//-------`(pm.)rStringRigidTermination`----------
// Bidirectional block implementing a right rigid string termination (no damping,
// just phase inversion).
//
// #### Usage
//
// ```
// chain(rStringRigidTermination : stringSegment : ...)
// ```
//----------------------------------
rStringRigidTermination = rTermination(basicBlock,*(-1));


//-------`(pm.)lStringRigidTermination`----------
// Bidirectional block implementing a left rigid string termination (no damping,
// just phase inversion).
//
// #### Usage
//
// ```
// chain(... : stringSegment : lStringRigidTermination)
// ```
//----------------------------------
lStringRigidTermination = lTermination(*(-1),basicBlock);


//-------`(pm.)elecGuitarBridge`----------
// Bidirectional block implementing a simple electric guitar bridge. This
// block is based on [`bridgeFilter`](#bridgeFilter). The bridge doesn't
// implement transmittance since it is not meant to be connected to a
// body (unlike acoustic guitar). It also partially sets the resonance
// duration of the string with the nuts used on the other side.
//
// #### Usage
//
// ```
// chain(... : stringSegment : elecGuitarBridge)
// ```
//----------------------------------
elecGuitarBridge = rTermination(basicBlock,-bridgeFilter(0.8,0.6));


//-------`(pm.)elecGuitarNuts`----------
// Bidirectional block implementing a simple electric guitar nuts. This
// block is based on [`bridgeFilter`](#bridgeFilter) and does essentially
// the same thing as [`elecGuitarBridge`](#elecguitarbridge), but on the
// other side of the chain. It also partially sets the resonance duration of
// the string with the bridge used on the other side.
//
// #### Usage
//
// ```
// chain(elecGuitarNuts : stringSegment : ...)
// ```
//----------------------------------
elecGuitarNuts = lTermination(-bridgeFilter(0.8,0.6),basicBlock);


//-------`(pm.)guitarBridge`----------
// Bidirectional block implementing a simple acoustic guitar bridge. This
// bridge damps more hight frequencies than
// [`elecGuitarBridge`](#elecguitarbridge) and implements a transmittance
// filter. It also partially sets the resonance duration of the string with
// the nuts used on the other side.
//
// #### Usage
//
// ```
// chain(... : stringSegment : guitarBridge)
// ```
//----------------------------------
guitarBridge = rTermination(basicBlock,reflectance) : _,transmittance,_
with {
    reflectance = -bridgeFilter(0.4,0.5);
    transmittance = _; // TODO
};


//-------`(pm.)guitarNuts`----------
// Bidirectional block implementing a simple acoustic guitar nuts. This
// nuts damps more hight frequencies than
// [`elecGuitarNuts`](#elecguitarnuts) and implements a transmittance
// filter. It also partially sets the resonance duration of the string with
// the bridge used on the other side.
//
// #### Usage
//
// ```
// chain(guitarNuts : stringSegment : ...)
// ```
//----------------------------------
guitarNuts = lTermination(-bridgeFilter(0.4,0.5),basicBlock);


//-------`(pm.)idealString`----------
// An "ideal" string with rigid terminations and where the plucking position
// and the pick-up position are the same. Since terminations are rigid, this
// string will ring forever.
//
// #### Usage
//
// ```
// 1-1' : idealString(length,reflexion,xPosition,excitation)
// ```
//
// With:
// * `length`: the length of the string in meters
// * `pluckPosition`: the plucking position (0.001-0.999)
// * `excitation`: the input signal for the excitation.
//----------------------------------------------------------
idealString(length,pluckPosition,excitation) = wg
with {
    maxStringLength = maxLength;
    lengthTuning = 0.08; // tuned "by hand"
    tunedLength = length-lengthTuning;
    nUp = tunedLength*pluckPosition; // upper string segment length
    nDown = tunedLength*(1-pluckPosition); // lower string segment length
    wg = chain(lStringRigidTermination : stringSegment(maxStringLength,nUp) :
    in(excitation) : out : stringSegment(maxStringLength,nDown) :
    rStringRigidTermination); // waveguide chain
};


//-------`(pm.)ks`----------
// A Karplus-Strong string (in that case, the string is implemented as a
// one dimension waveguide).
//
// #### Usage
//
// ```
// ks(length,damping,excitation) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `damping`: string damping (0-1)
// * `excitation`: excitation signal
//----------------------------------
ks(length,damping,excitation) = endChain(ksChain)
with {
    maxStringLength = maxLength;
    lengthTuning = 0.05; // tuned "by hand"
    tunedLength = length-lengthTuning;
    refCoef = (1-damping)*0.2+0.8;
    refFilter = ksReflexionFilter*refCoef;
    ksChain = terminations(_,chain(in(excitation) :
    stringSegment(maxStringLength,tunedLength) : out),refFilter);
};


//-------`(pm.)ks_ui_MIDI`----------
// Ready-to-use, MIDI-enabled Karplus-Strong string with buil-in UI.
//
// #### Usage
//
// ```
// ks_ui_MIDI : _
// ```
//----------------------------------
ks_ui_MIDI = gate : impulseExcitation*gain : ks( (freq : f2l), damping )
with {
    f = hslider("v:karplus/h:[0]params/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:karplus/h:[0]params/[1]bend[style:knob][hidden:1][midi:pitchwheel]"
    ,0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:karplus/h:[0]params/[2]gain[style:knob]",0.8,0,1,0.01);
    s = hslider("v:karplus/h:[0]params/[3]sustain[hidden:1][midi:ctrl 64][style:knob]"
    ,0,0,1,1);
    damping = hslider("v:karplus/h:[0]params/[1]damping[midi:ctrl 1][style:knob]"
    ,0.01,0,1,0.01) : si.smoo;
    t = button("v:karplus/[1]gate");

    gate = t+s : min(1);
    freq = f*bend;
};


//-------`(pm.)elecGuitarModel`----------
// A simple electric guitar model (without audio effects, of course) with
// selectable pluck position.
// This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function. Pitch is changed by
// changing the length of the string and not through a finger model.
//
// #### Usage
//
// ```
// elecGuitarModel(length,pluckPosition,mute,excitation) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `mute`: mute coefficient (1 for no mute and 0 for instant mute)
// * `excitation`: excitation signal
//----------------------------------
elecGuitarModel(length,pluckPosition,mute,excitation) = endChain(egChain)
with {
    maxStringLength = maxLength;
    lengthTuning = 0.11; // tuned "by hand"
    stringL = length-lengthTuning;
    muteBlock = *(mute),*(mute),_;
    egChain = chain(
        elecGuitarNuts :
        openStringPick(stringL,0.05,pluckPosition,excitation) :
        muteBlock :
        elecGuitarBridge);
};


//-------`(pm.)elecGuitar`----------
// A simple electric guitar model with steel strings (based on
// [`elecGuitarModel`](#elecguitarmodel)) implementing an excitation
// model.
// This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function.
//
// #### Usage
//
// ```
// elecGuitar(length,pluckPosition,trigger) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `mute`: mute coefficient (1 for no mute and 0 for instant mute)
// * `gain`: gain of the pluck (0-1)
// * `trigger`: trigger signal (1 for on, 0 for off)
//----------------------------------
elecGuitar(stringLength,pluckPosition,mute,gain,trigger) =
pluckString(stringLength,1,1,1,gain,trigger) :
elecGuitarModel(stringLength,pluckPosition,mute);


//-------`(pm.)elecGuitar_ui_MIDI`----------
// Ready-to-use MIDI-enabled electric guitar physical model with built-in UI.
//
// #### Usage
//
// ```
// elecGuitar_ui_MIDI : _
// ```
//----------------------------------
elecGuitar_ui_MIDI = elecGuitar(stringLength,pluckPosition,1,gain,gate)*outGain
with {
    f = hslider("v:elecGuitar/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:elecGuitar/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel][style:knob]"
    ,0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:elecGuitar/h:[0]midi/[2]gain[style:knob]",0.8,0,1,0.01);
    s = hslider("v:elecGuitar/h:[0]midi/[3]sustain[hidden:1]
    [midi:ctrl 64][style:knob]",0,0,1,1);
    pluckPosition = hslider("v:elecGuitar/[1]pluckPosition[midi:ctrl 1]",0.8,0,1,0.01) : si.smoo;
    outGain = hslider("v:elecGuitar/[2]outGain",0.5,0,1,0.01);
    t = button("v:elecGuitar/[3]gate");
    gate = t+s : min(1);
    freq = f*bend;
    stringLength = freq : f2l;
};


//-------`(pm.)guitarBody`----------
// WARNING: not implemented yet!
// Bidirectional block implementing a simple acoustic guitar body.
//
// #### Usage
//
// ```
// chain(... : guitarBody)
// ```
//----------------------------------
// TODO: not implemented yet
guitarBody = reflectance,transmittance,_
with {
    transmittance = _;
    reflectance = _;
};


//-------`(pm.)guitarModel`----------
// A simple acoustic guitar model with steel strings and selectable excitation
// position. This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function. Pitch is changed by
// changing the length of the string and not through a finger model.
// WARNING: this function doesn't currently implement a body (just strings and
// bridge).
//
// #### Usage
//
// ```
// guitarModel(length,pluckPosition,excitation) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `excitation`: excitation signal
//----------------------------------
guitarModel(length,pluckPosition,excitation) = endChain(egChain)
with {
    maxStringLength = maxLength;
    lengthTuning = 0.1; // tuned "by hand"
    stringL = length-lengthTuning;
    egChain = chain(guitarNuts : steelString(stringL,pluckPosition,excitation) :
    guitarBridge : guitarBody : out);
};


//-------`(pm.)guitar`----------
// A simple acoustic guitar model with steel strings (based on
// [`guitarModel`](#guitarmodel)) implementing an excitation model.
// This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function.
//
// #### Usage
//
// ```
// guitar(length,pluckPosition,trigger) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `gain`: gain of the excitation
// * `trigger`: trigger signal (1 for on, 0 for off)
//----------------------------------
guitar(stringLength,pluckPosition,gain,trigger) =
pluckString(stringLength,1,1.5,1,gain,trigger) : guitarModel(stringLength, pluckPosition);


//-------`(pm.)guitar_ui_MIDI`----------
// Ready-to-use MIDI-enabled steel strings acoustic guitar physical model with
// built-in UI.
//
// #### Usage
//
// ```
// guitar_ui_MIDI : _
// ```
//----------------------------------
guitar_ui_MIDI = guitar(stringLength,pluckPosition,gain,gate)*outGain
with {
    f = hslider("v:guitar/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:guitar/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:guitar/h:[0]midi/[2]gain[style:knob]",0.8,0,1,0.01);
    s = hslider("v:guitar/h:[0]midi/[3]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    pluckPosition = hslider("v:guitar/pluckPosition[midi:ctrl 1]"
    ,0.8,0,1,0.01) : si.smoo;
    outGain = hslider("v:guitar/outGain",0.5,0,1,0.01);
    t = button("v:guitar/[4]gate");
    gate = t+s : min(1);
    freq = f*bend;
    stringLength = freq : f2l;
};


//-------`(pm.)nylonGuitarModel`----------
// A simple acoustic guitar model with nylon strings and selectable excitation
// position. This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function. Pitch is changed by
// changing the length of the string and not through a finger model.
// WARNING: this function doesn't currently implement a body (just strings and
// bridge).
//
// #### Usage
//
// ```
// nylonGuitarModel(length,pluckPosition,excitation) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `excitation`: excitation signal
//----------------------------------
nylonGuitarModel(length,pluckPosition,excitation) = endChain(egChain)
with {
    maxStringLength = maxLength; // meters
    lengthTuning = 0.11;
    stringL = length-lengthTuning;
    egChain = chain(guitarNuts : nylonString(stringL,pluckPosition,excitation) :
    guitarBridge : guitarBody : out);
};


//-------`(pm.)nylonGuitar`----------
// A simple acoustic guitar model with nylon strings (based on
// [`nylonGuitarModel`](#nylonguitarmodel)) implementing an excitation model.
// This model implements a single string. Additional strings should be created
// by making a polyphonic application out of this function.
//
// #### Usage
//
// ```
// nylonGuitar(length,pluckPosition,trigger) : _
// ```
//
// Where:
//
// * `length`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `gain`: gain of the excitation (0-1)
// * `trigger`: trigger signal (1 for on, 0 for off)
//----------------------------------
nylonGuitar(stringLength,pluckPosition,gain,trigger) =
pluckString(stringLength,1,1.5,1,gain,trigger) : nylonGuitarModel(stringLength, pluckPosition);


//-------`(pm.)nylonGuitar_ui_MIDI`----------
// Ready-to-use MIDI-enabled nylon strings acoustic guitar physical model with
// built-in UI.
//
// #### Usage
//
// ```
// nylonGuitar_ui_MIDI : _
// ```
//----------------------------------
nylonGuitar_ui_MIDI = nylonGuitar(stringLength,pluckPosition,gain,gate)*outGain
with {
    f = hslider("v:nylonGuitar/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:nylonGuitar/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:nylonGuitar/h:[0]midi/[2]gain[style:knob]",0.8,0,1,0.01);
    s = hslider("v:nylonGuitar/h:[0]midi/[3]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    pluckPosition = hslider("v:nylonGuitar/pluckPosition[midi:ctrl 1]"
    ,0.8,0,1,0.01) : si.smoo;
    outGain = hslider("v:nylonGuitar/outGain",0.5,0,1,0.01);
    t = button("v:nylonGuitar/[4]gate");
    gate = t+s : min(1);
    freq = f*bend;
    stringLength = freq : f2l;
};


//-------`(pm.)modeInterpRes`----------
// Modular string instrument resonator based on IR measurements made on 3D 
// printed models. The 2D space allowing for the control of the shape and the
// scale of the model is enabled by interpolating between modes parameters.
// More information about this technique/project can be found here: 
// <https://ccrma.stanford.edu/~rmichon/3dPrintingModeling/>.
//
// #### Usage
//
// ```
// _ : modeInterpRes(nModes,x,y) : _
// ```
//
// Where:
//
// * `nModes`: number of modeled modes (40 max)
// * `x`: shape of the resonator (0: square, 1: square with rounded corners, 2: round)
// * `y`: scale of the resonator (0: small, 1: medium, 2: large)
//----------------------------------
modeInterpRes(nModes,x,y) = _ <: par(i,min(40,nModes),modeFilter(modeFreq(i),modeT60(i),modeGain(i))) :> /(nModes)
with {
    // modes parameters
    roundBigCenterFreq(i) = ba.take(i+1,(286.134404,476.384455,602.872587,831.823345,862.773611,1024.985313,1209.932428,1515.548598,1777.819591,2001.806241,2674.135433,3154.282144,3329.555657,3431.385589,3673.656995,4228.801853,4389.842060,4552.745406,4766.151594,4854.354993,4970.516137,5064.203325,5200.399916,5347.650831,5459.807495,5630.144253,5870.258034,6121.961406,6244.522341,6677.052569,7233.100711,7520.031251,7677.902736,7805.138717,9918.905945,10418.988301,11296.101954,11432.719867,13422.238365,13679.232778));
    roundBigCenterGain(i) = ba.take(i+1,(1.000000,0.609711,0.222375,0.025505,0.026108,0.026840,0.605767,0.368186,0.093526,0.016656,0.159600,0.013836,0.044000,0.023103,0.107971,0.013696,0.043615,0.039848,0.041903,0.043671,0.055018,0.041730,0.062832,0.030464,0.028049,0.012880,0.047154,0.084707,0.050532,0.032735,0.031938,0.023554,0.017534,0.017268,0.026202,0.013725,0.015611,0.021294,0.011486,0.011875));
    roundBigCenterT60(i) = ba.take(i+1,(0.283544,0.076723,0.287423,0.058391,0.075172,0.120794,0.131540,0.081173,0.088549,0.034991,0.060237,0.014153,0.074709,0.068893,0.033354,0.032013,0.019088,0.020992,0.026703,0.018040,0.015815,0.016671,0.028732,0.017929,0.021520,0.024622,0.022706,0.023138,0.027323,0.020556,0.024577,0.023669,0.024714,0.024420,0.017733,0.017881,0.018715,0.017737,0.015660,0.016464));
    roundMidCenterFreq(i) = ba.take(i+1,(380.910663,615.862831,933.855656,1128.769825,1231.353452,1404.948092,1764.641239,2095.586536,2451.220548,2566.435588,2781.467120,2946.757398,3522.047891,3748.375254,3892.352583,4340.896246,4454.851770,4629.302932,4806.417393,5024.545475,5604.684752,5671.112515,6533.147087,6751.949816,6862.781990,6951.226279,7180.863112,7901.631311,8074.662174,8979.395914,10075.249374,10521.350963,11964.859921,14227.213322,15501.455250,15587.852956,16290.138756,16418.172746,16663.266694,17991.614994));
    roundMidCenterGain(i) = ba.take(i+1,(1.000000,0.437495,0.110598,0.025858,0.010841,0.022400,0.245241,0.187990,0.029650,0.058174,0.007035,0.013188,0.009767,0.059284,0.160848,0.007463,0.016816,0.008689,0.031644,0.045608,0.019773,0.018546,0.013495,0.018482,0.020130,0.015417,0.031109,0.007624,0.008526,0.032539,0.009846,0.007605,0.019866,0.007583,0.008991,0.009194,0.010128,0.008557,0.009333,0.007494));
    roundMidCenterT60(i) = ba.take(i+1,(0.193654,0.076056,0.197517,0.019475,0.043087,0.063829,0.074308,0.059791,0.082616,0.081949,0.042855,0.036140,0.068564,0.041711,0.043799,0.025988,0.057494,0.028353,0.018586,0.018651,0.025380,0.021250,0.019217,0.016326,0.019469,0.020754,0.019487,0.019590,0.011348,0.015969,0.015674,0.016305,0.026772,0.010597,0.009926,0.007595,0.008676,0.008274,0.011471,0.012424));
    roundSmallCenterFreq(i) = ba.take(i+1,(511.119501,860.897453,977.675690,1245.454667,1408.823199,1604.478891,1825.857776,1960.985105,2340.204377,2594.337460,2700.041698,3277.766643,3388.510632,3505.100505,4151.457765,4257.308449,4425.319123,4526.878869,4629.688075,4819.012204,5493.628166,6323.383769,6514.513248,6709.979636,6823.104590,7009.748398,7198.916858,9019.248381,9492.252653,9706.451249,9858.889433,10080.422278,10252.000811,10389.864946,12327.892342,12465.475093,13491.187794,15072.291233,15176.746743,20143.302052));
    roundSmallCenterGain(i) = ba.take(i+1,(1.000000,0.686570,0.362616,0.178685,0.113981,0.158727,0.024067,0.023659,0.085113,0.192224,0.174412,0.160279,0.042275,0.042204,0.027219,0.040505,0.064866,0.069074,0.029473,0.033054,0.055311,0.073694,0.071341,0.068366,0.026735,0.037359,0.026008,0.025958,0.025327,0.052859,0.082911,0.038156,0.046079,0.034839,0.022738,0.024480,0.158056,0.026848,0.028009,0.031387));
    roundSmallCenterT60(i) = ba.take(i+1,(0.200041,0.030199,0.050855,0.020135,0.020101,0.026910,0.043249,0.040814,0.058128,0.047490,0.051189,0.069063,0.047344,0.049448,0.024720,0.027275,0.028407,0.030171,0.026606,0.030289,0.040772,0.023856,0.020889,0.032859,0.027543,0.026427,0.032154,0.020576,0.032990,0.019344,0.019827,0.018571,0.017213,0.019002,0.020440,0.021329,0.023027,0.018833,0.018344,0.018811));
    semBigCenterFreq(i) = ba.take(i+1,(318.812824,545.646717,727.103717,813.406232,938.026296,1067.472803,1239.311077,1478.409235,1567.366603,1840.215865,2162.299536,2302.660906,2394.006715,2499.397129,2724.662564,2876.713036,3051.447083,3219.116214,3378.667296,3507.736722,3687.907343,3827.904993,4039.858029,4208.742422,4385.651803,4574.818448,4744.136436,4866.893882,5090.532579,5200.283111,5448.339367,5751.719886,6040.059676,6500.665192,6780.883377,6878.732474,7127.463954,7701.025454,10337.623684,17411.019601));
    semBigCenterGain(i) = ba.take(i+1,(1.000000,0.386065,0.297923,0.235171,0.346593,0.219093,0.050655,0.310031,0.175387,0.180464,0.064716,0.070721,0.027399,0.076973,0.086367,0.024188,0.100553,0.028776,0.027974,0.040200,0.031680,0.021694,0.035362,0.031821,0.027726,0.037219,0.026308,0.016164,0.022938,0.046669,0.031128,0.027430,0.032404,0.018356,0.026887,0.041003,0.012779,0.028785,0.013346,0.017281));
    semBigCenterT60(i) = ba.take(i+1,(0.195061,0.269748,0.108389,0.132557,0.127125,0.113759,0.159514,0.058527,0.067333,0.071996,0.067001,0.054194,0.052818,0.050935,0.045373,0.051611,0.051410,0.061272,0.050844,0.029415,0.021966,0.023483,0.022866,0.023418,0.033572,0.026764,0.017334,0.021394,0.025984,0.026146,0.018801,0.018242,0.022440,0.020939,0.011937,0.011739,0.020314,0.013848,0.016504,0.085830));
    semMidCenterFreq(i) = ba.take(i+1,(316.511541,548.598445,682.927606,764.088220,1145.795222,1280.068992,1524.562996,1677.146625,2083.900770,2285.783756,2452.226740,2779.274140,3196.405966,3569.262928,3798.251779,4071.264178,4413.799330,4550.364404,4966.398087,5156.799448,5363.124121,5611.038360,5995.650074,6162.624225,6583.973014,7029.525719,7207.893991,7706.775019,7844.032020,10298.804956,10564.531850,11615.688447,12056.168362,12661.205157,12876.405687,13082.890228,13851.437211,14452.383106,15683.910228,16667.950883));
    semMidCenterGain(i) = ba.take(i+1,(0.999641,0.798181,0.200080,0.280473,0.042761,0.079955,0.286673,0.162652,1.000000,0.309945,0.097527,0.023964,0.124529,0.155314,0.075600,0.178980,0.089168,0.134962,0.120730,0.036885,0.024020,0.102628,0.085255,0.059439,0.023339,0.027443,0.141618,0.044053,0.024269,0.046682,0.034465,0.068930,0.070502,0.025912,0.076088,0.024972,0.038501,0.051135,0.027801,0.022469));
    semMidCenterT60(i) = ba.take(i+1,(0.210670,0.062329,0.097470,0.136802,0.077790,0.025877,0.054632,0.060485,0.061625,0.061352,0.060614,0.030867,0.039440,0.030054,0.037956,0.030948,0.046137,0.045345,0.025326,0.025515,0.022514,0.043142,0.018663,0.015057,0.052844,0.016478,0.025993,0.017880,0.014006,0.022971,0.017720,0.032173,0.015391,0.009217,0.028936,0.011895,0.012607,0.018141,0.018341,0.018180));
    semSmallCenterFreq(i) = ba.take(i+1,(440.578370,637.619026,828.182637,978.279267,1088.401911,1321.315700,1648.932125,1817.216790,1938.686740,2173.897496,2450.967026,2582.312643,2748.003128,2858.292897,2973.003533,3235.482621,3348.561906,3496.046199,3706.507771,3934.307821,4418.408376,5239.069796,5428.018005,5771.733971,6106.378549,6304.512786,6788.234913,6905.001294,7054.694449,7227.988097,7348.344488,8450.449049,9000.316477,10930.327581,11271.451255,11418.301881,13448.385700,13938.733921,14334.394664,17413.449273));
    semSmallCenterGain(i) = ba.take(i+1,(1.000000,0.297197,0.715097,0.178811,0.229363,0.183731,0.107685,0.120764,0.060436,0.503361,0.352264,0.056432,0.042926,0.090906,0.224894,0.128151,0.257120,0.076416,0.096971,0.078669,0.065573,0.057421,0.085378,0.085397,0.050578,0.046534,0.213480,0.090142,0.078632,0.089706,0.276299,0.061182,0.038941,0.045447,0.043447,0.051347,0.055647,0.071836,0.040813,0.066853));
    semSmallCenterT60(i) = ba.take(i+1,(0.136710,0.070656,0.049965,0.034104,0.023973,0.032990,0.049605,0.026798,0.024777,0.040671,0.039262,0.035240,0.019214,0.053001,0.045917,0.045864,0.048770,0.040981,0.020747,0.022565,0.031089,0.026788,0.013828,0.023917,0.022485,0.035319,0.017955,0.014831,0.018902,0.014661,0.032984,0.027729,0.025067,0.021211,0.016910,0.027488,0.013164,0.012325,0.019985,0.079428));
    squareBigCenterFreq(i) = ba.take(i+1,(222.231629,336.343611,555.041069,772.219662,911.261442,1096.774273,1203.102824,1404.011001,1540.420421,1628.944840,1994.925663,2197.106949,2533.354244,2672.470165,2821.185962,2967.436266,3342.189608,3479.723368,3631.705350,3806.882908,3910.078404,4173.779128,4292.901543,4519.288440,4631.053559,4855.405142,5079.828139,5182.392108,5328.291818,5636.817979,5983.629129,6104.013800,6184.547942,6602.614840,6718.247643,6835.844044,6882.254514,7220.774814,7283.263950,7407.672066));
    squareBigCenterGain(i) = ba.take(i+1,(0.581227,1.000000,0.167366,0.663184,0.127173,0.141424,0.232076,0.043454,0.168681,0.244908,0.193343,0.147543,0.062214,0.112585,0.383899,0.053139,0.104651,0.093521,0.029220,0.027646,0.020979,0.031887,0.052610,0.045000,0.023542,0.049426,0.029063,0.049062,0.015097,0.030925,0.032215,0.032146,0.050862,0.016861,0.019228,0.021428,0.021039,0.016988,0.017971,0.033065));
    squareBigCenterT60(i) = ba.take(i+1,(0.417113,0.247980,0.103791,0.131234,0.194764,0.130569,0.119615,0.103862,0.102785,0.079975,0.068871,0.120830,0.066401,0.064964,0.037272,0.057999,0.038435,0.036927,0.036970,0.069147,0.036096,0.049136,0.046268,0.029452,0.025066,0.020942,0.040527,0.027374,0.030986,0.019875,0.028632,0.032344,0.033709,0.024938,0.063833,0.032297,0.028473,0.025768,0.024877,0.026385));
    squareMidCenterFreq(i) = ba.take(i+1,(306.243391,426.381556,523.626577,637.813655,855.576390,992.200454,1098.863347,1226.922784,1545.950522,1704.776255,1839.163057,1983.832482,2190.814571,2445.409860,2672.907579,2782.642755,2983.850376,3268.313978,3458.038380,3590.850048,3783.143602,4123.002328,4250.310528,4408.442518,4531.916055,4865.211680,5248.660726,5387.011592,5780.593786,5880.954673,6095.291499,6432.605387,6596.403501,6820.917687,6907.055108,7030.829711,7337.883035,7602.739152,7926.836412,12248.346109));
    squareMidCenterGain(i) = ba.take(i+1,(1.000000,0.292592,0.554710,0.359595,0.041223,0.283249,0.182057,0.078666,0.898203,0.105021,0.037156,0.128493,0.302569,0.134081,0.090181,0.116704,0.041711,0.268207,0.079651,0.034128,0.296838,0.051959,0.124454,0.077586,0.084915,0.126202,0.019697,0.017643,0.077301,0.068367,0.023764,0.015924,0.022163,0.017954,0.029917,0.021449,0.019083,0.026076,0.017965,0.017018));
    squareMidCenterT60(i) = ba.take(i+1,(0.236061,0.179001,0.084440,0.165197,0.083836,0.055543,0.067029,0.054547,0.046195,0.055390,0.052760,0.079021,0.073765,0.059744,0.071780,0.050510,0.036326,0.036162,0.039268,0.031777,0.037990,0.072567,0.045824,0.026516,0.027684,0.034229,0.040867,0.037800,0.026894,0.023151,0.023381,0.021277,0.023057,0.020207,0.019672,0.016358,0.024757,0.017036,0.016600,0.010195));
    squareSmallCenterFreq(i) = ba.take(i+1,(375.008194,527.945584,775.002560,981.584262,1125.203834,1287.701398,1394.790058,1577.138666,1727.825366,1860.700270,1950.809059,2020.169063,2185.798879,2382.343983,2587.140759,2792.652735,2970.443881,3047.515682,3332.695227,3416.636119,3513.876416,3994.492272,4128.261485,4247.828681,4425.131019,4686.893838,4731.141754,5553.904972,5702.206389,6415.071752,6552.422900,6785.461423,6990.104471,7256.116515,7884.508530,8787.479467,10529.849261,10641.861836,12075.018006,17405.170359));
    squareSmallCenterGain(i) = ba.take(i+1,(0.991070,1.000000,0.515946,0.549911,0.115279,0.157813,0.463743,0.206985,0.131491,0.244054,0.160114,0.116274,0.161678,0.048839,0.071332,0.072482,0.097678,0.083415,0.057484,0.059259,0.034662,0.064479,0.061452,0.110522,0.045770,0.053318,0.051225,0.034998,0.051538,0.081127,0.039422,0.226011,0.053013,0.052620,0.041956,0.058147,0.061262,0.046118,0.038256,0.124151));
    squareSmallCenterT60(i) = ba.take(i+1,(0.102256,0.088069,0.026433,0.131998,0.050149,0.052897,0.066531,0.060621,0.045315,0.047413,0.058375,0.059098,0.029950,0.018537,0.024086,0.032984,0.036868,0.032205,0.023710,0.023430,0.015516,0.022773,0.025265,0.028872,0.019743,0.036565,0.032704,0.025234,0.028383,0.018356,0.018609,0.017534,0.018699,0.015762,0.022800,0.013151,0.012300,0.018114,0.013192,0.085670));
    
    // computing modes freq in function of x and y
    modeFreq(i) = zxd + (zxu-zxd)*select2(y<1,(y-1),y)
    with {
        zx0 = squareSmallCenterFreq(i) + (semSmallCenterFreq(i)-squareSmallCenterFreq(i))*x;
        zx1 = semSmallCenterFreq(i) + (roundSmallCenterFreq(i)-semSmallCenterFreq(i))*(x-1);
        zx2 = squareMidCenterFreq(i) + (semMidCenterFreq(i)-squareMidCenterFreq(i))*x;
        zx3 = semMidCenterFreq(i) + (roundMidCenterFreq(i)-semMidCenterFreq(i))*(x-1);
        zx4 = squareBigCenterFreq(i) + (semBigCenterFreq(i)-squareBigCenterFreq(i))*x;
        zx5 = semBigCenterFreq(i) + (roundBigCenterFreq(i)-semBigCenterFreq(i))*(x-1);
        zxd = select2(y<1,select2(x<1,zx3,zx2),select2(x<1,zx1,zx0));
        zxu = select2(y<1,select2(x<1,zx5,zx4),select2(x<1,zx3,zx2));
    };

    // computing modes gain in function of x and y
    modeGain(i) = zxd + (zxu-zxd)*select2(y<1,(y-1),y) 
    with {
        zx0 = squareSmallCenterGain(i) + (semSmallCenterGain(i)-squareSmallCenterGain(i))*x;
        zx1 = semSmallCenterGain(i) + (roundSmallCenterGain(i)-semSmallCenterGain(i))*(x-1);
        zx2 = squareMidCenterGain(i) + (semMidCenterGain(i)-squareMidCenterGain(i))*x;
        zx3 = semMidCenterGain(i) + (roundMidCenterGain(i)-semMidCenterGain(i))*(x-1);
        zx4 = squareBigCenterGain(i) + (semBigCenterGain(i)-squareBigCenterGain(i))*x;
        zx5 = semBigCenterGain(i) + (roundBigCenterGain(i)-semBigCenterGain(i))*(x-1);
        zxd = select2(y<1,select2(x<1,zx3,zx2),select2(x<1,zx1,zx0));
        zxu = select2(y<1,select2(x<1,zx5,zx4),select2(x<1,zx3,zx2));
    };

    // computing modes T60 in function of x and y
    modeT60(i) = zxd + (zxu-zxd)*select2(y<1,(y-1),y) //: min(0.1)
    with {
        zx0 = squareSmallCenterT60(i) + (semSmallCenterT60(i)-squareSmallCenterT60(i))*x;
        zx1 = semSmallCenterT60(i) + (roundSmallCenterT60(i)-semSmallCenterT60(i))*(x-1);
        zx2 = squareMidCenterT60(i) + (semMidCenterT60(i)-squareMidCenterT60(i))*x;
        zx3 = semMidCenterT60(i) + (roundMidCenterT60(i)-semMidCenterT60(i))*(x-1);
        zx4 = squareBigCenterT60(i) + (semBigCenterT60(i)-squareBigCenterT60(i))*x;
        zx5 = semBigCenterT60(i) + (roundBigCenterT60(i)-semBigCenterT60(i))*(x-1);
        zxd = select2(y<1,select2(x<1,zx3,zx2),select2(x<1,zx1,zx0));
        zxu = select2(y<1,select2(x<1,zx5,zx4),select2(x<1,zx3,zx2));
    };
};


//-------`(pm.)modularInterpBody`----------
// Bidirectional block implementing a modular string instrument resonator 
// (see [`modeInterpRes`](#pm.modeinterpres)).
//
// #### Usage
//
// ```
// chain(... : modularInterpBody(nModes,shape,scale) : ...)
// ```
//
// Where:
//
// * `nModes`: number of modeled modes (40 max)
// * `shape`: shape of the resonator (0: square, 1: square with rounded corners, 2: round)
// * `scale`: scale of the resonator (0: small, 1: medium, 2: large)
//----------------------------------
modularInterpBody(nModes,shape,scale) = _,modeInterpRes(nModes,shape,scale),_;


//-------`(pm.)modularInterpStringModel`----------
// String instrument model with a modular body (see 
// [`modeInterpRes`](#pm.modeinterpres) and 
// <https://ccrma.stanford.edu/~rmichon/3dPrintingModeling/>).
//
// #### Usage
//
// ```
// modularInterpStringModel(length,pluckPosition,shape,scale,bodyExcitation,stringExcitation) : _
// ```
//
// Where:
//
// * `stringLength`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `shape`: shape of the resonator (0: square, 1: square with rounded corners, 2: round)
// * `scale`: scale of the resonator (0: small, 1: medium, 2: large)
// * `bodyExcitation`: excitation signal for the body
// * `stringExcitation`: excitation signal for the string
//----------------------------------
modularInterpStringModel(stringLength,pluckPosition,shape,scale,bodyExcitation,stringExcitation) = endChain(egChain)*0.5
with {
    maxStringLength = maxLength;
    lengthTuning = 0.1; // tuned "by hand"
    stringL = stringLength-lengthTuning;
    nBodyModes = 40;
    egChain = chain(guitarNuts : steelString(stringL,pluckPosition,stringExcitation) :
    guitarBridge : inRightWave(bodyExcitation) : modularInterpBody(nBodyModes,shape,scale) : out);
};


//-------`(pm.)modularInterpInstr`----------
// String instrument with a modular body (see 
// [`modeInterpRes`](#pm.modeinterpres) and 
// <https://ccrma.stanford.edu/~rmichon/3dPrintingModeling/>).
//
// #### Usage
//
// ```
// modularInterpInstr(stringLength,pluckPosition,shape,scale,gain,tapBody,triggerString) : _
// ```
//
// Where:
//
// * `stringLength`: the length of the string in meters
// * `pluckPosition`: pluck position (0-1) (1 is on the bridge)
// * `shape`: shape of the resonator (0: square, 1: square with rounded corners, 2: round)
// * `scale`: scale of the resonator (0: small, 1: medium, 2: large)
// * `gain`: of the string excitation
// * `tapBody`: send an impulse in the body of the instrument where the string is connected (1 for on, 0 for off)
// * `triggerString`: trigger signal for the string (1 for on, 0 for off)
//----------------------------------
modularInterpInstr(stringLength,pluckPosition,shape,scale,gain,tapBody,triggerString) = 
(tapBody : ba.impulsify), pluckString(stringLength,1,1,1,gain,triggerString) : 
modularInterpStringModel(stringLength,pluckPosition,shape,scale);


//-------`(pm.)modularInterpInstr_ui_MIDI`----------
// Ready-to-use MIDI-enabled string instrument with a modular body (see 
// [`modeInterpRes`](#pm.modeinterpres) and 
// <https://ccrma.stanford.edu/~rmichon/3dPrintingModeling/>)
// with built-in UI.
//
// #### Usage
//
// ```
// modularInterpInstr_ui_MIDI : _
// ```
//----------------------------------
modularInterpInstr_ui_MIDI = modularInterpInstr(stringLength,pluckPosition,shape,scale,gain,tapBody,gate)*outGain
with {
    f = hslider("v:modularInterpInstr/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:modularInterpInstr/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:modularInterpInstr/h:[0]midi/[2]gain[style:knob]",0.8,0,1,0.01);
    s = hslider("v:modularInterpInstr/h:[0]midi/[3]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    shape = hslider("v:modularInterpInstr/h:[1]body/[0]shape[style:knob]",0,0,1,0.01)*2;
    scale = hslider("v:modularInterpInstr/h:[1]body/[1]scale[style:knob]",0,0,1,0.011)*2;
    tapBody = button("v:modularInterpInstr/h:[1]body/[2]tapBody");
    pluckPosition = hslider("v:modularInterpInstr/[2]pluckPosition[midi:ctrl 1]"
    ,0.8,0,1,0.01) : si.smoo;
    outGain = hslider("v:modularInterpInstr/[3]outGain",0.5,0,1,0.01);
    t = button("v:modularInterpInstr/[4]gate");
    gate = t+s : min(1);
    freq = f*bend;
    stringLength = freq : f2l;
};


//=============================Bowed String Instruments===================================
// Low and high level basic string instruments parts. Most of the elements in
// this section can be used in a bidirectional chain.
//========================================================================================

//-------`(pm.)bowTable`----------
// Extremely basic bow table that can be used to implement a wide range of
// bow types for many different bowed string instruments (violin, cello, etc.).
//
// #### Usage
//
// ```
// excitation : bowTable(offeset,slope) : _
// ```
//
// Where:
//
// * `excitation`: an excitation signal
// * `offset`: table offset
// * `slope`: table slope
//----------------------------------
bowTable(offset,slope) = pow(abs(sample) + 0.75, -4) : min(1)
with {
    sample = +(offset)*slope;
};


//-------`(pm.)violinBowTable`----------
// Violin bow table based on [`bowTable`](#bowtable).
//
// #### Usage
//
// ```
// bowVelocity : violinBowTable(bowPressure) : _
// ```
//
// Where:
//
// * `bowVelocity`: velocity of the bow/excitation signal (0-1)
// * `bowPressure`: bow pressure on the string (0-1)
//----------------------------------
violinBowTable(bowPressure) = bowTable(0,tableSlope)
with {
    tableSlope = 5 - (4*bowPressure);
};


//-------`(pm.)bowInteraction`----------
// Bidirectional block implementing the interaction of a bow in a
// [`chain`](#chain).
//
// #### Usage
//
// ```
// chain(... : stringSegment : bowInteraction(bowTable) : stringSegment : ...)
// ```
//
// Where:
//
// * `bowTable`: the bow table
//----------------------------------
bowInteraction(b) = (_,_ <: b,_,_ :> _,_),_;


//-------`(pm.)violinBow`----------
// Bidirectional block implementing a violin bow and its interaction with
// a string.
//
// #### Usage
//
// ```
// chain(... : stringSegment : violinBow(bowPressure,bowVelocity) : stringSegment : ...)
// ```
//
// Where:
//
// * `bowVelocity`: velocity of the bow / excitation signal (0-1)
// * `bowPressure`: bow pressure on the string (0-1)
//----------------------------------
violinBow(bowPressure,bowVelocity) = bowInteraction(bowSystem)
with {
    bowSystem = + : bowVelocity-_ <: *(violinBowTable(bowPressure)) <: _,_;
};


//-------`(pm.)violinBowedString`----------
// Violin bowed string bidirectional block with controllable bow position.
// Terminations are not implemented in this model.
//
// #### Usage
//
// ```
// chain(nuts : violinBowedString(stringLength,bowPressure,bowVelocity,bowPosition) : bridge)
// ```
//
// Where:
//
// * `stringLength`: the length of the string in meters
// * `bowVelocity`: velocity of the bow / excitation signal (0-1)
// * `bowPressure`: bow pressure on the string (0-1)
// * `bowPosition`: the position of the bow on the string (0-1)
//----------------------------------
violinBowedString(stringLength,bowPressure,bowVelocity,bowPosition) =
    chain(
        stringSegment(maxStringLength,ntbd) :
        violinBow(bowPressure,bowVelocity) :
        stringSegment(maxStringLength,btbd)
    )
    with {
        maxStringLength = maxLength;
        ntbd = stringLength*bowPosition; // upper portion of the string length
        btbd = stringLength*(1-bowPosition); // lower portion of the string length
    };


//-------`(pm.)violinNuts`----------
// Bidirectional block implementing simple violin nuts. This function is
// based on [`bridgeFilter`](#bridgefilter).
//
// #### Usage
//
// ```
// chain(violinNuts : stringSegment : ...)
// ```
//----------------------------------
violinNuts = lTermination(-bridgeFilter(0.6,0.1),basicBlock);


//-------`(pm.)violinBridge`----------
// Bidirectional block implementing a simple violin bridge. This function is
// based on [`bridgeFilter`](#bridgefilter).
//
// #### Usage
//
// ```
// chain(... : stringSegment : violinBridge
// ```
//----------------------------------
// TODO:
//     * reflectance is not implemented yet
violinBridge = rTermination(basicBlock,reflectance) : _,transmittance,_
with {
    reflectance = -bridgeFilter(0.2,0.9);
    transmittance = _;
};


//-------`(pm.)violinBody`----------
// Bidirectional block implementing a simple violin body (just a simple
// resonant lowpass filter).
//
// #### Usage
//
// ```
// chain(... : stringSegment : violinBridge : violinBody)
// ```
//----------------------------------
// TODO:
//     * reflectance is not implemented yet
violinBody = reflectance,transmittance,_
with {
    transmittance = fi.resonbp(500,2,1);
    reflectance = _;
};


//-------`(pm.)violinModel`----------
// Ready-to-use simple violin physical model. This model implements a single
// string. Additional strings should be created
// by making a polyphonic application out of this function. Pitch is changed
// by changing the length of the string (and not through a finger model).
//
// #### Usage
//
// ```
// violinModel(stringLength,bowPressure,bowVelocity,bridgeReflexion,
// bridgeAbsorption,bowPosition) : _
// ```
//
// Where:
//
// * `stringLength`: the length of the string in meters
// * `bowVelocity`: velocity of the bow / excitation signal (0-1)
// * `bowPressure`: bow pressure on the string (0-1))
// * `bowPosition`: the position of the bow on the string (0-1)
//----------------------------------
violinModel(stringLength,bowPressure,bowVelocity,bowPosition) =
    endChain(modelChain)
with {
    stringTuning = 0.08;
    stringL = stringLength-stringTuning;
    modelChain = chain(
        violinNuts :
        violinBowedString(stringL,bowPressure,bowVelocity,bowPosition) :
        violinBridge :
        violinBody :
        out
    );
};


//-------`(pm.)violin_ui`----------
// Ready-to-use violin physical model with built-in UI.
//
// #### Usage
//
// ```
// violinModel_ui : _
// ```
//----------------------------------
violin_ui = violinModel(stringLength,bowPress,bowVel,bowPos)*outGain
with {
    stringLength = hslider("v:violin/v:[0]string/[0]length",0.75,0,2,0.01) : si.smoo;
    bowVel = hslider("v:violin/v:[1]bow/[0]velocity",0,0,1,0.01) : si.smoo;
    bowPress = hslider("v:violin/v:[1]bow/[1]pressure",0.5,0,1,0.01) : si.smoo;
    bowPos = hslider("v:violin/v:[1]bow/[2]position",0.7,0,1,0.01) : si.smoo;
    outGain = hslider("v:violin/outGain",0.5,0,1,0.01);
};


//-------`(pm.)violin_ui_MIDI`----------
// Ready-to-use MIDI-enabled violin physical model with built-in UI.
//
// #### Usage
//
// ```
// violin_ui_MIDI : _
// ```
//----------------------------------
violin_ui_MIDI = violinModel(stringLength,bowPress,bowVel,bowPos)*outGain
with {
    f = hslider("v:violin/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:violin/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:violin/h:[0]midi/[2]gain[style:knob]
    ",0.6,0,1,0.01);
    envAttack = hslider("v:violin/h:[0]midi/[3]envAttack[style:knob]
    ",1,0,30,0.01)*0.001;
    s = hslider("v:violin/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    bowPress = hslider("v:violin/h:[1]bow/[0]pressure[style:knob]
    [midi:ctrl 1]",0.5,0,1,0.01) : si.smoo;
    bowPos = hslider("v:violin/h:[1]bow/[1]position[style:knob]
    ",0.7,0,1,0.01) : si.smoo;
    vibratoFreq = hslider("v:violin/h:[2]otherParams/[0]vibratoFrequency
    [style:knob]",6,1,10,0.01);
    vibratoGain = hslider("v:violin/h:[2]otherParams/[1]vibratoGain
    [style:knob]",0.5,0,1,0.01)*0.01;
    outGain = hslider("v:violin/h:[2]otherParams/[2]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:violin/[3]gate");

    gate = t+s : min(1);
    vibrato = 1+os.osc(vibratoFreq)*vibratoGain*envelope;
    freq = f*bend*vibrato;
    envelope = gate*gain : si.smooth(ba.tau2pole(envAttack));

    stringLength = freq : f2l;
    bowVel = envelope;
};


//=================================Wind Instruments=======================================
// Low and high level basic wind instruments parts. Most of the elements in
// this section can be used in a bidirectional chain.
//========================================================================================

//-------`(pm.)openTube`----------
// A tube segment without terminations (same as [`stringSegment`](#stringsegment)).
//
// #### Usage
//
// ```
// chain(A : openTube(maxLength,length) : B)
// ```
//
// Where:
//
// * `maxLength`: the maximum length of the tube in meters (should be static)
// * `length`: the length of the tube in meters
//----------------------------------
openTube = stringSegment;


//-------`(pm.)reedTable`----------
// Extremely basic reed table that can be used to implement a wide range of
// single reed types for many different instruments (saxophone, clarinet, etc.).
//
// #### Usage
//
// ```
// excitation : reedTable(offeset,slope) : _
// ```
//
// Where:
//
// * `excitation`: an excitation signal
// * `offset`: table offset
// * `slope`: table slope
//----------------------------------
reedTable(offset,slope) = reedTable : min(1) : max(-1)
with {
    reedTable = *(slope) + offset;
};


//-------`(pm.)fluteJetTable`----------
// Extremely basic flute jet table.
//
// #### Usage
//
// ```
// excitation : fluteJetTable : _
// ```
//
// Where:
//
// * `excitation`: an excitation signal
//----------------------------------
fluteJetTable = _ <: *(* : -(1)) : clipping
with {
    clipping = min(1) : max(-1);
};


//-------`(pm.)brassLipsTable`----------
// Simple brass lips/mouthpiece table. Since this implementation is very basic
// and that the lips and tube of the instrument are coupled to each other, the
// length of that tube must be provided here.
//
// #### Usage
//
// ```
// excitation : brassLipsTable(tubeLength,lipsTension) : _
// ```
//
// Where:
//
// * `excitation`: an excitation signal (can be DC)
// * `tubeLength`: length in meters of the tube connected to the mouthpiece
// * `lipsTension`: tension of the lips (0-1) (default: 0.5)
//----------------------------------
brassLipsTable(tubeLength,lipsTension) = *(0.03) : lipFilter <: * : clipping
with {
    clipping = min(1) : max(-1);
    freq = (tubeLength : l2f)*pow(4,(2*lipsTension)-1);
    filterR = 0.997;
    a1 = -2*filterR*cos(ma.PI*2*freq/ma.SR);
    lipFilter = fi.tf2(1,0,0,a1,pow(filterR,2)); // resonance with same freq as tube
};


//-------`(pm.)clarinetReed`----------
// Clarinet reed based on [`reedTable`](#reedtable) with controllable
// stiffness.
//
// #### Usage
//
// ```
// excitation : clarinetReed(stiffness) : _
// ```
//
// Where:
//
// * `excitation`: an excitation signal
// * `stiffness`: reed stiffness (0-1)
//----------------------------------
clarinetReed(stiffness) = reedTable(0.7,tableSlope)
with {
    tableSlope = -0.44 + 0.26*stiffness;
};


//-------`(pm.)clarinetMouthPiece`----------
// Bidirectional block implementing a clarinet mouthpiece as well as the various
// interactions happening with traveling waves. This element is ready to be
// plugged to a tube...
//
// #### Usage
//
// ```
// chain(clarinetMouthPiece(reedStiffness,pressure) : tube : etc.)
// ```
//
// Where:
//
// * `pressure`: the pressure of the air flow (DC) created by the virtual performer (0-1).
// This can also be any kind of signal that will directly injected in the mouthpiece
// (e.g., breath noise, etc.).
// * `reedStiffness`: reed stiffness (0-1)
//----------------------------------
clarinetMouthPiece(reedStiffness,pressure) = lTermination(reedInteraction,in(pressure))
with {
    reedInteraction = *(-1) <: *(clarinetReed(reedStiffness));
};


//-------`(pm.)brassLips`----------
// Bidirectional block implementing a brass mouthpiece as well as the various
// interactions happening with traveling waves. This element is ready to be
// plugged to a tube...
//
// #### Usage
//
// ```
// chain(brassLips(tubeLength,lipsTension,pressure) : tube : etc.)
// ```
//
// Where:
//
// * `tubeLength`: length in meters of the tube connected to the mouthpiece
// * `lipsTension`: tension of the lips (0-1) (default: 0.5)
// * `pressure`: the pressure of the air flow (DC) created by the virtual performer (0-1).
// This can also be any kind of signal that will directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
brassLips(tubeLength,lipsTension,pressure) = lTermination(mpInteraction,basicBlock)
with {
    absorption = *(0.85); // absorption coefficient
    p = pressure*0.3; // scaling pressure
    mpInteraction = absorption <:
    (p-_ : brassLipsTable(tubeLength,lipsTension) <: *(p),1-_),_ : _,* : + : fi.dcblocker;
};


//-------`(pm.)fluteEmbouchure`----------
// Bidirectional block implementing a flute embouchure as well as the various
// interactions happening with traveling waves. This element is ready to be
// plugged between tubes segments...
//
// #### Usage
//
// ```
// chain(... : tube : fluteEmbouchure(pressure) : tube : etc.)
// ```
//
// Where:
//
// * `pressure`: the pressure of the air flow (DC) created by the virtual
// performer (0-1).
// This can also be any kind of signal that will directly injected in the
// mouthpiece (e.g., breath noise, etc.).
//----------------------------------
fluteEmbouchure(pressure) =
(_ <: _,_),_,_ : _,*(0.5)+(pressure-_*0.5 : fluteJetTable),_;


//-------`(pm.)wBell`----------
// Generic wind instrument bell bidirectional block that should be placed at
// the end of a [`chain`](#chain).
//
// #### Usage
//
// ```
// chain(... : wBell(opening))
// ```
//
// Where:
//
// * `opening`: the "opening" of bell (0-1)
//----------------------------------
wBell(opening) = rTermination(basicBlock,si.smooth(opening));


//-------`(pm.)fluteHead`----------
// Simple flute head implementing waves reflexion.
//
// #### Usage
//
// ```
// chain(fluteHead : tube : ...)
// ```
//----------------------------------
fluteHead = lTermination(*(absorption),basicBlock)
with {
    absorption = 0.95; // same as for foot
};


//-------`(pm.)fluteFoot`----------
// Simple flute foot implementing waves reflexion and dispersion.
//
// #### Usage
//
// ```
// chain(... : tube : fluteFoot)
// ```
//----------------------------------
fluteFoot = rTermination(basicBlock,*(absorption) : dispersion)
with {
    dispersion = si.smooth(0.7); // just a simple lowpass
    absorption = 0.95;           // same as for head
};


//-------`(pm.)clarinetModel`----------
// A simple clarinet physical model without tone holes (pitch is changed by
// changing the length of the tube of the instrument).
//
// #### Usage
//
// ```
// clarinetModel(length,pressure,reedStiffness,bellOpening) : _
// ```
//
// Where:
//
// * `tubeLength`: the length of the tube in meters
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will directly injected in the mouthpiece
// (e.g., breath noise, etc.).
// * `reedStiffness`: reed stiffness (0-1)
// * `bellOpening`: the opening of bell (0-1)
//----------------------------------
clarinetModel(tubeLength,pressure,reedStiffness,bellOpening) = endChain(modelChain)
with {
    lengthTuning = 0.05; // empirical adjustment of the tuning of the tube
    maxTubeLength = maxLength;
    tunedLength = tubeLength/2-lengthTuning; // not really sure why we had to shift octave here
    modelChain =
        chain(
            clarinetMouthPiece(reedStiffness,pressure) :
            openTube(maxTubeLength,tunedLength) :
            wBell(bellOpening) : out
        );
};


//-------`(pm.)clarinetModel_ui`----------
// Same as [`clarinetModel`](#clarinetModel) but with a built-in UI. This function
// doesn't implement a virtual "blower", thus `pressure` remains an argument here.
//
// #### Usage
//
// ```
// clarinetModel_ui(pressure) : _
// ```
//
// Where:
//
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will be directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
clarinetModel_ui(pressure) = clarinetModel(tubeLength,pressure,reedStiffness,bellOpening)*outGain
with {
    tubeLength = hslider("v:clarinetModel/[0]tubeLength",0.8,0.01,3,0.01) : si.smoo;
    reedStiffness = hslider("v:clarinetModel/[1]reedStiffness",0.5,0,1,0.01);
    bellOpening = hslider("v:clarinetModel/[2]bellOpening",0.5,0,1,0.01);
    outGain = hslider("v:clarinetModel/[3]outGain",0.5,0,1,0.01);
};


//-------`(pm.)clarinet_ui`----------
// Ready-to-use clarinet physical model with built-in UI based on
// [`clarinetModel`](#clarinetmodel).
//
// #### Usage
//
// ```
// clarinet_ui : _
// ```
//----------------------------------
clarinet_ui = hgroup("clarinet",blower_ui : clarinetModel_ui);


//-------`(pm.)clarinet_ui_MIDI`----------
// Ready-to-use MIDI compliant clarinet physical model with built-in UI.
//
// #### Usage
//
// ```
// clarinet_ui_MIDI : _
// ```
//----------------------------------
clarinet_ui_MIDI =
clarinetModel(tubeLength,blow,reedStiffness,bellOpening)*outGain
with {
    f = hslider("v:clarinet/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:clarinet/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:clarinet/h:[0]midi/[2]gain[style:knob]
    ",0.6,0,1,0.01);
    envAttack = hslider("v:clarinet/h:[0]midi/[3]envAttack[style:knob]
    ",1,0,30,0.01)*0.001;
    s = hslider("v:clarinet/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    reedStiffness = hslider("v:clarinet/h:[1]otherParams/[0]reedStiffness
    [midi:ctrl 1][style:knob]",0.5,0,1,0.01);
    bellOpening = hslider("v:clarinet/h:[1]otherParams/[1]bellOpening
    [midi:ctrl 3][style:knob]",0.5,0,1,0.01);
    vibratoFreq = hslider("v:clarinet/h:[1]otherParams/[2]vibratoFreq
    [style:knob]",5,1,10,0.01);
    vibratoGain = hslider("v:clarinet/h:[1]otherParams/[3]vibratoGain
    [style:knob]",0.25,0,1,0.01)*0.01;
    outGain = hslider("v:clarinet/h:[1]otherParams/[4]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:clarinet/[2]gate");

    gate = t+s : min(1);
    vibrato = 1+os.osc(vibratoFreq)*vibratoGain*envelope;
    freq = f*bend*vibrato;
    envelope = gate*gain : si.smooth(ba.tau2pole(envAttack));

    tubeLength = freq : f2l;
    pressure = envelope; // TODO: double vibrato here!!
    blow = blower(pressure,0.05,2000,vibratoFreq,vibratoGain);
};


//-------`(pm.)brassModel`----------
// A simple generic brass instrument physical model without pistons
// (pitch is changed by changing the length of the tube of the instrument).
// This model is kind of hard to control and might not sound very good if
// bad parameters are given to it...
//
// #### Usage
//
// ```
// brassModel(tubeLength,lipsTension,mute,pressure) : _
// ```
//
// Where:
//
// * `tubeLength`: the length of the tube in meters
// * `lipsTension`: tension of the lips (0-1) (default: 0.5)
// * `mute`: mute opening at the end of the instrument (0-1) (default: 0.5)
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
brassModel(tubeLength,lipsTension,mute,pressure) = endChain(brassChain)
with {
    maxTubeLength = maxLength;
    lengthTuning = 0; // Not that important for that one because of lips tension
    tunedLength = tubeLength + lengthTuning;
    brassChain = chain(brassLips(tunedLength,lipsTension,pressure) : openTube(maxTubeLength,tunedLength) : wBell(mute) : out);
};


//-------`(pm.)brassModel_ui`----------
// Same as [`brassModel`](#brassModel) but with a built-in UI. This function
// doesn't implement a virtual "blower", thus `pressure` remains an argument here.
//
// #### Usage
//
// ```
// brassModel_ui(pressure) : _
// ```
//
// Where:
//
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will be directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
brassModel_ui(pressure) = brassModel(tubeLength,lipsTension,mute,pressure)
with {
    tubeLength = hslider("v:brassModel/[1]tubeLength",0.5,0.01,2.5,0.01) : si.smoo;
    lipsTension = hslider("v:brassModel/[2]lipsTension",0.5,0,1,0.01) : si.smoo;
    mute = hslider("v:brassModel/[3]mute",0.5,0,1,0.01) : si.smoo;
};


//-------`(pm.)brass_ui`----------
// Ready-to-use brass instrument physical model with built-in UI based on
// [`brassModel`](#brassmodel).
//
// #### Usage
//
// ```
// brass_ui : _
// ```
//----------------------------------
brass_ui = hgroup("brass",blower_ui : brassModel_ui);


//-------`(pm.)brass_ui_MIDI`----------
// Ready-to-use MIDI-controllable brass instrument physical model with built-in UI.
//
// #### Usage
//
// ```
// brass_ui_MIDI : _
// ```
//----------------------------------
brass_ui_MIDI = brassModel(tubeLength,lipsTension,mute,pressure)*outGain
with {
    f = hslider("v:brass/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:brass/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:brass/h:[0]midi/[2]gain[style:knob]
    ",0.5,0,1,0.01);
    envAttack = hslider("v:brass/h:[0]midi/[3]envAttack[style:knob]
    ",1,0,30,0.01)*0.001;
    s = hslider("v:brass/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    lipsTension = hslider("v:brass/h:[1]otherParams/[0]lipsTension[style:knob]
    [midi:ctrl 1]",0.5,0,1,0.01) : si.smoo;
    mute = hslider("v:brass/h:[1]otherParams/[1]mute[style:knob]
    ",0.5,0,1,0.01) : si.smoo;
    vibratoFreq = hslider("v:brass/h:[1]otherParams/[2]vibratoFreq[style:knob]
    ",5,1,10,0.01);
    vibratoGain = hslider("v:brass/h:[1]otherParams/[3]vibratoGain[style:knob]
    ",0.5,0,1,0.01)*0.04;
    outGain = hslider("v:brass/h:[1]otherParams/[4]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:brass/[2]gate");

    gate = t+s : min(1);
    vibrato = 1+os.osc(vibratoFreq)*vibratoGain*envelope;
    freq = f*bend;
    envelope = gate*gain : si.smooth(ba.tau2pole(envAttack));

    tubeLength = freq : f2l;
    pressure = envelope*vibrato;
};


//-------`(pm.)fluteModel`----------
// A simple generic flute instrument physical model without tone holes
// (pitch is changed by changing the length of the tube of the instrument).
//
// #### Usage
//
// ```
// fluteModel(tubeLength,mouthPosition,pressure) : _
// ```
//
// Where:
//
// * `tubeLength`: the length of the tube in meters
// * `mouthPosition`: position of the mouth on the embouchure (0-1) (default: 0.5)
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
// TODO: this model is out of tune and we're not really sure why
fluteModel(tubeLength,mouthPosition,pressure) = endChain(fluteChain) : fi.dcblocker
with {
    maxTubeLength = maxLength;
    tubeTuning = 0.27; // set "by hand"
    tLength = tubeLength+tubeTuning; // global tube length
    embouchurePos = 0.27 + (mouthPosition-0.5)*0.4; // position of the embouchure on the tube
    tted = tLength*embouchurePos; // head to embouchure distance
    eted = tLength*(1-embouchurePos); // embouchure to foot distance
    fluteChain = chain(fluteHead : openTube(maxTubeLength,tted) : fluteEmbouchure(pressure) : openTube(maxTubeLength,eted) : fluteFoot : out);
};


//-------`(pm.)fluteModel_ui`----------
// Same as [`fluteModel`](#fluteModel) but with a built-in UI. This function
// doesn't implement a virtual "blower", thus `pressure` remains an argument here.
//
// #### Usage
//
// ```
// fluteModel_ui(pressure) : _
// ```
//
// Where:
//
// * `pressure`: the pressure of the air flow created by the virtual performer (0-1).
// This can also be any kind of signal that will be directly injected in the mouthpiece
// (e.g., breath noise, etc.).
//----------------------------------
fluteModel_ui(pressure) =
fluteModel(tubeLength,mouthPosition,pressure)*outGain
with {
    tubeLength = hslider("v:fluteModel/[0]tubeLength",0.8,0.01,3,0.01) : si.smoo;
    mouthPosition = hslider("v:fluteModel/[1]mouthPosition",0.5,0,1,0.01) : si.smoo;
    outGain = hslider("v:fluteModel/[2]outGain",0.5,0,1,0.01);
};


//-------`(pm.)flute_ui`----------
// Ready-to-use flute physical model with built-in UI based on
// [`fluteModel`](#flutemodel).
//
// #### Usage
//
// ```
// flute_ui : _
// ```
//----------------------------------
flute_ui = hgroup("flute",blower_ui : fluteModel_ui);


//-------`(pm.)flute_ui_MIDI`----------
// Ready-to-use MIDI-controllable flute physical model with built-in UI.
//
// #### Usage
//
// ```
// flute_ui_MIDI : _
// ```
//----------------------------------
flute_ui_MIDI = fluteModel(tubeLength,mouthPosition,blow)*outGain
with {
    f = hslider("v:flute/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:flute/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:flute/h:[0]midi/[2]gain[style:knob]
    ",0.9,0,1,0.01);
    envAttack = hslider("v:flute/h:[0]midi/[3]envAttack[style:knob]
    ",1,0,30,0.01)*0.001;
    s = hslider("v:flute/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    mouthPosition = hslider("v:flute/h:[1]otherParams/[0]mouthPosition
    [style:knob][midi:ctrl 1]",0.5,0,1,0.01)
    : si.smoo;
    vibratoFreq = hslider("v:flute/h:[1]otherParams/[1]vibratoFreq[style:knob]
    ",5,1,10,0.01);
    vibratoGain = hslider("v:flute/h:[1]otherParams/[2]vibratoGain[style:knob]
    ",0.5,0,1,0.01)*0.04;
    outGain = hslider("v:flute/h:[1]otherParams/[3]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:flute/[2]gate");

    gate = t+s : min(1);
    freq = f*bend;
    envelope = gate*gain : si.smooth(ba.tau2pole(envAttack));

    tubeLength = freq : f2l;
    pressure = envelope;
    blow = blower(pressure,0.05,2000,vibratoFreq,vibratoGain);
};


//=====================================Exciters===========================================
// Various kind of excitation signal generators.
//========================================================================================

//-------`(pm.)impulseExcitation`--------------
// Creates an impulse excitation of one sample.
//
// #### Usage
//
// ```
// gate = button('gate');
// impulseExcitation(gate) : chain;
// ```
//
// Where:
//
// * `gate`: a gate button
//--------------------------------------
impulseExcitation(trigger) = trigger : ba.impulsify;


//-------`(pm.)strikeModel`--------------
// Creates a filtered noise excitation.
//
// #### Usage
//
// ```
// gate = button('gate');
// strikeModel(LPcutoff,HPcutoff,sharpness,gain,gate) : chain;
// ```
//
// Where:
//
// * `HPcutoff`: highpass cutoff frequency
// * `LPcutoff`: lowpass cutoff frequency
// * `sharpness`: sharpness of the attack and release (0-1)
// * `gain`: gain of the excitation
// * `gate`: a gate button/trigger signal (0/1)
//--------------------------------------
strikeModel(HPcutoff,LPcutoff,sharpness,gain,trigger) =
no.noise : fi.highpass(2,HPcutoff) : fi.lowpass(2,LPcutoff) :
*(en.ar(att,rel,trigger))*gain
with {
    att = 0.002*sharpness;
    rel = att;
};


//-------`(pm.)strike`--------------
// Strikes generator with controllable excitation position.
//
// #### Usage
//
// ```
// gate = button('gate');
// strike(exPos,sharpness,gain,gate) : chain;
// ```
//
// Where:
//
// * `exPos`: excitation position wiht 0: for max low freqs and 1: for max high
// freqs. So, on membrane for example, 0 would be the middle and 1 the edge
// * `sharpness`: sharpness of the attack and release (0-1)
// * `gain`: gain of the excitation
// * `gate`: a gate button/trigger signal (0/1)
//--------------------------------------
strike(exPos,sharpness,gain,trigger) = strikeModel(HPcutoff,LPcutoff,sharpness,gain,trigger)
with {
    HPcutoff = 40+exPos*500;
    LPcutoff = 500+exPos*15000;
};


//-------`(pm.)pluckString`--------------
// Creates a plucking excitation signal.
//
// #### Usage
//
// ```
// trigger = button('gate');
// pluckString(stringLength,cutoff,maxFreq,sharpness,trigger)
// ```
//
// Where:
//
// * `stringLength`: length of the string to pluck
// * `cutoff`: cutoff ratio (1 for default)
// * `maxFreq`: max frequency ratio (1 for default)
// * `sharpness`: sharpness of the attack and release (1 for default)
// * `gain`: gain of the excitation (0-1)
// * `trigger`: trigger signal (1 for on, 0 for off)
//--------------------------------------
pluckString(stringLength,cutoff,maxFreq,sharpness,gain,trigger) =
no.noise : fi.lowpass(2,cutoffreq) : *(en.ar(att,rel,trigger))*gain
with {
    freq = stringLength : l2f;
    maxF = 2000*maxFreq;
    att = 0.002*sharpness*pow((1 - freq/maxF),2);
    rel = att;
    cutoffreq = freq*5*cutoff;
};


//-------`(pm.)blower`--------------
// A virtual blower creating a DC signal with some breath noise in it.
//
// #### Usage
//
// ```
// blower(pressure,breathGain,breathCutoff) : _
// ```
//
// Where:
//
// * `pressure`: pressure (0-1)
// * `breathGain`: breath noise gain (0-1) (recommended: 0.005)
// * `breathCutoff`: breath cuttoff frequency (Hz) (recommended: 2000)
//--------------------------------------
blower(pressure,breathGain,breathCutoff,vibratoFreq,vibratoGain) = pressure + vibrato + breathNoise
with {
    vibrato = os.osc(vibratoFreq)*vibratoGain;
    breathNoise = no.noise : fi.lowpass(2,breathCutoff) : *(pressure*breathGain);
};


//-------`(pm.)blower_ui`--------------
// Same as [`blower`](#blower) but with a built-in UI.
//
// #### Usage
//
// ```
// blower : somethingToBeBlown
// ```
//--------------------------------------
blower_ui = blower(pressure,breathGain,breathCutoff,vibratoFreq,vibratoGain)
with {
    pressure = hslider("v:blower/[0]pressure",0,0,1,0.01) : si.smoo;
    breathGain = hslider("v:blower/[1]breathGain",0.1,0,1,0.01)*0.05;
    breathCutoff = hslider("v:blower/[2]breathCutoff",2000,20,20000,0.1);
    vibratoFreq = hslider("v:blower/[3]vibratoFreq",5,0.1,10,0.1);
    vibratoGain = hslider("v:blower/[4]vibratoGain",0.25,0,1,0.01)*0.03;
};


//============================Modal Percussions===========================================
// High and low level functions for modal synthesis of percussion instruments.
//========================================================================================

//-------`(pm.)djembeModel`----------
// Dirt-simple djembe modal physical model. Mode parameters are empirically
// calculated and don't correspond to any measurements or 3D model. They
// kind of sound good though :).
//
// #### Usage
//
// ```
// excitation : djembeModel(freq)
// ```
//
// Where:
//
// * `excitation`: excitation signal
// * `freq`: fundamental frequency of the bar
//----------------------------------
djembeModel(freq) = _ <: par(i,nModes,modeFilter(modeFreqs(i),modeT60s(i),modeGains(i))) :> /(nModes)
with {
    nModes = 20;
    theta = 0; // angle
    modeFreqs(i) = freq + 200*i;
    modeT60s(i) = (nModes-i)*0.03;
    modeGains(i) = cos((i+1)*theta)/float(i+1)*(1/(i+1));
};


//-------`(pm.)djembe`----------
// Dirt-simple djembe modal physical model. Mode parameters are empirically
// calculated and don't correspond to any measurements or 3D model. They
// kind of sound good though :).
//
// This model also implements a virtual "exciter".
//
// #### Usage
//
// ```
// djembe(freq,strikePosition,strikeSharpness,gain,trigger)
// ```
//
// Where:
//
// * `freq`: fundamental frequency of the model
// * `strikePosition`: strike position (0 for the middle of the membrane and
// 1 for the edge)
// * `strikeSharpness`: sharpness of the strike (0-1, default: 0.5)
// * `gain`: gain of the strike
// * `trigger`: trigger signal (0: off, 1: on)
//----------------------------------
djembe(freq,strikePosition,strikeSharpness,gain,trigger) =
strike(strikePosition,strikeSharpness,gain,trigger) : djembeModel(freq);

//-------`(pm.)djembe_ui_MIDI`----------
// Simple MIDI controllable djembe physical model with built-in UI.
//
// #### Usage
//
// ```
// djembe_ui_MIDI : _
// ```
//----------------------------------
djembe_ui_MIDI =
djembe(freq,strikePosition,strikeSharpness,gain,gate)*outGain
with {
    freq = hslider("v:djembe/h:[0]midi/[0]freq[style:knob]",60,50,100,0.01);
    gain = hslider("v:djembe/h:[0]midi/[2]gain[style:knob]",1,0,1,0.01);
    strikePosition = hslider("v:djembe/h:[1]otherParams/[0]strikePosition
    [midi:ctrl 1][style:knob]",0.5,0,1,0.01);
    strikeSharpness = hslider("v:djembe/h:[1]otherParams/[1]strikeSharpness[style:knob]",0.5,0.01,5,0.01);
    outGain = hslider("v:djembe/h:[1]otherParams/[2]outGain
    [style:knob]",1,0,1,0.01);
    gate = button("v:djembe/[3]gate");
};


//-------`(pm.)marimbaBarModel`----------
// Generic marimba tone bar modal model.
//
// This model was generated using
// `mesh2faust` from a 3D CAD model of a marimba tone bar
// (`libraries/modalmodels/marimbaBar`). The corresponding CAD model is that
// of a C2 tone bar (original fundamental frequency: ~65Hz). While
// `marimbaBarModel` allows to translate the harmonic content of the generated
// sound by providing a frequency (`freq`), mode transposition has limits and
// the model will sound less and less like a marimba tone bar as it
// diverges from C2. To make an accurate model of a marimba, we'd want to have
// an independent model for each bar...
//
// This model contains 5 excitation positions going linearly from the center
// bottom to the center top of the bar. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : marimbaBarModel(freq,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: excitation signal
// * `freq`: fundamental frequency of the bar
// * `exPos`: excitation position (0-4)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
marimbaBarModel(freq,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nModes = 50;
    nExPos = 5;
    modesFreqRatios(n) = ba.take(n+1,(1,3.31356,3.83469,8.06313,9.44778,14.1169,18.384,21.0102,26.1775,28.9944,37.0728,37.8703,40.0634,47.6439,51.019,52.43,58.286,63.5486,65.3628,66.9587,74.5301,78.692,80.8375,89.978,92.9661,95.1914,97.4807,110.62,112.069,113.826,119.356,127.045,129.982,132.259,133.477,144.549,149.438,152.033,153.166,155.597,158.183,168.105,171.863,174.464,178.937,181.482,185.398,190.369,192.19,195.505));
    modesFreqs(i) = freq*modesFreqRatios(i);
    modesGains(p,n) = waveform{1,0.776725,0.625723,0.855223,0.760159,0.698373,0.768011,0.641127,0.244034,0.707754,0.634013,0.247527,0.660849,0.450396,0.567783,0.106361,0.716814,0.66392,0.291208,0.310599,0.801495,0.635292,0.307435,0.874124,0.497668,0.487088,0.459115,0.733455,0.541818,0.441318,0.31392,0.40309,0.685353,0.60314,0.400552,0.453511,0.634386,0.291547,0.131605,0.368507,0.839907,0.60216,0.288296,0.57967,0.0242493,0.262746,0.368588,0.890284,0.408963,0.556072,0.884427,0.83211,0.612015,0.757176,0.919477,1,0.827963,0.89241,0.0357408,0.480789,0.752872,0.0546301,0.235937,0.362938,0.444472,0.101751,0.703418,0.453136,0.316629,0.490394,0.982508,0.551622,0.602009,0.666957,0.77683,0.905662,0.0987197,0.402968,0.829452,0.307645,0.64048,0.983971,0.584205,0.650365,0.334447,0.58357,0.540191,0.672534,0.245712,0.687298,0.883058,0.79295,0.600619,0.572682,0.122612,0.388248,0.290658,0.380255,0.290967,0.567819,0.0737721,0.42099,0.0786578,0.393995,0.268983,0.260614,0.494086,0.238026,0.0987824,0.277879,0.440563,0.0770212,0.450591,0.128137,0.0368275,0.128699,0.329605,0.374512,0.36359,0.272594,0.379052,0.305241,0.0741129,0.345728,0.29935,0.221284,0.0261391,0.293202,0.361885,0.11433,0.239005,0.434156,0.329583,0.21946,0.284175,0.198555,0.431976,0.302985,1,0.146221,0.140701,0.264243,0.185997,0.426322,0.30478,0.34399,0.19543,0.386955,0.1876,0.172812,0.0434115,0.303761,0.069454,0.453943,0.832451,0.317817,0.940601,1,0.180658,0.737921,0.832297,0.402352,0.126786,0.594398,0.485455,0.32447,0.365102,0.777922,0.588272,0.401353,0.610735,0.158693,0.0746072,0.825099,0.925459,0.65377,0.260792,0.719384,0.559908,0.37259,0.360035,0.622939,0.210271,0.444595,0.311286,0.464309,0.557231,0.52408,0.0701056,0.320749,0.19446,0.727609,0.522062,0.394004,0.235035,0.395646,0.494796,0.517318,0.109752,0.692849,0.00632009,0.0207583,0.00306107,0.0637191,0.081661,0.03511,0.127814,0.202294,0.0764145,0.263127,0.400199,0.267278,0.633385,1,0.739902,0.413763,0.41811,0.612715,0.672374,0.339674,0.21172,0.459645,0.1025,0.32589,0.148154,0.265442,0.0974305,0.286438,0.275213,0.109111,0.575089,0.370283,0.29411,0.259826,0.0648719,0.583418,0.282663,0.182004,0.117421,0.417727,0.16965,0.24853,0.122819,0.185486,0.0433618,0.373849,0.252768,0.195103,0.0927835,0.166543},int(p*nModes+n) : rdtable : select2(modesFreqs(n)<(ma.SR/2-1),0);
    modesT60s(i) = t60*pow(1-(modesFreqRatios(i)/195.955)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)marimbaResTube`----------
// Simple marimba resonance tube.
//
// #### Usage
//
// ```
// marimbaResTube(tubeLength,excitation)
// ```
//
// Where:
//
// * `tubeLength`: the length of the tube in meters
// * `excitation`: the excitation signal (audio in)
//----------------------------------
marimbaResTube(tubeLength,excitation) = endChain(tubeChain)
with {
    maxTubeLength = maxLength;
    lengthTuning = 0.04;
    tunedLength = tubeLength-lengthTuning;
    endTubeReflexion = si.smooth(0.95)*0.99;
    tubeChain =
        chain(
            in(excitation) :
            terminations(endTubeReflexion,
                openTube(maxTubeLength,tunedLength),
                endTubeReflexion) :
            out
        );
};


//-------`(pm.)marimbaModel`----------
// Simple marimba physical model implementing a single tone bar connected to
// tube. This model is scalable and can be adapted to any size of bar/tube
// (see [`marimbaBarModel`](#marimbabarmodel) to know more about the
// limitations of this type of system).
//
// #### Usage
//
// ```
// excitation : marimbaModel(freq,exPos) : _
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `freq`: the frequency of the bar/tube couple
// * `exPos`: excitation position (0-4)
//----------------------------------
marimbaModel(freq,exPos) =
marimbaBarModel(freq,exPos,maxT60,T60Decay,T60Slope) : marimbaResTube(resTubeLength)
with {
    resTubeLength = freq : f2l;
    maxT60 = 0.1;
    T60Decay = 1;
    T60Slope = 5;
};


//-------`(pm.)marimba`----------
// Simple marimba physical model implementing a single tone bar connected to
// tube. This model is scalable and can be adapted to any size of bar/tube
// (see [`marimbaBarModel`](#marimbabarmodel) to know more about the
// limitations of this type of system).
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// marimba(freq,strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `freq`: the frequency of the bar/tube couple
// * `strikePosition`: strike position (0-4)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
marimba(freq,strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
marimbaModel(freq,strikePosition);


//-------`(pm.)marimba_ui_MIDI`----------
// Simple MIDI controllable marimba physical model with built-in UI
// implementing a single tone bar connected to
// tube. This model is scalable and can be adapted to any size of bar/tube
// (see [`marimbaBarModel`](#marimbabarmodel) to know more about the
// limitations of this type of system).
//
// #### Usage
//
// ```
// marimba_ui_MIDI : _
// ```
//----------------------------------
marimba_ui_MIDI =
marimba(freq,strikePosition,strikeCutoff,strikeSharpness,gain,gate)*outGain
with {
    freq = hslider("v:marimba/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    gain = hslider("v:marimba/h:[0]midi/[2]gain[style:knob]",1,0,1,0.01);
    strikePosition = nentry("v:marimba/h:[1]otherParams/[0]strikePosition
    [midi:ctrl 1]",0,0,4,1);
    strikeCutoff = hslider("v:marimba/h:[1]otherParams/[1]strikeCutOff
    [midi:ctrl 1][style:knob]",6500,20,20000,1);
    strikeSharpness = hslider("v:marimba/h:[1]otherParams/[2]strikeSharpness
    [style:knob]",0.5,0.01,5,0.01);
    outGain = hslider("v:marimba/h:[1]otherParams/[2]outGain
    [style:knob]",0.8,0,1,0.01);
    gate = button("v:marimba/[3]gate");
};


//-------`(pm.)churchBellModel`----------
// Generic church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/churchBell`.
//
// Modeled after T. Rossing and R. Perrin, Vibrations of Bells, Applied
// Acoustics 2, 1987.
//
// Model height is 301 mm.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : churchBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
churchBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(451.918,455,864.643,871.402,1072.47,1073.98,1292.23,1292.48,1504.6,1532.41,1646.2,1647,1677.83,1678.13,1866.63,1882.08,1985.2,1989.87,2114.93,2356.81,2444.32,2446.22,2528.77,2530.92,2668.86,2669.63,2738.06,2749.03,2750.53,2753.11,2827.19,2842.32,2992.56,2996.84,3172.61,3330.79,3390.33,3403.1,3516.43,3538.54,3582.79,3583.82,3730.04,3739.35,3758.66,3903.78,3942.59,3981.74,3983.25,4033.97));
    modesGains(p,n) = waveform{0.525285,0.814174,0.483261,0.296745,0.975056,0.472244,0.409501,0.425364,0.687559,0.288381,0.309285,0.123054,0.286333,0.576706,0.908322,0.626974,0.0801852,0.309835,0.45143,0.132845,0.470635,0.417008,0.265112,0.0752802,0.46347,0.47181,0.275324,0.547027,0.512519,0.394078,0.595404,0.941306,0.392501,0.381435,0.391232,0.118924,0.339495,0.101421,0.241755,0.0873255,0.378944,0.637705,0.171946,0.149859,0.23329,0.54181,1,0.115554,0.244172,0.574329,0.606171,0.938397,0.392555,0.277359,0.86857,0.432489,0.408856,0.407932,0.299815,0.256659,0.549572,0.406347,0.312331,0.627578,0.670167,0.524648,0.406926,0.637524,0.555837,1,0.818979,0.705347,0.678141,0.427382,0.674404,0.636105,0.643635,0.699136,0.836201,0.613085,0.319019,0.725259,0.545519,0.479861,0.49836,0.488654,0.861672,0.314287,0.671052,0.531905,0.421781,0.815066,0.772032,0.488722,0.0896674,0.291286,0.65873,0.635632,0.679357,0.459497,0.36024,0.582289,0.650605,0.49095,0.38191,0.157261,0.479624,0.477491,0.174435,0.013094,0.879113,0.608069,0.268877,0.604479,0.24513,0.170507,0.292888,0.545849,0.476646,0.922316,0.669192,0.578094,0.578797,0.311396,0.60121,0.549955,1,0.66573,0.980115,0.537848,0.0740531,0.252472,0.25575,0.223974,0.0865103,0.138209,0.198623,0.0453034,0.432453,0.292407,0.39441,0.857659,0.271668,0.201545,0.583994,0.0602378,0.190618,0.849505,0.975542,0.17314,0.206472,0.344793,0.761011,0.558125,0.117245,0.0338485,0.337597,0.336646,0.174253,0.23017,0.934873,0.593647,0.393225,0.683704,0.0566093,0.0405012,0.148972,0.338722,0.283419,0.394007,0.237475,0.269964,0.428313,0.177499,0.462585,0.443963,0.981793,0.408239,0.676527,0.402865,0.0163303,0.0515114,0.34139,0.311135,0.613276,0.805884,0.95329,0.406091,0.578705,0.386785,0.434103,0.77526,1,0.635909,0.782052,0.0137183,0.0387725,0.618964,0.857071,0.131522,0.184988,0.299495,0.789212,0.603114,0.0704989,0.0129339,0.252481,0.254121,0.189206,0.357713,0.950308,0.552573,0.466454,0.77736,0.0307886,0.0251943,0.378886,0.740187,0.247637,0.235201,0.493045,0.51785,0.883954,0.429473,0.409433,0.415266,0.940198,0.282334,0.43789,0.375385,0.0157366,0.0171763,0.485555,0.461015,0.858958,0.907991,0.935191,0.37551,1,0.585493,0.269981,0.423053,0.666067,0.43509,0.790252,0.00889586,0.0208844,0.449735,0.790808,0.159856,0.089599,0.161546,0.528168,0.380642,0.0206276,0.00726426,0.0315352,0.0315841,0.197649,0.475057,0.517232,0.360922,0.421204,0.63134,0.00952139,0.016105,0.499615,0.922958,0.214983,0.0655141,0.50397,0.514848,1,0.483619,0.254027,0.228372,0.436105,0.233125,0.152242,0.279513,0.00906739,0.0132332,0.451257,0.388566,0.737226,0.479378,0.233037,0.103767,0.845609,0.644127,0.261359,0.371457,0.527229,0.381373,0.334492,0.00833749,0.00861982,0.255919,0.254197,0.0872333,0.0461512,0.113018,0.345986,0.236344,0.01078,0.00816506,0.405181,0.38718,0.343681,0.816493,0.259082,0.211906,0.432455,0.696886,0.00576224,0.013131,0.455969,0.811609,0.426544,0.128489,0.215937,0.233934,0.72307,0.351623,0.394231,0.323766,0.168803,0.276932,0.264684,0.227703,0.00680935,0.0170703,0.603018,0.476461,0.585925,0.71696,1,0.576527,0.475524,0.447322,0.356902,0.597573,0.697246,0.505333,0.285421,0.0147193,0.0141618,0.136188,0.0336537,0.216437},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/4035.44)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)churchBell`----------
// Generic church bell modal model.
//
// Modeled after T. Rossing and R. Perrin, Vibrations of Bells, Applied
// Acoustics 2, 1987.
//
// Model height is 301 mm.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// churchBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
churchBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
churchBellModel(50,strikePosition,30,1,2.5);


//-------`(pm.)churchBell_ui`----------
// Church bell physical model based on [`churchBell`](#pmchurchbell) with
// built-in UI.
//
// #### Usage
//
// ```
// churchBell_ui : _
// ```
//----------------------------------
churchBell_ui =
churchBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
  strikePosition = nentry("v:churchBell/[0]strikePosition",
    0,0,4,1);
    strikeCutoff = hslider("v:churchBell/[1]strikeCutOff",
    6500,20,20000,1);
    strikeSharpness = hslider("v:churchBell/[2]strikeSharpness",
    0.5,0.01,5,0.01);
    gain = hslider("v:churchBell/[3]gain",1,0,1,0.01);
    gate = button("v:churchBell/[4]gate");
};


//-------`(pm.)englishBellModel`----------
// English church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/englishBell`.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : englishBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
englishBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(259.429,261.855,493.331,494.472,603.048,604.742,724.566,724.777,786.985,931.893,932.846,943.429,946.533,947.109,952.782,1024.24,1034.86,1211.59,1212.21,1236.28,1340.21,1341.15,1452.63,1453.38,1453.77,1457.45,1459.23,1460.7,1482.81,1483.04,1519.55,1547.27,1565.95,1576.95,1697.47,1826.61,1892.56,1896.09,1899.86,1908.29,2003.95,2005.07,2005.89,2005.99,2067.42,2082.42,2086.81,2120.56,2129.79,2133.23));
    modesGains(p,n) = waveform{0.694274,0.471698,0.17239,0.585446,0.397986,0.919579,0.531947,0.100205,0.639469,0.672209,0.144345,0.416595,0.124108,0.380591,0.256578,0.646031,0.8522,0.0224376,0.382762,0.143925,0.368691,0.430556,0.32318,0.142956,0.274521,0.713824,0.442794,0.352473,0.247756,0.415152,1,0.401869,0.197981,0.27951,0.210249,0.36974,0.369227,0.155769,0.272368,0.335712,0.31645,0.714103,0.285781,0.22006,0.827704,0.206342,0.180177,0.311478,0.197607,0.575475,0.473311,0.587232,0.50288,0.337308,0.304514,0.429039,0.351522,0.341373,0.175081,0.561748,0.439172,0.323164,0.540518,0.536523,0.0743865,0.213417,0.358012,0.474494,0.310274,0.839413,0.241372,0.202343,0.480634,0.995685,0.37374,0.133998,0.520674,0.207514,1,0.101239,0.279536,0.185985,0.436293,0.62411,0.334519,0.283585,0.179317,0.353847,0.449545,0.574128,0.135172,0.538275,0.476424,0.832903,0.164198,0.188562,0.135978,0.390128,0.131045,0.312065,0.142139,0.0255901,0.266947,0.371607,0.0168435,0.0249468,0.508917,0.35441,0.283348,0.628155,0.292478,0.35835,0.342569,0.441237,0.886699,0.0174698,0.00641843,0.55532,0.880129,0.0306909,0.290081,0.248816,0.981736,0.324624,0.213676,0.432885,0.0981559,0.444149,0.395554,0.525069,0.0771308,0.0488804,0.591321,0.108044,0.443802,0.740318,0.599438,0.293093,1,0.141662,0.910031,0.226126,0.299702,0.341472,0.0568061,0.222494,0.918718,0.199478,0.21662,0.107759,0.692324,0.556336,0.281718,0.430832,0.341656,0.608095,0.342129,0.311312,0.229953,0.695087,0.0761489,0.349818,0.361706,0.577611,0.147797,0.327376,0.465715,0.342902,0.521381,0.836828,0.241112,0.284394,0.539316,0.143408,0.51702,1,0.236336,0.480333,0.676744,0.807582,0.468621,0.236953,0.411604,0.579251,0.425098,0.37599,0.461176,0.27653,0.462368,0.613004,0.666849,0.954715,0.161507,0.170433,0.290461,0.117303,0.365133,0.233794,0.194568,0.338874,0.523381,0.39835,0.447839,0.652891,0.15708,0.340354,0.44217,0.338764,0.184397,0.771607,0.160502,0.324487,0.477499,0.831519,0.0168764,0.154264,0.201865,0.237786,0.390066,0.880164,0.284234,0.162837,0.437557,0.227846,0.39934,1,0.158107,0.396903,0.513029,0.676457,0.204282,0.0895575,0.55541,0.732486,0.125062,0.171229,0.081646,0.0541394,0.39616,0.454288,0.466863,0.928842,0.155393,0.26285,0.113453,0.133163,0.302021,0.612616,0.228392,0.195617,0.287236,0.198762,0.499884,0.809139,0.00929925,0.0840299,0.286446,0.182112,0.186044,0.754138,0.279556,0.266948,0.494291,1,0.321284,0.0230981,0.0375537,0.262531,0.602204,0.489925,0.633077,0.407409,0.422256,0.0910641,0.357935,0.550179,0.106192,0.132366,0.376231,0.351151,0.0420302,0.0372183,0.696335,0.840821,0.855235,0.249248,0.860011,0.340481,0.28558,0.363039,0.324122,0.515699,0.228131,0.172592,0.0188723,0.168243,0.995105,0.741759,0.107093,0.070349,0.136636,0.0780455,0.315748,0.502201,0.0190422,0.033914,0.225724,0.160236,0.184101,0.564203,0.247317,0.284225,0.327153,0.651443,0.593471,0.0163899,0.0141048,0.52113,1,0.105109,0.530936,0.363724,0.924809,0.25041,0.69288,0.414122,0.0793658,0.347813,0.441731,0.476428,0.0808834,0.0581638,0.557239,0.556515,0.746084,0.582228,0.177231,0.158425,0.850903,0.755271,0.673614,0.31706,0.0515386,0.201898,0.0577938,0.232031,0.734988,0.141594,0.267062,0.145807},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/2137.04)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)englishBell`----------
// English church bell modal model.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// englishBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
englishBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
englishBellModel(50,strikePosition,30,1,3);


//-------`(pm.)englishBell_ui`----------
// English church bell physical model based on [`englishBell`](#pmenglishbell) with
// built-in UI.
//
// #### Usage
//
// ```
// englishBell_ui : _
// ```
//----------------------------------
englishBell_ui =
englishBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
    strikePosition = nentry("v:englishBell/[0]strikePosition",0,0,6,1);
    strikeCutoff = hslider("v:englishBell/[1]strikeCutOff",6500,20,20000,1);
    strikeSharpness = hslider("v:englishBell/[2]strikeSharpness",0.5,0.01,5,0.01);
    gain = hslider("v:englishBell/[3]gain",1,0,1,0.01);
    gate = button("v:englishBell/[4]gate");
};


//-------`(pm.)frenchBellModel`----------
// French church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/frenchBell`.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : frenchBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
frenchBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(439.077,440.305,606.452,611.776,709.617,709.877,856.966,899.051,1064.37,1064.75,1073.03,1074.43,1076.58,1145.42,1148.9,1197.09,1199.91,1387.7,1400.16,1402.37,1512.06,1513.06,1592.47,1593.35,1663.45,1666.47,1691.32,1693.77,1804.83,1808.84,1820.07,1847.14,1910.73,1927.7,2035.07,2038.73,2054.66,2055.67,2078.69,2091.76,2115.78,2127.01,2209.86,2214.05,2296.15,2298.08,2298.97,2307.62,2404.44,2408.99));
    modesGains(p,n) = waveform{0.97289,0.542393,0.495832,0.897966,0.552367,0.557895,0.614213,0.353694,0.436039,0.441024,0.674913,0.566754,0.755008,0.69903,0.164398,0.91004,0.628373,0.201718,0.429517,0.503715,0.871174,0.106886,0.761173,0.673602,0.291937,0.58859,0.31528,0.413081,0.274464,0.494062,0.696121,0.61201,0.382757,0.995113,0.228806,0.198449,0.595847,0.306263,0.252397,0.0883567,0.236086,1,0.245278,0.379388,0.198824,0.548892,0.492764,0.420871,0.794637,0.605634,1,0.604159,0.399841,0.799347,0.507187,0.50981,0.477383,0.310226,0.426976,0.437623,0.735712,0.630728,0.625785,0.651168,0.277738,0.850858,0.578079,0.645354,0.3837,0.370598,0.782542,0.181325,0.614391,0.740684,0.342441,0.586186,0.286909,0.405197,0.259215,0.566983,0.748219,0.655898,0.36826,0.940814,0.336157,0.413702,0.561557,0.402176,0.117698,0.329369,0.254571,0.870706,0.260981,0.274122,0.206247,0.6453,0.400758,0.363622,0.636834,0.584566,0.975603,0.616937,0.295401,0.650447,0.464221,0.465057,0.312467,0.238358,0.383695,0.399674,0.753912,0.662012,0.504268,0.599248,0.378665,0.725363,0.493214,1,0.332836,0.265889,0.674073,0.320401,0.440018,0.769782,0.316419,0.529968,0.380509,0.578676,0.249013,0.591384,0.761717,0.687057,0.324437,0.818864,0.505369,0.672485,0.461783,0.426198,0.0678875,0.435329,0.347955,0.708394,0.293322,0.328986,0.258254,0.893512,0.320131,0.433554,0.459302,0.542213,0.817241,0.51623,0.205302,0.467354,0.388683,0.388216,0.171262,0.150865,0.29172,0.311414,0.658876,0.570647,0.383619,0.502455,0.364114,0.532313,0.352989,1,0.261984,0.219591,0.544031,0.337199,0.279173,0.668303,0.208439,0.39923,0.418674,0.648618,0.234133,0.504729,0.645347,0.572851,0.232828,0.614292,0.485272,0.666264,0.31657,0.320355,0.191421,0.340131,0.342069,0.538371,0.281131,0.393115,0.251394,0.890725,0.310644,0.5037,0.29909,0.442478,0.733128,0.455217,0.199322,0.315699,0.375856,0.37649,0.029145,0.0200283,0.279578,0.3168,0.655957,0.546843,0.349666,0.470249,0.353765,0.286794,0.180185,1,0.210831,0.280133,0.535853,0.376488,0.15367,0.634745,0.0510449,0.485575,0.593111,0.917884,0.380477,0.422925,0.599373,0.311421,0.135654,0.359954,0.295748,0.474439,0.353339,0.116743,0.454313,0.112858,0.35931,0.483897,0.301608,0.577342,0.262663,0.794986,0.54958,0.808086,0.152511,0.439591,0.535941,0.308018,0.419837,0.579191,0.250428,0.25212,0.10286,0.288332,0.599679,0.665108,0.636285,0.495234,0.251613,0.208527,0.0939072,0.458415,0.318952,0.23535,0.215189,0.313412,0.0472787,0.0386893,0.0863359,0.222346,0.361511,0.997037,0.238509,0.38214,1,0.203554,0.472087,0.509015,0.206371,0.441013,0.794008,0.971623,0.796384,0.55043,0.68778,0.554298,0.0436545,0.0595577,0.214685,0.58147,0.27732,0.368466,0.121382,0.152739,0.0782244,0.190497,0.316269,0.180456,0.405196,0.518716,0.159365,0.165808,0.194372,0.614764,0.553415,0.644794,0.44178,0.44149,0.221588,1,0.579371,0.606339,0.529167,0.0214475,0.417046,0.663152,0.894841,0.475752,0.168289,0.46859,0.304604,0.666369,0.308337,0.295091,0.768945,0.350911,0.49068,0.516584,0.400117,0.820687,0.545486,0.709357,0.509759,0.472006,0.152455,0.981265,0.459134,0.698666,0.247154,0.623253,0.255508,0.830815,0.088709,0.126082,0.0770623,0.263328},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/2555.3)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)frenchBell`----------
// French church bell modal model.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
//frenchBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
frenchBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
frenchBellModel(50,strikePosition,30,1,3);


//-------`(pm.)frenchBell_ui`----------
// French church bell physical model based on [`frenchBell`](#pmfrenchbell) with
// built-in UI.
//
// #### Usage
//
// ```
// frenchBell_ui : _
// ```
//----------------------------------
frenchBell_ui =
frenchBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
    strikePosition = nentry("v:frenchBell/[0]strikePosition",0,0,6,1);
    strikeCutoff = hslider("v:frenchBell/[1]strikeCutOff",6500,20,20000,1);
    strikeSharpness = hslider("v:frenchBell/[2]strikeSharpness",0.5,0.01,5,0.01);
    gain = hslider("v:frenchBell/[3]gain",1,0,1,0.01);
    gate = button("v:frenchBell/[4]gate");
};


//-------`(pm.)germanBellModel`----------
// German church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/germanBell`.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : germanBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
germanBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(238.909,240.235,406.25,408.214,561.105,562.923,749.632,749.88,781.389,800.68,801.547,884.486,942.52,942.715,958.571,960.581,1016.87,1018.29,1022.29,1206.2,1206.87,1213.93,1222.77,1231.68,1345.65,1355.18,1356.55,1357.3,1389.29,1391.11,1420.17,1424.05,1424.31,1426.54,1582.3,1631.07,1638.15,1720.28,1726.7,1803.79,1827.08,1829.44,1922.87,1926.94,1940.96,1944.47,1948.2,1969.86,1992.24,2011.85));
    modesGains(p,n) = waveform{0.761411,0.797803,0.586057,0.225039,0.926354,0.813875,0.529871,0.158915,0.662516,0.30841,0.0454785,0.556082,0.784713,0.213035,0.956227,0.795165,0.313508,0.158531,0.301586,0.16832,0.281267,0.154237,0.217235,0.247239,0.522688,0.170303,0.320614,0.154857,0.414941,0.788364,0.713299,0.381592,0.707011,1,0.0106946,0.359763,0.0746156,0.431124,0.118053,0.52729,0.255303,0.646855,0.201145,0.919908,0.359389,0.253494,0.865695,0.0829263,0.222589,0.0415736,0.735774,0.769678,0.381416,0.197543,0.682389,0.614766,0.249525,0.289883,0.300319,0.443691,0.275272,0.370218,0.74879,0.161247,0.634353,0.498498,0.221988,0.350637,0.436817,0.436463,0.52508,0.842646,0.195324,0.224755,0.545681,0.353778,0.269044,0.327814,0.448952,0.852305,0.380503,1,0.458063,0.565058,0.354635,0.659529,0.449485,0.802014,0.283446,0.521563,0.374639,0.652112,0.181042,0.468394,0.430116,0.157207,0.414995,0.292737,0.487785,0.47768,0.481445,0.556988,0.561132,0.378265,0.291462,0.32325,0.482602,0.328248,0.196792,0.712254,0.389806,0.307751,0.792876,0.199098,0.288838,0.146811,0.0178444,0.394017,0.588212,0.489226,0.701486,0.940303,0.213191,0.0778845,0.474463,0.609972,0.16338,0.408376,0.330864,0.88999,0.0798101,1,0.167172,0.208727,0.227549,0.825903,0.531897,0.580946,0.215339,0.247319,0.268148,0.504088,0.160955,0.170316,0.230769,0.151191,0.0825031,0.141328,0.288079,0.499676,0.420355,0.487764,0.75371,0.492248,0.181364,0.20594,0.52497,0.325641,0.236827,0.931842,0.472981,0.312162,0.901032,0.270478,0.167112,0.0709698,0.166004,0.290745,0.425893,0.403633,0.581772,0.855694,0.0325587,0.0568359,0.241923,0.79931,0.181962,0.530283,0.41561,1,0.0291501,0.831155,0.119755,0.102188,0.132037,0.76603,0.442221,0.1749,0.142175,0.0140794,0.4375,0.85186,0.196836,0.0907522,0.551673,0.256528,0.0424377,0.490363,0.108178,0.503276,0.226584,0.312758,0.824079,0.48412,0.0347338,0.0619428,0.255097,0.145273,0.190359,0.995297,0.50342,0.217575,0.793975,0.357179,0.0409728,0.0187931,0.414458,0.15844,0.346233,0.0771673,0.175111,0.380567,0.497409,0.316164,0.488669,1,0.420657,0.442333,0.24706,0.796407,0.0104146,0.227079,0.10684,0.0289873,0.76757,0.0643122,0.0423098,0.715522,0.472117,0.392577,0.660433,0.803575,0.188653,0.0191653,0.873212,0.448719,0.0147128,0.618967,0.30718,0.345899,0.131394,0.193225,0.715283,0.40562,0.040637,0.0503336,0.0544331,0.0452023,0.152302,0.894549,0.443343,0.172071,0.647667,0.328993,0.0416014,0.0229488,0.551071,0.423544,0.862362,0.363526,0.769929,0.123452,0.710702,0.394895,0.556638,1,0.450925,0.422849,0.092187,0.413892,0.00995372,0.321146,0.160624,0.0300681,0.915385,0.4755,0.506951,0.387549,0.275591,0.468281,0.184776,0.186837,0.188195,0.0202785,0.697136,0.408862,0.0116935,0.349823,0.202301,0.461379,0.0559023,0.0582571,0.260608,0.181875,0.0271739,0.0249578,0.685089,0.470384,0.283194,0.33105,0.138349,0.338068,0.481992,0.178242,0.0155942,0.0110435,0.783771,0.442707,0.616478,0.381542,0.510892,0.045985,0.303119,0.0731909,0.547715,0.348941,0.149981,0.302158,0.284482,0.398177,0.00413049,0.180739,0.062839,0.0133459,0.347088,0.57637,0.240764,0.978481,0.452755,0.529742,0.340471,0.662282,0.444305,0.0429901,1,0.36194,0.0183372,0.626893,0.55285,0.384936},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/2016.94)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)germanBell`----------
// German church bell modal model.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 1 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// germanBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
germanBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
germanBellModel(50,strikePosition,30,1,2.5);


//-------`(pm.)germanBell_ui`----------
// German church bell physical model based on [`germanBell`](#pmgermanbell) with
// built-in UI.
//
// #### Usage
//
// ```
// germanBell_ui : _
// ```
//----------------------------------
germanBell_ui =
germanBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
    strikePosition = nentry("v:germanBell/[0]strikePosition",0,0,6,1);
    strikeCutoff = hslider("v:germanBell/[1]strikeCutOff",6500,20,20000,1);
    strikeSharpness = hslider("v:germanBell/[2]strikeSharpness",0.5,0.01,5,0.01);
    gain = hslider("v:germanBell/[3]gain",1,0,1,0.01);
    gate = button("v:germanBell/[4]gate");
};


//-------`(pm.)russianBellModel`----------
// Russian church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/russianBell`.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 2 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : russianBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
russianBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(136.491,136.992,258.958,260.419,316.489,318.411,393.852,393.949,454.677,467.493,503.594,503.637,530.46,531.046,541.16,546.158,578.335,579.863,660.026,708.34,716.915,717.23,775.839,776.066,783.356,783.464,788.826,789.463,800.348,806.758,876.788,880.414,887.297,888.222,916.815,919.711,937.686,984.606,1057.2,1058.45,1065.48,1082.27,1083.93,1086.55,1086.77,1108.2,1116.17,1149,1150.54,1199.54));
    modesGains(p,n) = waveform{0.925507,0.59752,0.0965671,0.45412,0.989773,0.593498,0.512541,0.124241,0.705411,0.292396,0.673399,0.302181,0.026234,0.286249,0.556267,1,0.250426,0.107711,0.427299,0.336295,0.616257,0.21442,0.0845294,0.231363,0.522724,0.559114,0.34847,0.854197,0.835576,0.735036,0.288494,0.117122,0.409686,0.363575,0.484943,0.170862,0.420531,0.164793,0.233847,0.861232,0.214037,0.283462,0.173153,0.876122,0.607809,0.294745,0.143142,0.332009,0.491878,0.626104,0.962027,0.584298,0.213653,0.420452,0.812329,0.545172,0.380744,0.331536,0.386801,0.248909,0.688756,0.313904,0.377894,0.595846,0.412274,0.739626,0.541148,0.393005,0.656637,0.847672,0.930653,0.288289,0.404938,0.657989,0.763949,0.574085,0.282241,1,0.429131,0.572049,0.734868,0.577477,0.520789,0.355593,0.890067,0.272391,0.448223,0.423969,0.392237,0.856091,0.0583794,0.784967,0.359527,0.576567,0.201513,0.642013,0.419308,0.340667,0.42319,0.860812,0.69402,0.423568,0.376987,0.568453,0.502302,0.280716,0.464041,0.395601,0.19985,0.0950398,0.64279,0.228326,0.484911,0.842353,0.161404,0.403432,0.655549,0.409098,0.699392,0.87605,1,0.185606,0.414255,0.695205,0.732612,0.478298,0.24947,0.927739,0.213135,0.227382,0.976352,0.642745,0.376311,0.260674,0.811034,0.264631,0.239979,0.261897,0.191642,0.78167,0.390679,0.382437,0.206714,0.22472,0.0676332,0.502611,0.301455,0.241029,0.224505,0.721193,0.436348,0.254062,0.480496,0.772371,0.210681,0.103415,0.485338,0.378334,0.228484,0.114877,0.68676,0.296942,0.50742,0.99747,0.0377103,0.132212,0.319547,0.192032,0.520574,0.585844,0.419362,0.0847317,0.134544,0.194762,0.616689,0.39736,0.298809,0.914746,0.0513371,0.0543569,0.989159,0.602892,0.145582,0.102149,0.129022,0.116485,0.279988,0.39427,0.290139,1,0.708655,0.780555,0.387526,0.0395217,0.00964067,0.24368,0.0740625,0.619315,0.447395,0.292497,0.295458,0.163748,0.521375,0.864533,0.0773408,0.0334231,0.345121,0.27886,0.223876,0.307756,0.763106,0.424707,0.487201,0.968962,0.00911747,0.0324653,0.334275,0.166787,0.312498,0.264262,0.35488,0.37559,0.261094,0.565006,0.474935,0.31352,0.251249,0.735352,0.0189072,0.0153634,0.786225,0.503299,0.369038,0.250765,0.673019,0.387573,0.742025,0.749056,0.261714,0.751868,0.771816,1,0.459484,0.0124402,0.0105394,0.964109,0.457052,0.532316,0.407128,0.697241,0.1522,0.0895893,0.406175,0.65104,0.0233951,0.010107,0.0722347,0.066634,0.195489,0.41674,0.654451,0.382782,0.305396,0.631501,0.00162802,0.0140906,0.762915,0.410245,0.189436,0.0604462,0.624941,0.439781,0.516273,0.896678,0.273298,0.202868,0.0996022,0.257657,0.0174508,0.0157859,0.429433,0.202184,0.443111,0.343811,0.447562,0.268694,0.753551,0.529426,0.0772973,0.097696,0.446414,0.261678,0.169035,0.0116219,0.0207399,1,0.488856,0.413029,0.252661,0.148369,0.0919644,0.0330634,0.268764,0.441849,0.0139873,0.0108584,0.657799,0.488248,0.375433,0.958179,0.761492,0.40949,0.151709,0.314931,0.0036118,0.013307,1,0.605343,0.550506,0.363516,0.255278,0.137537,0.448881,0.867615,0.483247,0.308622,0.348444,0.534835,0.0157716,0.0128965,0.147608,0.0762611,0.762224,0.511585,0.985863,0.540227,0.691691,0.905296,0.397521,0.69794,0.423289,0.924613,0.491559,0.00567911,0.0106002,0.647916,0.324182,0.579449,0.451936,0.877897},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/1201.01)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)russianBell`----------
// Russian church bell modal model.
//
// Modeled after D.Bartocha and Baron, Influence of Tin Bronze Melting and
// Pouring Parameters on Its Properties and Bell' Tone, Archives of Foundry
// Engineering, 2016.
//
// Model height is 2 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// russianBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
russianBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
russianBellModel(50,strikePosition,30,1,3);


//-------`(pm.)russianBell_ui`----------
// Russian church bell physical model based on [`russianBell`](#pmrussianbell) with
// built-in UI.
//
// #### Usage
//
// ```
// russianBell_ui : _
// ```
//----------------------------------
russianBell_ui =
russianBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
    strikePosition = nentry("v:russianBell/[0]strikePosition",0,0,6,1);
    strikeCutoff = hslider("v:russianBell/[1]strikeCutOff",6500,20,20000,1);
    strikeSharpness = hslider("v:russianBell/[2]strikeSharpness",0.5,0.01,5,0.01);
    gain = hslider("v:russianBell/[3]gain",1,0,1,0.01);
    gate = button("v:russianBell/[4]gate");
};

//-------`(pm.)standardBellModel`----------
// Standard church bell modal model generated by `mesh2faust` from
// `libraries/modalmodels/standardBell`.
//
// Modeled after T. Rossing and R. Perrin, Vibrations of Bells, Applied
// Acoustics 2, 1987.
//
// Model height is 1.8 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// #### Usage
//
// ```
// excitation : standardBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope)
// ```
//
// Where:
//
// * `excitation`: the excitation signal
// * `nModes`: number of synthesized modes (max: 50)
// * `exPos`: excitation position (0-6)
// * `t60`: T60 in seconds (recommended value: 0.1)
// * `t60DecayRatio`: T60 decay ratio (recommended value: 1)
// * `t60DecaySlope`: T60 decay slope (recommended value: 5)
//----------------------------------
standardBellModel(nModes,exPos,t60,t60DecayRatio,t60DecaySlope) = _ <: par(i,nModes,modeFilter(modesFreqs(i),modesT60s(i),modesGains(int(exPos),i))) :> /(nModes)
with {
    nExPos = 7;
    modesFreqs(n) = ba.take(n+1,(490.25,493.646,924.838,927.779,1181.21,1186.94,1348.84,1349.5,1560.33,1635.97,1706.73,1712.89,1745.05,1745.25,2005.51,2025.47,2053.88,2142.37,2151.4,2408.16,2534.11,2536.42,2623.3,2628.4,2711.57,2712.46,2823.23,2827.22,2863.42,2874.19,2923,2925.69,3032.52,3042.15,3208.57,3392.52,3485.92,3493.65,3539.8,3550.56,3678.71,3719.04,3722.59,3786.28,3789.38,3993.59,3998.43,4123.41,4164.83,4187.98));
    modesGains(p,n) = waveform{0.691911,0.622333,0.548651,0.463306,0.826946,0.749513,0.2242,0.642678,0.760442,0.326054,0.276463,0.359344,0.18258,0.686765,0.457159,0.839015,0.845338,0.372377,0.306417,0.147381,0.359707,0.653537,0.27553,0.401233,0.435417,0.251481,0.190062,0.773372,0.315014,0.228812,0.521512,0.411542,0.720762,1,0.286502,0.338938,0.119995,0.432289,0.409677,0.156272,0.298871,0.250786,0.640776,0.209431,0.17001,0.390014,0.301698,0.799413,0.980581,0.385,0.82544,0.818894,0.349616,0.235396,0.783164,0.821914,0.28411,0.430286,0.507671,0.326254,0.260488,0.273364,0.20518,0.714852,0.47995,0.803637,0.683943,0.355371,0.406924,0.656257,0.423025,0.413515,0.38636,0.384787,0.389448,0.813367,0.234988,1,0.311268,0.350245,0.403856,0.646143,0.500485,0.833553,0.431768,0.467064,0.298979,0.487413,0.514907,0.369383,0.106197,0.494224,0.816079,0.535807,0.379873,0.380201,0.606306,0.516117,0.748449,0.556948,0.587066,0.584423,0.394866,0.341121,0.433458,0.455987,0.361237,0.42939,0.122969,0.133175,0.505176,0.513985,0.0554619,0.604942,0.372074,0.381126,0.314354,0.499636,0.518711,0.923792,0.259544,0.576517,0.553915,0.585444,0.245369,1,0.117757,0.977318,0.652862,0.509314,0.14855,0.506402,0.180059,0.356005,0.38681,0.279354,0.205792,0.551055,0.689107,0.445724,0.306857,0.324747,0.603621,0.394466,0.288613,0.264697,0.60612,0.20274,0.267271,0.925656,0.439228,0.425884,0.626633,0.547204,0.230022,0.225654,0.392697,0.493474,0.149857,0.0604048,0.693889,0.740271,0.175485,0.704998,0.329732,0.153026,0.125744,0.286995,0.278878,0.812372,0.0562174,0.241479,0.294525,0.358834,0.171047,0.847604,0.17228,0.97521,0.892073,0.613987,0.0659213,0.301583,0.0610847,0.125438,0.145151,0.180086,0.124231,0.260161,0.337573,0.203743,0.655798,0.425893,0.902347,0.500686,0.311173,0.215561,0.349591,0.0854218,0.0805062,1,0.338652,0.295396,0.698314,0.664972,0.118983,0.0881905,0.31158,0.391136,0.151915,0.239504,0.685742,0.884332,0.288516,0.768688,0.274851,0.0490311,0.0357865,0.293303,0.249461,0.493771,0.340984,0.467623,0.216631,0.255235,0.0988695,0.46198,0.147247,0.640196,1,0.551938,0.0453732,0.189907,0.0197542,0.0309217,0.769837,0.360418,0.384041,0.867434,0.398948,0.171848,0.748652,0.301957,0.860611,0.958674,0.54903,0.272753,0.372753,0.0180728,0.0292353,0.8502,0.224583,0.214805,0.670319,0.586433,0.0435142,0.0388574,0.144811,0.157061,0.155569,0.418334,0.673656,0.749573,0.337354,0.747254,0.255997,0.0239656,0.0310719,0.721087,0.700616,0.199051,0.511844,0.849485,0.700682,0.778658,0.171289,0.261973,0.129228,0.328597,0.781821,0.583813,0.0806713,0.416876,0.0118202,0.00868563,1,0.461884,0.186882,0.641364,0.994705,0.501902,0.566449,0.0678845,0.139737,0.462582,0.318656,0.233947,0.495941,0.0314028,0.0146478,0.70432,0.124953,0.132549,0.457126,0.378636,0.0169362,0.0195494,0.204155,0.294401,0.271367,0.730857,0.459322,0.433078,0.325171,0.734536,0.416205,0.012873,0.0388489,0.821567,0.863683,0.0920531,0.393972,0.539544,0.832052,0.842732,0.241144,0.479558,0.283092,0.477845,0.385473,0.436587,0.144308,0.642395,0.0215791,0.00779029,0.563714,0.838279,0.410004,0.829086,1,0.630598,0.0233729,0.496217,0.711042,0.914266,0.695042,0.331894,0.898442,0.028568,0.0174966,0.482846},int(p*nModes+n) : rdtable;
    modesT60s(i) = t60*pow(1-(modesFreqs(i)/4191.95)*t60DecayRatio,t60DecaySlope);
};


//-------`(pm.)standardBell`----------
// Standard church bell modal model.
//
// Modeled after T. Rossing and R. Perrin, Vibrations of Bells, Applied
// Acoustics 2, 1987.
//
// Model height is 1.8 m.
//
// This model contains 7 excitation positions going linearly from the
// bottom to the top of the bell. Obviously, a model with more excitation
// position could be regenerated using `mesh2faust`.
//
// This function also implement a virtual exciter to drive the model.
//
// #### Usage
//
// ```
// standardBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) : _
// ```
//
// Where:
//
// * `strikePosition`: strike position (0-6)
// * `strikeCutoff`: cuttoff frequency of the strike genarator (recommended: ~7000Hz)
// * `strikeSharpness`: sharpness of the strike (recommended: ~0.25)
// * `gain`: gain of the strike (0-1)
// * `trigger` signal (0: off, 1: on)
//----------------------------------
standardBell(strikePosition,strikeCutoff,strikeSharpness,gain,trigger) =
strikeModel(10,strikeCutoff,strikeSharpness,gain,trigger) :
standardBellModel(50,strikePosition,30,1,2.5);


//-------`(pm.)standardBell_ui`----------
// Standard church bell physical model based on [`standardBell`](#pmstandardbell) with
// built-in UI.
//
// #### Usage
//
// ```
// standardBell_ui : _
// ```
//----------------------------------
standardBell_ui =
standardBell(strikePosition,strikeCutoff,strikeSharpness,gain,gate)
with {
    strikePosition = nentry("v:standardBell/[0]strikePosition",0,0,6,1);
    strikeCutoff = hslider("v:standardBell/[1]strikeCutOff",6500,20,20000,1);
    strikeSharpness = hslider("v:standardBell/[2]strikeSharpness",0.5,0.01,5,0.01);
    gain = hslider("v:standardBell/[3]gain",1,0,1,0.01);
    gate = button("v:standardBell/[4]gate");
};


//==============================Vocal Synthesis===========================================
// Vocal synthesizer functions (source/filter, fof, etc.).
//========================================================================================

//-------`(pm.)formantValues`----------
// Formant data values in an environment.
//
// The formant data used here come from the CSOUND manual
// <http://www.csounds.com/manual/html/>.
//
// #### Usage
//
// ```
// ba.take(j+1,formantValues.f(i)) : _
// ba.take(j+1,formantValues.g(i)) : _
// ba.take(j+1,formantValues.bw(i)) : _
// ```
//
// Where:
//
// * `i`: formant number
// * `j`: (voiceType*nFormants)+vowel
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3:
// soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
//--------------------------------------
formantValues = environment {
    f(0) = (800,400,350,450,325,600,400,250,400,350,660,440,270,430,370,800,
        350,270,450,325,650,400,290,400,350); // formant 0 freqs
    f(1) = (1150,1600,1700,800,700,1040,1620,1750,750,600,1120,1800,1850,820,630,
        1150,2000,2140,800,700,1080,1700,1870,800,600); // formant 1 freqs
    f(2) = (2800,2700,2700,2830,2530,2250,2400,2600,2400,2400,2750,2700,2900,2700,2750,
        2900,2800,2950,2830,2700,2650,2600,2800,2600,2700); // formant 2 freqs
    f(3) = (3500,3300,3700,3500,3500,2450,2800,3050,2600,2675,3000,3000,3350,3000,3000,
        3900,3600,3900,3800,3800,2900,3200,3250,2800,2900); // formant 3 freqs
    f(4) = (4950,4950,4950,4950,4950,2750,3100,3340,2900,2950,3350,3300,3590,3300,3400,
        4950,4950,4950,4950,4950,3250,3580,3540,3000,3300); // formant 4 freqs
    g(0) = (1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1); // formant 0 gains
    g(1) = (0.630957,0.063096,0.100000,0.354813,0.251189,0.446684,0.251189,0.031623,
        0.281838,0.100000,0.501187,0.199526,0.063096,0.316228,0.100000,
        0.501187,0.100000,0.251189,0.281838,0.158489,0.501187,0.199526,0.177828,
        0.316228,0.100000); // formant 1 gains
    g(2) = (0.100000,0.031623,0.031623,0.158489,0.031623,0.354813,0.354813,0.158489,
        0.089125,0.025119,0.070795,0.125893,0.063096,0.050119,0.070795,
        0.025119,0.177828,0.050119,0.079433,0.017783,0.446684,0.251189,0.125893,
        0.251189,0.141254); // formant 2 gains
    g(3) = (0.015849,0.017783,0.015849,0.039811,0.010000,0.354813,0.251189,0.079433,
        0.100000,0.039811,0.063096,0.100000,0.015849,0.079433,0.031623,
        0.100000,0.010000,0.050119,0.079433,0.010000,0.398107,0.199526,0.100000,
        0.251189,0.199526); // formant 3 gains
    g(4) = (0.001000,0.001000,0.001000,0.001778,0.000631,0.100000,0.125893,0.039811,
        0.010000,0.015849,0.012589,0.100000,0.015849,0.019953,0.019953,
        0.003162,0.001585,0.006310,0.003162,0.001000,0.079433,0.100000,0.031623,
        0.050119,0.050119); // formant 4 gains
    bw(0) = (80,60,50,70,50,60,40,60,40,40,80,70,40,40,40,80,60,60,40,50,
    50,70,40,70,40); // formant 0 bandwidths
    bw(1) = (90,80,100,80,60,70,80,90,80,80,90,80,90,80,60,90,100,90,80,60,
        90,80,90,80,60); // formant 1 bandwidths
    bw(2) = (120,120,120,100,170,110,100,100,100,100,120,100,100,100,100,
        120,120,100,100,170,120,100,100,100,100); // formant 2 bandwidths
    bw(3) = (130,150,150,130,180,120,120,120,120,120,130,120,120,120,120,
        130,150,120,120,180,130,120,120,130,120); // formant 3 bandwidths
    bw(4) = (140,200,200,135,200,130,120,120,120,120,140,120,120,120,120,
        140,200,120,120,200,140,120,120,135,120); // formant 4 bandwidths
};


// array of values used to multiply BWs by to get attack Bws for FOF version.
// min/max values per vowel (AEIOU) and per gender (M/F). Index by:
// gender*5 + vowel;
// values were chosen based on informal listening tests
bwMultMins = (1.0, 1.25, 1.25, 1.0, 1.5, 2.0, 3.0, 3.0, 2.0, 2.0);
bwMultMaxes = (10.0, 2.5, 2.5, 10.0, 4.0, 15.0, 12.0, 12.0, 12.0, 12.0);


// minimum/maximum frequency values per gender (M/F) used in the calculation
// of the attack Bws from the release Bws in the FOF version
// values are based on arbitrary maximum/minimum singing values
// in Hz for male/female voices
minGenderFreq = (82.41,174.61);
maxGenderFreq = (523.25,1046.5);


//--------------`(pm.)voiceGender`-----------------
// Calculate the gender for the provided `voiceType` value. (0: male, 1: female)
//
// #### Usage
//
// ```
// voiceGender(voiceType) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
//---------------------------------------------
declare voiceGender author "Mike Olsen";

voiceGender(voiceType) = ba.if(voiceType == 0,1,ba.if(voiceType == 3,1,0));


//-----------`(pm.)skirtWidthMultiplier`------------
// Calculates value to multiply bandwidth to obtain `skirtwidth`
// for a Fof filter.
//
// #### Usage
//
// ```
// skirtWidthMultiplier(vowel,freq,gender) : _
// ```
//
// Where:
//
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `freq`: the fundamental frequency of the excitation signal
// * `gender`: gender of the voice used in the fof filter (0: male, 1: female)
//---------------------------------------------
declare skirtWidthMultiplier author "Mike Olsen";

skirtWidthMultiplier(vowel,freq,gender) = (multMax-multMin)*skirtParam+multMin
with {
    nVowels = 5;
    index = gender*nVowels + vowel;
    multMin = bwMultMins : ba.selectn(10,index);
    multMax = bwMultMaxes : ba.selectn(10,index);
    freqMin = minGenderFreq : ba.selectn(2,gender);
    freqMax = maxGenderFreq : ba.selectn(2,gender);
    skirtParam = ba.if(freq <= freqMin,0.0,ba.if(freq >= freqMax,1.0,
                 (1.0/(freqMax-freqMin))*(freq-freqMin)));
};


//--------------`(pm.)autobendFreq`-----------------
// Autobends the center frequencies of formants 1 and 2 based on
// the fundamental frequency of the excitation signal and leaves
// all other formant frequencies unchanged. Ported from `chant-lib`.
//
// #### Reference
//
// <https://ccrma.stanford.edu/~rmichon/chantLib/>.
//
// #### Usage
//
// ```
// _ : autobendFreq(n,freq,voiceType) : _
// ```
//
// Where:
//
// * `n`: formant index
// * `freq`: the fundamental frequency of the excitation signal
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * input is the center frequency of the corresponding formant
//---------------------------------------------
declare autobendFreq author "Mike Olsen";

autobendFreq(n,freq,voiceType) = autobend(n)
with {
    autobend(0) = _ <: ba.if(_ <= freq,freq,_);
    autobend(1) = _ <: ba.if(voiceType != 2,
                _ <: ba.if((_ >= 1300)&(freq >= 200),
                _ -(freq-200)*(_-1300)/1050,
                ba.if(_ <= (30 + 2*freq),30 + 2*freq,_)), _);
    autobend(n) = _;
};


//--------------`(pm.)vocalEffort`-----------------
// Changes the gains of the formants based on the fundamental
// frequency of the excitation signal. Higher formants are
// reinforced for higher fundamental frequencies.
// Ported from `chant-lib`.
//
// #### Reference
//
// <https://ccrma.stanford.edu/~rmichon/chantLib/>.
//
// #### Usage
//
// ```
// _ : vocalEffort(freq,gender) : _
// ```
//
// Where:
//
// * `freq`: the fundamental frequency of the excitation signal
// * `gender`: the gender of the voice type (0: male, 1: female)
// * input is the linear amplitude of the formant
//---------------------------------------------
declare vocalEffort author "Mike Olsen";

vocalEffort(freq,gender) = _ <: ba.if(gender == 0,*(3+1.1*(400-freq)/300),*(0.8+1.05*(1000-freq)/1250));


//-------------------------`(pm.)fof`--------------------------
// Function to generate a single Formant-Wave-Function.
// 
// #### Reference
//
// <https://ccrma.stanford.edu/~mjolsen/pdfs/smc2016_MOlsenFOF.pdf>.
//
// #### Usage
//
// ```
// _ : fof(fc,bw,a,g) : _
// ```
//
// Where:
//
// * `fc`: formant center frequency,
// * `bw`: formant bandwidth (Hz),
// * `sw`: formant skirtwidth (Hz)
// * `g`: linear scale factor (g=1 gives 0dB amplitude response at fc)
// * input is an impulse signal to excite filter
//---------------------------------------------------------
declare fof author "Mike Olsen";

fof(fc,bw,sw,g) = _ <: (_',_) : (f * s)
with {
    T = 1/ma.SR;               // sample period
    pi = ma.PI;             // pi
    u1 = exp(-sw*pi*T);     // exponential controlling rise
    u2 = exp(-bw*pi*T);     // exponential controlling decay
    a1 = -1*(u1+u2);        // a1 filter coefficient
    a2 = u1*u2;             // a2 filter coefficient
    G0 = 1/(1+a1+a2);       // magnitude at DC
    b0 = g/G0;              // filter gain
    s  = os.hs_oscsin(fc);     // hardsyncing wavetable oscillator
    f  = fi.tf2(b0,0,0,a1,a2); // biquad filter
};


//-------------------------`(pm.)fofSH`-------------------------
// FOF with sample and hold used on `bw` and a parameter
// used in the filter-cycling FOF function `fofCycle`.
// 
// #### Reference
//
// <https://ccrma.stanford.edu/~mjolsen/pdfs/smc2016_MOlsenFOF.pdf>.
//
// #### Usage
//
// ```
// _ : fofSH(fc,bw,a,g) : _
// ```
//
// Where: all parameters same as for [`fof`](#fof)
//---------------------------------------------------------
declare fofSH author "Mike Olsen";

fofSH(fc,bw,a,g) = _ <: (((_,bw):ba.sAndH),((_,a):ba.sAndH),_) : (fc,_,_,g,_') : fof;


//----------------------`(pm.)fofCycle`-------------------------
// FOF implementation where time-varying filter parameter noise is
// mitigated by using a cycle of `n` sample and hold FOF filters.
// 
// #### Reference
//
// <https://ccrma.stanford.edu/~mjolsen/pdfs/smc2016_MOlsenFOF.pdf>.
//
// #### Usage
//
// ```
// _ : fofCycle(fc,bw,a,g,n) : _
// ```
//
// Where:
//
// * `n`: the number of FOF filters to cycle through
// * all other parameters are same as for [`fof`](#fof)
//---------------------------------------------------------
declare fofCycle author "Mike Olsen";

fofCycle(fc,bw,a,g,n) = _ : ba.cycle(n) : par(i,n,fofSH(fc,bw,a,g)) :> _;


//----------------------`(pm.)fofSmooth`-------------------------
// FOF implementation where time-varying filter parameter
// noise is mitigated by lowpass filtering the filter
// parameters `bw` and `a` with [smooth](#smooth).
//
// #### Usage
//
// ```
// _ : fofSmooth(fc,bw,sw,g,tau) : _
// ```
//
// Where:
//
// * `tau`: the desired smoothing time constant in seconds
// * all other parameters are same as for [`fof`](#fof)
//---------------------------------------------------------
declare fofSmooth author "Mike Olsen";

fofSmooth(fc,bw,sw,g,tau) = fof(fc,bw2,sw2,g)
with {
    sw2 = sw : si.smooth(ba.tau2pole(tau));
    bw2 = bw : si.smooth(ba.tau2pole(tau));
};


//-------`(pm.)formantFilterFofCycle`--------------
// Formant filter based on a single FOF filter.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. A cycle of `n` fof filters with sample-and-hold is
// used so that the fof filter parameters can be varied in realtime.
// This technique is more robust but more computationally expensive than
// [`formantFilterFofSmooth`](#formantFilterFofSmooth).Voice type can be
// selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterFofCycle(voiceType,vowel,nFormants,i,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor,
//    3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `nFormants`: number of formant regions in frequency domain, typically 5
// * `i`: formant number (i.e. 0 - 4) used to index formant data value arrays
// * `freq`: fundamental frequency of excitation signal. Used to calculate
//         rise time of envelope
//--------------------------------------
declare formantFilterFofCycle author "Mike Olsen";

formantFilterFofCycle(voiceType,vowel,nFormants,i,freq) =
    fofCycle(formantFreq(i),formantBw(i),formantSw(i),formantGain(i),n)
with {
    n = 3; // number of fof filters to cycle between
    index = (voiceType*nFormants)+vowel; // index of formant values
    // formant center frequency using autobend correction
    formantFreq(i) = ba.listInterp(formantValues.f(i),index) : autobendFreq(i,freq,voiceType);
    // formant amplitude using vocal effort correction
    formantGain(i) = ba.listInterp(formantValues.g(i),index) : vocalEffort(freq,gender);
    formantBw(i) = ba.listInterp(formantValues.bw(i),index); // formant bandwidth
    // formant skirtwidth
    formantSw(i) = skirtWidthMultiplier(vowel,freq,gender)*formantBw(i);
    gender = voiceGender(voiceType); // gender of voice
};


//-------`(pm.)formantFilterFofSmooth`--------------
// Formant filter based on a single FOF filter.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Fof filter parameters are lowpass filtered
// to mitigate possible noise from varying them in realtime.
// Voice type can be selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterFofSmooth(voiceType,vowel,nFormants,i,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor,
//           3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `nFormants`: number of formant regions in frequency domain, typically 5
// * `i`: formant number (i.e. 1 - 5) used to index formant data value arrays
// * `freq`: fundamental frequency of excitation signal. Used to calculate
//         rise time of envelope
//--------------------------------------
declare formantFilterFofSmooth author "Mike Olsen";

formantFilterFofSmooth(voiceType,vowel,nFormants,i,freq) =
    fofSmooth(formantFreq(i),formantBw(i),formantSw(i),formantGain(i),tau)
with {
    tau = 0.001;
    index = (voiceType*nFormants)+vowel; // index of formant values
    // formant center frequency using autobend correction
    formantFreq(i) = ba.listInterp(formantValues.f(i),index) : autobendFreq(i,freq,voiceType);
    // formant amplitude using vocal effort correction
    formantGain(i) = ba.listInterp(formantValues.g(i),index) : vocalEffort(freq,gender);
    formantBw(i) = ba.listInterp(formantValues.bw(i),index); // formant bandwidth
    // formant skirtwidth
    formantSw(i) = skirtWidthMultiplier(vowel,freq,gender)*formantBw(i);
    gender = voiceGender(voiceType); // gender of voice
};


//-------`(pm.)formantFilterBP`--------------
// Formant filter based on a single resonant bandpass filter.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterBP(voiceType,vowel,nFormants,i,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `nFormants`: number of formant regions in frequency domain, typically 5
// * `i`: formant index used to index formant data value arrays
// * `freq`: fundamental frequency of excitation signal.
//--------------------------------------
formantFilterBP(voiceType,vowel,nFormants,i,freq) =
    fi.resonbp(formantFreq(i),formantFreq(i)/formantBw(i),
           formantGain(i))
with {
    index = (voiceType*nFormants)+vowel; // index of formant values
    // formant center frequency using autobend correction
    formantFreq(i) = ba.listInterp(formantValues.f(i),index) : autobendFreq(i,freq,voiceType);
    // formant amplitude using vocal effort correction
    formantGain(i) = ba.listInterp(formantValues.g(i),index) : vocalEffort(freq,gender);
    formantBw(i) = ba.listInterp(formantValues.bw(i),index); // formant bandwidth
    gender = voiceGender(voiceType); // gender of voice
};


//-------`(pm.)formantFilterbank`--------------
// Formant filterbank which can use different types of filterbank
// functions and different excitation signals. Formant parameters are
// linearly interpolated allowing to go smoothly from one vowel to another.
// Voice type can be selected but must correspond to the frequency range
// of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterbank(voiceType,vowel,formantGen,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `formantGen`: the specific formant filterbank function
//  (i.e. FormantFilterbankBP, FormantFilterbankFof,...)
// * `freq`: fundamental frequency of excitation signal. Needed for FOF
//  version to calculate rise time of envelope
//--------------------------------------
declare formantFilterbank author "Mike Olsen";

formantFilterbank(voiceType,vowel,formantGen,freq) =
    _ <: par(i,nFormants,formantGen(voiceType,vowel,nFormants,i,freq)) :> _
with {
    nFormants = 5;
};


//-----`(pm.)formantFilterbankFofCycle`-----
// Formant filterbank based on a bank of fof filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterbankFofCycle(voiceType,vowel,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `freq`: the fundamental frequency of the excitation signal. Needed to calculate the skirtwidth 
// of the FOF envelopes and for the autobendFreq and vocalEffort functions
//-------------------------------------
declare formantFilterbankFofCycle author "Mike Olsen";

formantFilterbankFofCycle(voiceType,vowel,freq) =
formantFilterbank(voiceType,vowel,formantFilterFofCycle,freq);


//-----`(pm.)formantFilterbankFofSmooth`----
// Formant filterbank based on a bank of fof filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterbankFofSmooth(voiceType,vowel,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `freq`: the fundamental frequency of the excitation signal. Needed to
// calculate the skirtwidth of the FOF envelopes and for the
// autobendFreq and vocalEffort functions
//-------------------------------------
declare formantFilterbankFofSmooth author "Mike Olsen";

formantFilterbankFofSmooth(voiceType,vowel,freq) =
formantFilterbank(voiceType,vowel,formantFilterFofSmooth,freq);


//-------`(pm.)formantFilterbankBP`--------------
// Formant filterbank based on a bank of resonant bandpass filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the provided source to be realistic.
//
// #### Usage
//
// ```
// _ : formantFilterbankBP(voiceType,vowel,freq) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u)
// * `freq`: the fundamental frequency of the excitation signal. Needed for the autobendFreq and vocalEffort functions
//--------------------------------------
formantFilterbankBP(voiceType,vowel,freq) =
formantFilterbank(voiceType,vowel,formantFilterBP,freq);


//-------`(pm.)SFFormantModel`--------------
// Simple formant/vocal synthesizer based on a source/filter model. The `source`
// and `filterbank` must be specified by the user. `filterbank` must take the same
// input parameters as [`formantFilterbank`](#formantFilterbank) (`BP`/`FofCycle`
// /`FofSmooth`).
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the synthesized voice to be realistic.
//
// #### Usage
//
// ```
// SFFormantModel(voiceType,vowel,exType,freq,gain,source,filterbank,isFof) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u
// * `exType`: voice vs. fricative sound ratio (0-1 where 1 is 100% fricative)
// * `freq`: the fundamental frequency of the source signal
// * `gain`: linear gain multiplier to multiply the source by
// * `isFof`: whether model is FOF based (0: no, 1: yes)
//--------------------------------------
declare SFFormantModel author "Mike Olsen";

SFFormantModel(voiceType,vowel,exType,freq,gain,source,filterbank,isFof) =
    excitation : resonance
with {
    breath = no.noise;
    excitation = ba.if(isFof,source,source*(1-exType) + breath*exType :
            *(gain));
    resonance = filterbank(voiceType,vowel,freq) <: ba.if(isFof,*(gain),_);
};


//-------`(pm.)SFFormantModelFofCycle`-------
// Simple formant/vocal synthesizer based on a source/filter model. The source
// is just a periodic impulse and the "filter" is a bank of FOF filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the synthesized voice to be realistic. This model
// does not work with noise in the source signal so exType has been removed
// and model does not depend on SFFormantModel function.
//
// #### Usage
//
// ```
// SFFormantModelFofCycle(voiceType,vowel,freq,gain) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u
// * `freq`: the fundamental frequency of the source signal
// * `gain`: linear gain multiplier to multiply the source by
//---------------------------------------
declare SFFormantModelFofCycle author "Mike Olsen";

SFFormantModelFofCycle(voiceType,vowel,freq,gain) =
SFFormantModel(voiceType,vowel,0,freq,gain,os.lf_imptrain(freq),
formantFilterbankFofCycle,1);


//-------`(pm.)SFFormantModelFofSmooth`-------
// Simple formant/vocal synthesizer based on a source/filter model. The source
// is just a periodic impulse and the "filter" is a bank of FOF filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the synthesized voice to be realistic.
//
// #### Usage
//
// ```
// SFFormantModelFofSmooth(voiceType,vowel,freq,gain) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u
// * `freq`: the fundamental frequency of the source signal
// * `gain`: linear gain multiplier to multiply the source by
//---------------------------------------
declare SFFormantModelFofSmooth author "Mike Olsen";

SFFormantModelFofSmooth(voiceType,vowel,freq,gain) =
SFFormantModel(voiceType,vowel,0,freq,gain,os.lf_imptrain(freq),
formantFilterbankFofSmooth,1);


//-------`(pm.)SFFormantModelBP`--------------
// Simple formant/vocal synthesizer based on a source/filter model. The source
// is just a sawtooth wave and the "filter" is a bank of resonant bandpass filters.
// Formant parameters are linearly interpolated allowing to go smoothly from
// one vowel to another. Voice type can be selected but must correspond to
// the frequency range of the synthesized voice to be realistic.
//
// The formant data used here come from the CSOUND manual
// <http://www.csounds.com/manual/html/>.
//
// #### Usage
//
// ```
// SFFormantModelBP(voiceType,vowel,exType,freq,gain) : _
// ```
//
// Where:
//
// * `voiceType`: the voice type (0: alto, 1: bass, 2: countertenor, 3: soprano, 4: tenor)
// * `vowel`: the vowel (0: a, 1: e, 2: i, 3: o, 4: u
// * `exType`: voice vs. fricative sound ratio (0-1 where 1 is 100% fricative)
// * `freq`: the fundamental frequency of the source signal
// * `gain`: linear gain multiplier to multiply the source by
//---------------------------------------
SFFormantModelBP(voiceType,vowel,exType,freq,gain) =
SFFormantModel(voiceType,vowel,exType,freq,gain,os.sawtooth(freq),
formantFilterbankBP,0);


//-------`(pm.)SFFormantModelFofCycle_ui`----------
// Ready-to-use source-filter vocal synthesizer with built-in user interface.
//
// #### Usage
//
// ```
// SFFormantModelFofCycle_ui : _
// ```
//----------------------------------
declare SFFormantModelFofCycle_ui author "Mike Olsen";

SFFormantModelFofCycle_ui = SFFormantModelFofCycle(voiceType,vowel,freq2,gain*corrFactor)
with {
    freq1 = hslider("v:vocal/[0]freq",440,50,1000,0.01);
    gain = hslider("v:vocal/[1]gain",0.9,0,1,0.01);
    corrFactor = 75.0;
    voiceType = hslider("v:vocal/[2]voiceType",0,0,4,1);
    vowel = hslider("v:vocal/[3]vowel",0,0,4,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/[5]vibratoFreq",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/[6]vibratoGain",0.5,0,1,0.01)*0.1;
    freq2 = freq1*(os.osc(vibratoFreq)*vibratoGain+1);
};


//-------`(pm.)SFFormantModelFofSmooth_ui`----------
// Ready-to-use source-filter vocal synthesizer with built-in user interface.
//
// #### Usage
//
// ```
// SFFormantModelFofSmooth_ui : _
// ```
//----------------------------------
declare SFFormantModelFofSmooth_ui author "Mike Olsen";

SFFormantModelFofSmooth_ui = SFFormantModelFofSmooth(voiceType,vowel,freq2,gain*corrFactor)
with {
    freq1 = hslider("v:vocal/[0]freq",440,50,1000,0.01);
    gain = hslider("v:vocal/[1]gain",0.9,0,1,0.01);
    corrFactor = 25.0;
    voiceType = hslider("v:vocal/[2]voiceType",0,0,4,1);
    vowel = hslider("v:vocal/[3]vowel",0,0,4,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/[5]vibratoFreq",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/[6]vibratoGain",0.5,0,1,0.01)*0.1;
    freq2 = freq1*(os.osc(vibratoFreq)*vibratoGain+1);
};


//-------`(pm.)SFFormantModelBP_ui`----------
// Ready-to-use source-filter vocal synthesizer with built-in user interface.
//
// #### Usage
//
// ```
// SFFormantModelBP_ui : _
// ```
//----------------------------------
SFFormantModelBP_ui = SFFormantModelBP(voiceType,vowel,fricative,freq2,gain)
with {
    freq1 = hslider("v:vocal/[0]freq",440,50,1000,0.01);
    gain = hslider("v:vocal/[1]gain",0.9,0,1,0.01);
    voiceType = hslider("v:vocal/[2]voiceType",0,0,4,1);
    vowel = hslider("v:vocal/[3]vowel",0,0,4,0.01) : si.smoo;
    fricative = hslider("v:vocal/[4]fricative",0,0,1,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/[5]vibratoFreq",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/[6]vibratoGain",0.5,0,1,0.01)*0.1;
    freq2 = freq1*(os.osc(vibratoFreq)*vibratoGain+1);
};


//-------`(pm.)SFFormantModelFofCycle_ui_MIDI`----------
// Ready-to-use MIDI-controllable source-filter vocal synthesizer.
//
// #### Usage
//
// ```
// SFFormantModelFofCycle_ui_MIDI : _
// ```
//----------------------------------
SFFormantModelFofCycle_ui_MIDI = SFFormantModelFofCycle(voiceType,vowel,freq2,envelope)*outGain
with {
    freq1 = hslider("v:vocal/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:vocal/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:vocal/h:[0]midi/[2]gain[style:knob]
    ",0.9,0,1,0.01);
    corrFactor = 75.0;
    envAttack = hslider("v:vocal/h:[0]midi/[3]envAttack[style:knob]
    ",10,0,30,0.01)*0.001;
    s = hslider("v:vocal/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    voiceType = hslider("v:vocal/h:[1]otherParams/[0]voiceType
    [style:knob]",0,0,4,1);
    vowel = hslider("v:vocal/h:[1]otherParams/[1]vowel
    [style:knob][midi:ctrl 1]",0,0,4,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/h:[1]otherParams/[3]vibratoFreq[style:knob]
    ",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/h:[1]otherParams/[4]vibratoGain[style:knob]
    ",0.5,0,1,0.01)*0.1;
    outGain = hslider("v:vocal/h:[1]otherParams/[5]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:vocal/[2]gate");

    gate = t+s : min(1);
    freq2 = freq1*bend*(os.osc(vibratoFreq)*vibratoGain+1);
    envelope = gate*gain*corrFactor : si.smooth(ba.tau2pole(envAttack));
};


//-------`(pm.)SFFormantModelFofSmooth_ui_MIDI`----------
// Ready-to-use MIDI-controllable source-filter vocal synthesizer.
//
// #### Usage
//
// ```
// SFFormantModelFofSmooth_ui_MIDI : _
// ```
//----------------------------------
SFFormantModelFofSmooth_ui_MIDI = SFFormantModelFofSmooth(voiceType,vowel,freq2,envelope)*outGain
with {
    freq1 = hslider("v:vocal/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:vocal/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:vocal/h:[0]midi/[2]gain[style:knob]
    ",0.9,0,1,0.01);
    corrFactor = 25.0;
    envAttack = hslider("v:vocal/h:[0]midi/[3]envAttack[style:knob]
    ",10,0,30,0.01)*0.001;
    s = hslider("v:vocal/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    voiceType = hslider("v:vocal/h:[1]otherParams/[0]voiceType
    [style:knob]",0,0,4,1);
    vowel = hslider("v:vocal/h:[1]otherParams/[1]vowel
    [style:knob][midi:ctrl 1]",0,0,4,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/h:[1]otherParams/[3]vibratoFreq[style:knob]
    ",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/h:[1]otherParams/[4]vibratoGain[style:knob]
    ",0.5,0,1,0.01)*0.1;
    outGain = hslider("v:vocal/h:[1]otherParams/[5]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:vocal/[2]gate");

    gate = t+s : min(1);
    freq2 = freq1*bend*(os.osc(vibratoFreq)*vibratoGain+1);
    envelope = gate*gain*corrFactor : si.smooth(ba.tau2pole(envAttack));
};


//-------`(pm.)SFFormantModelBP_ui_MIDI`----------
// Ready-to-use MIDI-controllable source-filter vocal synthesizer.
//
// #### Usage
//
// ```
// SFFormantModelBP_ui_MIDI : _
// ```
//----------------------------------
SFFormantModelBP_ui_MIDI = SFFormantModelBP(voiceType,vowel,fricative,freq2,envelope)*outGain
with {
    freq1 = hslider("v:vocal/h:[0]midi/[0]freq[style:knob]",440,50,1000,0.01);
    bend = ba.semi2ratio(hslider("v:vocal/h:[0]midi/[1]bend[hidden:1][midi:pitchwheel]
    [style:knob]",0,-2,2,0.01)) : si.polySmooth(gate,0.999,1);
    gain = hslider("v:vocal/h:[0]midi/[2]gain[style:knob]
    ",0.9,0,1,0.01);
    envAttack = hslider("v:vocal/h:[0]midi/[3]envAttack[style:knob]
    ",10,0,30,0.01)*0.001;
    s = hslider("v:vocal/h:[0]midi/[4]sustain[hidden:1][midi:ctrl 64]
    [style:knob]",0,0,1,1);
    voiceType = hslider("v:vocal/h:[1]otherParams/[0]voiceType
    [style:knob]",0,0,4,1);
    vowel = hslider("v:vocal/h:[1]otherParams/[1]vowel
    [style:knob][midi:ctrl 1]",0,0,4,0.01) : si.smoo;
    fricative = hslider("v:vocal/h:[1]otherParams/[2]fricative
    [style:knob]",0,0,1,0.01) : si.smoo;
    vibratoFreq = hslider("v:vocal/h:[1]otherParams/[3]vibratoFreq[style:knob]
    ",6,1,10,0.01);
    vibratoGain = hslider("v:vocal/h:[1]otherParams/[4]vibratoGain[style:knob]
    ",0.5,0,1,0.01)*0.1;
    outGain = hslider("v:vocal/h:[1]otherParams/[5]outGain[style:knob]
    ",0.5,0,1,0.01);
    t = button("v:vocal/[2]gate");

     gate = t+s : min(1);
    freq2 = freq1*bend*(os.osc(vibratoFreq)*vibratoGain+1);
    envelope = gate*gain : si.smooth(ba.tau2pole(envAttack));
};


//=============================== Misc Functions =========================================
// Various miscellaneous functions.
//========================================================================================

//-------`(pm.)allpassNL`--------------
// Bidirectional block adding nonlinearities in both directions in a chain.
// Nonlinearities are created by modulating the coefficients of a passive
// allpass filter by the signal it is processing.
//
// #### Usage
//
// ```
// chain(... : allpassNL(nonlinearity) : ...)
// ```
//
// Where:
//
// * `nonlinearity`: amount of nonlinearity to be added (0-1)
//--------------------------------------
allpassNL(nonlinearity) = par(i,2,nlf),_
with {
    nlf = _ <: fi.allpassn(2,par(i,2,*(nonlinearity)*ma.PI));
};


//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// UNCATEGORIZED FUNCTIONS (TODO)
//%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

//-------`(pm.)modalModel`--------------
//
// Implement multiple resonance modes using resonant bandpass filters.
//
// #### Usage
//
// ```
// _ : modalModel(n, freqs, t60s, gains) : _
// ```
//
// Where:
//
// * `n`: number of given modes
// * `freqs` : list of filter center freqencies
// * `t60s` : list of mode resonance durations (in seconds)
// * `gains` : list of mode gains (0-1)
//
// For example, to generate a model with 2 modes (440 Hz and 660 Hz, a
// fifth) where the higher one decays faster and is attenuated:
//
// ```
// os.impulse : modalModel(2, (440, 660),
//                            (0.5, 0.25),
//                            (ba.db2linear(-1), ba.db2linear(-6)) : _
// ```
//
// Further reading: [Grumiaux et. al., 2017:
// Impulse-Response and CAD-Model-Based Physical Modeling in
// Faust](https://raw.githubusercontent.com/grame-cncm/faust/master-dev/tools/physicalModeling/ir2dsp/lacPaper2017.pdf)
//
//--------------------------------------
modalModel(n,modeFreqs,modeRes,modeGains) = _ <: par(i,n,modeFilter(freqs(i),res(i),gain(i))) :> _
with {
    freqs(i) = ba.take(i+1,modeFreqs);
    res(i) = ba.take(i+1,modeRes);
    gain(i) = ba.take(i+1,modeGains);
};

//-----------------------------`(pm.)rk_solve`----------------------------
// Solves the system of ordinary differential equations of any order using
// the explicit Runge-Kutta methods.
//
// #### Usage
//
// ```
// rk_solve(ts,ks, ni,h, eq,iv) : si.bus(outputs(eq))
// ```
//
// Where:
//
// * `ts,ks` : the Butcher tableau (see below)
// * `ni` : number of iterations at each tick, compile time constant
//          ni > 1 can improve accuracy but will degrade performance
// * `h`  : time step, run time constant, e.g. 1/ma.SR
// * `eq` : list of derivative functions
// * `iv` : list of initial values
//
// `rk_solve()` with the "standard" 1-4 tableaux and ni = 1:
// ```
// rk_solve_1 = rk_solve((0), (1), 1);
// rk_solve_2 = rk_solve((0,1/2), (1/2, 0,1), 1);
// rk_solve_3 = rk_solve((0,1/2,1), (1/2,-1,2, 1/6,2/3,1/6), 1);
// rk_solve_4 = rk_solve((0,1/2,1/2,1), (1/2,0,1/2,0,0,1, 1/6,1/3,1/3,1/6), 1);
// ```
//
// #### Example test program
//
// Suppose we have a system of differential equations:
// ```
// dx/dt = dx_dt(t,x,y,z)
// dy/dt = dy_dt(t,x,y,z)
// dz/dt = dz_dt(t,x,y,z)
// ```
// with initial conditions:
// ```
//    x(0) = x0
//    y(0) = y0
//    z(0) = z0
// ```
// and we want to solve it using this Butcher tableau:
// ```
//  0 |
// c2 | a21
// c3 | a31 a32
// c4 | a41 a42 a43
// -------------------
//    | b1  b1  b3  b4
// ```
//
// ```
// EQ(t,x,y,z) = dx_dt(t,x,y,z),
//               dy_dt(t,x,y,z),
//               dz_dt(t,x,y,z);
//
// IV = x0, y0, z0;
//
// TS = 0, c2, c3, c4;
// KS = a21,
//      a31, a32,
//      a41, a42, a43,
//      b1,  b2,  b3,  b4;
//
// process = rk_solve(TS,KS, 1,1/ma.SR, EQ,IV);
// ```
// Less abstract example which can actually be compiled/tested:
//
// ```
// // Lotka-Volterra equations parameterized by a,b,c,d:
// LV(a,b,c,d, t,x,y) =
//     a*x - b*x*y,
//     c*x*y - d*y;
//
// // Solved using the "standard" fourth-order method:
// process = rk_solve_4(
//     0.01,                  // time step
//     LV(0.1,0.02,0.03,0.4), // LV() with random parameters
//     (3,4)                  // initial values
// );
// ```
//
// #### References
//
// * <https://wikipedia.org/wiki/Runge%E2%80%93Kutta_methods>
//------------------------------------------------------------------------
rk_solve(ts,ks, ni,h, eq,iv) = doit
with {
    T = ba.time * h;
    O = outputs(ts);
    N = outputs(eq);

    // solution for the current tick, calculated at the previous
    // tick or `iv` if T == 0 (first tick).
    curr = par(i,N, ck_iv(ba.pick(iv,i))) with {
        ck_iv = case {
            (0) => _; (0.0) => _; // optimization
            (x) => select2(ba.time, x);
        };
    };

    // Why not simply `doit = tick ~ curr` ? This would output the
    // solution for T+h time. And we can't just add `: par(i,N,mem)`,
    // this won't work at T == 0 time if iv != 0, see ck_iv above.
    doit = (si.bus(N) <: tick, si.bus(N)) ~ curr :
            si.block(N), si.bus(N);

    // repeat `step` `ni` times using the `h/ni` time step, outputs
    // the solution for T+h time (next tick).
    tick = seq(i,ni, step(ts,ks, T+i*h/ni, h/ni, eq));

    // main function, outputs the solution for t+h time.
    step(ts,ks, t,h, eq) =
        seq(i,O, push_k(i+1, vec(i), eqt(i))) :
        ro.interleave(N,O+1) : par(i,N, vec(O))
    with {
        push_k(K, vec, eqt) = si.bus(N*K) <: si.bus(N*K), push
        with {
            push = ro.interleave(N,K) : par(i,N, vec) : eqt;
        };

        vec(K) = _, par(i,K, *(ba.pick(ks, K*(K-1)/2 + i))) :> _;
        eqt(K) = eq(t + h * ba.pick(ts, K)) : par(i,N, *(h));
    };
};

// rk_solve() with the "standard" 1-4 tableaux and ni = 1
rk_solve_1 = rk_solve((0), (1), 1);
rk_solve_2 = rk_solve((0,1/2), (1/2, 0,1), 1);
rk_solve_3 = rk_solve((0,1/2,1), (1/2,-1,2, 1/6,2/3,1/6), 1);
rk_solve_4 = rk_solve((0,1/2,1/2,1), (1/2,0,1/2,0,0,1, 1/6,1/3,1/3,1/6), 1);