相位误差补偿表LUT-创建

LUT

1.相位误差补偿表：就是建立一个查找表对由gamma效应造成的投影仪光栅非正弦性进行一个补偿
2.相关资料可以看李中伟博士的论文 或者其他资料
3.步骤
3.1 获得理想相位图 （用多步相移代替理想相位）
3.2 获得实际相位图
3.3 求取误差补偿表 并将一个周期的区域划分为256个周期
3.4 对实际相位图进行补偿
4. 这里只是一个最基础的流程 存在很多问题 可以在此基础上进行优化。

%%  第一部分 读取图片并对图片进行预处理

W = 140;
H = 1000;
C1 = cell(1,8);  % 存储原始的八张图片
for i = 1:8
    str = ['D:\halcon-image\lut-test\tttt\11\Pic_',num2str(i+33),'.bmp'];   % 拼接地址
    img = imread(str);  % 读取图片
    img = img(1:W,1:H);
    img = mat2gray(img);  % 归一化
    % 对图像进行平滑处理 中值滤波
    img=medfilt2(img);
    C1{1,i} =  img;
end
% figure(1);
% imshow(C1{1,1});
C2 = cell(1,4);  % 存储四步相移图片
j = 1;
for i = 1:4
  C2{1,i} = C1{1,j} ;
  j = j+2;
end

%% 第二部分 用八步相移  代替理想相位图

% 首先测试  八步相移解包裹相位
% 1.计算包裹相位
pha_st = atan2(-sqrt(2)*(C1{1,2} - C1{1,6}) + sqrt(2)*(C1{1,8} - C1{1,4}) + sqrt(2)*(C1{1,7} - C1{1,3}) ...
, sqrt(2)*(C1{1,2} - C1{1,6}) + sqrt(2)*(C1{1,8} - C1{1,4}) + sqrt(2)*(C1{1,1} - C1{1,5}));
% 映射到0-2pi
pha_st(pha_st < 0) = pha_st(pha_st < 0) + 2*pi;
% pha_st = mat2gray(pha_st);
figure(2);
mesh(pha_st);
title("八步相移包裹相位");


%% 第三部分  未补偿前的解包裹相位
% 四步相移的解包裹相位
pha_t = atan2(C2{1,4}-C2{1,2}, C1{1,1}- C1{1,3});
% 映射到0-2pi
pha_t(pha_t < 0) = pha_t(pha_t < 0) + 2*pi;
% pha_t = mat2gray(pha_t);
figure(3);
mesh(pha_t);
title("四步相移包裹相位");

%% 计算理想相位和实际相位之间的方差
pre_re = sum(sum((pha_st - pha_t).^2));
%% 第四部分  计算相位差  
L = pha_st - pha_t;   % 相位误差
L = L(100,:);
figure(5);
plot(L);
title("第五百行的相位-误差");

% 90-174列为一个周期
figure(6);
hold on;plot(pha_st(100,90:174));plot(pha_t(100,90:174));hold off;
% figure(7);
% plot(L(119:204));
% 理想相位
yt = pha_st(100,90:174);
% 绘制lut
x = pha_t(100,90:174);
y = L(90:174); % 相位误差


% 对一个周期进行划分
[fitresult, gof] = createFit2(x,y);
       a1 =      fitresult.a1;
       b1 =      fitresult.b1;
       c1 =      fitresult.c1;
       a2 =      fitresult.a2;
       b2 =      fitresult.b2;
       c2 =      fitresult.c2;
       a3 =      fitresult.a3;
       b3 =      fitresult.b3;
       c3 =      fitresult.c3;
       a4 =      fitresult.a4;
       b4 =      fitresult.b4;
       c4 =      fitresult.c4;
       a5 =      fitresult.a5;
       b5 =      fitresult.b5;
       c5 =      fitresult.c5;
       a6 =      fitresult.a6;
       b6 =      fitresult.b6;
       c6 =      fitresult.c6;
       a7 =      fitresult.a7;
       b7 =      fitresult.b7;
       c7 =      fitresult.c7;

% 划分区域  绘制函数
X_temp = [0:(2*pi)/256:2*pi];
Y = zeros(1,256);
for  i = 1:256
    z = X_temp(i);
    Y(i) = a1*sin(b1*z+c1) + a2*sin(b2*z+c2) + a3*sin(b3*z+c3) + a4*sin(b4*z+c4) + ...
    a5*sin(b5*z+c5) + a6*sin(b6*z+c6) + a7*sin(b7*z+c7);
end
figure(20);
plot(Y);

% 对数据进行补偿  存在的问题 ？？？ 将包裹相位都映射到了0-1
all = pha_t(:,:);
figure(10);
mesh(all);
for i = 1:W
    for j = 1:H
        k = round(256*all(i,j)/(2*pi) + 0.5 );      % 计算某一点的相位补偿
        all(i,j) = all(i,j) + Y(k);
    end
end
figure(11);
hold on; plot(all(100,:));plot(pha_st(100,:)); hold off;

%% 计算补偿后的包裹相位 误差
dl_re = sum(sum((pha_st - all).^2));
figure(12);
mesh(all); 
% for i = 1:W
%      hang(i,:) = unwrap(hang(i,:));
%  end
%  for j = 1:H
%      hang(:,j) = unwrap(hang(:,j));
%  end

function [fitresult, gof] = createFit2(x, y)
%CREATEFIT2(X,Y)
%  Create a fit.
%
%  Data for 'untitled fit 1' fit:
%      X Input : x
%      Y Output: y
%  Output:
%      fitresult : a fit object representing the fit.
%      gof : structure with goodness-of fit info.
%
%  另请参阅 FIT, CFIT, SFIT.

%  由 MATLAB 于 01-May-2022 01:42:53 自动生成


%% Fit: 'untitled fit 1'.
[xData, yData] = prepareCurveData( x, y );

% Set up fittype and options.
ft = fittype( 'sin7' );
opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
opts.Display = 'Off';
opts.Lower = [-Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf -Inf 0 -Inf];
opts.StartPoint = [0.545195364682761 0.509584466101639 0.183654103781012 0.481755378886401 2.03833786440656 -2.33206741722098 0.292919653173777 3.05750679660983 1.19340952787291 0.163023007114948 1.01916893220328 2.01497389987387 0.126085400824674 10.1916893220328 -2.54982977296841 0.126084344183955 7.13418252542295 2.6876591999933 0.0641616429693231 4.07667572881311 1.96605815925331];

% Fit model to data.
[fitresult, gof] = fit( xData, yData, ft, opts );

% Plot fit with data.
figure( 'Name', 'untitled fit 1' );
h = plot( fitresult, xData, yData );
legend( h, 'y vs. x', 'untitled fit 1', 'Location', 'NorthEast', 'Interpreter', 'none' );
% Label axes
xlabel( 'x', 'Interpreter', 'none' );
ylabel( 'y', 'Interpreter', 'none' );
grid on