音频 matlab 挑战：混响、FDN、滤波
Audio matlab challenge: Reverberation, FDN, Filtering
为了提供一些上下文，我试图在 matlab 中从头开始实现以下信号图，即反馈延迟网络 (FDN)。图片：FDN
使用合适的矩阵，延迟长度无关紧要，当馈入狄拉克脉冲时，几乎会产生白噪声。
我已经设法在代码中做到了这一点，但我的目标是另一个，因此我的问题。我想在每个延迟线 z^-m 之后应用过滤器 h(z)。图片：h(z)
更具体地说，我想在每条延迟线之后应用一个三分之一倍频程级联图形均衡器。目的是在整个结构上产生依赖于频率的衰减，并因此产生依赖于延迟的衰减。我已经成功设计了 SOS 形式的过滤器，但我的问题是：如何在结构中应用它？我假设在某个地方使用 sosfilt()，但我不确定。
我并没有为了目的而降低系统的顺序。阶数为16（16x16矩阵，16条延迟线，31x16双二阶滤波器）
第一个代码指的是无损 FDN，可安全运行，会产生白噪声。我评论了我在循环中引入过滤的失败尝试：% Filtering
不幸的是，我不能 post 所有 GEQ 条目，但我会在最后留下 8 个对应于前 8 个延迟。
所以，问题是我如何编写代码来过滤白噪声，在整个 FDN 结构中实现与频率相关的衰减。此外，尽管它可能在计算上效率低下，但我更愿意在没有更高级别功能的情况下应用它并基于我已经拥有的功能，即：适用于 GNU Octave
编辑：假设您必须使用差分方程逐个样本应用带通二阶滤波，您将如何递归地对 31 个串联波段执行此操作？一个显示在第二个代码部分。
% Feedback Delay Network - Mono
fs = 44100;
in = [ 1; 0 ]; % Dirac Impulse
in = [in; zeros(3*fs,1)]; % Space for reverb
% Householder matrix N=16
A = eye(4) - ( 2/4 * ones(4,1) * ones(1,4) );
MTX = .5 * [A -A -A -A; -A A -A -A; -A -A A -A; -A -A -A A];
N = size(MTX,1); % Matrix order
delays = [1500 1547 1602 1668 1745 1838 1947 2078 2232,...
2415 2623 2890 3196 3559 3989 4500]; % N=16 delays in samples
load('GEQ.mat'); % Third octave graphic equalizer calculated based
% on an atenuation-per-sample and scaled by delay.
% To be applied before or after each delayline
% Initialize signals
delayBuffer = max(delays) + max(delays)/10;
bufferN = zeros(delayBuffer,N); % Delay buffers
FB = zeros(1,N); % Feedback matrix output
in_dl = zeros(1,N); % Delay input
out_dl = zeros(1,N); % Delay output
nSample = length(in); % Number of samples
out = zeros(nSample,1); % Output
% FDN Computation
for sample = 1:nSample % each sample
for n = 1:N % each delayline
in_dl(n) = in(sample,1) + FB(n); % Input and Feedback matrix sum
[out_dl(n),bufferN(:,n)] = funcDelay( in_dl(n), bufferN(:,n),...
sample, delays(n) ); % Delaying
% Filtering
%out_dl(n) = sosfilt( GEQ(:,:,n), out_dl(n) );
end
out(sample,1) = 1/sqrt(2) * sum(out_dl); % Delay output sum
FB = out_dl * MTX; % Feedback matrix output recalculation
end
% Used function
function [out,buffer] = funcDelay(in,buffer,n,delay)
% Circular buffer indices
len = length(buffer);
indexC = mod(n-1, len) + 1; % Current
indexD = mod(n-delay-1, len) + 1; % Delay
out = buffer(indexD,1);
% Stores output on appropriate index
buffer(indexC,1) = in;
end
%sound(out,fs,16)
第二个代码段：将滤波器差分方程应用于信号。关于递归地为 31 个过滤器编码的建议？
in = (rand(1,100)*2)-1; % Example noise 100 samples
samples = length(in);
out = zeros(1,samples);
b0 = GEQ(1,1,1); % Coeffs extracted from actual GEQ
b1 = GEQ(1,2,1); a1 = GEQ(1,5,1);
b2 = GEQ(1,3,1); a2 = GEQ(1,6,1);
out(1) = b0 * in(1);
out(2) = b0 * in(2) + b1 * in(1) - a1 * out(1);
for n = 3:samples
out(n) = b0*in(n) + b1*in(n-1) + b2*in(n-2) - a1*out(n-1) - a2*out(n-2);
end
谢谢！
GEQ(:,:,1) = [0.6444 -1.2882 0.6438 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9958 0.9960 1.0000 -1.9957 0.9959;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9802 1.0000 -1.9757 0.9789;
1.0000 -1.9702 0.9752 1.0000 -1.9685 0.9735;
1.0000 -1.9609 0.9688 1.0000 -1.9587 0.9667;
1.0000 -1.9483 0.9608 1.0000 -1.9458 0.9584;
1.0000 -1.9311 0.9508 1.0000 -1.9281 0.9478;
1.0000 -1.9070 0.9381 1.0000 -1.9039 0.9350;
1.0000 -1.8738 0.9228 1.0000 -1.8698 0.9187;
1.0000 -1.8275 0.9043 1.0000 -1.8215 0.8980;
1.0000 -1.7608 0.8807 1.0000 -1.7538 0.8732;
1.0000 -1.6659 0.8520 1.0000 -1.6580 0.8432;
1.0000 -1.5308 0.8178 1.0000 -1.5209 0.8059;
1.0000 -1.3382 0.7756 1.0000 -1.3278 0.7616;
1.0000 -1.0671 0.7226 1.0000 -1.0607 0.7118;
1.0000 -0.7061 0.6745 1.0000 -0.6929 0.6388;
1.0000 -0.2324 0.6083 1.0000 -0.2311 0.5703;
1.0000 0.3354 0.5587 1.0000 0.3047 0.4869;
1.0000 0.9408 0.5246 1.0000 0.8392 0.4163;
1.0000 1.5310 0.6212 1.0000 1.2251 0.3584];
GEQ(:,:,2) = [0.6356 -1.2706 0.6350 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9789;
1.0000 -1.9702 0.9752 1.0000 -1.9684 0.9734;
1.0000 -1.9609 0.9689 1.0000 -1.9587 0.9666;
1.0000 -1.9483 0.9608 1.0000 -1.9458 0.9583;
1.0000 -1.9311 0.9509 1.0000 -1.9280 0.9478;
1.0000 -1.9070 0.9382 1.0000 -1.9039 0.9350;
1.0000 -1.8739 0.9228 1.0000 -1.8697 0.9186;
1.0000 -1.8276 0.9044 1.0000 -1.8214 0.8979;
1.0000 -1.7609 0.8808 1.0000 -1.7537 0.8731;
1.0000 -1.6660 0.8522 1.0000 -1.6579 0.8431;
1.0000 -1.5310 0.8180 1.0000 -1.5208 0.8057;
1.0000 -1.3384 0.7758 1.0000 -1.3276 0.7614;
1.0000 -1.0672 0.7227 1.0000 -1.0606 0.7116;
1.0000 -0.7063 0.6751 1.0000 -0.6927 0.6382;
1.0000 -0.2324 0.6088 1.0000 -0.2311 0.5697;
1.0000 0.3359 0.5598 1.0000 0.3042 0.4858;
1.0000 0.9423 0.5263 1.0000 0.8375 0.4146;
1.0000 1.5349 0.6247 1.0000 1.2195 0.3537];
GEQ(:,:,3) = [0.6255 -1.2504 0.6249 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9788;
1.0000 -1.9702 0.9752 1.0000 -1.9684 0.9734;
1.0000 -1.9610 0.9689 1.0000 -1.9587 0.9666;
1.0000 -1.9483 0.9609 1.0000 -1.9458 0.9583;
1.0000 -1.9312 0.9509 1.0000 -1.9280 0.9477;
1.0000 -1.9071 0.9382 1.0000 -1.9038 0.9349;
1.0000 -1.8739 0.9229 1.0000 -1.8697 0.9185;
1.0000 -1.8277 0.9045 1.0000 -1.8213 0.8978;
1.0000 -1.7610 0.8810 1.0000 -1.7536 0.8730;
1.0000 -1.6662 0.8523 1.0000 -1.6577 0.8429;
1.0000 -1.5312 0.8182 1.0000 -1.5206 0.8055;
1.0000 -1.3385 0.7761 1.0000 -1.3274 0.7612;
1.0000 -1.0674 0.7229 1.0000 -1.0605 0.7115;
1.0000 -0.7066 0.6757 1.0000 -0.6925 0.6376;
1.0000 -0.2324 0.6095 1.0000 -0.2311 0.5690;
1.0000 0.3363 0.5610 1.0000 0.3037 0.4845;
1.0000 0.9440 0.5282 1.0000 0.8355 0.4126;
1.0000 1.5395 0.6288 1.0000 1.2128 0.3482];
GEQ(:,:,4) = [0.6136 -1.2265 0.6130 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9939 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9809 0.9829;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9788;
1.0000 -1.9702 0.9753 1.0000 -1.9684 0.9734;
1.0000 -1.9610 0.9689 1.0000 -1.9586 0.9665;
1.0000 -1.9484 0.9609 1.0000 -1.9457 0.9582;
1.0000 -1.9312 0.9510 1.0000 -1.9279 0.9476;
1.0000 -1.9072 0.9383 1.0000 -1.9037 0.9348;
1.0000 -1.8740 0.9230 1.0000 -1.8696 0.9184;
1.0000 -1.8278 0.9046 1.0000 -1.8211 0.8977;
1.0000 -1.7612 0.8811 1.0000 -1.7534 0.8728;
1.0000 -1.6663 0.8525 1.0000 -1.6575 0.8427;
1.0000 -1.5314 0.8184 1.0000 -1.5204 0.8053;
1.0000 -1.3388 0.7764 1.0000 -1.3272 0.7609;
1.0000 -1.0675 0.7232 1.0000 -1.0604 0.7112;
1.0000 -0.7069 0.6764 1.0000 -0.6922 0.6368;
1.0000 -0.2325 0.6103 1.0000 -0.2310 0.5681;
1.0000 0.3369 0.5624 1.0000 0.3030 0.4830;
1.0000 0.9460 0.5304 1.0000 0.8332 0.4102;
1.0000 1.5449 0.6336 1.0000 1.2047 0.3416];
GEQ(:,:,5) = [0.5999 -1.1993 0.5993 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9939 1.0000 -1.9928 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9826 0.9846 1.0000 -1.9809 0.9829;
1.0000 -1.9772 0.9803 1.0000 -1.9756 0.9788;
1.0000 -1.9703 0.9753 1.0000 -1.9683 0.9733;
1.0000 -1.9611 0.9690 1.0000 -1.9586 0.9665;
1.0000 -1.9485 0.9610 1.0000 -1.9456 0.9582;
1.0000 -1.9313 0.9511 1.0000 -1.9278 0.9475;
1.0000 -1.9072 0.9384 1.0000 -1.9037 0.9348;
1.0000 -1.8741 0.9231 1.0000 -1.8695 0.9183;
1.0000 -1.8280 0.9048 1.0000 -1.8210 0.8975;
1.0000 -1.7614 0.8813 1.0000 -1.7532 0.8726;
1.0000 -1.6665 0.8527 1.0000 -1.6573 0.8425;
1.0000 -1.5316 0.8187 1.0000 -1.5201 0.8050;
1.0000 -1.3390 0.7767 1.0000 -1.3270 0.7605;
1.0000 -1.0677 0.7234 1.0000 -1.0602 0.7109;
1.0000 -0.7072 0.6773 1.0000 -0.6918 0.6359;
1.0000 -0.2325 0.6112 1.0000 -0.2310 0.5672;
1.0000 0.3376 0.5640 1.0000 0.3022 0.4813;
1.0000 0.9484 0.5331 1.0000 0.8304 0.4073;
1.0000 1.5511 0.6393 1.0000 1.1953 0.3338];
GEQ(:,:,6) = [0.5839 -1.1672 0.5833 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9886 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9826 0.9846 1.0000 -1.9808 0.9828;
1.0000 -1.9772 0.9804 1.0000 -1.9756 0.9787;
1.0000 -1.9703 0.9753 1.0000 -1.9683 0.9733;
1.0000 -1.9611 0.9691 1.0000 -1.9585 0.9664;
1.0000 -1.9485 0.9611 1.0000 -1.9456 0.9581;
1.0000 -1.9314 0.9512 1.0000 -1.9277 0.9475;
1.0000 -1.9073 0.9385 1.0000 -1.9036 0.9347;
1.0000 -1.8742 0.9232 1.0000 -1.8693 0.9182;
1.0000 -1.8281 0.9050 1.0000 -1.8208 0.8973;
1.0000 -1.7616 0.8815 1.0000 -1.7530 0.8724;
1.0000 -1.6668 0.8530 1.0000 -1.6571 0.8422;
1.0000 -1.5319 0.8191 1.0000 -1.5198 0.8046;
1.0000 -1.3393 0.7771 1.0000 -1.3266 0.7601;
1.0000 -1.0678 0.7237 1.0000 -1.0600 0.7106;
1.0000 -0.7076 0.6783 1.0000 -0.6914 0.6348;
1.0000 -0.2325 0.6123 1.0000 -0.2310 0.5660;
1.0000 0.3384 0.5659 1.0000 0.3013 0.4792;
1.0000 0.9513 0.5363 1.0000 0.8271 0.4039;
1.0000 1.5586 0.6460 1.0000 1.1838 0.3244];
GEQ(:,:,7) = [0.5656 -1.1306 0.5651 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9978;
1.0000 -1.9975 0.9976 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9893 0.9900 1.0000 -1.9886 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9827 0.9847 1.0000 -1.9808 0.9828;
1.0000 -1.9773 0.9804 1.0000 -1.9755 0.9787;
1.0000 -1.9704 0.9754 1.0000 -1.9682 0.9732;
1.0000 -1.9612 0.9691 1.0000 -1.9584 0.9663;
1.0000 -1.9486 0.9611 1.0000 -1.9455 0.9580;
1.0000 -1.9315 0.9513 1.0000 -1.9276 0.9473;
1.0000 -1.9074 0.9386 1.0000 -1.9035 0.9345;
1.0000 -1.8744 0.9234 1.0000 -1.8692 0.9180;
1.0000 -1.8283 0.9052 1.0000 -1.8206 0.8971;
1.0000 -1.7618 0.8818 1.0000 -1.7527 0.8721;
1.0000 -1.6670 0.8533 1.0000 -1.6568 0.8419;
1.0000 -1.5323 0.8195 1.0000 -1.5194 0.8041;
1.0000 -1.3397 0.7776 1.0000 -1.3262 0.7595;
1.0000 -1.0681 0.7241 1.0000 -1.0598 0.7102;
1.0000 -0.7080 0.6795 1.0000 -0.6910 0.6335;
1.0000 -0.2326 0.6136 1.0000 -0.2309 0.5646;
1.0000 0.3393 0.5682 1.0000 0.3002 0.4768;
1.0000 0.9546 0.5400 1.0000 0.8231 0.3999;
1.0000 1.5672 0.6537 1.0000 1.1702 0.3133];
GEQ(:,:,8) = [0.5444 -1.0883 0.5439 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9978;
1.0000 -1.9975 0.9976 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9911 0.9916;
1.0000 -1.9893 0.9901 1.0000 -1.9886 0.9894;
1.0000 -1.9862 0.9874 1.0000 -1.9855 0.9867;
1.0000 -1.9827 0.9847 1.0000 -1.9807 0.9827;
1.0000 -1.9773 0.9805 1.0000 -1.9755 0.9786;
1.0000 -1.9705 0.9755 1.0000 -1.9681 0.9731;
1.0000 -1.9613 0.9692 1.0000 -1.9583 0.9662;
1.0000 -1.9487 0.9612 1.0000 -1.9454 0.9579;
1.0000 -1.9316 0.9514 1.0000 -1.9275 0.9472;
1.0000 -1.9076 0.9387 1.0000 -1.9033 0.9344;
1.0000 -1.8745 0.9235 1.0000 -1.8690 0.9179;
1.0000 -1.8286 0.9054 1.0000 -1.8203 0.8968;
1.0000 -1.7621 0.8821 1.0000 -1.7524 0.8718;
1.0000 -1.6674 0.8537 1.0000 -1.6564 0.8415;
1.0000 -1.5327 0.8200 1.0000 -1.5190 0.8036;
1.0000 -1.3401 0.7782 1.0000 -1.3258 0.7589;
1.0000 -1.0683 0.7246 1.0000 -1.0595 0.7098;
1.0000 -0.7085 0.6810 1.0000 -0.6904 0.6319;
1.0000 -0.2327 0.6152 1.0000 -0.2309 0.5630;
1.0000 0.3403 0.5708 1.0000 0.2990 0.4740;
1.0000 0.9586 0.5445 1.0000 0.8184 0.3951;
1.0000 1.5774 0.6630 1.0000 1.1537 0.2999];
因此，一个样本一个样本地从 FDN 的狄拉克脉冲中滤除噪声的效率相当低的方法是再添加 2 个缓冲区和计算 31 个级联双二阶滤波器的差分方程的方法
（任何提高计算速度的建议都在下面评论）
% Feedback Delay Network - Mono
% Creates impulse response based on reverberation times
fs = 44100;
in = [ 1; 0 ]; % Dirac Impulse
in = [in; zeros(3*fs,1)]; % Space for reverb
% Householder matrix N=16
A = eye(4) - ( 2/4 * ones(4,1) * ones(1,4) );
MTX = .5 * [A -A -A -A; -A A -A -A; -A -A A -A; -A -A -A A];
N = size(MTX,1); % Matrix order
delays = randi([dmin dmax],N); % N delays
load('GEQ.mat'); % Third octave graphic equalizer calculated based
% on an atenuation-per-sample and scaled by delays.
% SOS Form Size 31x6xN
% Initialize signals
delayBuffer = max(delays) + max(delays)/10;
bufferN = zeros(delayBuffer,N); % Delay buffers
buffFin = zeros(31,3,N); % New buffers for filter outputs
buffFout = zeros(31,3,N);
FB = zeros(1,N); % Feedback matrix output
in_dl = zeros(1,N); % Delay input
out_dl = zeros(1,N); % Delay output
out_fdl = zeros(1,N); % Filtered delay output
nSample = length(in); % Number of samples
out = zeros(nSample,1); % Output
% FDN Computation
for sample = 1:nSample % each sample
for n = 1:N % each delayline
in_dl(n) = in(sample,1) + FB(n); % Input and Feedback matrix sum
% Delaying
[out_dl(n),bufferN(:,n)] = funcDelay( in_dl(n), bufferN(:,n),...
sample, delays(n) );
% Filtering
[out_fdl(n),buffFin(:,:,n),buffFout(:,:,n)] = funcFilterGEQ(...
GEQ(:,:,n), out_dl(n), buffFin(:,:,n), buffFout(:,:,n));
end
out(sample,1) = sum(out_fdl); % Filtered delay output sum
FB = out_fdl * MTX; % Feedback matrix output recalculation
end
% Used functions
function [out,buffer] = funcDelay( in, buffer, n, delay)
% Circular buffer indices
len = length(buffer);
indexC = mod(n-1, len) + 1; % Current
indexD = mod(n-delay-1, len) + 1; % Delay
out = buffer(indexD,1);
% Stores output on appropriate index
buffer(indexC,1) = in;
end
function [outfilt,buffin,buffout] = funcFilterGEQ( GEQ, in, buffin, buffout)
% Sample based filter cascading
nBands = size(GEQ,1);
out = zeros(1,nBands+1);
out(1) = in; % Stores value pre-filtered
for b = 1:nBands % Performs series bandpass filtering
[out(b+1),buffin(b,:),buffout(b,:)] = funcBandpassFilt( out(b),...
GEQ(b,1:3), GEQ(b,5:6), buffin(b,:), buffout(b,:) );
end
outfilt = out(end); % Outputs final value post-filtered
end
function [out,buffin,buffout] = funcBandpassFilt( in, bs, as, buffin, buffout)
% Sample based filtering
buffin(3) = buffin(2); % Valid for biquad filters
buffin(2) = buffin(1);
buffin(1) = in; % Sequential indexing
buffout(1) = bs(1)*buffin(1) + bs(2)*buffin(2) + bs(3)*buffin(3)...
- as(1)*buffout(2) - as(2)*buffout(3);
out = buffout(1);
% Outputs calculation based on 3 latest values
buffout(3) = buffout(2);
buffout(2) = buffout(1);
buffout(1) = 0;
end
为了提供一些上下文，我试图在 matlab 中从头开始实现以下信号图，即反馈延迟网络 (FDN)。图片：FDN
使用合适的矩阵，延迟长度无关紧要，当馈入狄拉克脉冲时，几乎会产生白噪声。
我已经设法在代码中做到了这一点，但我的目标是另一个，因此我的问题。我想在每个延迟线 z^-m 之后应用过滤器 h(z)。图片：h(z)
更具体地说，我想在每条延迟线之后应用一个三分之一倍频程级联图形均衡器。目的是在整个结构上产生依赖于频率的衰减，并因此产生依赖于延迟的衰减。我已经成功设计了 SOS 形式的过滤器，但我的问题是：如何在结构中应用它？我假设在某个地方使用 sosfilt()，但我不确定。
我并没有为了目的而降低系统的顺序。阶数为16（16x16矩阵，16条延迟线，31x16双二阶滤波器）
第一个代码指的是无损 FDN，可安全运行，会产生白噪声。我评论了我在循环中引入过滤的失败尝试：% Filtering
不幸的是，我不能 post 所有 GEQ 条目，但我会在最后留下 8 个对应于前 8 个延迟。
所以，问题是我如何编写代码来过滤白噪声，在整个 FDN 结构中实现与频率相关的衰减。此外，尽管它可能在计算上效率低下，但我更愿意在没有更高级别功能的情况下应用它并基于我已经拥有的功能，即：适用于 GNU Octave
编辑：假设您必须使用差分方程逐个样本应用带通二阶滤波，您将如何递归地对 31 个串联波段执行此操作？一个显示在第二个代码部分。
% Feedback Delay Network - Mono
fs = 44100;
in = [ 1; 0 ]; % Dirac Impulse
in = [in; zeros(3*fs,1)]; % Space for reverb
% Householder matrix N=16
A = eye(4) - ( 2/4 * ones(4,1) * ones(1,4) );
MTX = .5 * [A -A -A -A; -A A -A -A; -A -A A -A; -A -A -A A];
N = size(MTX,1); % Matrix order
delays = [1500 1547 1602 1668 1745 1838 1947 2078 2232,...
2415 2623 2890 3196 3559 3989 4500]; % N=16 delays in samples
load('GEQ.mat'); % Third octave graphic equalizer calculated based
% on an atenuation-per-sample and scaled by delay.
% To be applied before or after each delayline
% Initialize signals
delayBuffer = max(delays) + max(delays)/10;
bufferN = zeros(delayBuffer,N); % Delay buffers
FB = zeros(1,N); % Feedback matrix output
in_dl = zeros(1,N); % Delay input
out_dl = zeros(1,N); % Delay output
nSample = length(in); % Number of samples
out = zeros(nSample,1); % Output
% FDN Computation
for sample = 1:nSample % each sample
for n = 1:N % each delayline
in_dl(n) = in(sample,1) + FB(n); % Input and Feedback matrix sum
[out_dl(n),bufferN(:,n)] = funcDelay( in_dl(n), bufferN(:,n),...
sample, delays(n) ); % Delaying
% Filtering
%out_dl(n) = sosfilt( GEQ(:,:,n), out_dl(n) );
end
out(sample,1) = 1/sqrt(2) * sum(out_dl); % Delay output sum
FB = out_dl * MTX; % Feedback matrix output recalculation
end
% Used function
function [out,buffer] = funcDelay(in,buffer,n,delay)
% Circular buffer indices
len = length(buffer);
indexC = mod(n-1, len) + 1; % Current
indexD = mod(n-delay-1, len) + 1; % Delay
out = buffer(indexD,1);
% Stores output on appropriate index
buffer(indexC,1) = in;
end
%sound(out,fs,16)
第二个代码段：将滤波器差分方程应用于信号。关于递归地为 31 个过滤器编码的建议？
in = (rand(1,100)*2)-1; % Example noise 100 samples
samples = length(in);
out = zeros(1,samples);
b0 = GEQ(1,1,1); % Coeffs extracted from actual GEQ
b1 = GEQ(1,2,1); a1 = GEQ(1,5,1);
b2 = GEQ(1,3,1); a2 = GEQ(1,6,1);
out(1) = b0 * in(1);
out(2) = b0 * in(2) + b1 * in(1) - a1 * out(1);
for n = 3:samples
out(n) = b0*in(n) + b1*in(n-1) + b2*in(n-2) - a1*out(n-1) - a2*out(n-2);
end
谢谢！
GEQ(:,:,1) = [0.6444 -1.2882 0.6438 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9958 0.9960 1.0000 -1.9957 0.9959;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9802 1.0000 -1.9757 0.9789;
1.0000 -1.9702 0.9752 1.0000 -1.9685 0.9735;
1.0000 -1.9609 0.9688 1.0000 -1.9587 0.9667;
1.0000 -1.9483 0.9608 1.0000 -1.9458 0.9584;
1.0000 -1.9311 0.9508 1.0000 -1.9281 0.9478;
1.0000 -1.9070 0.9381 1.0000 -1.9039 0.9350;
1.0000 -1.8738 0.9228 1.0000 -1.8698 0.9187;
1.0000 -1.8275 0.9043 1.0000 -1.8215 0.8980;
1.0000 -1.7608 0.8807 1.0000 -1.7538 0.8732;
1.0000 -1.6659 0.8520 1.0000 -1.6580 0.8432;
1.0000 -1.5308 0.8178 1.0000 -1.5209 0.8059;
1.0000 -1.3382 0.7756 1.0000 -1.3278 0.7616;
1.0000 -1.0671 0.7226 1.0000 -1.0607 0.7118;
1.0000 -0.7061 0.6745 1.0000 -0.6929 0.6388;
1.0000 -0.2324 0.6083 1.0000 -0.2311 0.5703;
1.0000 0.3354 0.5587 1.0000 0.3047 0.4869;
1.0000 0.9408 0.5246 1.0000 0.8392 0.4163;
1.0000 1.5310 0.6212 1.0000 1.2251 0.3584];
GEQ(:,:,2) = [0.6356 -1.2706 0.6350 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9789;
1.0000 -1.9702 0.9752 1.0000 -1.9684 0.9734;
1.0000 -1.9609 0.9689 1.0000 -1.9587 0.9666;
1.0000 -1.9483 0.9608 1.0000 -1.9458 0.9583;
1.0000 -1.9311 0.9509 1.0000 -1.9280 0.9478;
1.0000 -1.9070 0.9382 1.0000 -1.9039 0.9350;
1.0000 -1.8739 0.9228 1.0000 -1.8697 0.9186;
1.0000 -1.8276 0.9044 1.0000 -1.8214 0.8979;
1.0000 -1.7609 0.8808 1.0000 -1.7537 0.8731;
1.0000 -1.6660 0.8522 1.0000 -1.6579 0.8431;
1.0000 -1.5310 0.8180 1.0000 -1.5208 0.8057;
1.0000 -1.3384 0.7758 1.0000 -1.3276 0.7614;
1.0000 -1.0672 0.7227 1.0000 -1.0606 0.7116;
1.0000 -0.7063 0.6751 1.0000 -0.6927 0.6382;
1.0000 -0.2324 0.6088 1.0000 -0.2311 0.5697;
1.0000 0.3359 0.5598 1.0000 0.3042 0.4858;
1.0000 0.9423 0.5263 1.0000 0.8375 0.4146;
1.0000 1.5349 0.6247 1.0000 1.2195 0.3537];
GEQ(:,:,3) = [0.6255 -1.2504 0.6249 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9938 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9873 1.0000 -1.9856 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9810 0.9830;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9788;
1.0000 -1.9702 0.9752 1.0000 -1.9684 0.9734;
1.0000 -1.9610 0.9689 1.0000 -1.9587 0.9666;
1.0000 -1.9483 0.9609 1.0000 -1.9458 0.9583;
1.0000 -1.9312 0.9509 1.0000 -1.9280 0.9477;
1.0000 -1.9071 0.9382 1.0000 -1.9038 0.9349;
1.0000 -1.8739 0.9229 1.0000 -1.8697 0.9185;
1.0000 -1.8277 0.9045 1.0000 -1.8213 0.8978;
1.0000 -1.7610 0.8810 1.0000 -1.7536 0.8730;
1.0000 -1.6662 0.8523 1.0000 -1.6577 0.8429;
1.0000 -1.5312 0.8182 1.0000 -1.5206 0.8055;
1.0000 -1.3385 0.7761 1.0000 -1.3274 0.7612;
1.0000 -1.0674 0.7229 1.0000 -1.0605 0.7115;
1.0000 -0.7066 0.6757 1.0000 -0.6925 0.6376;
1.0000 -0.2324 0.6095 1.0000 -0.2311 0.5690;
1.0000 0.3363 0.5610 1.0000 0.3037 0.4845;
1.0000 0.9440 0.5282 1.0000 0.8355 0.4126;
1.0000 1.5395 0.6288 1.0000 1.2128 0.3482];
GEQ(:,:,4) = [0.6136 -1.2265 0.6130 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9984 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9939 1.0000 -1.9929 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9895;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9825 0.9845 1.0000 -1.9809 0.9829;
1.0000 -1.9771 0.9803 1.0000 -1.9757 0.9788;
1.0000 -1.9702 0.9753 1.0000 -1.9684 0.9734;
1.0000 -1.9610 0.9689 1.0000 -1.9586 0.9665;
1.0000 -1.9484 0.9609 1.0000 -1.9457 0.9582;
1.0000 -1.9312 0.9510 1.0000 -1.9279 0.9476;
1.0000 -1.9072 0.9383 1.0000 -1.9037 0.9348;
1.0000 -1.8740 0.9230 1.0000 -1.8696 0.9184;
1.0000 -1.8278 0.9046 1.0000 -1.8211 0.8977;
1.0000 -1.7612 0.8811 1.0000 -1.7534 0.8728;
1.0000 -1.6663 0.8525 1.0000 -1.6575 0.8427;
1.0000 -1.5314 0.8184 1.0000 -1.5204 0.8053;
1.0000 -1.3388 0.7764 1.0000 -1.3272 0.7609;
1.0000 -1.0675 0.7232 1.0000 -1.0604 0.7112;
1.0000 -0.7069 0.6764 1.0000 -0.6922 0.6368;
1.0000 -0.2325 0.6103 1.0000 -0.2310 0.5681;
1.0000 0.3369 0.5624 1.0000 0.3030 0.4830;
1.0000 0.9460 0.5304 1.0000 0.8332 0.4102;
1.0000 1.5449 0.6336 1.0000 1.2047 0.3416];
GEQ(:,:,5) = [0.5999 -1.1993 0.5993 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9980 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9973 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9935 0.9939 1.0000 -1.9928 0.9932;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9887 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9826 0.9846 1.0000 -1.9809 0.9829;
1.0000 -1.9772 0.9803 1.0000 -1.9756 0.9788;
1.0000 -1.9703 0.9753 1.0000 -1.9683 0.9733;
1.0000 -1.9611 0.9690 1.0000 -1.9586 0.9665;
1.0000 -1.9485 0.9610 1.0000 -1.9456 0.9582;
1.0000 -1.9313 0.9511 1.0000 -1.9278 0.9475;
1.0000 -1.9072 0.9384 1.0000 -1.9037 0.9348;
1.0000 -1.8741 0.9231 1.0000 -1.8695 0.9183;
1.0000 -1.8280 0.9048 1.0000 -1.8210 0.8975;
1.0000 -1.7614 0.8813 1.0000 -1.7532 0.8726;
1.0000 -1.6665 0.8527 1.0000 -1.6573 0.8425;
1.0000 -1.5316 0.8187 1.0000 -1.5201 0.8050;
1.0000 -1.3390 0.7767 1.0000 -1.3270 0.7605;
1.0000 -1.0677 0.7234 1.0000 -1.0602 0.7109;
1.0000 -0.7072 0.6773 1.0000 -0.6918 0.6359;
1.0000 -0.2325 0.6112 1.0000 -0.2310 0.5672;
1.0000 0.3376 0.5640 1.0000 0.3022 0.4813;
1.0000 0.9484 0.5331 1.0000 0.8304 0.4073;
1.0000 1.5511 0.6393 1.0000 1.1953 0.3338];
GEQ(:,:,6) = [0.5839 -1.1672 0.5833 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9979;
1.0000 -1.9975 0.9975 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9892 0.9900 1.0000 -1.9886 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9826 0.9846 1.0000 -1.9808 0.9828;
1.0000 -1.9772 0.9804 1.0000 -1.9756 0.9787;
1.0000 -1.9703 0.9753 1.0000 -1.9683 0.9733;
1.0000 -1.9611 0.9691 1.0000 -1.9585 0.9664;
1.0000 -1.9485 0.9611 1.0000 -1.9456 0.9581;
1.0000 -1.9314 0.9512 1.0000 -1.9277 0.9475;
1.0000 -1.9073 0.9385 1.0000 -1.9036 0.9347;
1.0000 -1.8742 0.9232 1.0000 -1.8693 0.9182;
1.0000 -1.8281 0.9050 1.0000 -1.8208 0.8973;
1.0000 -1.7616 0.8815 1.0000 -1.7530 0.8724;
1.0000 -1.6668 0.8530 1.0000 -1.6571 0.8422;
1.0000 -1.5319 0.8191 1.0000 -1.5198 0.8046;
1.0000 -1.3393 0.7771 1.0000 -1.3266 0.7601;
1.0000 -1.0678 0.7237 1.0000 -1.0600 0.7106;
1.0000 -0.7076 0.6783 1.0000 -0.6914 0.6348;
1.0000 -0.2325 0.6123 1.0000 -0.2310 0.5660;
1.0000 0.3384 0.5659 1.0000 0.3013 0.4792;
1.0000 0.9513 0.5363 1.0000 0.8271 0.4039;
1.0000 1.5586 0.6460 1.0000 1.1838 0.3244];
GEQ(:,:,7) = [0.5656 -1.1306 0.5651 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9978;
1.0000 -1.9975 0.9976 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9912 0.9917;
1.0000 -1.9893 0.9900 1.0000 -1.9886 0.9894;
1.0000 -1.9861 0.9874 1.0000 -1.9855 0.9868;
1.0000 -1.9827 0.9847 1.0000 -1.9808 0.9828;
1.0000 -1.9773 0.9804 1.0000 -1.9755 0.9787;
1.0000 -1.9704 0.9754 1.0000 -1.9682 0.9732;
1.0000 -1.9612 0.9691 1.0000 -1.9584 0.9663;
1.0000 -1.9486 0.9611 1.0000 -1.9455 0.9580;
1.0000 -1.9315 0.9513 1.0000 -1.9276 0.9473;
1.0000 -1.9074 0.9386 1.0000 -1.9035 0.9345;
1.0000 -1.8744 0.9234 1.0000 -1.8692 0.9180;
1.0000 -1.8283 0.9052 1.0000 -1.8206 0.8971;
1.0000 -1.7618 0.8818 1.0000 -1.7527 0.8721;
1.0000 -1.6670 0.8533 1.0000 -1.6568 0.8419;
1.0000 -1.5323 0.8195 1.0000 -1.5194 0.8041;
1.0000 -1.3397 0.7776 1.0000 -1.3262 0.7595;
1.0000 -1.0681 0.7241 1.0000 -1.0598 0.7102;
1.0000 -0.7080 0.6795 1.0000 -0.6910 0.6335;
1.0000 -0.2326 0.6136 1.0000 -0.2309 0.5646;
1.0000 0.3393 0.5682 1.0000 0.3002 0.4768;
1.0000 0.9546 0.5400 1.0000 0.8231 0.3999;
1.0000 1.5672 0.6537 1.0000 1.1702 0.3133];
GEQ(:,:,8) = [0.5444 -1.0883 0.5439 1.0000 -1.9989 0.9989;
1.0000 -1.9987 0.9987 1.0000 -1.9987 0.9987;
1.0000 -1.9984 0.9985 1.0000 -1.9983 0.9983;
1.0000 -1.9980 0.9981 1.0000 -1.9978 0.9978;
1.0000 -1.9975 0.9976 1.0000 -1.9972 0.9973;
1.0000 -1.9967 0.9968 1.0000 -1.9966 0.9967;
1.0000 -1.9959 0.9960 1.0000 -1.9957 0.9958;
1.0000 -1.9949 0.9951 1.0000 -1.9944 0.9946;
1.0000 -1.9936 0.9939 1.0000 -1.9928 0.9931;
1.0000 -1.9915 0.9920 1.0000 -1.9911 0.9916;
1.0000 -1.9893 0.9901 1.0000 -1.9886 0.9894;
1.0000 -1.9862 0.9874 1.0000 -1.9855 0.9867;
1.0000 -1.9827 0.9847 1.0000 -1.9807 0.9827;
1.0000 -1.9773 0.9805 1.0000 -1.9755 0.9786;
1.0000 -1.9705 0.9755 1.0000 -1.9681 0.9731;
1.0000 -1.9613 0.9692 1.0000 -1.9583 0.9662;
1.0000 -1.9487 0.9612 1.0000 -1.9454 0.9579;
1.0000 -1.9316 0.9514 1.0000 -1.9275 0.9472;
1.0000 -1.9076 0.9387 1.0000 -1.9033 0.9344;
1.0000 -1.8745 0.9235 1.0000 -1.8690 0.9179;
1.0000 -1.8286 0.9054 1.0000 -1.8203 0.8968;
1.0000 -1.7621 0.8821 1.0000 -1.7524 0.8718;
1.0000 -1.6674 0.8537 1.0000 -1.6564 0.8415;
1.0000 -1.5327 0.8200 1.0000 -1.5190 0.8036;
1.0000 -1.3401 0.7782 1.0000 -1.3258 0.7589;
1.0000 -1.0683 0.7246 1.0000 -1.0595 0.7098;
1.0000 -0.7085 0.6810 1.0000 -0.6904 0.6319;
1.0000 -0.2327 0.6152 1.0000 -0.2309 0.5630;
1.0000 0.3403 0.5708 1.0000 0.2990 0.4740;
1.0000 0.9586 0.5445 1.0000 0.8184 0.3951;
1.0000 1.5774 0.6630 1.0000 1.1537 0.2999];
因此，一个样本一个样本地从 FDN 的狄拉克脉冲中滤除噪声的效率相当低的方法是再添加 2 个缓冲区和计算 31 个级联双二阶滤波器的差分方程的方法 （任何提高计算速度的建议都在下面评论）
% Feedback Delay Network - Mono
% Creates impulse response based on reverberation times
fs = 44100;
in = [ 1; 0 ]; % Dirac Impulse
in = [in; zeros(3*fs,1)]; % Space for reverb
% Householder matrix N=16
A = eye(4) - ( 2/4 * ones(4,1) * ones(1,4) );
MTX = .5 * [A -A -A -A; -A A -A -A; -A -A A -A; -A -A -A A];
N = size(MTX,1); % Matrix order
delays = randi([dmin dmax],N); % N delays
load('GEQ.mat'); % Third octave graphic equalizer calculated based
% on an atenuation-per-sample and scaled by delays.
% SOS Form Size 31x6xN
% Initialize signals
delayBuffer = max(delays) + max(delays)/10;
bufferN = zeros(delayBuffer,N); % Delay buffers
buffFin = zeros(31,3,N); % New buffers for filter outputs
buffFout = zeros(31,3,N);
FB = zeros(1,N); % Feedback matrix output
in_dl = zeros(1,N); % Delay input
out_dl = zeros(1,N); % Delay output
out_fdl = zeros(1,N); % Filtered delay output
nSample = length(in); % Number of samples
out = zeros(nSample,1); % Output
% FDN Computation
for sample = 1:nSample % each sample
for n = 1:N % each delayline
in_dl(n) = in(sample,1) + FB(n); % Input and Feedback matrix sum
% Delaying
[out_dl(n),bufferN(:,n)] = funcDelay( in_dl(n), bufferN(:,n),...
sample, delays(n) );
% Filtering
[out_fdl(n),buffFin(:,:,n),buffFout(:,:,n)] = funcFilterGEQ(...
GEQ(:,:,n), out_dl(n), buffFin(:,:,n), buffFout(:,:,n));
end
out(sample,1) = sum(out_fdl); % Filtered delay output sum
FB = out_fdl * MTX; % Feedback matrix output recalculation
end
% Used functions
function [out,buffer] = funcDelay( in, buffer, n, delay)
% Circular buffer indices
len = length(buffer);
indexC = mod(n-1, len) + 1; % Current
indexD = mod(n-delay-1, len) + 1; % Delay
out = buffer(indexD,1);
% Stores output on appropriate index
buffer(indexC,1) = in;
end
function [outfilt,buffin,buffout] = funcFilterGEQ( GEQ, in, buffin, buffout)
% Sample based filter cascading
nBands = size(GEQ,1);
out = zeros(1,nBands+1);
out(1) = in; % Stores value pre-filtered
for b = 1:nBands % Performs series bandpass filtering
[out(b+1),buffin(b,:),buffout(b,:)] = funcBandpassFilt( out(b),...
GEQ(b,1:3), GEQ(b,5:6), buffin(b,:), buffout(b,:) );
end
outfilt = out(end); % Outputs final value post-filtered
end
function [out,buffin,buffout] = funcBandpassFilt( in, bs, as, buffin, buffout)
% Sample based filtering
buffin(3) = buffin(2); % Valid for biquad filters
buffin(2) = buffin(1);
buffin(1) = in; % Sequential indexing
buffout(1) = bs(1)*buffin(1) + bs(2)*buffin(2) + bs(3)*buffin(3)...
- as(1)*buffout(2) - as(2)*buffout(3);
out = buffout(1);
% Outputs calculation based on 3 latest values
buffout(3) = buffout(2);
buffout(2) = buffout(1);
buffout(1) = 0;
end