% LETTURA, SCRITTURA E CONVERSIONE DI IMMAGINI
% ESEMPIO DI FILTRO DI SOBEL PER EDGE-DETECTION

% "edgedemo" apre un dimostrativo Matlab con vari metodi (tra cui Sobel)

clear all;

% load trees;   % carica l'immagine nella matrice "X" e la mappa colore in "map"
% im1= ind2gray(X,map);	% converte in scala di grigi

% lettura immagine (pout.tif ? di libreria Matlab)
im1= imread('pout.tif');	% im1 ? una matrice (immagine B/N)
% vedere 'help imread' x formati immagine supportati e sintassi x immagini a colori
% display dell'immagine
figure(1); clf;
imshow(im1);

% matrice di Sobel
s= 0.25*[1 2 1; 0 0 0; -1 -2 -1];	% ho aggiunto *0.25 per avere sum(s(k,:))=1

im1h= conv2(im1,s);	% bordi orizzontali
figure(101); clf; imshow(im1h);
im1v= conv2(im1,s');	% bordi verticali
figure(102); clf; imshow(im1v);
% im2= abs(im1hB)+abs(im1vB);	% calcola dove l'una o l'altra sono grandi
im2= sqrt(im1h.^2+im1v.^2);	% stessa cosa...: teo. di Pitagora

figure(2); clf;
imshow(im2);
title('bordi difficili da visualizzare per immagini in scala di grigi');
% alcuni punti "troppo luminosi" impediscono di visualizzare i bordi

k=  0.05; % fattore di attenuazione luminosit?... va cercato iterativamente
im2att= k*im2;		% immagine "attenuata"
figure(3); clf;
imshow(im2att);
title('bordi enfatizzati smorzando i punti pi? luminosi');

% per vedere ancora meglio i bordi, applico una soglia all'immagine filtrata e attenuata
sogliaim2att= graythresh(im2att);	% trova la soglia che discrimina B e N
im2attthr= im2bw(im2att,sogliaim2att);	% converte in B/N (matrice di tipo 'logical')
figure(31); clf;
imshow(im2attthr);
title('qui si ? applicata una soglia all''immagine filtrata e attenuata');



% PROVIAMO A LAVORARE SU UN'IMMAGINA CONVERTITA IN B/N

% conversione dell'immagine da scala di grigi a B/N
sogliaim1= graythresh(im1);	% trova la soglia che discrimina B e N
im1bn= im2bw(im1,sogliaim1);	% converte in B/N (matrice di tipo 'logical'
figure(4); clf;
imshow(im1bn);

im1hbn= conv2(im1bn,s);	% bordi orizzontali
figure(401); clf; imshow(im1hbn);
im1vbn= conv2(im1bn,s');	% bordi verticali
figure(402); clf; imshow(im1vbn);
% im2bn= abs(im1hbnB)+abs(im1vbnB);	% calcola dove l'una o l'altra sono grandi
im2bn= sqrt(im1hbn.^2+im1vbn.^2);

figure(5); clf;
imshow(im2bn);
title('bordi molto pi? chiari per immagini B/N');



% visualizziamo la risposta in ampiezza del filtro di Sobel
figure(10); clf;
freqz2(s);
title('frequenze normalizzate alla ''freq. di Nyquist'' unitaria');
% e del suo trasposto
figure(11); clf;
freqz2(s');
title('matrice s'': ruolo delle due dimensioni scambiato');
% NOTA: a parte il fattore scalare 0.25, la matrice s può essere ottenuta come prodotto colonna*riga di due vettori:
% s= sc*sr; con sc= [1; 0; -1] e sr= [1 2 1];
% sc e sr sono due filtri FIR con M+1 elementi (M=2) che agiscono in vericale (sc) e orizzontale (sr)
% sc è un filtro FIR a fase lineare di tipo III(antisimmetrico con M pari) con zeri in z=+-1 (frequenze f=0 e f=+-fNyquist)
% (un FIR tipo III ha SEMPRE degli zeri in z=+-1)
% sr è un filtro FIR a fase lineare di tipo I (simmetrico con M pari) con due zeri coincidenti in z=-1 (frequenza f=+-fNyquist)
% Tali proprietà sono ben evidenti nella risposta in ampiezza del filtro, sulla frequenza orizzontale Fx e verticale Fy


% scrittura immagine
imwrite(im2att,'pout_bordi.jpg');	% l'estensione specifica il formato
% richiesta informazioni sull'immagine: nell'intestazione del formato .gif
imfinfo('pout_bordi.jpg')		% niente ; in modo da visualizzare il risultato