function varargout = datarequest(request,m,data,range,idata,loglikopt)
% DATAREQUEST  [Not a public function] Request data from database or datapack.
%
% Backend IRIS function.
% No help provided.

% -IRIS Toolbox.
% -Copyright (c) 2007-2012 Jaromir Benes.

if ~exist('idata','var')
    idata = ':';
end

if ~exist('loglikopt','var')
    loglikopt = [];
end

%**************************************************************************

nx = size(m.solution{1},1);
nb = size(m.solution{7},1);
nf = nx - nb;
nalt = size(m.Assign,3);

if isempty(data)
    data = struct();
end

dpackmean = [];
dpackmse = [];
dbasemean = [];
dbasemse = [];

if iscell(data)
    % Plain datapack.
    dpackmean = data;
elseif isstruct(data) ...
        && (isfield(data,'mean')  || isfield(data,'mean_'))
    if isfield(data,'mean') && isstruct(data.mean)
        dbasemean = data.mean;
        if isfield(data,'mse') && isa(data.mse,'tseries')
            dbasemse = data.mse;
        end
    end
    if isempty(dbasemean)
        if isfield(data,'mean_') && iscell(data.mean_)
            dpackmean = data.mean_;
        end
    end
    if isfield(data,'mse_') && iscell(data.mse_)
        dpackmse = data.mse_;
    end
elseif isstruct(data)
    % Plain database.
    dbasemean = data;
else
    error('iris:model','Unknown type of input data.');
end

