/************************************************************************
 ************************************************************************
  	FAUST library file
	Copyright (C) 2003-2011 GRAME, Centre National de Creation Musicale
    ---------------------------------------------------------------------
    This program is free software; you can redistribute it and/or modify
    it under the terms of the GNU Lesser General Public License as 
	published by the Free Software Foundation; either version 2.1 of the 
	License, or (at your option) any later version.

    This program is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU Lesser General Public License for more details.

    You should have received a copy of the GNU Lesser General Public
 	License along with the GNU C Library; if not, write to the Free
  	Software Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA
  	02111-1307 USA. 
 ************************************************************************
 ************************************************************************/

declare name "MaxMSP compatibility Library";
declare author "GRAME";
declare copyright "GRAME";
declare version "1.0";
declare license "LGPL"; 

import("music.lib");


atodb = db2lin;

//-------------------------------------------------------------------------
// 
// Implementation of MaxMSP filtercoeff 
// 
//		from : Cookbook formulae for audio EQ biquad filter coefficients
//        by : Robert Bristow-Johnson  <rbj@audioimagination.com>
//       URL : http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt
// 
//-------------------------------------------------------------------------

filtercoeff(f0, dBgain, Q) = environment 
{
	//----------------------------------------
	// biquad coeffs for various filters 
	// usage : filtercoeff(f0, dBgain, Q).LPF 
	//----------------------------------------
	
		LPF = rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
				    b0 =  (1 - cos(w0))/2;
				    b1 =   1 - cos(w0);
				    b2 =  (1 - cos(w0))/2;
				    a0 =   1 + alpha;
				    a1 =  -2*cos(w0);
				    a2 =   1 - alpha;
 				};
  				
		HPF = rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
				    b0 =  (1 + cos(w0))/2;
				    b1 = -1 - cos(w0);
				    b2 =  (1 + cos(w0))/2;
				    a0 =   1 + alpha;
				    a1 =  -2*cos(w0);
				    a2 =   1 - alpha;
  				};
  				
		BPF = rbjcoef( a0, a1, a2, b0, b1, b2 ) // constant 0 dB peak gain
  				with {
				    b0 =   alpha;
				    b1 =   0;
				    b2 =  -alpha;
				    a0 =   1 + alpha;
				    a1 =  -2*cos(w0);
				    a2 =   1 - alpha;
 				};
  				
  	  notch = rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
				    b0 =   1;
				    b1 =  -2*cos(w0);
				    b2 =   1;
				    a0 =   1 + alpha;
				    a1 =  -2*cos(w0);
				    a2 =   1 - alpha;
  				};
  				
		APF = rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
				    b0 =   1 - alpha;
				    b1 =  -2*cos(w0);
				    b2 =   1 + alpha;
				    a0 =   1 + alpha;
				    a1 =  -2*cos(w0);
				    a2 =   1 - alpha;
 				};
  				
  peakingEQ	= rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
            		b0 =   1 + alpha*A;
            		b1 =  -2*cos(w0);
            		b2 =   1 - alpha*A;
            		a0 =   1 + alpha/A;
            		a1 =  -2*cos(w0);
            		a2 =   1 - alpha/A;
  				};
  				
  peakNotch	= rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
            		b0 =   1 + alpha*G;
            		b1 =  -2*cos(w0);
            		b2 =   1 - alpha*G;
            		a0 =   1 + alpha/G;
            		a1 =  -2*cos(w0);
            		a2 =   1 - alpha/G;
  				};
  				
   lowShelf = rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
             		b0 =    A*( (A+1) - (A-1)*cos(w0) + 2*sqrt(A)*alpha );
            		b1 =  2*A*( (A-1) - (A+1)*cos(w0)                   );
            		b2 =    A*( (A+1) - (A-1)*cos(w0) - 2*sqrt(A)*alpha );
            		a0 =        (A+1) + (A-1)*cos(w0) + 2*sqrt(A)*alpha;
            		a1 =   -2*( (A-1) + (A+1)*cos(w0)                   );
            		a2 =        (A+1) + (A-1)*cos(w0) - 2*sqrt(A)*alpha;
  				};
  				
  highShelf	= rbjcoef( a0, a1, a2, b0, b1, b2 ) 
  				with {
				    b0 =    A*( (A+1) + (A-1)*cos(w0) + 2*sqrt(A)*alpha );
				    b1 = -2*A*( (A-1) + (A+1)*cos(w0)                   );
				    b2 =    A*( (A+1) + (A-1)*cos(w0) - 2*sqrt(A)*alpha );
				    a0 =        (A+1) - (A-1)*cos(w0) + 2*sqrt(A)*alpha;
				    a1 =    2*( (A-1) - (A+1)*cos(w0)                   );
				    a2 =        (A+1) - (A-1)*cos(w0) - 2*sqrt(A)*alpha;
				};

	// --------------------- implementation ------------------------------
	
	// convert rbj coeffs to biquad coeffs
	rbjcoef(a0,a1,a2,b0,b1,b2) = (b0/a0, b1/a0, b2/a0,-a1/a0,-a2/a0);
	
	// common values
