schroeder3.m

MATLAB实现REVERB音效算法。经典！

function [y,b,a]=schroeder3(x,tn,tg,td,cg,cd,ag,ad,k)

%This is a reverberator based on Schroeder's design which consists of a tapped delay line followed by a  
%reverberator with 4 parallel comb filters in series with 2 cascaded allpass filters
%
%The structure is:  [y,b,a] = schroeder3(x,tn,tg,td,cg,cd,ag,ad,k) 
%
%where x = the input signal
%      tn = the number of taps of the tapped delay line
%      tg = a vector which contains the gain of each tap of the tapped delay line
%      td = a vector which contains the delay length of each tap of the tapped delay line
%      cg = a vector of length 4 which contains the gain of each of the comb filters (should be less than 1 for stability)	
%      cd = a vector of length 4 which contains the delay of each of the comb filters 
%      ag = the gain of the allpass filters (this should be less than 1 for stability)
%      ad = a vector of length 2 which contains the delay of each of the allpass filters 
%      k = the gain factor of the direct signal added to the cascaded allpass filters
%      y = the output signal
%      b = the numerator coefficients of the transfer function
%      a = the denominator coefficients of the transfer function
%
%
% Gautham J. Mysore - [e-mail]
%


%send the input signal througth the tapped delay line
[ty,b1,a1] = tdl(x,tn,tg,td);

%increase the delay lengths so that the total delays of the allpass filters are greater than that of the tapped delay line
ad = ad + td(tn);

%send the input signal through the reverberator
[ay,b2,a2]=schroeder2(x,cg,cd,ag,ad,k);

%find the total output
y = ty + ay;

%make the number of coefficients of the filters equal
if length(b1) > length(b2)
    b2=[b2 zeros(1,length(b1)-length(b2))];
else    
    b1=[b1 zeros(1,length(b2)-length(b1))];
end

if length(a1) > length (a2)
    a2=[a2 zeros(1,length(a1)-length(a2))];
else    
    a1=[a1 zeros(1,length(a2)-length(a1))];
end

%find the composite filter coefficients
[b,a]=parallelcoefficients(b1,a1,b2,a2);

% normalize the output signal
y = y/max(y);