function modelpred = doPLSDApred(Xt,Yt,modpls)
    % Preparation of the class matrices and vectors
    if ~isempty(Yt)
        if size(Yt,2) == 1
            if size(Yt,1)==1
                Yt=Yt';
            end
            yclt=Yt;
            nclass=max(max(yclt), max(modpls.classification.TrueClass));
            Yt=zeros(length(yclt), nclass);
        for i=1:length(yclt)
            Yt(i,yclt(i))=1;
        end
        Yt=Yt(:,1:max(modpls.classification.TrueClass));
        else
            [~,yclt]=max(Yt,[],2);
        end            
    end
    
    % ===========================
    % Preprocessing
    % Apply to the unknown X the preprocessing parameters of the calibration set.
    Xp=apprepfn(Xt,modpls.options.prepr{1},modpls.preprpars{1,1}, modpls.preprpars{1,2});
    T=Xp*modpls.Rweights;
    Yhat=Xp*modpls.regcoef; 
    U=Yhat*modpls.loadings{2};
    modelpred.scores={T U}; 
    modelpred.loadings=modpls.loadings; 
    Ypred=unprepfn(Yhat,modpls.options.prepr{1},modpls.preprpars{2,1}, modpls.preprpars{2,2});
    modelpred.preds=Ypred; 

    if ~isempty(Yt)
        modelpred.YResiduals=Yt-Ypred; 
        modelpred.bias=mean(modelpred.YResiduals); 
        modelpred.RMSEP=sqrt(sum(modelpred.YResiduals.^2)/size(Ypred,1)); 
        RSS=sum(sum(modelpred.YResiduals.^2));
        TSS=sum(sum((Yt-repmat(mean(Yt), size(Yt,1), 1)).^2)); 
        modelpred.R2=1-(RSS./TSS(1));
        PRESS = sum(sum((Yt-Ypred).^2));
        modelpred.Q2=1-PRESS./TSS(1);    

        % class evaluation
        [~,class_true] = max(Yt');
        resthr = plsdafindthr(Ypred,class_true');
        modelpred.threshold = resthr;


    end 

    if strcmp(modpls.options.detailedoutput,'on')
         q=sum((Xp-T*modpls.loadings{1}').^2,2);
        t2=diag((size(modpls.scores{1},1)-1)*T/(modpls.scores{1}'*modpls.scores{1})*T'); %tsquare            
        modelpred.tsq=t2;
        modelpred.qres=q;
        modelpred.tsqlim=modpls.tsqlim;
        modelpred.qlim=modpls.qlim;
        modelpred.tsqr=t2/modpls.tsqlim; 
        modelpred.qr=q/modpls.qlim;        
    end

    % ===========================
    
    
    modelpred.classification.TrueClass=[];
    switch modpls.options.classcriterion
        case 'maxY'
            [~,YPredClass]=max(modelpred.preds,[],2);
            modelpred.classification.PredClass=YPredClass;
        case 'lda'
            if strcmp(modpls.options.ldax, 'scores')
                xx=modpls.scores{1};
                xxt=modelpred.scores{1};
            elseif strcmp(modpls.options.ldax, 'predy')
                xx=modpls.preds;
                xx=xx(:,1:end-1);
                xxt=modelpred.preds;
                xxt=xxt(:,1:end-1);
                
            end
            
            mlda= ldaf(xx,xxt,modpls.classification.TrueClass, modpls.options.prior);
            modelpred.classification.PredClass=mlda.PredClass;
            modelpred.classification.probability=mlda.PosteriorProb;
    end
    
    if ~isempty(Yt)
        modelpred.classification.TrueClass=yclt;
        cmatrix=cmatrixcalc(modelpred.classification.PredClass, yclt, size(modelpred.preds,2));
        modelpred.classification.ConfusionMatrix=cmatrix;
        ccrate(1:size(cmatrix,2))=100*diag(cmatrix)./sum(cmatrix,2);
        modelpred.classification.CorrClassRate=ccrate;
        modelpred.classification.CorrClassTot=100*sum(diag(cmatrix))/length(yclt);
        [~,class_true] = max(Yt');
        assigned_class = plsdafindclass(Ypred,modpls.threshold.class_thr);
        modelpred.class_param = calc_class_param(assigned_class',class_true');
    end    
    
    function cmatrix=cmatrixcalc(Ypred, Ytrue, nc)
        ncts=max(Ytrue);
    
        if ncts>nc
            cmatrix=zeros(ncts,nc);
        else
            cmatrix=zeros(nc,nc);
        end
    
        for i=1:length(Ytrue)
            cmatrix(Ytrue(i), Ypred(i))=cmatrix(Ytrue(i), Ypred(i))+1;
        end
    end
end