What's wrong with my STFT implementation?

Question

I'm studying the linear convolution for achieve more complex results in Audio FX.

My main goal is to implement a realtime convolution between two stream of sound in live perfomance (for example Sound 1 from a synthesizer and Sound2 from a guitar o microphone) .

The length of sounds are not determinate and for testing i recorded two part of these streams. You can think the file used in Matlab as a part of infinite streams.

I have two audio files (x1 and x2) with same fixed length (len1 and len2 = 673792 samples) and from DSP theory i know that the length of convolution array is lenOut = (len1+len2-1).

You can download the input audio files at these links

http://www.sun-art.org/download/matlab/test_01.aif

http://www.sun-art.org/download/matlab/test_02.aif

I wrote two OCTAVE procedures

1. with direct convolution using conv(n,m) function
2. with FFT product of input signals with zero padding

Both procedures gives the same results. You can download the output audio files (30 seconds long each) at these links:

http://www.sun-art.org/download/matlab/conv.wav

http://www.sun-art.org/download/matlab/convFFT.wav

The next step is to make a realtime convolution using STFT but before going in C++ i would like to improve it and test it in Matlab/Octave. (my reference are: "SPECTRAL AUDIO SIGNAL PROCESSING" by Julius O. Smith III and "DAFX: Digital Audio Effects" by Udo Zölzer).

I noticed that the result of STFT procedure is quite different from direnct convolution/FFT product procedure. You can download the output file at this link http://www.sun-art.org/download/matlab/STFTAudioConv.wav

The file is 30 seconds long, but from 15 seconds onwards the audio is zero. It'possible to achieve the result of function Matlab conv() with STFT procedure ? What's wrong with my STFT ?

Below you will find the three procedures.

1) LINEAR CONVOLUTION

    [x1, FS] = audioread('conv_test\test_01.aif');
    [x2, FS] = audioread('conv_test\test_02.aif');

    y        = conv(x1,x2);
    outFile  = y / max(abs(y));

    audiowrite ('conv_test\conv.wav', outFile, FS );  %Octave

2) FFT PRODUCT

    % ============= LOAD audio files x1 and x2==============
    [x1, FS] = audioread('conv_test\test_01.aif');
    [x2, FS] = audioread('conv_test\test_02.aif');

    % x1 and x2 are files with the same length / length(x1) = length(x2) = 673792 samples
    len1   = length(x1);
    len2   = length(x2);  
    lenOut = len1 + len2 - 1;

    N      = 2^nextpow2(lenOut);   %for FFT and padding

    % ============= PADDING _ & scaling =============
    inx1   = [x1;  zeros(N-len1,1)] / max( abs(x1) );
    inx2   = [x2;  zeros(N-len2,1)] / max( abs(x2) );

    %============= FFT of signal x1 =============

    f1     = fft(inx1); 
    r1     = abs(f1);
    phi1   = angle(f1);
    fft1   = (r1.* exp(i*phi1));

    %============= FFT of signal x2 =============

    f2     = fft(inx2);
    r2     = abs(f2);
    phi2   = angle(f2);
    fft2   = (r2.* exp(i*phi2));

    %============= FFT PRODUCT  and iFFT =============

    Yfft   = fft1 .* fft2;
    yout   = real( ifft(Yfft) );

    %============= Write and Save iFFT (linear convolution) =============

    outFile = yout(1:lenOut) / max( abs(yout) );
    audiowrite ('conv_test\convFFT.wav', outFile, FS);  %Octave

3) STFT

For Overlapp-Add i define an Hanning window with L = 1024 samples, an Hopsize of L/2 = 512 samples and an FFT size N = 2^nextpow2(L+1) = 2048 samples. In the while loop

1. take 1024 (L) samples of inputs with windowing
2. zero padding with (N-L) = 1024 zeros. The buffer length in the FFT is now 2048 samples long.
3. FFT of zero padded inputs
4. FFT product
5. iFFT
6. overlapp add the output
7. HopSize advance

