Real-time implementation of cascaded all-pass filters from given transfer function

I am working on real-time implementation of spring reverb based on scientific paper called the Parametric Spring Reverberation Effect by Välimäki, Vesa; Parker, Julian; Abel, Jonathan S. One block of this effect is a cascade of 100 stretched all-pass filters. Single stretched all-pass filter is defined by:

I calculated the difference equation of this filter and used it for implementing the cascade of M = 100 stretched all-pass filters in Matlab.

for i = 1:M
    for n = 1+K1+1:(length(x))
        y(n) = a1 * x(n) + a1*a2 * x(n-1) + a2 * x(n-K1) + ...
            x(n-1-K1) - a2 * y(n-1) - a1*a2 * y(n-K1) - a1 * y(n-1-K1);         
    end
    x = y;
end

However, the problem seems to be that I need to add another filters to the structure which are not cascaded (as you can see in the following block diagram), therefore they cannot be "repeated" like the all-pass filter.

How can I calculate the output of my all-pass filter cascade sample-by-sample?

I can imagine the multiplication of two or three transfer functions in order to obtain the transfer function of the whole cascade, but multiplying 100 transfer functions seems to be a crazy and wrong idea. There is probably a simple solution, but I just couldn't figure anything out. That's why I am asking you for help.

Answer 1

[Darshan Ramakant Bhat]

I think you can use the cascade function available in the signal processing toolbox. Please try to use this function in a loop to cascade the filter 100 times:

https://www.mathworks.com/help/signal/ref/dfilt.cascade.html