%% LABORATORIO 1 - METODI STATISTICI PER LA BIOINGEGNERIA
%% Canale B - A.A. 2024/2025
%% Docente: Enrico Longato

%% Esercizio 3

%% Direttive per la gestione del workspace e delle figure
clc        
clear all  
close all  

%% Creazione matrici A e B date

A = [-8     2     7     3
     -8    -9   -10     5
      9    -6   -10     3
     10    -3    -7    -1];

B = [ 2     5     6     1
     -3    -5    -7     4
      6    -1     0     2
     12    -2     0     4];


%% Calcolare il rango
rank_A = rank(A); % Righe/colonne ortogonali
rank_B = rank(B); % Terza e quarta riga combinazione lineare

%% Calcolare gli autovalori
eigenvalues_A = eig(A); % Sono complessi (poco male...)
eigenvalues_B = eig(B); % Uno è zero (e infatti il rango non è pieno)

% Verifichiamo che B non è invertibile
B_not_really_inv = inv(B);

%% Calcolare la matrice D = "A per inversa di A"
% Calcolo con le operazioni
inv_A = inv(A);
D = A*inv_A; % Non ci stupisce che sia l'identità, ma è "un po' sporca"

disp('Rispetto a quanto atteso, la differenza è:')
disp(D - eye(size(A))) % eye è l'identità

%% Utilizzare gli operatori / e \ in sostituzione al calcolo diretto
% MATLAB stesso ci suggerisce che A*inv(A) sia meglio scriverlo come A/A
D1 = A/A; % Infatti, il risultato è pulito.

disp(' ')
disp('Rispetto a quanto atteso, la differenza è:')
disp(D1 - eye(size(A))) % Nulla!

%% Calcolare il prodotto E = "inversa di A, per B" nei due modi
% Direttamente, il risultato è "sporco" come prima, ma è difficile vederlo
E = inv(A)*B;

% MATLAB ci suggerisce di fare come segue
E1 = A\B; % Notare che ora lo slash è il backslash \

% Permane una differenza simile a prima
disp(E1 - E)

%% Calcolare il prodotto elemento per elemento
F = A .* B; % Ricordarsi il punto

%% Verificare che è diverso dal prodotto matriciale
disp(F - A * B) % Tutt'altra cosa!

%% Trovare indici riga e colonna di ogni elemento < 0
i_neg = find(F < 0); % Lineare
[r_neg, c_neg] = find(F < 0); % Riga/colonna
 
% Mostriamo i valori e gli indici
for k = 1:length(i_neg)
    i_curr = i_neg(k); % Indice lineare corrente
    r_curr = r_neg(k); % Indice riga corrente
    c_curr = c_neg(k); % Indice colonna corrente

    disp(['Elemento in posizione lineare ', num2str(i_curr), ...
        ' e posizione riga/colonna (' num2str(r_curr), ...
        ', ' num2str(c_curr), ') uguale a ' num2str(F(r_curr, c_curr))])
end

%% Trovare minimo e massimo della matrice A e i loro indici
[min_A, i_min_A] = min(A(:)); % Sono il min e l'argmin
[r_min_A, c_min_A] = ind2sub(size(A), i_min_A); % Solita conversione

[max_A, i_max_A] = max(A(:)); % Sono il max e l'argmax
[r_max_A, c_max_A] = ind2sub(size(A), i_max_A); % Solita conversione


%% Ripetere l'operazione senza vettorizzare
[max_A_nonvec, i_max_A_nonvec] = max(A); % Notare che A non è "srotolata"

r_max_A_nonvec = i_max_A_nonvec; % Indice riga
c_max_A_nonvec = 1:size(A, 2); % Le colonne sono semplicemente 4

% Di default, le operazioni vettoriali sono svolte per colonne; in questo
% caso si poteva usare il parametro "dim", cioè il terzo, per ovviare e
% trovare il massimo di ciascuna riga.
[max_A_nonvec_by_rows, i_max_A_nonvec_by_rows] = max(A, [], 2);

r_max_A_nonvec_by_rows = 1:size(A, 1); % Le righe sono  4
c_max_A_nonvec_by_rows = i_max_A_nonvec_by_rows; % Indice colonna