X = tensor(activity_tensor);

%%

opts = struct();
opts.maxiters =1000;     % aumenta il numero massimo di iterazioni
opts.tol = 1e-16;         % richiedi maggiore accuratezza prima di fermarti
opts.printitn = 100;      % stampa log ogni 10 iterazioni

% Rank della decomposizione

S = {};
% Decomposizione CP con ALS usando opzioni personalizzate
max_R = 4;
fits = zeros(max_R, 1);


for R = 4:4
    fprintf('\n--- Rank %d ---\n', R);
    rng(13);
    [M, ~] = cp_nmu(X, R, opts);
    X_hat = full(M);
    fits(R) = 1 - norm(X - X_hat) / norm(X);
    

end



% Plot
figure;
plot(1:max_R, fits, '-o', 'LineWidth', 2);
xlabel('Numero componenti (R)');
ylabel('Fit');
title('Fit in funzione del numero di componenti - cp\_nmu');
grid on;


%Ricostruzione e errore
X_hat = full(M);
rel_error = norm(X - X_hat) / norm(X);
fprintf('Errore relativo di ricostruzione: %.6f\n', rel_error);

%%
visualize_components_full(M, labels);

%%
visualize_trials_grid_with_trends(M, labels, 1, 'robust');

%%
visualize_time_grid(M, 1)