光束远场强度整型的二元光学器件设计基于G-S算法

本程序实现了对远场光束的强度整形的二元光学器件的设计，程序采用窗口文件读取模式，从一张待设计的光强分布图程序将利用G-S算法设计出所需要的二元光学器件。

% This program uses G-S algorithm to realize the design of beam shaping
% device
%%                          Initiation
clear;
clc;
N=512;
M=512;
f=500;
lamda=560e-6;
L0=15;
r=5;
dx0=L0/N;
dy0=L0/M;
L=lamda*f*N/L0;
dx=L/N;
dy=L/M;
% - N and M are the number of horizontal and vertical sampling points
% respectively
% - f is the focal length and its nuit is micrometre
% - lamda is the wave length and its unit is micrometre
% - L0 is the size of binary optical element and its unit is micrometre
% - r is the clear radius of the binary optical element and its unit is 
% micrometre
% - dx0 and dy0 are the horizontal and vertical sample interval in binary
% optical element respectively and their unit are micrometre
% - L is the size of image plane and its unit is micrometre
% - dx and dy are the horizontal and vertical sample interval in image 
% plane
%%     Creating the orginal BOE palne and the ideal image plane 
[mn]=meshgrid（linspace（-N/2N/2-1N）linspace（-M/2M/2-1M））;
I0=zeros（NM）;
R=sqrt（（dx0*m）.^2+（dy0*n）.^2）;
index1=find（R<=r）;
I0（index1）=1;
energy0=sum（sum（I0））;
chemin=‘.\image\‘;
[nomchemin]=uigetfile（[chemin‘*.*‘][‘Input the picture‘]100100）;
[IMAP]=imread（[cheminnom]）;
I=rgb2gray（I）;
I=imresize（I[1024 1024]）;
I=im2double（~imbinarize（I））;
% The image sampler
I_i=zeros（NM）;
for i=1:1024/N:1024
    for j=1:1024/M:1024
        I_i（（i-1）*N/1024+1（j-1）*M/1024+1）=I（ij）;
    end
end
% energy conservation
energy=sum（sum（I_i））;
I_i=I_i*energy0/energy;
%%                          G-S algorithm
[I0_phaseerr]=GSA_Q（I0I_i8100）;
%%      Actual result of imaging and calculating the rms and ESEU
I1=abs（fftshift（fft2（sqrt（I0）.*exp（1i*I0_phase）））/sqrt（M*N））.^2;
I1_u=I1/max（max（I1））;
index2=find（I_i>0.8）;
ave_I1=sum（I1（index2））/（M*N）;
rms=1/（length（index2）-1）*sqrt（sum（（（I1（index2）-ave_I1）/ave_I1）.^2））;
eta=sum（I1（index2））/sum（sum（I1））;
%%                       Ploting the result
figure（1）;
mesh（I_i）;
title（‘The intensity distribution of ideal image plane‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（2）;
mesh（I1）;
title（‘The intensity distribution of actual image plane‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（3）;
bar3（I0_phase）;
title（‘The phase distribution‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（4）;
imshow（I_i）;
title（‘The ideal image plane‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（5）
imshow（I1_u）;
title（‘The actual image plane‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（6）
imshow（I0_phase）;
title（‘The phase of the binanry optical element‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
figure（7）
stairs（I0_phase（N/2:））;
title（‘The phase of the binanry optical element （Central section）‘‘FontName‘‘Times New Roman‘‘FontSize‘32）
xlabel（‘\itPeriods‘‘FontName‘‘Times New Roman‘‘FontSize‘20）;
ylabel（‘\itPhase/rad‘‘FontName‘‘Times New Roman‘‘FontSize‘20）;
axis（[1 M -4 4]）;
figure（8）
plot（err‘markersi

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        3293  2018-09-06 23:31  BOE_GSA\BOE_GSA.m
     文件         756  2018-06-08 11:07  BOE_GSA\GSA_Q.m