function  obj = train(obj,displayMode)
% TRAIN is used to train the coefficients of a global polynom. The usage is
%
% obj = train(obj)
%
% OUTPUTS:
%   obj     a global polynom object
%

% Tobias Ebert, 11-Oct-2011
% Institut für Mess- und Regelungstechnik, Universität Siegen, Deutschland
% Institute of Measurement and Control, University of Siegen, Germany
% Copyright (c) 2011 by Tobias Ebert

% handle the display mode
if exist('displayMode','var') && islogical(displayMode)
    obj.history.displayMode = displayMode;
end

obj.numberOfInputs = size(obj.input,2);
obj.numberOfOutputs = size(obj.output,2);

if isempty(obj.maxNumberOfParameters)
    error('globalPolynom:train','maxNumberOfParameters is empty')
end

if obj.maxNumberOfParameters > size(obj.input,1)/2
    warning('globalPolynom:train','maxNumberOfParameters will be reduzed to size(obj.input,1)/2')
    obj.maxNumberOfParameters = size(obj.input,1)/2;
end

if isempty(obj.xRegressorExponentMatrix)
    % create temporary exponent Matrix
    if obj.maxNumberOfParameters == 1
        %polyOrder = 0;
        exponentMatrix = 0;
    else
        polyOrder = 1;
        exponentMatrix = obj.buildExponentMatrix('polynomial',polyOrder,size(obj.input,2));
        %exponentMatrix = obj.buildPolynomialExponentMatrix(polyOrder,size(obj.input,2));
        while size(exponentMatrix,1)+1 < obj.maxNumberOfParameters
            polyOrder = polyOrder+1;
            exponentMatrix = obj.buildExponentMatrix('polynomial',polyOrder,size(obj.input,2));
            
            % delete higher orders of exponents
            idx2Delete = any(exponentMatrix > obj.xRegressorMaxPower,2);
            exponentMatrix = exponentMatrix(~idx2Delete,:);
            
        end
    end
    obj.xRegressorDegree = polyOrder;
    
    % sort the exponent Matrix, according to lowest order to highest
    % calculate the sum of all exponents
    sumExponents = sum(exponentMatrix,2);
    % get highest power for all expoents of each regressor
    maxExponents = max(exponentMatrix,[],2);
    % initialise vector for the new indice sequence to use
    idxOrder = [];
    for k = unique(sumExponents') % loop over sum of exponents
        % get the indices of the regressor to sort for, start with highest
        idx2Sort = (sumExponents == k);
        for kk = unique(maxExponents(idx2Sort))' % loop over max of exponents
            idxOrder = [idxOrder; find((maxExponents == kk) & idx2Sort)];
        end
    end
    % use the idxOrder vector to shape expoentMatrix in the new Order
    exponentMatrix = exponentMatrix(idxOrder,:);
    
    % save exponent matrix in obj
    obj.xRegressorExponentMatrix = exponentMatrix;
end
% create regressor
obj.xRegressor = obj.data2xRegressor(obj.input);



%% Initialisation, Iteration 1
tStart = tic;

% calculate parameters
obj.localModels(1).parameter = obj.estimateParametersLocal(obj.xRegressor, obj.output, [], 1, obj.dataWeighting);

% set last model as active
obj.leafModels(1) = true;

% set history
obj.history = writeHistory(obj.history,obj,tStart);

if obj.history.displayMode
    fprintf('\n\n%d. Iteration. Number of x regressors = %d ...', ...
        find(obj.leafModels), size(obj.xRegressor,2));
end
if obj.history.displayMode
    fprintf('\n   Subset selection used %d parameters. J = %f.',...
        obj.history.currentNumberOfParameters(obj.history.iteration),obj.history.globalLossFunction(obj.history.iteration))
end


%% loop
while obj.history.checkTerminationCriterions(obj) && ...
        obj.history.iteration <= size(obj.xRegressor,2)
    
    % calculate parameters
    obj.localModels(end+1).parameter = obj.estimateParametersLocal(obj.xRegressor, obj.output, [], obj.history.iteration+1, obj.dataWeighting);
    
    % set last model as active
    obj.leafModels = false(1,length(obj.localModels));
    obj.leafModels(end) = true;
    
    % set history
    obj.history = writeHistory(obj.history,obj,tStart);
    
    if obj.history.displayMode
        fprintf('\n\n%d. Iteration. Number of x regressors = %d ...', ...
            find(obj.leafModels), size(obj.xRegressor,2));
    end
    if obj.history.displayMode
        fprintf('\n   Subset selection used %d parameters. J = %f.',...
            obj.history.currentNumberOfParameters(obj.history.iteration),obj.history.globalLossFunction(obj.history.iteration))
    end
    
end

[~, obj.suggestedNet] = min(obj.history.penaltyLossFunction);

if obj.history.displayMode
    fprintf('\n\nNet %d with %d parameters is suggested as the model with the best complexity trade-off.\n\n', ...
        obj.suggestedNet, obj.history.currentNumberOfParameters(obj.suggestedNet));
end

end