% Example of how to use NCGP for deblurring Poisson-corrupted images   
% with Tikhonov regularization

clear all
close all

% Phantom test problem
load('galaxyNGC1288');

% fft of the psf image
TF = fft2(fftshift(psf));

% Set the regularization parameter value
eta = 0.005;

% initial iterate
rec0 = ones(size(gn))*sum(gn(:))/numel(gn);

% Save data in structure array 
data.obj = obj; % exact image
data.gn = gn; % blurred and noisy image
data.bg = bg; % background
data.TF = TF; % fft of the psf image
data.rec0 = rec0; % initial guess

% Set objective and constraint functions and save in structure fun_params
fun_params.obj = 'KLTKb_obj'; % objective function
fun_params.eta = eta; % regularization parameter 
fun_params.Amult = 'Amult'; % matrix-vector product with the Hessian 

% Set up parameters
params.max_fft = 1e3; % maximum number of FFTs
params.kmax = 20; % maximum number of iterations
params.tol_gp = 1e-4; % projected gradient tolerance
params.tol_f = 1e-2; % objective function tolerance  
params.pcg_max = 30; % maximum number of CG iterations
params.pcg_tol = 0.1; % CG relative tolerance
params.precond = 0; % if true=> diagonal preconditioner for CG  

% Newton-Cg Projection method
[x,histout] = NCGP(data,fun_params,params);

% Print numerical results and visualize images
NPoutput(x,histout,data)










