% Grayscale BM3D denoising demo file, based on 
% Y. Mäkinen, L. Azzari, A. Foi, 2020,
% "Collaborative Filtering of Correlated Noise: Exact Transform-Domain Variance for Improved Shrinkage and Patch Matching",
% in IEEE Transactions on Image Processing, vol. 29, pp. 8339-8354.
addpath '/Users/marketatkadlecova/Desktop/Tkadlecová/obrazky/'
f=@(x) 3.504*10^(-5)*(x+0.06)-1.874*10^(-6);
I_shep=im2double(imread('obrazky/Shep_orig.tif'));
A_shep=im2double(imread('noisy_shep.tif'));
I_parrots=im2double(imread('obrazky/Parrots_prog.tif'));
A_parrots=im2double(imread('noisy_parrots.tif'));
I_JIMA=im2double(imread('obrazky/JIMA.tif'));
A_JIMA=im2double(imread('noisy_JIMA.tif'));
load('local_var_JIMA.mat')
load('local_var_phantom.mat')
load('local_var_parrots.mat')

% ---
% The location of the BM3D files -- this folder only contains demo data
addpath('bm3d');

% Experiment specifications   
%imagename = 'cameraman256.png';

% Load noise-free image
%y = im2double(imread(imagename));

% % Possible noise types to be generated 'gw', 'g1', 'g2', 'g3', 'g4', 'g1w',
% % 'g2w', 'g3w', 'g4w'.
% noise_type =  'g3';
% noise_var = 0.02; % Noise variance
% seed = 0; % seed for pseudorandom noise realization
% % Generate noise with given PSD
% [noise, PSD, kernel] = getExperimentNoise(noise_type, noise_var, seed, size(y));
% % N.B.: For the sake of simulating a more realistic acquisition scenario,
% % the generated noise is *not* circulant. Therefore there is a slight
% % discrepancy between PSD and the actual PSD computed from infinitely many
% % realizations of this noise with different seeds.
% 
% % Generate noisy image corrupted by additive spatially correlated noise
% % with noise power spectrum PSD
% z = y + noise;

% Call BM3D With the default settings.
y_est = BM3D(A_shep, 0.85*max(sqrt(local_var1(:))),'refilter');
[SSIM_shep, MSE_shep, PSNR_shep]=FR_metrics(y_est,I_shep)
%0.9

% Call BM3D With the default settings.
y_parrots = BM3D(A_parrots, 0.85*max(sqrt(local_var2(:))),'refilter');
[SSIM_parrots, MSE_parrots, PSNR_parrots]=FR_metrics(y_parrots,I_parrots)


% Call BM3D With the default settings.
y_JIMA = BM3D(A_JIMA, 0.85*max(sqrt(local_var3(:))),'refilter');
[SSIM_JIMA, MSE_JIMA, PSNR_JIMA]=FR_metrics(y_JIMA,I_JIMA)



% For other settings, use BM3DProfile.
% profile = BM3DProfile(); % equivalent to profile = BM3DProfile('np');
% profile.gamma = 6;  % redefine value of gamma parameter
% y_est = BM3D(z, PSD, profile);

% Note: For white noise, you may instead of the PSD 
% also pass a standard deviation 
% y_est = BM3D(z, sqrt(noise_var));


