clear
close all
clc

% 'Dammann_numerical.mat' contains the transition points data in [Zhou1995]
% [Zhou1995] C. Zhou and L. Liu, "Numerical study of Dammann array 
% illuminators," Applied Optics, vol. 34, p. 5961, September 1995. 
load('Dammann_Numerical_Data.mat');
Dammann_numerical=Dammann_Numerical_Data.Dammann_numerical;
fprintf('Numerical data loaded.\n');

%% BEGIN of Dammann Grating Generation
fprintf('Generating Dammann Grating...\n');
Edge_A_Array_Size=6;
Edge_B_Array_Size=6;
Edge_C_Array_Size=6;

Angle=60;

Filename_Prefix=sprintf('Dammann_%d-by-%d-by-%d_%d_',Edge_A_Array_Size,Edge_B_Array_Size,Edge_C_Array_Size,Angle);

Edge_A_Transition_Points=Dammann_numerical{table2array(cell2table(Dammann_numerical(:,1)))==Edge_A_Array_Size,2};
Size_A_Transition_Points=size(Edge_A_Transition_Points,2);
Edge_B_Transition_Points=Dammann_numerical{table2array(cell2table(Dammann_numerical(:,1)))==Edge_B_Array_Size,2};
Size_B_Transition_Points=size(Edge_B_Transition_Points,2);
Edge_C_Transition_Points=Dammann_numerical{table2array(cell2table(Dammann_numerical(:,1)))==Edge_C_Array_Size,2};
Size_C_Transition_Points=size(Edge_C_Transition_Points,2);

Dammann_Pixel_A_Size=2048;
Dammann_Pixel_B_Size=2048;
Dammann_Pixel_C_Size=2048;
Dammann_Pixel_FINAL_X_Size=1024;
Dammann_Pixel_FINAL_Y_Size=1024;
Dammann_Gap=128;

Light_X_Size=ceil(1*Dammann_Pixel_FINAL_X_Size);
Light_Y_Size=ceil(1*Dammann_Pixel_FINAL_Y_Size);
Light_Sigma=0.5;

Dammann_1D_A=zeros(Dammann_Pixel_A_Size,Dammann_Pixel_A_Size);
Dammann_1D_A_FINAL=zeros(Dammann_Pixel_A_Size,Dammann_Pixel_A_Size);
Dammann_1D_B=zeros(Dammann_Pixel_B_Size,Dammann_Pixel_B_Size);
Dammann_1D_B_FINAL=zeros(Dammann_Pixel_B_Size,Dammann_Pixel_B_Size);
Dammann_1D_C=zeros(Dammann_Pixel_C_Size,Dammann_Pixel_C_Size);
Dammann_1D_C_FINAL=zeros(Dammann_Pixel_C_Size,Dammann_Pixel_C_Size);

for i=1:1:Dammann_Pixel_A_Size
    Dammann_1D_A(:,i)=mod(i-1,Dammann_Gap);
end
for i=1:1:Dammann_Pixel_B_Size
    Dammann_1D_B(i,:)=mod(i-1,Dammann_Gap);
end
for i=1:1:Dammann_Pixel_B_Size
    Dammann_1D_C(i,:)=mod(i-1,Dammann_Gap);
end
Dammann_1D_A=uint32(Dammann_1D_A*4294967296/Dammann_Gap);
Dammann_1D_B=uint32(Dammann_1D_B*4294967296/Dammann_Gap);
Dammann_1D_C=uint32(Dammann_1D_C*4294967296/Dammann_Gap);


for ia=1:1:Size_A_Transition_Points
    bin_A(:,:,ia)=imbinarize(Dammann_1D_A,Edge_A_Transition_Points(ia));
    if mod(ia,2)
        Dammann_1D_A_FINAL=logical(Dammann_1D_A_FINAL+bin_A(:,:,ia));
    else
        Dammann_1D_A_FINAL=logical(Dammann_1D_A_FINAL-bin_A(:,:,ia));
    end
end

for ib=1:1:Size_B_Transition_Points
    bin_B(:,:,ib)=imbinarize(Dammann_1D_B,Edge_B_Transition_Points(ib));
    if mod(ib,2)
        Dammann_1D_B_FINAL=logical(Dammann_1D_B_FINAL+bin_B(:,:,ib));
    else
        Dammann_1D_B_FINAL=logical(Dammann_1D_B_FINAL-bin_B(:,:,ib));
    end
end

for ic=1:1:Size_C_Transition_Points
    bin_C(:,:,ic)=imbinarize(Dammann_1D_C,Edge_C_Transition_Points(ic));
    if mod(ic,2)
        Dammann_1D_C_FINAL=logical(Dammann_1D_C_FINAL+bin_C(:,:,ic));
    else
        Dammann_1D_C_FINAL=logical(Dammann_1D_C_FINAL-bin_C(:,:,ic));
    end
end
crop=[Dammann_Pixel_FINAL_X_Size/2 Dammann_Pixel_FINAL_Y_Size/2 Dammann_Pixel_FINAL_X_Size-1 Dammann_Pixel_FINAL_Y_Size-1];
Dammann_1D_A_FINAL=imcrop(imrotate(Dammann_1D_A_FINAL,90),crop);
Dammann_1D_B_FINAL=imcrop(imrotate(Dammann_1D_B_FINAL,Angle,'crop'),crop);
Dammann_1D_C_FINAL=imcrop(imrotate(Dammann_1D_C_FINAL,180-Angle,'crop'),crop);

% figure
% imshow(Dammann_1D_A_FINAL);
% figure
% imshow(Dammann_1D_B_FINAL);
% figure
% imshow(Dammann_1D_C_FINAL);

% Dammann_2D=xor(xor(Dammann_1D_A_FINAL,Dammann_1D_B_FINAL),Dammann_1D_C_FINAL);
Dammann_2D=xor(Dammann_1D_B_FINAL,Dammann_1D_C_FINAL);

figure
imshow(Dammann_2D);
fprintf('\tDammann Grating Generated.\n');
% END of Dammann Grating Generation


%% START Image Export
fprintf('Exporting images...\n');
%  Dammann Grating
Dammann_2D_Image=uint16(Dammann_2D*65535);
Dammann_2D_Image_Filename=sprintf('%sPattern_Original_PX%d_PY%d_G%d',Filename_Prefix,Dammann_Pixel_FINAL_X_Size,Dammann_Pixel_FINAL_Y_Size,Dammann_Gap);
imwrite(Dammann_2D_Image,[pwd '\Dammann Grating\' Dammann_2D_Image_Filename '.png']);
fprintf('\tExporting finished.\n');

%% START Image Show
fprintf('Image Show...\n');
figure('Name',Dammann_2D_Image_Filename,'NumberTitle','off');
imshow(Dammann_2D_Image);
fprintf('\t-----Complete-----\n');
