相移法+多频外差相位解包

最新推荐文章于 2023-09-14 15:22:17 发布

马少爷

最新推荐文章于 2023-09-14 15:22:17 发布

阅读量2.1k

点赞数

分类专栏：工程三维视觉文章标签：计算机视觉三维视觉

本文链接：https://blog.csdn.net/qq_27353621/article/details/124236395

版权

三维视觉同时被 2 个专栏收录

72 篇文章 90 订阅

订阅专栏

工程

55 篇文章 6 订阅

订阅专栏

本文介绍了三频外差光栅的生成过程，使用Matlab代码详细展示了如何创建不同频率的光栅，并通过三频外差方法进行相位解包，从而获取绝对相位信息。高频条纹用于获取包裹相位，中频和低频条纹则用于确定阶次，最终实现相位展开。

摘要由CSDN通过智能技术生成

三频外差光栅生成
三组频率分别为76、57、37，76为高频条纹的频率，57为中频条纹的频率，37为低频条纹的频率。其中高频的条纹设置为4步相移，中频和低频的分别设置为3步相移，具体的Matlab生成代码如下：


R = 1140;  C = 912;   % 光栅尺寸1140*912

P1 = 1140/76;        % 频率76的光栅宽度
P2 = 1140/57;        % 频率57的光栅宽度
P3 = 1140/37;        % 频率37的光栅宽度

a = 126;           % 背景光强
b = 126;           % 调制光强

I0 = zeros(1140,912);   %定义图像尺寸大小
H1 = I0; H2 = I0; H3 = I0; H4 = I0;
H5 = I0; H6 = I0; H7 = I0; 
H8 = I0; H9 = I0; H10 = I0;

for i = 1:1140
    for j = 1:912
        H1(i,j) = a+b*cos(2*pi*(j-1)/P1);                %4步相移余弦光栅 --高频
        H2(i,j) = a+b*cos(2*pi*(j-1)/P1+1*2*pi/4);
        H3(i,j) = a+b*cos(2*pi*(j-1)/P1+2*2*pi/4);
        H4(i,j) = a+b*cos(2*pi*(j-1)/P1+3*2*pi/4);
        
        H5(i,j) = a+b*cos(2*pi*(j-1)/P2-2*pi/3);         %3步相移余弦光栅 --中频
        H6(i,j) = a+b*cos(2*pi*(j-1)/P2);
        H7(i,j) = a+b*cos(2*pi*(j-1)/P2+2*pi/3);
        
        H8(i,j) = a+b*cos(2*pi*(j-1)/P3-2*pi/3);         %3步相移余弦光栅 --低频
        H9(i,j) = a+b*cos(2*pi*(j-1)/P3);
        H10(i,j) = a+b*cos(2*pi*(j-1)/P3+2*pi/3);
 
    end
end 

%% ----------------------【光栅写出】--------------------------

for k = 1:10
    eval(['H','=H',num2str(k),';']);
    H = uint8(H);
    img_name = ['gratings/',num2str(k),'.bmp'];
    imwrite(H,Himg_name);
end

三频外差相位解包

每个频率的光栅求出的相位都是包裹相位，其形式如同下式的ϕ ( x , y ) 所示，这样的相位图不能直接用于双目匹配，为了获取绝对相位Φ ( x , y ) ，我们需要获取一个条纹阶次k ( x , y )，然后通过下式将包裹相位展开为绝对相位：

Φ ( x , y ) = ϕ ( x , y ) + 2 π × k ( x , y )

在这里插入图片描述
三频外差的高频条纹就是为了获取包裹相位，这里有人可能会问，为什么其他频率的不能当作包裹相位ϕ ( x , y ) 呢，因为条纹频率越高，信噪比就越大，因此高频的相位精度更好。中频和低频的光栅就是为了获取k ( x , y ) 的，其具体获取方式如下（matlab代码）：

% ----------------变量预定义
[R C] = size(H1);
phi1 = zeros(R,C);  phi2 = phi1;  phi3 = phi1;
k1 = phi1;  k2 =phi1;  k = phi1;
phi_eq2 = phi1; phi_eq1 = phi1;  phi_eq = phi1;
phase = phi1;
% ---------------高频包裹相位phi1
phi1 = atan2(H4-H2,H1-H3);
phi1(phi1<0) = phi1(phi1<0)+2*pi;

% ---------------中频包裹相位phi2
phi2 = atan2(sqrt(3)*(H5-H7),2*H6-H5-H7);
phi2(phi2<0) = phi2(phi2<0)+2*pi;

% ---------------低频包裹相位phi3
phi3 = atan2(sqrt(3)*(H8-H10),2*H9-H8-H10);
phi3(phi3<0) = phi3(phi3<0)+2*pi;