switch request
    case 'init'
        if nargout < 4
            [xinitmean,naninit] = dodata2xinit();
            xinitmse = [];
            ainitmean = doxinit2ainit();
            varargout{1} = ainitmean;
            varargout{2} = xinitmean;
            varargout{3} = naninit;
        else
            [xinitmean,naninit,xinitmse] = dodata2xinit();
            [ainitmean,ainitmse] = doxinit2ainit();
            varargout{1} = ainitmean;
            varargout{2} = xinitmean;
            varargout{3} = naninit;
            varargout{4} = ainitmse;
            varargout{5} = xinitmse;
        end
    case 'xinit'
        [varargout{1:nargout}] = dodata2xinit();
    case {'y','kalman'}
        y = dodata2y();
        y = y(:,:,idata);
        if ~isempty(loglikopt) && isstruct(loglikopt) ...
                && isfield(loglikopt,'domain') ...
                && strncmpi(loglikopt.domain,'f',1)
            y = fft(y.').';
        end
        varargout{1} = y;
    case 'e'
        varargout{1} = dodata2e();
    case 'x'
        varargout{1} = dodata2x();
    case 'y,x,e'
        data = {dodata2y(),dodata2x(),dodata2e()};
        ndata = max([size(data{1},3),size(data{2},3),size(data{3},3)]);
        % Make the size of all data arrays equal in 3rd dimension.
        if size(data{1},3) < ndata
            data{1} = cat(3,data{1}, ...
                data{1}(:,:,end*ones(1,ndata-size(data{1},3))));
        end
        if size(data{2},3) < ndata
            data{2} = cat(3,data{2}, ...
                data{2}(:,:,end*ones(1,ndata-size(data{2},3))));
        end
        if size(data{3},3) < ndata
            data{3} = cat(3,data{3}, ...
                data{3}(:,:,end*ones(1,ndata-size(data{3},3))));
        end
        varargout = data;
    case 'alpha'
        varargout{1} = dodata2alpha();
end

% Nested functions.

%**************************************************************************
    function [xinitmean,naninit,xinitmse] = dodata2xinit()
        xinitmean = nan(nb,1,nalt);
        xinitmse = [];
        if ~isempty(dpackmean)
            xinitmean = xxresize(dpackmean{2}(nf+1:end,:,:), ...
                dpackmean{4},range(1)-1,'mean');
        elseif ~isempty(dbasemean)
            realid = real(m.solutionid{2}(nf+1:end));
            imagid = imag(m.solutionid{2}(nf+1:end));
            xinitmean = db2array(dbasemean,m.name(realid),range(1)-1, ...
                imagid,m.log(realid),false,true,true,true);
            xinitmean = permute(xinitmean,[2,1,3]);
        end
        if nargout >= 3
            if ~isempty(dbasemse)
                xinitmse = dbasemse(range(1)-1);
                xinitmse = ipermute(xinitmse,[3,2,1,4]);
            elseif ~isempty(dpackmse)
                xinitmse = xxresize(dpackmse{2}(nf+1:end,nf+1:end,:,:), ...
                    dpackmse{4},range(1)-1,'mse');
            end
        end
        realid = real(m.solutionid{2}(nf+1:end));
        ndata = size(xinitmean,3);
        naninit = false(nb,1);
        for iloop = 1 : ndata
            if iloop <= nalt
                required = m.icondix(1,:,iloop);
                required = required(:);
            end
            naninit = naninit | (isnan(xinitmean(:,1,iloop)) & required);
        end
        naninit = m.name(unique(realid(naninit)));
    end
% dodata2xinit().

%***********************************************************************
% Get initial conditions for xb and alpha.
% Those that are not required are set to `NaN` in xinitmean, and
% to 0 when computing ainitmean.
    function [ainitmean,ainitmse] = doxinit2ainit()
        % Transform mean x to alpha.
        ndata = size(xinitmean,3);
        if ndata < nalt
            xinitmean(:,1,end+1:nalt) = ...
                xinitmean(:,1,end*ones(1,nalt-ndata));
            ndata = nalt;
        end
        ainitmean = xinitmean;
        for iloop = 1 : ndata
            if iloop <= nalt
                U = m.solution{7}(:,:,iloop);
                notrequired = ~m.icondix(1,:,iloop);
            end
            index = isnan(xinitmean(:,1,iloop)) & notrequired(:);
            ainitmean(index,1,iloop) = 0;
            ainitmean(:,1,iloop) = U\ainitmean(:,1,iloop);
        end
        % Transform MSE x to alpha.
        if nargout < 2 || isempty(xinitmse)
            ainitmse = xinitmse;
            return
        end
        ndata = size(xinitmse,4);
        if ndata < nalt
            xinitmse(:,:,1,end+1:nalt) = ...
                xinitmse(:,:,1,end*ones(1,nalt-ndata));
            ndata = nalt;
        end
        ainitmse = xinitmse;
        for iloop = 1 : ndata
            if iloop <= nalt
                U = m.solution{7}(:,:,iloop);
                Ut = U.';
                %required = m.icondix(1,:,iloop);
            end
            %xinitmse(~required,:,1,iloop) = NaN;
            %xinitmse(:,~required,1,iloop) = NaN;
            %ainitmse(~required,:,1,iloop) = 0;
            %ainitmse(:,~required,1,iloop) = 0;
            ainitmse(:,:,1,iloop) = U\ainitmse(:,:,1,iloop);
            ainitmse(:,:,1,iloop) = ainitmse(:,:,1,iloop)/Ut;
        end
    end
% doxinit2ainit().

%***********************************************************************
% Get observables.
    function y = dodata2y()
        if ~isempty(dbasemean)
            realid = real(m.solutionid{1});
            imagid = imag(m.solutionid{1});
            tmplog = m.log(realid);
            y = db2array(dbasemean,m.name(realid),range, ...
                imagid,tmplog,false,true,true,true);
            y = permute(y,[2,1,3]);
        else
            y = xxresize(dpackmean{1},dpackmean{4},range,'mean');
        end
    end
% dodata2y().

%***********************************************************************
% Get residuals. Set `NaN` residuals to zero.
    function e = dodata2e()
        if ~isempty(dbasemean)
            realid = real(m.solutionid{3});
            imagid = imag(m.solutionid{3});
            e = db2array(dbasemean,m.name(realid),range, ...
                imagid,m.log(realid),false,true,true,true);
            e = permute(e,[2,1,3]);
        else
            e = xxresize(dpackmean{3},dpackmean{4},range,'mean');
        end
        ereal = real(e);
        eimag = imag(e);
        ereal(isnan(ereal)) = 0;
        eimag(isnan(eimag)) = 0;
        e = ereal + 1i*eimag;
    end
% dodata2e().

%***********************************************************************
% Get current dates of transition variables.
% Set lags and leads to NaN.
    function x = dodata2x()
        realid = real(m.solutionid{2});
        imagid = imag(m.solutionid{2});
        index = imagid == 0;
        if ~isempty(dbasemean)
            realid = realid(index);
            imagid = imagid(index);
            tmplog = m.log(realid);
            tmpdata = db2array(dbasemean,m.name(realid),range, ...
                imagid,tmplog,false,true,true,true);
            tmpdata = permute(tmpdata,[2,1,3]);
            x = nan(length(index),size(tmpdata,2),size(tmpdata,3));
            x(index,:,:) = tmpdata;
        else
            x = xxresize(dpackmean{2},dpackmean{4},range,'mean');
            % Set lags and leads to NaN.
            x(~index,:,:) = NaN;
        end
    end
% dodata2x().

%***********************************************************************
    function a = dodata2alpha()
        if ~isempty(dbasemean)
            realid = real(m.solutionid{2});
            imagid = imag(m.solutionid{2});
            realid = realid(nf+1:end);
            imagid = imagid(nf+1:end);
            a = db2array(dbasemean,m.name(realid),range, ...
                imagid,m.log(realid),false,true,true,true);
            a = permute(a,[2,1,3]);
        else
            a = xxresize(dpackmean{2}(nf+1:end,:,:),dpackmean{4},range,'mean');
        end
        ndata = size(a,3);
        if ndata < nalt
            a(:,:,end+1:nalt) = a(:,:,end*ones(1,nalt-ndata));
            ndata = nalt;
        end
        for iloop = 1 : ndata
            if iloop <= nalt
                U = m.solution{7}(:,:,iloop);
            end
            a(:,:,iloop) = U\a(:,:,iloop);
        end
    end
% dodata2alpha().

end

%**************************************************************************
function x = xxresize(x,oldrange,newrange,type)
ismean = strcmp(type,'mean');
% Same input range and output range.
if round(oldrange(1) - newrange(1)) == 0 ...
        && round(oldrange(end) - newrange(end)) == 0
    return
end

% Resolve [-Inf,date] or [date,Inf].
if length(newrange) == 2
    if isinf(newrange(1))
        newrange(1) = oldrange(1);
    end
    if isinf(newrange(2))
        newrange(2) = oldrange(end);
    end
end
newrange = newrange(1) : newrange(end);

noldrange = length(oldrange);
index = round(newrange - oldrange(1)) + 1;
npre = sum(index < 1);
npost = sum(index > noldrange);
index(index < 1 | index > noldrange) = [];

if ismean
    tmpsize = [ ...
        size(x,1), ...
        size(x,2), ...
        size(x,3), ...
        ];
    x = [ ...
        nan(tmpsize(1),npre,tmpsize(3)), ...
        x(:,index,:), ...
        nan(tmpsize(1),npost,tmpsize(3)) ];
else
    tmpsize = [ ...
        size(x,1), ...
        size(x,2), ...
        size(x,3), ...
        size(x,4), ...
        ];
    x = cat(3, ...
        nan(tmpsize(1:2),npre,tmpsize(4)), ...
        x(:,:,index,:), ...
        nan(tmpsize(1:2),npost,tmpsize(4)) );
end
end
% xxresize().