function y = event_mix(x, a, delay)
    // Mixes in a single stereo signal a number 
    // of signals (events), each one of arbitrary
    // length and arbitrary initial sample, with
    // corresponding gains and delays. 
    //
    // Usage:
    //       y = event_mix(x, a, delay)
    //
    // where
    //       x:     cell array of row stereo signals
    //       a:     cell array of corresponding gains 
    //       delay: cell array of corresponding delays 
    //              expressed as a number of samples, 
    //              in ascending order 
    //       y:     output signal
    //
    // Notes: 1) All the cell arrays must have the
    //           same number of cells
    //        2) The cell array may be 2-dimensional 
    //           (and possibly N-dimensional) allowing 
    //           several layers (each one corresponds
    //           to a voice in a musical composition
    //           context or to a track in an audio 
    //           mixing context). Since not necessarily
    //           each voice or track has the same
    //           number of events, it is allowed that 
    //           some cells be empty. The overlapping
    //           nature of this algorithm would allow
    //           a one dimensional distribution of the
    //           events, since the delay of several 
    //           events could be the same. However, it 
    //           is better if the structure of the data
    //           represents the structure of the 
    //           composition.
    //        3) The cells of a may be scalars, two-row
    //           column vectors allowing different gains
    //           for each stereo channel or arrays of the 
    //           same size as the corresponding signal
    //           in cell array x, allowing amplitude
    //           modulation.
    //        4) Signals may also be monophonic (all of
    //           them)
    //        5) In this first version, no check is done
    //           of argument structure, it is the user's 
    //           own responsability to provide arguments
    //           with the required structure.
    //
    // Se also: wavmix()
    // 
    // -------------------------------------
    //
    // Author: Federico Miyara
    // Date:   2015-05-31
    //         2021-06-21
    
    // For testing purposes only  --DELETE--
    if 1==2
        x = {[1 2 3 4; 5 6 7 8], [9 10 11; 12 13 14], [15 16 ; 17 18]; ...
             [1 2 3; 4 5 6], [7 8 9; 10 11 12], []}
        a = {[2; 2], [1; 1], [3; 2]; ...
             [1; 2], [2; 1], []}
        delay = {0, 3, 7; ...
                 1, 2, []}     
    end
    
    // size of signal container
    Ny = size(x);

    // Detect whether there is at least 
    // one stereo signal
    channels = 1;
    for k=1:prod(Ny)
        if size(x{k})(1)==2
            channels = 2;
            break;
        end
    end
    
    // Adapt size of gain factors to the
    // sizes of the signals
    for k=1:prod(Ny)
        // If a{k} has two components it
        // means that is a stereo gain
        // (for instance the result of
        // panning). In that case it cannot
        // be just multiplied by x{k}. It 
        // is first necessary to replicate
        // to match the size of x{k} 
        if and(size(a{k})==[2,1])
            a{k} = repmat(a{k},[1,size(x{k})(2)]);
        end
    end
    
    // Initialize cell of last sample indexes
    last = cell(size(x));
    // Initialize maximum last sample index 
    maxlast = 1;
    // Compute the last sample index for each 
    // non-empty signal after applying the delay
    for n=1:Ny(1)
        for m=1:Ny(2)
            if ~isempty(x{n,m})
                // If there is no delay specified for
                // a non-empty signal set it equal to 0
                if isempty(delay{n,m})
                    delay{n,m} = 0;
                end
                // Get last sample index
                last{n,m} = delay{n,m} + size(x{n,m})(2);
                // Update the maximum last sample index
                if last{n,m}>maxlast
                    maxlast = last{n,m};
                end
            end
        end
    end

    // Initialize output signal. Length is computed
    // so that there is enough room for the longest 
    // combination of signal and delay
    y = zeros(channels, maxlast);     
                
    for k=1:prod(Ny)
        // The k-th signal is added at the correct
        // location with the intended gain
        if ~isempty(x{k})
            y(:, delay{k}+1:last{k}) = ...
                 y(:, delay{k}+1:last{k}) + a{k}.*x{k};
        end
    end

endfunction

