%%
N = 256;

% Bild
f = zeros(1,N);
f(7) = 0.1; f(8) = 0.5; f(9) = 0.9; f(10) = 1; f(11) = 0.9; f(12) = 0.5; f(13) = 0.1;
f(60:68) = 0.8;
f(121:129) = sin([1:9]/10*pi);
f(200:209) = 0.5;

% PSF (point spread function)
h = zeros(1,N);
h(1) = 1/5;
h(2) = 1/5;h(3) = 1/5;
h(N) = 1/5;h(N-1) = 1/5;

fh = fft(f);
hh = fft(h);
dh = fh.*hh;

% Daten
d = ifft(dh);

figure(1)
subplot(3,1,1),plot(d), ylim([-0.1,1]), title('d')
subplot(3,1,2),plot(f), ylim([-0.1,1]), title('f')
subplot(3,1,3),plot(h), ylim([-0.1,1]), title('h')

%%

d = d + 0.03*randn(1,N); %  Messfehler
dh = fft(d);
fh = dh./hh;
f  = ifft(fh); % Rekonstruiertes Bild

figure(2)
subplot(3,1,1),plot(d), ylim([-0.1,1]), title('d')
subplot(3,1,2),plot(f), ylim([-0.1,1]), title('f')
subplot(3,1,3),plot(h), ylim([-0.1,1]), title('h')