-
Notifications
You must be signed in to change notification settings - Fork 6
/
function_FunObj_dualbinary.m
42 lines (38 loc) · 1.39 KB
/
function_FunObj_dualbinary.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
function [loss, df ] = function_FunObj_dualbinary( phase, source, z, Nx, Ny, thresholdh, thresholdl, maskFun, kickmaskFun,fresnelKernelFun, useGPU)
%
% This function is called by Matlab's fmincon library.
% Computes loss and gradient or BINARY NOVOCGH with a Dual threshold-based cost
% function that enforces either Dark or Bright voxels, or relaxes constraints elsewhere !
if useGPU
df = zeros(Nx, Ny, 'gpuArray');
phase = gpuArray(phase);
source = gpuArray(source);
else
df = zeros(Nx, Ny);
end
loss = 0;
phase = reshape(phase, [Nx, Ny]);
objectField = source.*exp(1i * phase);
for i = 1 : numel(z)
HStack = fresnelKernelFun(i,i);
mask = maskFun(i,i);
kickmask = kickmaskFun(i,i);
imagez = fftshift(fft2(objectField .* HStack));
imageInten = abs(imagez.^2);
maskh = mask .* (imageInten < thresholdh);
maskl = kickmask .* (imageInten > thresholdl);
diffh = maskh .* (imageInten - thresholdh);
diffl = maskl .* (imageInten - thresholdl);
temph = imagez.*diffh;
temph = conj(HStack).*(Nx*Ny*ifft2(ifftshift(temph)));
templ = imagez.*diffl;
templ = conj(HStack).*(Nx*Ny*ifft2(ifftshift(templ)));
%Compute losses
loss = loss + sum(sum(diffh.^2 + diffl.^2));
%Compute gradient
df = df + temph + templ;
end
dfphase = source.*(- real(df).*sin(phase) + imag(df) .* cos(phase));
df = gather(real(dfphase(:)));
loss = gather(real(loss));
end