同轴重建（无需相移）2.zip_全息_全息 重建_全息图 恢复_相位恢复_相位重建迭代

%%%%%同轴全息迭代获得相位%%%%% %%%%方法1%%%% %%%%%席特立，2017年10月210日编写%%%%% clear all close all clc %% Step1：定义文件参数 Folder='D:\1毕业设计\资料\'; %%设置所需要处理系列全息图所在文件夹位置 Folder_result='D:\1毕业设计\资料\'; %%设置输出图片所在文件夹位置 Hologram_name1='11'; %%%%背景文件名，此时无干涉条纹，为纯背景必须设置为数字 Hologram_name2='22'; %%全息图文件名，必须设置为数字 %% Step2：建立归一化全息图（normalized hologram） %%读入全息图 File_string=[Folder num2str(Hologram_name1)]; Hologram1=double(imread([File_string,'.bmp'])); %% MATLAB读入图像 % Hologram1=double(rgb2gray(imread([File_string,'.bmp']))); %% MATLAB读入图像 % figure (1),imagesc(Hologram1), colormap(gray) %%显示读入的图像 d = find(Hologram1==0); Hologram1(d(1:end)) =1; % [row, col] = find( Hologram1 = 0 ); % row，col是所有非零元素的横坐标和纵坐标 % num = size(row, 1); % 有多少个非零元素 % for i = 1:num % R(row(i), col(i)) = 1; % 非零元素置一 % end File_string=[Folder num2str(Hologram_name2)]; Hologram2=double(imread([File_string,'.bmp'])); %% MATLAB读入图像 % Hologram2=double(rgb2gray(imread([File_string,'.bmp']))); %% MATLAB读入图像 % figure (2),imagesc(Hologram2), colormap(gray) %%显示读入的图像 %%得到标准化全息图（normalized hologram） N_Hologram=(Hologram2)./(Hologram1); [xx,yy]=size(N_Hologram); %%获得图像的尺寸 % figure (3),imagesc(N_Hologram), colormap(gray) axis off %%显示标准化全息图 %%%Tukey窗 L=xx; %window width v=yy; r=0.2; %taper ratio of tukey window %produce 2d tukey window w=tukeywin(L,r).*tukeywin(v,r)'; % %%%Hamming窗 % w=hamming(xx).*hamming(yy)'; N_Hologram=N_Hologram.*w; %% Step2：定义实验参数 %%定义激光波长、CCD等参数 format long lambda=9550e-8; %%激光波长，单位米 delta=80e-6; %%CCD像素尺寸 CCD_Sizex=delta*xx; %%CCD横向尺寸 CCD_Sizey=delta*yy; %%CCD纵向尺寸 Depth_range=[116]; %%手动定义全息图重建距离范围，单位mm Depth_steps=1; if Depth_steps==1 Slice_depths=Depth_range(1); else Slice_depths=linspace(Depth_range(1),Depth_range(2), Depth_steps); end Magnification=1; %%Magnificaton of the system [nn,mm]=meshgrid(1:yy,1:xx); %%设定渐变矩阵 水平方向为mm 竖直方向为nn %% Step3：建立迭代关系 %%像平面复分布 N_Intensity=sqrt(N_Hologram); % figure, imagesc(N_Intensity),colormap(gray),axis image E=1; %%初始像平面指数项（初始假定） for n=1:2; Imaging_field_1=N_Intensity.*E; %%全息图平面复分布 for Reconstruction_distance=(Slice_depths/1000) %%循环数值重建 H=exp(i*2*pi*-Reconstruction_distance*sqrt((1/lambda).^2-((mm-xx/2)/CCD_Sizex).^2-((nn-yy/2)/CCD_Sizey).^2)); %%设定传递函数 Reconstruction_field_1=ifft2(fftshift(fftshift(fft2(Imaging_field_1)).*H)); % % Intensity=(abs(Reconstruction_field_1).^2); %重建物光波的强度图 % Ifigure=figure; imagesc(abs(Reconstruction_field_1));colormap(gray) % % % Ifigure=figure; imagesc(Intensity,[min(Intensity(:)),0.8*max(Intensity(:))]);colormap(gray) % axis image; axis off; %%画出强度图 % hold on; % title (['背景全息图数值重建像的强度图（重建距离d= ' num2str(Reconstruction_distance*1000) 'mm）']) % hold off % % % Path=[Folder_result,'1-Rec_distance_',num2str(Reconstruction_distance),'.jpg']; % print(Ifigure,'-djpeg',Path); %%以截屏的方式存储上述强度图，路径为Path end % % %%赋值法1 % % % exp1=Reconstruction_field_1./abs(Reconstruction_field_1); % A=(abs(Reconstruction_field_1)); % B=Reconstruction_field_1./A; % % % if all(A(:))<1 % % % break % % % else % A(abs(A)>=1)=1; %%判定吸收率不能大于1小于0 % Reconstruction_field_2=A.*B; %%赋值法2 Reconstruction_field_2= Reconstruction_field_1; Reconstruction_field_2(abs(Reconstruction_field_2)>=1)=1; H=exp(i*2*pi*Reconstruction_distance*sqrt((1/lambda).^2-((mm-xx/2)/CCD_Sizex).^2-((nn-yy/2)/CCD_Sizey).^2)); %%设定传递函数 Reconstruction_field_3=ifft2(fftshift(fftshift(fft2(Reconstruction_field_2)).*H)); E=Reconstruction_field_3./abs(Reconstruction_field_3); % end end Ifigure=figure; imagesc(abs(Reconstruction_field_1(34:630,70:680)));colormap(gray) % % Ifigure=figure; imagesc(Intensity,[min(Intensity(:)),0.8*max(Intensity(:))]);colormap(gray) axis image; axis off; %%画出强度图 hold on; title (['背景全息图数值重建像的强度图（重建距离d= ' num2str(Reconstruction_distance*1000) 'mm）']) hold off % Path=[Folder_result,'迭代15-Rec_distance_',num2str(Reconstruction_distance),'.jpg']; % print(Ifigure,'-djpeg',Path); %%以截屏的方式存储上述强度图，路径为Path Phase_wrap=angle(Reconstruction_field_1(34:630,70:680)); %重建物光波的相位图 %2D 包裹相位图 Phase_wrap_figure=figure; % Phase_wrap=Phase_wrap(1:2:end, 1:2:end); % imshow(Phase_wrap,[ ]); colormap(jet);colorbar imagesc(Phase_wrap); colormap(gray);colorbar;axis off set(gca, 'ydir', 'normal'); set(gca,'FontSize',18,'FontName','Times Nes Roman','FontWeight','bold') % Path=[Folder_result,'迭代15-Wrapped Phase image','.jpg']; % print(Phase_wrap_figure,'-djpeg',Path); %%以截屏的方式存储上述强度图，路径为Path