function ivecs = fast_hybrid_ivec(T, ubm_wts, zo_stats, fo_stats)
% FAST_HYBRID_IVEC = Fast method to  estimate i-vectors (approximately)
%   ivecs = quick(T, ubm_wts, zo_stats, fo_stats) returns a matrix
%   of i-vectors given the matrix T, weights of mixtures of the ubm
%   (ubm_wts), zeroth order statistics (zo_stats) and first order 
%   statistics (fo_stats)
%
%   Input:
%       T           :   MD x R matrix
%       ubm_wts     :   M x 1 column vector of mixture weights
%                       that sum up to 1
%       zo_stats    :   M x n matrix with each column representing the
%                       zeroth order statistics of the recording for 
%                       which i-vector is being estimated
%       fo_stats    :   MD x n matrix of first order statistics. Each
%                       column represents the first order stats corresponding
%                       to the column in zo_stats
%
%   Output:
%       ivecs       :   R x n matrix where each column is an i-vector
%
%
%  
% Copyright © 2012 Indian Institute of Technology Madras. All Rights Reserved.
% 
% Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:

% 1. Redistributions of source code must retain the above copyright notice,
% this list of conditions and the following disclaimer.
%
% 2. Redistributions in binary form must reproduce the above copyright notice,
% this list of conditions and the following disclaimer in the documentation
% and/or other materials provided with the distribution.
%
% 3. The name of the author may not be used to endorse or promote products
% derived from this software without specific prior written permission.
%
% THIS SOFTWARE IS PROVIDED BY Indian Institute of Technology Madras "AS IS"
% AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
% IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
% DISCLAIMED. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
% INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
% LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
% PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
% LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE
% OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
% ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 


    % get dimensions
    M = max(size(ubm_wts));
    [MD R] = size(T);
    D = MD / M;

    % some sanity check
    if MD < R
        warning('supervector is shorted than i-vector');
    elseif M > MD
        warning('dimension mismatch - more mixtures than expected');
    elseif min(size(ubm_wts)) > 1
        warning('ubm weights is expected to be a column/row vector');
    elseif size(fo_stats,1) ~= MD
        warning(['dimension mismatch - fostats should have the same  ' ...
                ,'no. of rows as T'])
    elseif size(fo_stats, 2) ~= size(zo_stats, 2)
        warning (['dimension mismatch - fo_stats and zo_stats have ' ...
                 ,'different no. of columns']);   
    elseif size(zo_stats, 1) ~= M
        warning('zo_stats does not have M rows');
    end   

    % pre-compute matrices
    L = zeros(R, R);
    l = ([1:M]-1)*D+1;
    r = [1:M]*D;
    for m = 1:M
        L = L + ubm_wts(m) * T(l(m):r(m),:)' * T(l(m):r(m),:);
    end
    N = sum(zo_stats);
    size(N)

    [U S V] = svd(L);
    premat = U' * T';
    S = diag(S);
    unz = ones(R, 1);
    n = size(fo_stats,2);
    ivecs = zeros(R, n);

    % compute i-vectors
    for i = 1:n
        ivecs(:, i) = bsxfun(@rdivide, U, (N(i) * S + unz)) * (premat * fo_stats(:, i));        
    end

%   Written by: Srikanth Madikeri
%   Last revision: 01-03-2014
