function x=GMRF_sampler(theta,G,N)
% x=GMRF_sampler(theta,G,N)
%
% Gauss Markov Random Field (GMRF) sampler, generates synthetic texture
% based on the AR-model theta
%
% theta:	2D array of filter parameters
% G:		Number of gray levels in output image
% N:		Size of output image is NxN

% Prepare the parameter matrix A
waitbarhandle=waitbar(0,'Please wait while sampling GMRF',1);
A=zeros(N);
tc=(size(theta)+[1 1])/2;
ts=size(theta);
A(1:ts(1),1:ts(2))=-theta;
A=[A(tc(1):N,:) 
    A(1:tc(1)-1,:)];
A=[A(:,tc(2):N) A(:,1:tc(2)-1)];
A(1,1)=-A(1,1);

waitbar(0.0411);

eta=randn(N);
waitbar(0.0918);

eta=fft2(eta);
waitbar(0.3066);

A=fft2(A);
waitbar(0.4796);

x = eta./sqrt(A);
waitbar(0.5625);

x=ifft2(x);
waitbar(0.8732);

% Make both real and imaginary part positive
x=x - min(min(real(x)) + min(imag(x))) * (1+i);
x=round(scale(abs(x),0,G-1));

waitbar(1.0000);
close(waitbarhandle);
