function y = wavmix(x1, x2, a1, a2, delay, Fs)
    // Flexible mix of two audio signals
    // This function mixes two audio signals, either 
    // of them mono or stereo, column or row, even
    // if their lengths are different and the second
    // one has a delay.
    //
    // Usage:
    //       y = wavmix(x1, x2, a1, a2, delay, Fs)
    //
    // where
    //       x1:    first signal
    //       x2:    second signal
    //       a1:    first signal gain factor   
    //       a2:    second signal gain factor   
    //       delay: time delay of x2 respect to x1 in s
    //       Fs:    sampling rate in Hz
    //       y:     mixer output
    //
    // Arguments Fs, delay, a2 and a1 are optional. Their
    // default values are 44100, 0, 1 and 1. x1 and x2 may
    // be any combination of row, column, mono or stereo 
    // signals. However, the output signal is row if x1
    // is row and stereo if either x1 or x2 are stereo,
    // according to the following table:
    // 
    // x1                 x2                   y
    // -------------------------------------------------------
    // mono    row        mono    row          mono    row
    // mono    row        mono    column       mono    row
    // mono    row        stereo  row          stereo  row
    // mono    row        stereo  column       stereo  row
    // mono    column     mono    row          mono    column
    // mono    column     mono    column       mono    column
    // mono    column     stereo  row          stereo  column 
    // mono    column     stereo  column       stereo  column
    // stereo  row        mono    row          stereo  row
    // stereo  row        mono    column       stereo  row
    // stereo  row        stereo  row          stereo  row
    // stereo  row        stereo  column       stereo  row
    // stereo  column     mono    row          stereo  column
    // stereo  column     mono    column       stereo  column
    // stereo  column     stereo  row          stereo  column 
    // stereo  column     stereo  column       stereo  column
    //
    // In other words, stereo dominates over mono and the
    // first signal orientation dominates over the second
    // signal.
    // 
    // a1 and a2 may be scalars, two-component vectors, or
    // arrays of the same length as x1 and x2. In the first
    // case, the corresponding whole signal (either mono or
    // stereo) is scaled by the scalar. In the second case,
    // a1(1) affects the first channel of x1 and a1(2)
    // affects the second channel of x1 and similarly for 
    // x2. If x1 is mono (single channel), a1(1) and a1(2) 
    // act as panning factors (and similarly for x2). The
    // third case allows an amplitude modulation of the
    // signal. If a1 is a vector, it is multipled component
    // by component by the signal's samples. If it is a
    // two-column or two-row array, it acts on a component
    // by component basis as a time-variable panning factor
    // or as time-variable channel gains according to the 
    // dimension of the signal. 
    //
    // NOTE 1: If either x1, x2, a1 or a2 are stereo, the 
    //         result is also stereo. 
    // NOTE 2: If both x1 and a1 are mono but y is stereo,
    //         a1.*x1 is fed to both channels yielding
    //         pseudo stereo to comply with the stereo
    //         format. No attempt is done to delay, pan
    //         or decorrelate both channels, which is 
    //         responsibility of the user. 
    //
    // The function performs the following action:
    //
    // y(t) = a1 * x1(t) + a2 * x2(t-delay) 
    // 
    // or
    // 
    // y(t) = a1(t) * x1(t) + a2(t-delay) * x2(t-delay)
    // 
    // if either a1, a2 or both are variable. 
    //
    // NOTE 3: x2 is asumed to be 0 for t < 0.
    //
    // ----------------------------
    // Author: Federico Miyara 
    // Date:   2010-07-06
    //         2010-07-14
    //         2021-06-14
    
    rhs = argn(2);

    // Default values if not enough arguments
    if rhs<6
        Fs = 44100;
    end
    if rhs<5
        delay = 0;
    end
    if rhs<4
        a2 = 1;
    end
    if rhs<3
        a1 = 1;
    end

    // -----------
    // BUREAUCRACY
    // -----------
    //
    // (would be unnecessary if only everyone used only column signals...)

    // Size of signals and coefficients
    Sx1 = size(x1);
    Sx2 = size(x2);
    Sa1 = size(a1);
    Sa2 = size(a2);

    // Delay expressed in samples
    offset = round(delay*Fs);

    // Classification of signal types
    if Sx1(1)==1
        // x1 is mono row 
        L1 = Sx1(2);
        x1_wasrow = 1;
        x1_mono = 1;
        // Convert to mono column
        x1 = x1.';
    elseif Sx1(2)==1
        // x1 is mono column 
        L1 = Sx1(1);
        x1_wasrow = 0;
        x1_mono = 1;
    end

    if Sx1(1)==2
        // x1 is stereo row 
        L1 = Sx1(2);
        x1_wasrow = 1;
        x1_mono = 0;
        // Convert to stereo column
        x1 = x1.';
    elseif Sx1(2)==2
        // x1 is stereo column 
        L1 = Sx1(1);
        x1_wasrow = 0;
        x1_mono = 0;
    end

    if Sx2(1)==1
        // x2 is mono row 
        L2 = Sx2(2);
        x2_wasrow = 1;
        x2_mono = 1;
        // Convert to mono column
        x2 = x2.';
    elseif Sx2(2)==1
        // x2 is mono column 
        L2 = Sx2(1);
        x2_wasrow = 0;
        x2_mono = 1;
    end

    if Sx2(1)==2
        // x2 is stereo row 
        L2 = Sx2(2);
        x2_wasrow = 1;
        x2_mono = 0;
        // Convert to stereo column
        x2 = x2.';
    elseif Sx2(2)==2
        // x2 is stereo column 
        L2 = Sx2(1);
        x2_wasrow = 0;
        x2_mono = 0;
    end

    // Recalculate Size of signals and coefficients
    Sx1 = size(x1);
    Sx2 = size(x2);

    // Now we have all input signals as columns. Next we must deal
    // with the coefficients (scalars, panning or modulating signals).

    select prod(Sa1)
    case 1
        // If a1 is scalar then use it as a general factor for x1
        x1 = a1 * x1;
    case 2
        // If a1 is a pair of scalars, one per channel, then use
        // it as a panning pair or a gain pair
        if Sx1(2)==1
            // If x1 is a mono signal, then it must be panned according
            // to a1 
            x1 = [a1(1)*x1, a1(2)*x1];
        else
            // If x1 is a stereo signal, then both components of a1 act 
            // as different weighting factors for each channel  
            x1 = [a1(1)*x1(:,1), a1(2)*x1(:,2)];
        end
    else
        // a1 is a vector
        if (Sa1(1)==1) | (Sa1(2)==1)
            if Sa1(1)==1
                // If a1 is a row vector, convert to column, otherwise
                // do nothing since it is already a column
                a1 = a1.';
            end
            // Multiply by x1
            if x1_mono==1
                // If x1 is also mono, then just multiply a1 by x1 element
                // by element
                x1 = a1 .* x1;
            else
                // If x1 is stereo, first clone (¿kron?) a1 into a second 
                // column and then multiply element by element
                x1 = kron(a1, [1,1]) .* x1;
            end
        end
        // a1 is a two-vector array
        if (Sa1(1)==2) | (Sa1(2)==2)
            if Sa1(1)==2
                // If a1 is a two-row array, transpose to have a two-column
                // array, otherwise do nothing since it has already two
                // columns
                a1 = a1.';
            end
            // Multiply by x1
            if x1_mono==1
                // If x1 is mono, clone it to be stereo and then multiply 
                // by a1 element by element
                x1 = a1 .* kron(x1, [1 1]);
            else
                // If x1 is already stereo, then just multiply by a1 element
                // by element
                x1 = a1 .* x1;
            end
        end
    end

    select prod(Sa2)
    case 1
        // If a2 is a scalar use it as a general factor for x2
        x2 = a2 * x2;
    case 2
        // If a2 is a pair of scalars, one per channel, then use
        // it as a panning pair or a gain pair
        if Sx2(2)==1
            // If x2 is a mono signal, then it must be panned according
            // to a2 
            x2 = [a2(1)*x2, a2(2)*x2];
        else
            // If x1 is a stereo signal, then both components of a1 act 
            // as different weighting factors for each channel  
            x2 = [a2(1)*x2(:,1), a1(2)*x2(:,2)];
        end
    else
        // a2 is a vector
        if (Sa2(1)==1) | (Sa2(2)==1)
            if Sa2(1)==1
                // If a2 is a row vector, convert to column, otherwise
                // do nothing since it is already a column
                a2 = a2.';
            end
            // Multiply by x2
            if x2_mono==1
                // If x2 is also mono, then just multiply a2 by x2 element
                // by element
                x2 = a2 .* x2;
            else
                // If x2 is stereo, first clone (¿kron?) a2 into a second 
                // column and then multiply element by element
                x2 = kron(a2, [1,1]) .* x2;
            end
        end
        // a2 is a two-vector array
        if (Sa2(1)==2) | (Sa2(2)==2)
            if Sa2(1)==2
                // If a2 is a two-row array, transpose to have a two-column
                // array, otherwise do nothing since it has already two
                // columns
                a2 = a2.';
            end
            // Multiply by x2
            if x2_mono==1
                // If x2 is mono, clone it to be stereo and then multiply 
                // by a2
                x2 = a2 .* kron(x2, [1 1]);
            else
                // If x2 is already stereo, then just multiplyby a2 element 
                // by element
                x2 = a2 .* x2;
            end
        end
    end


    // --------------------------------------
    // MIXING ALGORITHM
    // --------------------------------------

    // Length of combined signal
    L = max(L1, offset + L2);

    // Initialisation of y
    if (size(x1,2)==1) & (size(x2,2)==1)
        // If both x1 and x2 are mono then y is mono
        y = zeros(L, 1);
    else
        // ... otherwise, y is stereo
        y = zeros(L, 2);
    end 

    // Convert x1 to stereo if necessary
    if size(x1, 2)==1 & size(y, 2)==2
        x1 = kron(x1, [1,1]);
    end

    // Convert x2 to stereo if necessary
    if size(x2, 2)==1 & size(y, 2)==2
        x2 = kron(x2, [1,1]);
    end

    // Mix of signals x1 and x2
    y(1:L1, :) = x1;
    y(offset+1:offset+L2, :) = y(offset+1:offset+L2, :) + x2;


    // If x1 was row mono or row stereo, transpose y 
    if x1_wasrow==1
        y = y.';
    end

endfunction
