Multi_frequency_heterodyne_phase_solution.rar

利用matlab实现三频外差法解相位，通过生成三组相应频率的四步相移灰度图，共12张原始图片，计算出其中的每一坐标的相位。

%模拟波浪图像
%clear all
%tic;
%L=2000;d=200;f0=1/16;w=2*pi*f0;
%x=1:512;y=x;
%[xy]=meshgrid（xy）;
%z=3*peaks（512）;
%figure（1）;mesh（z）;axis on
%xlabel（‘x（pixel）‘）;ylabel（‘y（pixel）‘）;zlabel（‘h（mm）‘）

%基准光栅
%J0=128+127*cos（w*x）;
%J0=mat2gray（J0）;
%figure（2）;imshow（J0）;axis on

%变形光栅
%J1=128+127*cos（w*（x+z*d./（L-z）））;
%J1=mat2gray（J1）;
%J1=imnoise（J1‘gaussian‘00.001）;
%figure（3）;imshow（J1）;axis on

%参考李中伟博士论文实现相移法加三频外差解相位仿真
% 程序开始

clc;

close all;

clear;

 

% 图片的初始化

width = 1024;  

heigth = 768;

 

% 三频率

% 这个可以参见李中伟的博士论文

freq = [70 64 59];  %像素单位为个数，可以看做频率；正弦函数为周期含义
f12 =abs（（freq（1）*freq（2））/（freq（1）-freq（2）））;
f23 = abs（（freq（2）*freq（3））/（freq（2）-freq（3）））;
f13 = abs（（freq（1）*freq（3））/（freq（1）-freq（3）））;
f123 = abs（（f12*f23）/f12-f23）;

% 利用分块矩阵C存储3组共计12张图

% 三种频率，四组相位

C = cell（34）;  

for i=1:3

    for j=1:4

        C{ij} = zeros（heigthwidth）;

    end

end

 

% 利用余弦函数计算12张图的灰度值

% 图像的生成

% 三种频率，四组相位

for i = 1:3 % 对应三种不同的频率

    for  j = 0:3 % 对应四种相位

        for k = 1:width

            C{ij+1}（:k） = 128+127*sin（2*pi*k*freq（i）/width+j*pi/2）;

        end

    end

end

 

% 对灰度值进行归一化处理

for i = 1:3

    for j = 1:4

        C{ij} = mat2gray（C{ij}）;

    end

end

 

% 显示12张图

 for i = 1:3

   for j = 1:4

       n = 4*（i-1）+j;

         h = figure（n）;

        imshow（C{ij}）;

     end

end

 

% 初始化三组处理后的图片灰度矩阵

% phi也是分块矩阵

% 存储相位主值图像

phi = cell（31）;

for i = 1:3

    phi{i1} = zeros（heigthwidth）;

end

 

% 求取相位差

% 计算每种频率对应的相位主值

% 输出三种频率的相位主值，用于相差计算

for i = 1:3 % 对于3组中的每一组图片，每一组相同频率的有四张图片

     I1 = C{i1};

     I2 = C{i2};

     I3 = C{i3};

     I4 = C{i4};

      for g = 1:heigth

          for k = 1:width          

            if     I4（gk）==I2（gk）&&I1（gk）>I3（gk） %四个特殊位置

                    phi{i1}（gk）=0;

            elseif I4（gk）==I2（gk）&&I1（gk）
                    phi{i1}（gk）=pi; 

            elseif I1（gk）==I3（gk）&&I4（gk）>I2（gk） %四个特殊位置

                    phi{i1}（gk）=pi/2;

            elseif I1（gk）==I3（gk）&&I4（gk）
                    phi{i1}（gk）=3*pi/2;

            elseif I1（gk）
                    phi{i1}（gk）=atan（（I4（gk）-I2（gk））./（I1（gk）-I3（gk）））+pi;

            elseif I1（gk）>I3（gk）&&I4（gk）>I2（gk） %第一象限

                    phi{i1}（gk）=atan（（I4（gk）-I2（gk））./（I1（gk）-I3（gk）））;

            elseif I1（gk）>I3（gk）&&I4（gk）
                    phi{i1}（gk）=atan（（I4（gk）-I2（gk））./（I1（gk）-I3（gk）））+2*pi;  

            end

          end             

     end

end

  

% 计算相差

% 保存矩阵，用于多频相差的计算

PH1 = phi{11};   %频率1

PH2 = phi{21};   %频率2

PH3 = phi{31};   %频率3

 

% 初始化相差变量

% 多频相差

PH12 = zeros（heigthwidth）;

PH23 = zeros（heigthwidth）;

PH13 = zeros（heigthwidth）;

PH123 = zeros（

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       6449  2019-09-29 13:03  Multi_frequency_heterodyne_phase_solution.m

----------- ---------  ---------- -----  ----

                 6449                    1