Mean shift The algorithm is a semi-automatic tracking method in the start-up tracking frame by manually determining the search window to select the motion target calculation kernel function weighted by the histogram distribution of the search window with the same method to calculate the histogram distribution of the corresponding window of the current frame with the principle of maximum similarity of two distributions to move the target's real position in the direction of maximum density increase 。
Weighted histogram
The traditional histogram only counts the number of pixels falling into the histogram interval, while the weighted histogram further considers the distance between the pixel and the target center, and the pixels away from the target center contribute less to the histogram.
The idea of weighted histogram with spatial position information is that when calculating histograms, each point is given a certain weight, and the weight is adjusted according to its distance from the center point, which can be modulated by Gauss kernel function or Epanechnikov kernal. The commonly used kernel density functions are epanechnikovkernal:
In the statistics of the distribution of color levels, for each pixel, according to its distance from the center of the window, the size of a given it a certain weight. The closer it is to the center point, the greater its weight, and the smaller the inverse. This also increases the robustness of the object description, as the closer to the perimeter of the object, the less the point is, in the case of occlusion or complex background effects.
T_w and T_h respectively for the width and height of the detection window H = POW ((double) t_w)/2,2) + POW (((double) t_h)/2,2); Bandwidth//initialization weights matrix and target histogram for (i = 0;i < t_w*t_h;i++) {m_wei[i] = 0.0; } for (i=0;i<4096;i++) {hist1[i] = 0.0; }//Generate Weights matrix for (i = 0;i < T_h; i++) {for (j = 0;j < T_w; J + +) {dist = Pow (i-(double) t_h/2,2) + POW (J-(double) t_w/2,2); M_wei[i * T_w + j] = 1-dist/h; C + = M_wei[i * T_w + j]; }}//Compute the target weights straight for (i = t_y;i < t_y + T_h; i++) {for (j = t_x;j < t_x + T_w; j + +) {//rgb color space quantization to 16*16*16 bins Q_r = ((U_char) current->imagedata[i * current->widthste P + J * 3 + 2])/16; Q_g = ((U_char) current->imagedata[i * current->widthstep + J * 3 + 1])/16; Q_b = ((U_char) current->imagedata[i * current->widthstep + J * 3 + 0])/16; Q_temp = Q_r * + q_g * + q_b; Hist1[q_temp] = Hist1[q_temp] + m_wei[(i-t_y) * t_w + (j-t_x)]; }}//Normalized histogram for (i=0;i<4096;i++) {hist1[i] = hist1[i]/C; }
color model Similarity measurementTarget model of the initial frameA feature space with a value of pixel color value is divided into multiple eigenvalues at a certain interval of color value then the probability of the first U eigenvalue in the initial frame containing the target's search window is
The x0 is the center pixel coordinate of the search window (n pixels), and Xi is the coordinates of the first pixel; K (| | x| | 2) is a kernel function, which is used to weighting the probability density of color distribution, why should it be weighted? by using the kernel function estimation method, we can gradually converge to any density function in the case of sufficient sampling, that is, it is possible to estimate the density of the data subjected to any distribution. This can ensure the convergence of the meanshift algorithm, in order to use the Meanshift vector to iterate. H represents the bandwidth of the kernel function, which is generally equal to half the width of the window; function B and δ determine whether the color value at XI is a eigenvalues u, and if Xi belongs to the eigenvalues U,the value of δ is 1, otherwise the value of Δ takes 0;C is a normalized constant coefficient that makes all eigenvalue probabilities and is 1. c=1/The calculation of weights is regulated by the kernel function K, which is also the introduction of the kernel function K, so that the histogram matching can be combined with the mean drift algorithm.
the model of the current frameThe probability of calculating the eigenvalues U of the Search window in the current frame (nth frame) is
y0 is the center pixel coordinate of the current frame Search window, Xi is the coordinates of the I-pixel, and ch corresponds to the C in the above formula.
Bhattacharrya distance of similarity function the PAP coefficient is used to measure the similarity between the target histogram and the candidate histogram, the larger the PAP coefficient, the more similar.
The similarity function describes the initial frame target model and the similarity metric of the current frame model defined as
therefore, in order to find the position of the next frame target in the field, it is to find a candidate window in the next frame image, so that the candidate histogram is the most similar to the target histogram, that is, maximizing the PAP coefficient. how can search within the domain maximize the PAP coefficients with the fewest iterations? This will use the Meanshift vector. Meanshift VectorTo make the maximum position of the search window in the current frame as the position of the current Frame Search window, the center of the window is y0, in order to find the local optimal target position in the Y0 neighborhood y1 the Taylor expansion similarity function can be approximated as
because the first item in P is not related to Y, you only need to maximize the second item. The second is exactly the probability density estimate weighted by the kernel function, so that, according to the mean drift theory, we can use Meanshift to iterate to find the peak of the probability density, that is, the maximum value of the second item. The mean shift vector can be deduced by finding the maximum value of the similarity function, that is, the y is biased and the partial derivative is 0.
so
where G () is the derivative of the weighted function K ()
Make y0 = y1, iterative calculation with the above, until | | y1-y0| | is less than a certain threshold or the maximum number of iterations has been reached. The Mean shift algorithm iterates over the last iteration to get the optimal position of the target in the current frame Y.
function [] = SELECT () close all; Clear all; %%%%%%%%%%%%%%%%%% the tracking target%%%%%%%%%%%%%%%%%%%%%%% i=imread (' result72.jpg ') according to a fully visible image of the target; Figure (1); Imshow (I); [Temp,rect]=imcrop (I); [A,b,c]=size (temp); %a:row,b:col%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% calculates the weight matrix of the target image%%%%%%%%%%%%%%%%%%%%%%% y (1) =a/2; Y (2) =b/2; Tic_x=rect (1) +rect (3)/2; Tic_y=rect (2) +rect (4)/2; M_wei=zeros (A, b);% Weight Matrix h=y (1) ^2+y (2) ^2;% bandwidth for i=1:a for J=1:b dist= (i-y (1)) ^2+ (J-y (2)) ^2; M_wei (i,j) =1-dist/h; %epanechnikov profile end End C=1/sum (sum (m_wei)),% normalization factor% calculates the target weight histogram qu%hist1=c*wei_hist (temp,m_wei,a,b);%targ ET model Hist1=zeros (1,4096); For the i=1:a for j=1:b%rgb color space quantization to 16*16*16 bins Q_r=fix (Double (temp (i,j,1))/16); %fix is the nearly 0 rounding function q_g=fix (double (temp (i,j,2))/16); Q_b=fix (Double (temp (i,j,3))/16); Q_temp=q_r*256+q_g*16+q_b; % set the percentage of red, green, and blue components per pixel hist1 (q_temp+1) = Hist1 (q_temp+1) +m_wei (I,J); % meterCalculates the weight of each pixel in the histogram statistic end end hist1=hist1*c; Rect (3) =ceil (Rect (3)); Rect (4) =ceil (Rect (4)); %%%%%%%%%%%%%%%%%%%%%%%%% Read sequence image myfile=dir (' D:\matlab7\work\mean shift\image\*.jpg '); Lengthfile=length (myfile); For L=1:lengthfile Im=imread (MyFile (L). name); num=0; y=[2,2]; %%%%%%%mean Shift Iteration while ((Y (1) ^2+y (2) ^2>0.5) &num<20)% iteration condition num=num+1; Temp1=imcrop (Im,rect); % calculation candidate Area histogram%hist2=c*wei_hist (temp1,m_wei,a,b);%target candidates pu Hist2=zeros (1,4096); For i=1:a for J=1:b Q_r=fix (double (Temp1 (i,j,1))/16); Q_g=fix (Double (Temp1 (i,j,2))/16); Q_b=fix (Double (Temp1 (i,j,3))/16); Q_TEMP1 (i,j) =q_r*256+q_g*16+q_b; Hist2 (Q_TEMP1 (i,j) +1) = Hist2 (Q_TEMP1 (i,j) +1) +m_wei (I,J); End End Hist2=hist2*c; Figure (2); Subplot (1,2,1); Plot (Hist2); Hold on; W=zeros (1,4096); For i=1:4096 if (Hist2 (i) ~=0)% is not equal to W (i) =sqrt (Hist1 (i)/hist2 (i)); else W (i) = 0; End end% variable initialization sum_w=0; xw=[0,0]; For i=1:a; For J=1:b sum_w=sum_w+w (UInt32 (Q_TEMP1 (i,j)) +1); Xw=xw+w (UInt32 (Q_TEMP1 (i,j)) +1) *[i-y (1) -0.5,j-y (2)-0.5]; End End Y=xw/sum_w; % Center Point Location Update rect (1) =rect (1) +y (2); Rect (2) =rect (2) +y (1); End%%% Tracking Track Matrix%%% Tic_x=[tic_x;rect (1) +rect (3)/2]; Tic_y=[tic_y;rect (2) +rect (4)/2]; V1=rect (1); V2=rect (2); V3=rect (3); V4=rect (4); %%% Display trace Results%%% subplot (1,2,2); Imshow (Uint8 (Im)); Title (' Target tracking results and their trajectories '); Hold on; Plot ([V1,v1+v3],[v2,v2],[v1,v1],[v2,v2+v4],[v1,v1+v3],[v2+v4,v2+v4],[v1+v3,v1+v3], [v2,v2+v4], ' linewidth ', 2, ' Color ', ' R '); Plot (tic_x,tic_y, ' linewidth ', 2, ' Color ', ' B '); End
A target tracking algorithm based on mean shift