//	alpha 	= sin(w0)/(2*Q);
//	w0 		= 2*PI*f0/Fs;
	alpha 	= sin(w0)/(2*max(0.001,Q));
	w0 		= 2*PI*max(0,f0)/Fs;
	Fs 		= SR;
	A  		= 10^(dBgain/40);     			// (for peaking and shelving EQ filters only)
	G  		= sqrt(max(0.00001, dBgain));   // When gain is a linear values (i.e. not in dB)
};


//-------------------------------------------------------------------------
// Implementation of MaxMSP biquad~
// y[n] = a0 * x[n] + a1 * x[n-1] + a2 * x[n-2] + b1 * y[n-1] + b2 * y[n-2] 
//-------------------------------------------------------------------------

biquad(x,a0,a1,a2,b1,b2) = x : conv3(a0, a1, a2) : + ~ conv2(b1, b2)
	with {
		conv2(c0,c1,x) = c0*x+c1*x';
		conv3(c0,c1,c2,x) = c0*x+c1*x'+c2*x'';
	};


//-------------------------------------------------------------------------
// 
// Filters using filtercoeff and biquad
//
//-------------------------------------------------------------------------


// Low Pass Filter
LPF(x, f0, gain, Q) 		= x , filtercoeff(f0,gain,Q).LPF : biquad;

// High Pass Filter
HPF(x, f0, gain, Q) 		= x , filtercoeff(f0,gain,Q).HPF : biquad;

// Band Pass Filter
BPF(x, f0, gain, Q) 		= x , filtercoeff(f0,gain,Q).BPF : biquad;

// notch Filter
notch(x, f0, gain, Q) 		= x , filtercoeff(f0,gain,Q).notch : biquad;
	
// All Pass Filter
APF(x, f0, gain, Q) 		= x , filtercoeff(f0,gain,Q).APF : biquad;
	
// ????
peakingEQ(x, f0, gain, Q) 	= x , filtercoeff(f0,gain,Q).peakingEQ : biquad;
	
// Max peakNotch is like peakingEQ but with a linear gain
peakNotch(x, f0, gain, Q) 	= x , filtercoeff(f0,gain,Q).peakNotch : biquad;
	
// ????
lowShelf(x, f0, gain, Q) 	= x , filtercoeff(f0,gain,Q).lowShelf : biquad;
	
// ????
highShelf(x, f0, gain, Q) 	= x , filtercoeff(f0,gain,Q).highShelf : biquad;




//-------------------------------------------------------------------------
// Implementation of Max/MSP line~. Generate signal ramp or envelope 
// 
// USAGE : line(value, time)
// 	value : the desired output value
//	time  : the interpolation time to reach this value (in milliseconds)
//
// NOTE : the interpolation process is restarted every time the desired
// output value changes. The interpolation time is sampled only then.
//-------------------------------------------------------------------------

line (value, time) = state ~ ( _ , _ ) : ! , _ 
	with {
		state (t , c) = nt , if( nt <= 0 , value , c + (value - c) / nt)
		with {
			nt = if( value != value' , samples, t - 1) ;
			samples = time * SR / 1000.0 ;
		} ;
	} ;