function [u, v,im1_hr,im1_hr_im2_hr,im1_hr_im2_orig] = optic_flow_brox(img1, img2)
% im1_hr 是王添加导出参数%%
%后面两个参数是我自己添加的，为了能在workspace中显示
% Copyright (c)2007 Visesh Chari <visesh [at] research.iiit.net>
% Centre for Visual Information Technology
% International Institute of Information Technology
% http://cvit.iiit.ac.in/
% http://students.iiit.ac.in/~ukvisesh
%
% This program is free software; you can redistribute it and/or modify
% it under the terms of the GNU General Public License as published by
% the Free Software Foundation; either version 2 of the License, or
% (at your option) any later version.
%
% The Software is provided "as is", without warranty of any kind.
alpha = 30.0 ; % Global smoothness variable. 平滑项的权重30.0
gamma = 80.0 ; % Global weight for derivatives.梯度不变的权重80.0
[ht, wt, dt] = size( img1 ) ;
num_levels = 40 ; % default40。 for face; Number of pyramid levels. Can also be calculated automatically
im1_hr = gaussianRescaling( img1, power( 0.95, num_levels ) ) ; % Rescaling images 先高斯平滑，再缩放图片
im2_hr = gaussianRescaling( img2, power( 0.95, num_levels ) ) ; % to the top of the
% laplacian pyramid.
% % % u = zeros(size(rgb2gray(im1_hr))); % Initialization.
% % % v = zeros(size(rgb2gray(im1_hr)));
u = zeros(size((im1_hr))); % Initialization.
v = zeros(size((im1_hr)));
miter = 0;
for i = num_levels-1 : -1 : 0
% % % I1 = rgb2gray(im1_hr) ;
% % % I2 = rgb2gray(im2_hr) ;
I1 = (im1_hr) ;
I2 = (im2_hr) ;
% Computing derivatives.
[Ikx Iky] = imgGrad( I2 ) ;
[Ikx2 Iky2] = imgGrad( I1 ) ;
Ikz = double(I2) - double(I1) ;
Ixz = double(Ikx) - double(Ikx2) ;
Iyz = double(Iky) - double(Iky2) ;
% Calling the processing for a particular resolution.
% Last two arguments are the outer and inner iterations, respectively.
% 1.8 is the omega value for the SOR iteration.
[du, dv] = resolutionProcess_brox( Ikz, Ikx, Iky, Ixz, Iyz, alpha, gamma, 1.8, u, v, 3, 500 ) ;
%{
%%%%%%%%%%保存矩阵为.mat格式%%%%%%%%%%%%%
miter = miter + 1;
path = sprintf('%su%d.mat','C:\Users\asus\Desktop\光流图\',miter);
save (path,'u');
path = sprintf('%sv%d.mat','C:\Users\asus\Desktop\光流图\',miter);
save (path,'v');
%%%%%%%%%%保存矩阵为.mat格式%%%%%%%%%%%%%
%}
% Adding up the optical flow.
u = u + du ;
v = v + dv ;
im1_hr = gaussianRescaling( img1, power( 0.95, i ) ) ;
im2_hr = gaussianRescaling( img2, power( 0.95, i ) ) ;
im2_orig = im2_hr ; % Original image without warping for comparison.
% Resize optical flow to current resolution.
u = imresize( u, [size(im1_hr, 1), size(im1_hr, 2)], 'bilinear' ) ;
v = imresize( v, [size(im1_hr, 1), size(im1_hr, 2)], 'bilinear' ) ;
im2_hr = uint8( mywarp_rgb( double( im2_hr ), u, v ) ) ; % taking im1_hr closer to im2_hr. warp校正后的图2
%im2_hr = mywarp_rgb( double( im2_hr ), u, v ) ;%如果把uint8删掉了，要报错
%%%%%%%%%%我自己加的%%%%%%%%%%%
%figure(3);
%subplot(3,1,1);imshow(im1_hr) ;
%subplot(3,1,2);imshow(im2_orig) ;
%subplot(3,1,3);imshow(im2_hr) ;
%%%%%%%%%%我自己加的%%%%%%%%%%%
% Displaying relevant figures.
%王修改mark
% figure(2);
% subplot(3, 3, 1); imshow(im1_hr,[]) ;title('图一');%%%加了【】
% %subplot(3, 3, 2); imshow(rgb2gray(im2_hr)) ;
% %subplot(3, 3, 3); imshow(rgb2gray(im2_orig)) ;
% subplot(3, 3, 2); imshow((im2_hr)) ;title('warp校正后的图');
% subplot(3, 3, 3); imshow((im2_orig),[]) ;title('图二');%%%加了【】
% subplot(3, 3, 4); imshow(uint8(double(im1_hr)-double(im2_hr))) ;title('图一和校正后的图求差');
% subplot(3, 3, 5); imshow(uint8(double(im1_hr)-double(im2_orig))) ;title('图一和图二求差');
% subplot(3, 3, 6); imshow(uint8(double(im2_hr)-double(im2_orig))) ;title('图二和校正后的图求差');
% subplot(3, 3, 7); imagesc(u) ;%上面直接用uint8直接转换，不够严谨吧，负的值直接置为0了
% subplot(3, 3, 8); imagesc(v) ;
% subplot(3, 3, 9); imagesc(sqrt(u.^2 + v.^2)) ;
%王修改mark
% pause;
im1_hr_im2_hr = double(im1_hr)-double(im2_hr);
im1_hr_im2_orig = double(im1_hr)-double(im2_orig);
%%%%%%%%序号的不同是因为上面做了金字塔的分层%%%%%%%%
% path = sprintf('%sim1_hr_im2_hr%d.mat','C:\Users\asus\Desktop\光流图\',miter);
% save (path,'im1_hr_im2_hr');
% path = sprintf('%sim1_hr_im2_orig%d.mat','C:\Users\asus\Desktop\光流图\',miter);
% save (path,'im1_hr_im2_orig');
%
% path = sprintf('%sim1_hr%d.mat','C:\Users\asus\Desktop\光流图\',miter);
% save (path,'im1_hr');%原图
% path = sprintf('%sim2_hr%d.mat','C:\Users\asus\Desktop\光流图\',miter);
% save (path,'im2_hr');%校正后图
% path = sprintf('%sim2_orig%d.mat','C:\Users\asus\Desktop\光流图\',miter);
% save (path,'im2_orig');%校前的后图
%
% ee = uint8(im1_hr_im2_hr)
% figure(3);imshow(rgb2gray(uint8(double(im1_hr)-double(im2_hr)))) ;title('图一和校正后的图求差');
% figure(4);imshow(rgb2gray(uint8(double(im1_hr)-double(im2_orig)))) ;title('图一和图二求差');
% figure(5); imagesc(u) ;
% figure(6); imagesc(v) ;
% figure(7); imagesc(sqrt(u.^2 + v.^2)) ;
end
评论0