在粒子滤波器框架下，选取区域协方差矩阵特征对运动独立的目标进行描述.

%***********************************************************************% % tracking_cm.m % % Description: This is a test programme using region covariance matrix % % for tracking the man under the complex environment. % % Author:Wu_Pian_Pian % % Date:2010-4-10 % %***********************************************************************% clear; % read the avi file and initialize avi_info = aviinfo('22.avi'); % read the avi file into the workspace frame_number = avi_info.NumFrames; % get the number of the frames avi_data = aviread('22'); % get the information of the avi file above % mov=aviread('22.avi'); %读入avi % fnum=size(mov,2); %读取电影的祯数，mov为1*temp % for i=1:fnum % % if (size(mov(i).cdata,3)==1) % % im=mov(i).cdata; % % else % im=mov(i).cdata; % im=rgb2gray(im); % end im=edge(im,'sobel',[],'both'); imshow(im); %end tic; % load the covariance matrix of the two target model load('C_guanche'); rectangle_a = round(20+2*randn(1,70)); rectangle_b = round(62+2*randn(1,70)); rectangle_x(1) = 188; % initialize the x-coordinate of the rectangle' center (left) rectangle_y(1) = 161; % initialize the y-coordinate of the rectangle' center %****************************************************** % Using particle filtering for tracking the target %****************************************************** for frame_th = 231:300 image_data = avi_data(frame_th); % read the current frame image_cdata = image_data.cdata; % convert the information of the current frame to a matrix imshow(frame_th) ; aviinfo(aviFile) %*************************************************** % Prediction ( Particle computation ) %*************************************************** % initiate N = 100; % Particle Number % represent the dynamic systems and predict the state of position if frame_th == 231 rectangle_x_temp = rectangle_x(frame_th-230) + 0.25*randn(1,N); rectangle_y_temp = rectangle_y(frame_th-230) + 0.01*randn(1,N); elseif frame_th == 232 rectangle_x_temp = rectangle_x(frame_th-231) + 0.25*randn(1,N); rectangle_y_temp = rectangle_y(frame_th-231) + 0.01*randn(1,N); elseif frame_th < 270 rectangle_x_temp = round(2*rectangle_x(frame_th-231) - rectangle_x(frame_th-232) + 2*randn(1,N)); rectangle_y_temp = round(2*rectangle_y(frame_th-231) - rectangle_y(frame_th-232) + 4*randn(1,N)); else rectangle_x_temp = round(2*rectangle_x(frame_th-231) - rectangle_x(frame_th-232) + 2.5*randn(1,N)); rectangle_y_temp = round(2*rectangle_y(frame_th-231) - rectangle_y(frame_th-232) + 6*randn(1,N)); end %*************************************** % Importance sampling %*************************************** % compute the covariance matrix of every region var_Gaussian = 0.5; %%%%%？？？ weight_sum = 0; % initialize the weight to be zero rectangle_x_temp = round(rectangle_x_temp); % 取整，但只限于这步 rectangle_y_temp = round(rectangle_y_temp); for particle_th = 1:N frame_covariance(:,:,particle_th) = covariance(image_cdata,... rectangle_x_temp(1,particle_th)-rectangle_a(frame_th-230),... rectangle_x_temp(1,particle_th)+rectangle_a(frame_th-230),... rectangle_y_temp(1,particle_th)-rectangle_b(frame_th-230),... rectangle_y_temp(1,particle_th)+rectangle_b(frame_th-230)); % Weight computation(协方差矩阵距离度量） distance_covariance(particle_th) = distance(C_guanche_third(:,:,frame_th-230),... frame_covariance(:,:,particle_th),5); invert_difference(particle_th) = (1/((sqrt(2*pi))*var_Gaussian))*... exp(-distance_covariance(particle_th)^2/(2*var_Gaussian^2)); weight_sum = weight_sum + invert_difference(particle_th); end % Weight normalization importance_weight = invert_difference./weight_sum; %********************************************************* % Resampling ( update and propagate the particles ) %********************************************************* % Classical Resampling Algorithm u = rand(N+1,1); t = -log(u); x_temp = rectangle_x_temp'; y_temp = rectangle_y_temp'; x_sample = 192.*ones(size(x_temp)); y_sample = 144.*ones(size(y_temp)); T = cumsum(t); Q = cumsum(importance_weight'); i = 1; j = 1; while j <= N, if (Q(j,1)*T(N,1)) > T(i,1) x_sample(i,1) = x_temp(j,1); y_sample(i,1) = y_temp(j,1); i = i+1; else j = j+1; end; end; %************************************************************** % Computation of the estimates ( find the optimal consequence ) %************************************************************** rectangle_x(frame_th-230) = mean(x_sample); rectangle_y(frame_th-230) = mean(y_sample); [row] = find(invert_difference == max((max(invert_difference)))); C_third(:,:,frame_th-230) = frame_covariance(:,:,min(row)); % update the template x_vec = rectangle_x_temp(min(row)); y_vec = rectangle_y_temp(min(row)); xx(frame_th-230) = x_vec(1); % the optimal x-coordinate of the rectangle' center yy(frame_th-230) = y_vec(1); % the optimal y-coordinate of the rectangle' center fprintf('this frame is OK: %d\n',frame_th); toc; end %********************************************************** % generate the output avi file %********************************************************** fprintf('generate the output avi file'); tic; fig = figure(1); aviobj = avifile('tracking_cm.avi'); for frame_th = 231:300 image_data = avi_data(frame_th); % read the current frame image_cdata = image_data.cdata; % convert the information of the current frame to a matrix point1_x = xx(frame_th-230) + rectangle_a(frame_th-230); % point1 is the right upper point of the rectangle point2_x = xx(frame_th-230) - rectangle_a(frame_th-230); % point2 is the left upper point of the rectangle point3_x = xx(frame_th-230) - rectangle_a(frame_th-230); % point3 is the left bottom point of the rectangle point4_x = xx(frame_th-230) + rectangle_a(frame_th-230); % point4 is the right bottom point of the rectangle point1_y = yy(frame_th-230) + rectangle_b(frame_th-230); % point2_y = yy(frame_th-230) + rectangle_b(frame_th-230); % point3_y = yy(frame_th-230) - rectangle_b(frame_th-230); % point4_y = yy(frame_th-230) - rectangle_b(frame_th-230); % imshow(image_cdata); % show the current frame hold on; line([point1_x point2_x],[point1_y point2_y],'Color','y','LineWidth',2); line([point1_x point4_x],[point1_y point4_y],'Color','y','LineWidth',2); line([point2_x point3_x],[point2_y point3_y],'Color','y','LineWidth',2); line([point3_x point4_x],[point3_y point4_y],'Color','y','LineWidth',2); hold off; FFF = getframe; aviobj = addframe(aviobj,FFF); end aviobj = close(aviobj); toc;