# Filtrace obrazu
function cviceni10

#DFT
obr=imread("obrcb.jpg");
figure(1)
imshow(obr)
spektrum=ifftshift(fft2(fftshift(obr)));
hustota=abs(spektrum);
hustota=hustota/max(max(hustota))*255;
figure(2)
imshow(hustota);



#souradnice obrazu
N=size(obr)(1);
idx=-N/2:N/2-1;
osa=idx/(N/2);
[xx,yy]=meshgrid(osa);

#rekonstrukce
filtr=ones(N);
fspektrum=filtr.*spektrum;
recon=abs(ifftshift(ifft2(fftshift(fspektrum))));
recon=recon/max(max(recon))*255;
imshow(uint8(recon))

#dolni propust
filtr=exp(-80*(xx.^2+yy.^2));
fspektrum=filtr.*spektrum;
recon=abs(ifftshift(ifft2(fftshift(fspektrum))));
recon=recon/max(max(recon))*255;
imshow(uint8(recon))



#horni propust
filtr=1-exp(-50*(xx.^2+yy.^2));
fspektrum=filtr.*spektrum;
recon=abs(ifftshift(ifft2(fftshift(fspektrum))));
recon=recon/max(max(recon))*255;
imshow(uint8(recon))


#pasmova propust
filtr=exp(-200*(sqrt(xx.^2+yy.^2)-0.1).^2);
#imshow(filtr);break
fspektrum=filtr.*spektrum;
recon=abs(ifftshift(ifft2(fftshift(fspektrum))));
recon=recon/max(max(recon))*255;
imshow(uint8(recon))



#filtrace periodickych horizontalnich struktur
filtr=zeros(N);
for i=1:4
   filtr+=exp(-180*(xx.^2+(yy-i*0.21).^2));
   filtr+=exp(-180*(xx.^2+(yy+i*0.21).^2));
endfor
filtr=1-filtr;
imshow(filtr.*hustota);

fspektrum=filtr.*spektrum;
recon=abs(ifftshift(ifft2(fftshift(fspektrum))));
recon=recon/max(max(recon))*255;
imshow(uint8(recon))
break

clear recon
clear obr
clear spektrum
#barevny
obr=imread("obr.jpg");
for i=1:3
  spektrum{i}=ifftshift(fft2(fftshift(obr(:,:,i))));
  fspektrum=filtr.*spektrum{i};
  tem=abs(ifftshift(ifft2(fftshift(fspektrum))));
  recon(:,:,i)=uint8(tem/max(max(tem))*255);
endfor
imshow(recon)

clear filtr
#filtry aplikovany na barvy
filtr{1}=ones(N);
filtr{2}=ones(N);
filtr{3}=1-exp(-80000*(xx.^2+yy.^2));
for i=1:3
  fspektrum=filtr{i}.*spektrum{i};
  tem=abs(ifftshift(ifft2(fftshift(fspektrum))));
  recon(:,:,i)=uint8(tem/max(max(tem))*255);
endfor
imshow(recon)
break

endfunction