% ---------------高频-中频外差相位phi_eq1
phi_eq1 = phi1-phi2;
phi_eq1(phi_eq1<0) = phi_eq1(phi_eq1<0)+2*pi;
% ---------------中频-低频外差相位phi_eq2
phi_eq2 = phi2-phi3;
phi_eq2(phi_eq2<0) = phi_eq2(phi_eq2<0)+2*pi;
% ---------------高频-中频-低频外差相位phi_eq（三个频率的综合外差）
phi_eq =  phi_eq2-phi_eq1;
phi_eq(phi_eq<0) = phi_eq(phi_eq<0)+2*pi;

L_12 = P1*P2/(P2-P1);          % 外差相位phi_eq1波长
L_23 = P2*P3/(P3-P2);          % 外差相位phi_eq2波长
L_eq = L_23*L_12/(L_12-L_23);  % 外差相位phi_eq波长

k1 = round(((L_12/P1)*phi_eq1-phi1)/(2*pi));  
k2 = round(((L_eq/L_12)*phi_eq-phi_eq1)/(2*pi));
k = k2*(L_12./P1) + k1;

代码完整版如下：


R = 1140;  C = 912;   % 光栅尺寸1140*912

P1 = 1140/76;        % 频率76的光栅宽度
P2 = 1140/57;        % 频率57的光栅宽度
P3 = 1140/37;        % 频率37的光栅宽度

a = 126;           % 背景光强
b = 126;           % 调制光强

I0 = zeros(1140,912);   %定义图像尺寸大小
H1 = I0; H2 = I0; H3 = I0; H4 = I0;
H5 = I0; H6 = I0; H7 = I0; 
H8 = I0; H9 = I0; H10 = I0;

for i = 1:1140
    for j = 1:912
        H1(i,j) = a+b*cos(2*pi*(j-1)/P1);                %4步相移余弦光栅 --高频
        H2(i,j) = a+b*cos(2*pi*(j-1)/P1+1*2*pi/4);
        H3(i,j) = a+b*cos(2*pi*(j-1)/P1+2*2*pi/4);
        H4(i,j) = a+b*cos(2*pi*(j-1)/P1+3*2*pi/4);
        
        H5(i,j) = a+b*cos(2*pi*(j-1)/P2-2*pi/3);         %3步相移余弦光栅 --中频
        H6(i,j) = a+b*cos(2*pi*(j-1)/P2);
        H7(i,j) = a+b*cos(2*pi*(j-1)/P2+2*pi/3);
        
        H8(i,j) = a+b*cos(2*pi*(j-1)/P3-2*pi/3);         %3步相移余弦光栅 --低频
        H9(i,j) = a+b*cos(2*pi*(j-1)/P3);
        H10(i,j) = a+b*cos(2*pi*(j-1)/P3+2*pi/3);
 
    end
end 

%% ----------------------【光栅写出】--------------------------

for k = 1:10
    eval(['H','=H',num2str(k),';']);
    H = uint8(H);
    img_name = ['',num2str(k),'.bmp'];
    imwrite(H,img_name);
end

% ----------------变量预定义
[R C] = size(H1);
phi1 = zeros(R,C);  phi2 = phi1;  phi3 = phi1;
k1 = phi1;  k2 =phi1;  k = phi1;
phi_eq2 = phi1; phi_eq1 = phi1;  phi_eq = phi1;
phase = phi1;
% ---------------高频包裹相位phi1
phi1 = atan2(H4-H2,H1-H3);
phi1(phi1<0) = phi1(phi1<0)+2*pi;

% ---------------中频包裹相位phi2
phi2 = atan2(sqrt(3)*(H5-H7),2*H6-H5-H7);
phi2(phi2<0) = phi2(phi2<0)+2*pi;

% ---------------低频包裹相位phi3
phi3 = atan2(sqrt(3)*(H8-H10),2*H9-H8-H10);
phi3(phi3<0) = phi3(phi3<0)+2*pi;

% ---------------高频-中频外差相位phi_eq1
phi_eq1 = phi1-phi2;
phi_eq1(phi_eq1<0) = phi_eq1(phi_eq1<0)+2*pi;
% ---------------中频-低频外差相位phi_eq2
phi_eq2 = phi2-phi3;
phi_eq2(phi_eq2<0) = phi_eq2(phi_eq2<0)+2*pi;
% ---------------高频-中频-低频外差相位phi_eq（三个频率的综合外差）
phi_eq =  phi_eq2-phi_eq1;
phi_eq(phi_eq<0) = phi_eq(phi_eq<0)+2*pi;

L_12 = P1*P2/(P2-P1);          % 外差相位phi_eq1波长
L_23 = P2*P3/(P3-P2);          % 外差相位phi_eq2波长
L_eq = L_23*L_12/(L_12-L_23);  % 外差相位phi_eq波长

k1 = round(((L_12/P1)*phi_eq1-phi1)/(2*pi));  
k2 = round(((L_eq/L_12)*phi_eq-phi_eq1)/(2*pi));
k = k2*(L_12./P1) + k1;

% -----------------------相位展开，得绝对相位phase---------------------------
phase = phi1 + k*2*pi;
figure,imshow(mat2gray(phase));
title('绝对相位phase');

绝对相位为：
在这里插入图片描述