Fixing the windows lenght L = 1024 implicitly I also fixed the signal (x1) and filter (x2) size L and FFT size = 2*L .

My inpulse response (x2) as the other signal is varying-time. My processor has two input buffer (left and right channel) for example 16 samples long. I collect L = 1024 sample in new two auxilaty buffers (in1 and in2 in Script 3) and after windowing i extend them with zero-pad until the FFT size ( = 2*L)

    % ============= LOAD audio files x1 and x2==============
    [x1, FS] = audioread('conv_test\test_01.aif');
    [x2, FS] = audioread('conv_test\test_02.aif');

    % x1 and x2 are files with the same length / length(x1) = length(x2) = 673792 samples
    len1     = length(x1);
    len2     = length(x2);  
    lenOut   = len1 + len2 - 1;

    L         = 1024;               % length of analisys windows
    w1        = hanningz(L);        % w = .5*(1 - cos(2*pi*(0:L-1)'/(L)));
    hop       = L/2;
    N         = 2^nextpow2(L+1);    % FFT length

    out       = zeros(lenOut,1);    %init of output
    start     = 0;                  %index of while loop
    pout      = 0;
    pend      = (length(x1) - L);   %end of while loop

    % ================== OVERLAPP ADD ==================

    while start<pend

      in1    = x1(start+1:start+L) .* w1;       %windowing
      in2    = x2(start+1:start+L) .* w1;

      max1   = max(abs(in1));
      max2   = max(abs(in2));

      buf1   = [in1; zeros(N-L,1)] / max1;  % zero padding
      buf2   = [in2; zeros(N-L,1)] / max2;

      f1     = fft(buf1);                   % FFT input 1
      r1     = abs(f1);                     
      phi1   = angle(f1);
      ft1    = (r1.* exp(i*phi1));

      f2     = fft(buf2);                   % FFT input 2
      r2     = abs(f2);
      phi2   = angle(f2);
      ft2    = (r2.* exp(i*phi2));  

      YConv  = ft1 .* ft2;                  % Convolution - FFT Product

      y      = real( ifft(YConv) );         % iFFT

      out(pout+1:pout+N) = out(pout+1:pout+N) + y;  % Overlapp-Ddd

      start  = start + hop;
      pout   = pout + hop;

    end

    % ================== Write and Save ==================

    out = out / max(abs(out));        
    audiowrite ('conv_test\STFTAudioConv.wav', out, FS ); %Octave

Answer 1

A few points

1. It is very unusual to convolve two audio files. You typically convolve one audio file with the impulse response of a filter. Thinking about "long file" and "short-ish" impulse response makes it easier to understand
2. For overlap add the FFT size is determined by the length of the filter. For example, if your file length is 670k samples and your impulse response is 16k samples, you would chose an FFT size of twice the short one, i.e. 32k. You multiply the FFT of each input frame with the FFT of the entire filter.
3. You don't window. Or to be precise: you use a rectangular window.
4. If your impulse response is too long for a single FFT (for latency, memory or CPU reasonse), you can chop this up as well, but that's more complicated and you should not do this until you have Overall Add figured out.
5. Overlap add is NOT the same as STFT. (see jithin's comment)
6. Your code has a lot of unnecessary stuff in it. You extract amplitude and phase from the FFT of the block and put them together again. Why? $ft1 = f1$
7. You normalize the input to the maxium of a windowed frame. That's just plain wrong

Update based on added info

My main goal is to implement a realime convolution between two stream of sound in live perfomance(for example sound 1 from a synthesizer and sound2 from a guitar o microphone) So the length of sounds are non determinate. For testing i recorded two part of these streams. You can think the file used in Matlab as a part of infinite streams.

Sorry, that won't work. You can't convolve two infinite sequences with a finite amount of latency. In essence one of the sequences must be time reversed, and you can't reverse an infinite sequence without waiting infinitely long.