From 30814aa76a51dc017ba74b24feb87cfd6075a96c Mon Sep 17 00:00:00 2001
From: Jeffsan <wang.chen@zoho.com>
Date: Wed, 22 Nov 2017 09:39:43 +0800
Subject: [PATCH] add tracking sources

---
 tracking/computeHistogram.m                  |  24 ++
 tracking/cropFilterResponse.m                |  26 ++
 tracking/display_statistics.m                |  14 +
 tracking/fhog.m                              |  76 +++++
 tracking/freezeColors.m                      | 275 ++++++++++++++++++
 tracking/gaussianResponse.m                  |  15 +
 tracking/gaussian_correlation.m              |  26 ++
 tracking/gaussian_correlation_scale_single.m |  10 +
 tracking/getCenterLikelihood.m               |  30 ++
 tracking/getColourMap.m                      |  17 ++
 tracking/getDesiredResponse.m                |  15 +
 tracking/getFeatureMap.m                     |  31 +++
 tracking/getScalePatch.m                     |  18 ++
 tracking/getScaleSubwindow.m                 |  38 +++
 tracking/getScaleSubwindow_v1.m              |  25 ++
 tracking/getSubwindow.m                      |  42 +++
 tracking/get_axis_aligned_BB.m               |  14 +
 tracking/get_feature_map.m                   |  61 ++++
 tracking/initializeAllAreas.m                |  40 +++
 tracking/isToolboxAvailable.m                |  41 +++
 tracking/mergeResponses.m                    |  21 ++
 tracking/params.txt                          |  36 +++
 tracking/precision_plot.m                    |  46 ++++
 tracking/readParams.m                        |  26 ++
 tracking/run_KCC.m                           |  41 +++
 tracking/trackerMain_otb_wangchen.m          | 276 +++++++++++++++++++
 tracking/unfreezeColors.m                    | 106 +++++++
 tracking/updateHistModel.m                   |  32 +++
 tracking/utility/affparam2geom.m             |  39 +++
 tracking/utility/affparam2mat.m              |  30 ++
 tracking/utility/affparaminv.m               |  17 ++
 tracking/utility/affwarpimg.m                |  24 ++
 tracking/utility/box_to_p.m                  |  39 +++
 tracking/utility/interp2.cpp                 | 110 ++++++++
 tracking/utility/interp2.dll                 | Bin 0 -> 24576 bytes
 tracking/utility/interp2.mexa64              | Bin 0 -> 12982 bytes
 tracking/utility/interp2.mexglx              | Bin 0 -> 10719 bytes
 tracking/utility/interp2.mexw64              | Bin 0 -> 10752 bytes
 tracking/utility/p_to_box.m                  |  38 +++
 tracking/utility/warpimg.m                   |  29 ++
 tracking/wrappers.hpp                        |  42 +++
 41 files changed, 1790 insertions(+)
 create mode 100644 tracking/computeHistogram.m
 create mode 100644 tracking/cropFilterResponse.m
 create mode 100644 tracking/display_statistics.m
 create mode 100644 tracking/fhog.m
 create mode 100644 tracking/freezeColors.m
 create mode 100644 tracking/gaussianResponse.m
 create mode 100644 tracking/gaussian_correlation.m
 create mode 100644 tracking/gaussian_correlation_scale_single.m
 create mode 100644 tracking/getCenterLikelihood.m
 create mode 100644 tracking/getColourMap.m
 create mode 100644 tracking/getDesiredResponse.m
 create mode 100644 tracking/getFeatureMap.m
 create mode 100644 tracking/getScalePatch.m
 create mode 100644 tracking/getScaleSubwindow.m
 create mode 100644 tracking/getScaleSubwindow_v1.m
 create mode 100644 tracking/getSubwindow.m
 create mode 100644 tracking/get_axis_aligned_BB.m
 create mode 100644 tracking/get_feature_map.m
 create mode 100644 tracking/initializeAllAreas.m
 create mode 100644 tracking/isToolboxAvailable.m
 create mode 100644 tracking/mergeResponses.m
 create mode 100644 tracking/params.txt
 create mode 100644 tracking/precision_plot.m
 create mode 100644 tracking/readParams.m
 create mode 100644 tracking/run_KCC.m
 create mode 100644 tracking/trackerMain_otb_wangchen.m
 create mode 100644 tracking/unfreezeColors.m
 create mode 100644 tracking/updateHistModel.m
 create mode 100644 tracking/utility/affparam2geom.m
 create mode 100644 tracking/utility/affparam2mat.m
 create mode 100644 tracking/utility/affparaminv.m
 create mode 100644 tracking/utility/affwarpimg.m
 create mode 100644 tracking/utility/box_to_p.m
 create mode 100644 tracking/utility/interp2.cpp
 create mode 100644 tracking/utility/interp2.dll
 create mode 100644 tracking/utility/interp2.mexa64
 create mode 100644 tracking/utility/interp2.mexglx
 create mode 100644 tracking/utility/interp2.mexw64
 create mode 100644 tracking/utility/p_to_box.m
 create mode 100644 tracking/utility/warpimg.m
 create mode 100644 tracking/wrappers.hpp

diff --git a/tracking/computeHistogram.m b/tracking/computeHistogram.m
new file mode 100644
index 0000000..9207538
--- /dev/null
+++ b/tracking/computeHistogram.m
@@ -0,0 +1,24 @@
+function histogram = computeHistogram(patch, mask, n_bins, grayscale_sequence)
+%COMPUTEHISTOGRAM creates a colour (or grayscale) histogram of an image patch
+% MASK has the same size as the image patch and selects what should
+% be used when computing the histogram (i.e. out-of-frame regions are ignored)
+
+	[h, w, d] = size(patch);
+
+	assert(all([h w]==size(mask)) == 1, 'mask and image are not the same size');
+
+	bin_width = 256/n_bins;
+
+	% convert image to 1d array with same n channels of img patch
+	patch_array = reshape(double(patch), w*h, d);
+	% compute to which bin each pixel (for all 3 channels) belongs to
+	bin_indices = floor(patch_array/bin_width) + 1;
+
+	if grayscale_sequence
+		histogram = accumarray(bin_indices, mask(:), [n_bins 1])/sum(mask(:));
+	else
+		% the histogram is a cube of side n_bins
+		histogram = accumarray(bin_indices, mask(:), [n_bins n_bins n_bins])/sum(mask(:));
+	end
+
+end
diff --git a/tracking/cropFilterResponse.m b/tracking/cropFilterResponse.m
new file mode 100644
index 0000000..0d2e734
--- /dev/null
+++ b/tracking/cropFilterResponse.m
@@ -0,0 +1,26 @@
+function new_response = cropFilterResponse(response_cf, response_size)
+%CROPFILTERRESPONSE makes RESPONSE_CF of size RESPONSE_SIZE (i.e. same size of colour response)
+
+    [h,w] = size(response_cf);
+    b = response_size(1);
+    a = response_size(2);
+
+    % a and b must be odd, as we want an exact center
+    if ~all_odd([a, b])
+        error('dimensions must be odd');
+    end
+    half_width = floor(a/2);
+    half_height = floor(b/2);
+
+    new_response = response_cf(...
+        mod_one(-half_height:half_height, h), ...
+        mod_one(-half_width:half_width, w));
+end
+
+function y = mod_one(a, b)
+    y = mod(a-1, b)+1;
+end
+
+function y = all_odd(x)
+    y = all(mod(x, 2) == 1);
+end
diff --git a/tracking/display_statistics.m b/tracking/display_statistics.m
new file mode 100644
index 0000000..abfcd0f
--- /dev/null
+++ b/tracking/display_statistics.m
@@ -0,0 +1,14 @@
+function display_statistics(name, precision, reinitializations, fps)
+    %print statistics to command window
+    str = sprintf('Overlap:% 1.3f', precision);
+
+    if ~isnan(reinitializations),
+        if reinitializations == round(reinitializations),  %prettier for integer values
+            str = [str, sprintf(', reinitializations:% 2i', reinitializations)];
+        else
+            str = [str, sprintf(', reinitializations:% 1.2f', reinitializations)];
+        end
+    end
+
+    fprintf('%12s - %s, FPS:% 4.2f\n', name, str, fps)
+end
diff --git a/tracking/fhog.m b/tracking/fhog.m
new file mode 100644
index 0000000..f96fc0f
--- /dev/null
+++ b/tracking/fhog.m
@@ -0,0 +1,76 @@
+function H = fhog( I, binSize, nOrients, clip, crop )
+% Efficiently compute Felzenszwalb's HOG (FHOG) features.
+%
+% A fast implementation of the HOG variant used by Felzenszwalb et al.
+% in their work on discriminatively trained deformable part models.
+%  http://www.cs.berkeley.edu/~rbg/latent/index.html
+% Gives nearly identical results to features.cc in code release version 5
+% but runs 4x faster (over 125 fps on VGA color images).
+%
+% The computed HOG features are 3*nOrients+5 dimensional. There are
+% 2*nOrients contrast sensitive orientation channels, nOrients contrast
+% insensitive orientation channels, 4 texture channels and 1 all zeros
+% channel (used as a 'truncation' feature). Using the standard value of
+% nOrients=9 gives a 32 dimensional feature vector at each cell. This
+% variant of HOG, refered to as FHOG, has been shown to achieve superior
+% performance to the original HOG features. For details please refer to
+% work by Felzenszwalb et al. (see link above).
+%
+% This function is essentially a wrapper for calls to gradientMag()
+% and gradientHist(). Specifically, it is equivalent to the following:
+%  [M,O] = gradientMag( I,0,0,0,1 ); softBin = -1; useHog = 2;
+%  H = gradientHist(M,O,binSize,nOrients,softBin,useHog,clip);
+% See gradientHist() for more general usage.
+%
+% This code requires SSE2 to compile and run (most modern Intel and AMD
+% processors support SSE2). Please see: http://en.wikipedia.org/wiki/SSE2.
+%
+% USAGE
+%  H = fhog( I, [binSize], [nOrients], [clip], [crop] )
+%
+% INPUTS
+%  I        - [hxw] color or grayscale input image (must have type single)
+%  binSize  - [8] spatial bin size
+%  nOrients - [9] number of orientation bins
+%  clip     - [.2] value at which to clip histogram bins
+%  crop     - [0] if true crop boundaries
+%
+% OUTPUTS
+%  H        - [h/binSize w/binSize nOrients*3+5] computed hog features
+%
+% EXAMPLE
+%  I=imResample(single(imread('peppers.png'))/255,[480 640]);
+%  tic, for i=1:100, H=fhog(I,8,9); end; disp(100/toc) % >125 fps
+%  figure(1); im(I); V=hogDraw(H,25,1); figure(2); im(V)
+%
+% EXAMPLE
+%  % comparison to features.cc (requires DPM code release version 5)
+%  I=imResample(single(imread('peppers.png'))/255,[480 640]); Id=double(I);
+%  tic, for i=1:100, H1=features(Id,8); end; disp(100/toc)
+%  tic, for i=1:100, H2=fhog(I,8,9,.2,1); end; disp(100/toc)
+%  figure(1); montage2(H1); figure(2); montage2(H2);
+%  D=abs(H1-H2); mean(D(:))
+%
+% See also hog, hogDraw, gradientHist
+%
+% Piotr's Image&Video Toolbox      Version 3.23
+% Copyright 2013 Piotr Dollar.  [pdollar-at-caltech.edu]
+% Please email me if you find bugs, or have suggestions or questions!
+% Licensed under the Simplified BSD License [see external/bsd.txt]
+
+%Note: modified to be more self-contained
+
+if( nargin<2 ), binSize=8; end
+if( nargin<3 ), nOrients=9; end
+if( nargin<4 ), clip=.2; end
+if( nargin<5 ), crop=0; end
+
+softBin = -1; useHog = 2; b = binSize;
+
+[M,O]=gradientMex('gradientMag',I,0,1);
+
+H = gradientMex('gradientHist',M,O,binSize,nOrients,softBin,useHog,clip);
+
+if( crop ), e=mod(size(I),b)<b/2; H=H(2:end-e(1),2:end-e(2),:); end
+
+end
diff --git a/tracking/freezeColors.m b/tracking/freezeColors.m
new file mode 100644
index 0000000..6c7b465
--- /dev/null
+++ b/tracking/freezeColors.m
@@ -0,0 +1,275 @@
+function freezeColors(varargin)
+% freezeColors  Lock colors of plot, enabling multiple colormaps per figure. (v2.3)
+%
+%   Problem: There is only one colormap per figure. This function provides
+%       an easy solution when plots using different colomaps are desired 
+%       in the same figure.
+%
+%   freezeColors freezes the colors of graphics objects in the current axis so 
+%       that subsequent changes to the colormap (or caxis) will not change the
+%       colors of these objects. freezeColors works on any graphics object 
+%       with CData in indexed-color mode: surfaces, images, scattergroups, 
+%       bargroups, patches, etc. It works by converting CData to true-color rgb
+%       based on the colormap active at the time freezeColors is called.
+%
+%   The original indexed color data is saved, and can be restored using
+%       unfreezeColors, making the plot once again subject to the colormap and
+%       caxis.
+%
+%
+%   Usage:
+%       freezeColors        applies to all objects in current axis (gca),
+%       freezeColors(axh)   same, but works on axis axh.
+%
+%   Example:
+%       subplot(2,1,1); imagesc(X); colormap hot; freezeColors
+%       subplot(2,1,2); imagesc(Y); colormap hsv; freezeColors etc...
+%
+%       Note: colorbars must also be frozen. Due to Matlab 'improvements' this can
+%				no longer be done with freezeColors. Instead, please
+%				use the function CBFREEZE by Carlos Adrian Vargas Aguilera
+%				that can be downloaded from the MATLAB File Exchange
+%				(http://www.mathworks.com/matlabcentral/fileexchange/24371)
+%
+%       h=colorbar; cbfreeze(h), or simply cbfreeze(colorbar)
+%
+%       For additional examples, see test/test_main.m
+%
+%   Side effect on render mode: freezeColors does not work with the painters
+%       renderer, because Matlab doesn't support rgb color data in
+%       painters mode. If the current renderer is painters, freezeColors
+%       changes it to zbuffer. This may have unexpected effects on other aspects
+%	      of your plots.
+%
+%       See also unfreezeColors, freezeColors_pub.html, cbfreeze.
+%
+%
+%   John Iversen (iversen@nsi.edu) 3/23/05
+%
+
+%   Changes:
+%   JRI (iversen@nsi.edu) 4/19/06   Correctly handles scaled integer cdata
+%   JRI 9/1/06   should now handle all objects with cdata: images, surfaces, 
+%                scatterplots. (v 2.1)
+%   JRI 11/11/06 Preserves NaN colors. Hidden option (v 2.2, not uploaded)
+%   JRI 3/17/07  Preserve caxis after freezing--maintains colorbar scale (v 2.3)
+%   JRI 4/12/07  Check for painters mode as Matlab doesn't support rgb in it.
+%   JRI 4/9/08   Fix preserving caxis for objects within hggroups (e.g. contourf)
+%   JRI 4/7/10   Change documentation for colorbars
+
+% Hidden option for NaN colors:
+%   Missing data are often represented by NaN in the indexed color
+%   data, which renders transparently. This transparency will be preserved
+%   when freezing colors. If instead you wish such gaps to be filled with 
+%   a real color, add 'nancolor',[r g b] to the end of the arguments. E.g. 
+%   freezeColors('nancolor',[r g b]) or freezeColors(axh,'nancolor',[r g b]),
+%   where [r g b] is a color vector. This works on images & pcolor, but not on
+%   surfaces.
+%   Thanks to Fabiano Busdraghi and Jody Klymak for the suggestions. Bugfixes 
+%   attributed in the code.
+
+% Free for all uses, but please retain the following:
+%   Original Author:
+%   John Iversen, 2005-10
+%   john_iversen@post.harvard.edu
+
+appdatacode = 'JRI__freezeColorsData';
+
+[h, nancolor] = checkArgs(varargin);
+
+%gather all children with scaled or indexed CData
+cdatah = getCDataHandles(h);
+
+%current colormap
+cmap = colormap;
+nColors = size(cmap,1);
+cax = caxis;
+
+% convert object color indexes into colormap to true-color data using 
+%  current colormap
+for hh = cdatah',
+    g = get(hh);
+    
+    %preserve parent axis clim
+    parentAx = getParentAxes(hh);
+    originalClim = get(parentAx, 'clim');    
+   
+    %   Note: Special handling of patches: For some reason, setting
+    %   cdata on patches created by bar() yields an error,
+    %   so instead we'll set facevertexcdata instead for patches.
+    if ~strcmp(g.Type,'patch'),
+        cdata = g.CData;
+    else
+        cdata = g.FaceVertexCData; 
+    end
+    
+    %get cdata mapping (most objects (except scattergroup) have it)
+    if isfield(g,'CDataMapping'),
+        scalemode = g.CDataMapping;
+    else
+        scalemode = 'scaled';
+    end
+    
+    %save original indexed data for use with unfreezeColors
+    siz = size(cdata);
+    setappdata(hh, appdatacode, {cdata scalemode});
+
+    %convert cdata to indexes into colormap
+    if strcmp(scalemode,'scaled'),
+        %4/19/06 JRI, Accommodate scaled display of integer cdata:
+        %       in MATLAB, uint * double = uint, so must coerce cdata to double
+        %       Thanks to O Yamashita for pointing this need out
+        idx = ceil( (double(cdata) - cax(1)) / (cax(2)-cax(1)) * nColors);
+    else %direct mapping
+        idx = cdata;
+        %10/8/09 in case direct data is non-int (e.g. image;freezeColors)
+        % (Floor mimics how matlab converts data into colormap index.)
+        % Thanks to D Armyr for the catch
+        idx = floor(idx);
+    end
+    
+    %clamp to [1, nColors]
+    idx(idx<1) = 1;
+    idx(idx>nColors) = nColors;
+
+    %handle nans in idx
+    nanmask = isnan(idx);
+    idx(nanmask)=1; %temporarily replace w/ a valid colormap index
+
+    %make true-color data--using current colormap
+    realcolor = zeros(siz);
+    for i = 1:3,
+        c = cmap(idx,i);
+        c = reshape(c,siz);
+        c(nanmask) = nancolor(i); %restore Nan (or nancolor if specified)
+        realcolor(:,:,i) = c;
+    end
+    
+    %apply new true-color color data
+    
+    %true-color is not supported in painters renderer, so switch out of that
+    if strcmp(get(gcf,'renderer'), 'painters'),
+        set(gcf,'renderer','zbuffer');
+    end
+    
+    %replace original CData with true-color data
+    if ~strcmp(g.Type,'patch'),
+        set(hh,'CData',realcolor);
+    else
+        set(hh,'faceVertexCData',permute(realcolor,[1 3 2]))
+    end
+    
+    %restore clim (so colorbar will show correct limits)
+    if ~isempty(parentAx),
+        set(parentAx,'clim',originalClim)
+    end
+    
+end %loop on indexed-color objects
+
+
+% ============================================================================ %
+% Local functions
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% getCDataHandles -- get handles of all descendents with indexed CData
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+function hout = getCDataHandles(h)
+% getCDataHandles  Find all objects with indexed CData
+
+%recursively descend object tree, finding objects with indexed CData
+% An exception: don't include children of objects that themselves have CData:
+%   for example, scattergroups are non-standard hggroups, with CData. Changing
+%   such a group's CData automatically changes the CData of its children, 
+%   (as well as the children's handles), so there's no need to act on them.
+
+error(nargchk(1,1,nargin,'struct'))
+
+hout = [];
+if isempty(h),return;end
+
+ch = get(h,'children');
+for hh = ch'
+    g = get(hh);
+    if isfield(g,'CData'),     %does object have CData?
+        %is it indexed/scaled?
+        if ~isempty(g.CData) && isnumeric(g.CData) && size(g.CData,3)==1, 
+            hout = [hout; hh]; %#ok<AGROW> %yes, add to list
+        end
+    else %no CData, see if object has any interesting children
+            hout = [hout; getCDataHandles(hh)]; %#ok<AGROW>
+    end
+end
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% getParentAxes -- return handle of axes object to which a given object belongs
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+function hAx = getParentAxes(h)
+% getParentAxes  Return enclosing axes of a given object (could be self)
+
+error(nargchk(1,1,nargin,'struct'))
+%object itself may be an axis
+if strcmp(get(h,'type'),'axes'),
+    hAx = h;
+    return
+end
+
+parent = get(h,'parent');
+if (strcmp(get(parent,'type'), 'axes')),
+    hAx = parent;
+else
+    hAx = getParentAxes(parent);
+end
+
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+%% checkArgs -- Validate input arguments
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+function [h, nancolor] = checkArgs(args)
+% checkArgs  Validate input arguments to freezeColors
+
+nargs = length(args);
+error(nargchk(0,3,nargs,'struct'))
+
+%grab handle from first argument if we have an odd number of arguments
+if mod(nargs,2),
+    h = args{1};
+    if ~ishandle(h),
+        error('JRI:freezeColors:checkArgs:invalidHandle',...
+            'The first argument must be a valid graphics handle (to an axis)')
+    end
+    % 4/2010 check if object to be frozen is a colorbar
+    if strcmp(get(h,'Tag'),'Colorbar'),
+      if ~exist('cbfreeze.m'),
+        warning('JRI:freezeColors:checkArgs:cannotFreezeColorbar',...
+            ['You seem to be attempting to freeze a colorbar. This no longer'...
+            'works. Please read the help for freezeColors for the solution.'])
+      else
+        cbfreeze(h);
+        return
+      end
+    end
+    args{1} = [];
+    nargs = nargs-1;
+else
+    h = gca;
+end
+
+%set nancolor if that option was specified
+nancolor = [nan nan nan];
+if nargs == 2,
+    if strcmpi(args{end-1},'nancolor'),
+        nancolor = args{end};
+        if ~all(size(nancolor)==[1 3]),
+            error('JRI:freezeColors:checkArgs:badColorArgument',...
+                'nancolor must be [r g b] vector');
+        end
+        nancolor(nancolor>1) = 1; nancolor(nancolor<0) = 0;
+    else
+        error('JRI:freezeColors:checkArgs:unrecognizedOption',...
+            'Unrecognized option (%s). Only ''nancolor'' is valid.',args{end-1})
+    end
+end
+
+
diff --git a/tracking/gaussianResponse.m b/tracking/gaussianResponse.m
new file mode 100644
index 0000000..b26204d
--- /dev/null
+++ b/tracking/gaussianResponse.m
@@ -0,0 +1,15 @@
+function y = gaussianResponse(rect_size, sigma)
+%GAUSSIANRESPONSE create the (fixed) target response of the correlation filter response
+    half = floor((rect_size-1) / 2);
+    i_range = -half(1):half(1);
+    j_range = -half(2):half(2);
+    [i, j] = ndgrid(i_range, j_range);
+    i_mod_range = mod_one(i_range, rect_size(1));
+    j_mod_range = mod_one(j_range, rect_size(2));
+    y = zeros(rect_size);
+    y(i_mod_range, j_mod_range) = exp(-(i.^2 + j.^2) / (2 * sigma^2));
+end
+
+function y = mod_one(a, b)
+    y = mod(a-1, b)+1;
+end
diff --git a/tracking/gaussian_correlation.m b/tracking/gaussian_correlation.m
new file mode 100644
index 0000000..9c7411b
--- /dev/null
+++ b/tracking/gaussian_correlation.m
@@ -0,0 +1,26 @@
+function kf = gaussian_correlation(xf, yf, sigma)
+%GAUSSIAN_CORRELATION Gaussian Kernel at all shifts, i.e. kernel correlation.
+%   Evaluates a Gaussian kernel with bandwidth SIGMA for all relative
+%   shifts between input images X and Y, which must both be MxN. They must 
+%   also be periodic (ie., pre-processed with a cosine window). The result
+%   is an MxN map of responses.
+%
+%   Inputs and output are all in the Fourier domain.
+%
+%   Joao F. Henriques, 2014
+%   http://www.isr.uc.pt/~henriques/
+	
+	N = size(xf,1) * size(xf,2);
+	xx = xf(:)' * xf(:) / N;  %squared norm of x
+	yy = yf(:)' * yf(:) / N;  %squared norm of y
+	
+	%cross-correlation term in Fourier domain
+	xyf = xf .* conj(yf);
+	xy = sum(real(ifft2(xyf)), 3);  %to spatial domain
+	
+	%calculate gaussian response for all positions, then go back to the
+	%Fourier domain
+	kf = fft2(exp(-1 / sigma^2 * max(0, (xx + yy - 2 * xy) / numel(xf))));
+
+end
+
diff --git a/tracking/gaussian_correlation_scale_single.m b/tracking/gaussian_correlation_scale_single.m
new file mode 100644
index 0000000..9ee189d
--- /dev/null
+++ b/tracking/gaussian_correlation_scale_single.m
@@ -0,0 +1,10 @@
+function kf = gaussian_correlation_scale_single(base, y, sigma)
+	
+    N = size(y,1) * size(y,2);%*size(y,3);
+    for i =1:size(base,2)
+        n = sum(sum(sum((base(:,i) - y).^2)))/N^2;
+        k(i) = exp(-n / sigma^2 ) ;
+    end
+    
+    kf = fft(k');
+end
diff --git a/tracking/getCenterLikelihood.m b/tracking/getCenterLikelihood.m
new file mode 100644
index 0000000..e273027
--- /dev/null
+++ b/tracking/getCenterLikelihood.m
@@ -0,0 +1,30 @@
+function center_likelihood = getCenterLikelihood(object_likelihood, m)
+%GETCENTERLIKELIHOOD computes the sum over rectangles of size M.
+% CENTER_LIKELIHOOD is the 'colour response'
+    [h,w] = size(object_likelihood);
+    n1 = h - m(1) + 1;
+    n2 = w - m(2) + 1;
+
+    %% integral images
+    % compute summed area table
+%     SAT = zeros(h, w);
+%     for y=1:h
+%         for x=1:w
+%             if x>1, SAT_left = SAT(y,x-1);
+%             else    SAT_left = 0; end
+%             if y>1, SAT_up = SAT(y-1,x);
+%             else    SAT_up = 0; end
+%             if y>1 && x>1, SAT_left_up = SAT(y-1,x-1);
+%             else    SAT_left_up = 0; end
+%
+%             SAT(y,x) = object_likelihood(y,x) + SAT_left + SAT_up - SAT_left_up;
+%         end
+%     end
+%     SAT = padarray(SAT,[1 1], 'pre');
+
+%% equivalent MATLAB function
+    SAT = integralImage(object_likelihood);
+    i = 1:n1;
+    j = 1:n2;
+    center_likelihood = (SAT(i,j) + SAT(i+m(1), j+m(2)) - SAT(i+m(1), j) - SAT(i, j+m(2))) / prod(m);
+end
diff --git a/tracking/getColourMap.m b/tracking/getColourMap.m
new file mode 100644
index 0000000..7380881
--- /dev/null
+++ b/tracking/getColourMap.m
@@ -0,0 +1,17 @@
+function [P_O] = getColourMap(patch, bg_hist, fg_hist, n_bins, grayscale_sequence)
+%% GETCOLOURMAP computes pixel-wise probabilities (PwP) given PATCH and models BG_HIST and FG_HIST
+    % check whether the patch has 3 channels
+    [h, w, d] = size(patch);
+    % figure out which bin each pixel falls into
+    bin_width = 256/n_bins;
+    % convert image to d channels array
+    patch_array = reshape(double(patch), w*h, d);
+    % to which bin each pixel (for all d channels) belongs to
+    bin_indices = floor(patch_array/bin_width) + 1;
+    % Get pixel-wise posteriors (PwP)
+    P_bg = getP(bg_hist, h, w, bin_indices, grayscale_sequence);
+    P_fg = getP(fg_hist, h, w, bin_indices, grayscale_sequence);
+
+    % Object-likelihood map
+    P_O = P_fg ./ (P_fg + P_bg);
+end
diff --git a/tracking/getDesiredResponse.m b/tracking/getDesiredResponse.m
new file mode 100644
index 0000000..94a5f02
--- /dev/null
+++ b/tracking/getDesiredResponse.m
@@ -0,0 +1,15 @@
+function [desiredResponse, x0, y0] = getDesiredResponse(saliency_map, bg_area, target_sz, p)
+% GETDESIREDRESPONSE is a wrapper for getGaussian, that computes the (object-centered) desired response
+
+% If the saliency patch is empty, use a patch centered desired response only for one frame
+if(sum(sum(saliency_map))==0)
+    p.gaussianType = 'bboxCentered';
+    [desiredResponse, ~, ~] = getGaussian(saliency_map, bg_area, target_sz, p);
+    x0 = 0;
+    y0 = 0;
+    p.gaussianType = 'objectCentered';
+else
+    [desiredResponse, x0, y0] = getGaussian(saliency_map, bg_area, target_sz, p);
+end
+
+end
\ No newline at end of file
diff --git a/tracking/getFeatureMap.m b/tracking/getFeatureMap.m
new file mode 100644
index 0000000..7cf9927
--- /dev/null
+++ b/tracking/getFeatureMap.m
@@ -0,0 +1,31 @@
+function out = getFeatureMap(im_patch, feature_type, cf_response_size, hog_cell_size)
+
+% code from DSST
+
+% allocate space
+switch feature_type
+    case 'fhog'
+        temp = fhog(single(im_patch), hog_cell_size);
+        h = cf_response_size(1);
+        w = cf_response_size(2);
+        out = zeros(h, w, 28, 'single');
+        out(:,:,2:28) = temp(:,:,1:27);
+        if hog_cell_size > 1
+            im_patch = mexResize(im_patch, [h, w] ,'auto');
+        end
+        % if color image
+        if size(im_patch, 3) > 1
+            im_patch = rgb2gray(im_patch);
+        end
+        out(:,:,1) = single(im_patch)/255 - 0.5;
+    case 'gray'
+        if hog_cell_size > 1, im_patch = mexResize(im_patch,cf_response_size,'auto');   end
+        if size(im_patch, 3) == 1
+            out = single(im_patch)/255 - 0.5;
+        else
+            out = single(rgb2gray(im_patch))/255 - 0.5;
+        end        
+end
+        
+end
+
diff --git a/tracking/getScalePatch.m b/tracking/getScalePatch.m
new file mode 100644
index 0000000..aec42a2
--- /dev/null
+++ b/tracking/getScalePatch.m
@@ -0,0 +1,18 @@
+function out = getScalePatch(img, pos, target_size, search_size, rotation, scale_model_sz, hog_scale_cell_size)
+
+% code from DSST
+padding=0;
+if size(img,3)>1
+img= single(rgb2gray(img))/255.0;
+end
+
+tmp_sz = floor((target_size * (1 + padding))*search_size);
+param0 = [pos(2), pos(1), tmp_sz(2)/target_size(2), rotation/180*pi,...
+            tmp_sz(1)/target_size(2)/(target_size(1)/target_size(2)),0];
+param0 = affparam2mat(param0); 
+patch = warpimg(double(img), param0, target_size);
+patch = fhog(single(patch), hog_scale_cell_size);
+patch = patch(:,:,1:31);
+out = patch(:);
+
+
diff --git a/tracking/getScaleSubwindow.m b/tracking/getScaleSubwindow.m
new file mode 100644
index 0000000..3113105
--- /dev/null
+++ b/tracking/getScaleSubwindow.m
@@ -0,0 +1,38 @@
+function out = getScaleSubwindow(im, pos, base_target_sz, scale_factors, scale_window, scale_model_sz, hog_scale_cell_size)
+
+% code from DSST
+
+num_scales = length(scale_factors);
+
+for s = 1:num_scales
+    patch_sz = floor(base_target_sz * scale_factors(s));
+    %make sure the size is not to small
+    patch_sz = max(patch_sz, 2);
+ 
+    xs = floor(pos(2)) + (1:patch_sz(2)) - floor(patch_sz(2)/2);
+    ys = floor(pos(1)) + (1:patch_sz(1)) - floor(patch_sz(1)/2);
+    
+    %check for out-of-bounds coordinates, and set them to the values at
+    %the borders
+    xs(xs < 1) = 1;
+    ys(ys < 1) = 1;
+    xs(xs > size(im,2)) = size(im,2);
+    ys(ys > size(im,1)) = size(im,1);
+    
+    %extract image
+    im_patch = im(ys, xs, :);
+    
+    % resize image to model size
+    im_patch_resized = mexResize(im_patch, scale_model_sz, 'auto');
+    
+    % extract scale features
+    temp_hog = fhog(single(im_patch_resized), hog_scale_cell_size);
+    temp = temp_hog(:,:,1:31);
+    
+    if s == 1
+        out = zeros(numel(temp), num_scales, 'single');
+    end
+    
+    % window
+    out(:,s) = temp(:) * scale_window(s); 
+end
diff --git a/tracking/getScaleSubwindow_v1.m b/tracking/getScaleSubwindow_v1.m
new file mode 100644
index 0000000..2153a7c
--- /dev/null
+++ b/tracking/getScaleSubwindow_v1.m
@@ -0,0 +1,25 @@
+function out = getScaleSubwindow_v1(img, pos, target_size, search_size, rotation, scale_model_sz, hog_scale_cell_size)
+
+% code from DSST
+padding=0;
+num_scales = length(search_size);
+if size(img,3)>1
+img= single(rgb2gray(img))/255.0;
+end
+for s = 1:num_scales
+    tmp_sz = floor((target_size * (1 + padding))*search_size(s));
+    param0 = [pos(2), pos(1), tmp_sz(2)/target_size(2), rotation/180*pi,...
+                tmp_sz(1)/target_size(2)/(target_size(1)/target_size(2)),0];
+    param0 = affparam2mat(param0); 
+    patch = warpimg(double(img), param0, target_size);
+    patch = fhog(single(patch), hog_scale_cell_size);
+    patch = patch(:,:,1:31);
+    
+    if s == 1
+        out = zeros(numel(patch), num_scales, 'single');
+    end
+    
+    % window
+    out(:,s) = patch(:);%* scale_window(s);
+end
+
diff --git a/tracking/getSubwindow.m b/tracking/getSubwindow.m
new file mode 100644
index 0000000..2fef312
--- /dev/null
+++ b/tracking/getSubwindow.m
@@ -0,0 +1,42 @@
+function im_patch = getSubwindow(im, pos, model_sz, scaled_sz)
+%GET_SUBWINDOW Obtain image sub-window, padding is done by replicating border values.
+%   Returns sub-window of image IM centered at POS ([y, x] coordinates),
+%   with size MODEL_SZ ([height, width]). If any pixels are outside of the image,
+%   they will replicate the values at the borders
+
+% with 3 input, no scale. With 4 params, scale adaptation
+if nargin < 4, sz = model_sz;   
+else, sz = scaled_sz;    
+end
+
+%make sure the size is not to small
+sz = max(sz, 2);
+%if sz(1) < 1, sz(1) = 2; end;
+%if sz(2) < 1, sz(2) = 2; end;
+
+%xs = floor(pos(2)) + (1:sz(2)) - floor(sz(2)/2);
+%ys = floor(pos(1)) + (1:sz(1)) - floor(sz(1)/2);
+xs = round(pos(2) + (1:sz(2)) - sz(2)/2);
+ys = round(pos(1) + (1:sz(1)) - sz(1)/2);
+
+%check for out-of-bounds coordinates, and set them to the values at
+%the borders
+xs(xs < 1) = 1;
+ys(ys < 1) = 1;
+xs(xs > size(im,2)) = size(im,2);
+ys(ys > size(im,1)) = size(im,1);
+
+%extract image
+im_patch_original = im(ys, xs, :);
+
+% (if rescaling is introduced) resize image to model size
+% im_patch = imresize(im_patch, model_sz, 'bilinear');
+if nargin>=4
+%     im_patch = mexResize(im_patch_original, model_sz, 'auto');
+    im_patch = mexResize(im_patch_original, model_sz, 'auto');
+else
+    im_patch = im_patch_original;
+end
+
+end
+
diff --git a/tracking/get_axis_aligned_BB.m b/tracking/get_axis_aligned_BB.m
new file mode 100644
index 0000000..7cc45a9
--- /dev/null
+++ b/tracking/get_axis_aligned_BB.m
@@ -0,0 +1,14 @@
+function [cx, cy, w, h] = get_axis_aligned_BB(region)
+% GETAXISALIGNEDBB extracts an axis aligned bbox from the ground truth REGION with same area as the rotated one
+    cx = mean(region(1:2:end));
+    cy = mean(region(2:2:end));
+    x1 = min(region(1:2:end));
+    x2 = max(region(1:2:end));
+    y1 = min(region(2:2:end));
+    y2 = max(region(2:2:end));
+    A1 = norm(region(1:2) - region(3:4)) * norm(region(3:4) - region(5:6));
+    A2 = (x2 - x1) * (y2 - y1);
+    s = sqrt(A1/A2);
+    w = s * (x2 - x1) + 1;
+    h = s * (y2 - y1) + 1;
+end
diff --git a/tracking/get_feature_map.m b/tracking/get_feature_map.m
new file mode 100644
index 0000000..baadc66
--- /dev/null
+++ b/tracking/get_feature_map.m
@@ -0,0 +1,61 @@
+function out = get_feature_map(im_patch, features, w2c)
+
+% out = get_feature_map(im_patch, features, w2c)
+%
+% Extracts the given features from the image patch. w2c is the
+% Color Names matrix, if used.
+
+if nargin < 3
+    w2c = [];
+end
+
+% the names of the features that can be used
+valid_features = {'gray', 'cn'};
+
+% the dimension of the valid features
+feature_levels = [1 10]';
+
+num_valid_features = length(valid_features);
+used_features = false(num_valid_features, 1);
+
+% get the used features
+for i = 1:num_valid_features
+    used_features(i) = any(strcmpi(valid_features{i}, features));
+end
+
+% total number of used feature levels
+num_feature_levels = sum(feature_levels .* used_features);
+
+level = 0;
+
+% If grayscale image
+if size(im_patch, 3) == 1
+    % Features that are available for grayscale sequances
+    
+    % Grayscale values (image intensity)
+    out = single(im_patch)/255 - 0.5;
+else
+    % Features that are available for color sequances
+    
+    % allocate space (for speed)
+    out = zeros(size(im_patch, 1), size(im_patch, 2), num_feature_levels, 'single');
+    
+    % Grayscale values (image intensity)
+    if used_features(1)
+        out(:,:,level+1) = single(rgb2gray(im_patch))/255 - 0.5;
+        level = level + feature_levels(1);
+    end
+    
+    % Color Names
+    if used_features(2)
+        if isempty(w2c)
+            % load the RGB to color name matrix if not in input
+            temp = load('w2crs');
+            w2c = temp.w2crs;
+        end
+        
+        % extract color descriptor
+        out(:,:,level+(1:10)) = im2c(single(im_patch), w2c, -2);
+        level = level + feature_levels(2);
+    end
+end
\ No newline at end of file
diff --git a/tracking/initializeAllAreas.m b/tracking/initializeAllAreas.m
new file mode 100644
index 0000000..95ecd6b
--- /dev/null
+++ b/tracking/initializeAllAreas.m
@@ -0,0 +1,40 @@
+function [params, bg_area, fg_area, area_resize_factor] = initializeAllAreas(im, params)
+
+	% we want a regular frame surrounding the object
+	avg_dim = sum(params.target_sz)/2;
+	% size from which we extract features
+	bg_area = round(params.target_sz + avg_dim);
+	% pick a "safe" region smaller than bbox to avoid mislabeling
+	fg_area = round(params.target_sz - avg_dim * params.inner_padding);
+	% saturate to image size
+	if(bg_area(2)>size(im,2)), bg_area(2)=size(im,2)-1; end
+	if(bg_area(1)>size(im,1)), bg_area(1)=size(im,1)-1; end
+	% make sure the differences are a multiple of 2 (makes things easier later in color histograms)
+	bg_area = bg_area - mod(bg_area - params.target_sz, 2);
+	fg_area = fg_area + mod(bg_area - fg_area, 2);
+
+	% Compute the rectangle with (or close to) params.fixedArea and
+	% same aspect ratio as the target bbox
+	area_resize_factor = sqrt(params.fixed_area/prod(bg_area));
+	params.norm_bg_area = round(bg_area * area_resize_factor);
+	% Correlation Filter (HOG) feature space
+	% It smaller that the norm bg area if HOG cell size is > 1
+	params.cf_response_size = floor(params.norm_bg_area / params.hog_cell_size);
+	% given the norm BG area, which is the corresponding target w and h?
+ 	norm_target_sz_w = 0.75*params.norm_bg_area(2) - 0.25*params.norm_bg_area(1);
+ 	norm_target_sz_h = 0.75*params.norm_bg_area(1) - 0.25*params.norm_bg_area(2);
+%    norm_target_sz_w = params.target_sz(2) * params.norm_bg_area(2) / bg_area(2);
+%	norm_target_sz_h = params.target_sz(1) * params.norm_bg_area(1) / bg_area(1);
+    params.norm_target_sz = round([norm_target_sz_h norm_target_sz_w]);
+	% distance (on one side) between target and bg area
+	norm_pad = floor((params.norm_bg_area - params.norm_target_sz) / 2);
+	radius = min(norm_pad);
+	% norm_delta_area is the number of rectangles that are considered.
+	% it is the "sampling space" and the dimension of the final merged resposne
+	% it is squared to not privilege any particular direction
+	params.norm_delta_area = (2*radius+1) * [1, 1];
+	% Rectangle in which the integral images are computed.
+	% Grid of rectangles ( each of size norm_target_sz) has size norm_delta_area.
+	params.norm_pwp_search_area = params.norm_target_sz + params.norm_delta_area - 1;
+
+end
diff --git a/tracking/isToolboxAvailable.m b/tracking/isToolboxAvailable.m
new file mode 100644
index 0000000..3c1913e
--- /dev/null
+++ b/tracking/isToolboxAvailable.m
@@ -0,0 +1,41 @@
+function result = isToolboxAvailable(toolboxName,action)
+%% check if an specific toolbox exists
+%
+% This small piece of code checks to find whether a specific toolbox is
+% installed on MATLAB.
+%
+% The input is the name of the toolbox as a string
+%
+% The output is 1 if the toolbox is found and 0 otherwise
+%
+% Example:
+%
+%   result = isToolboxAvailable('Image Processing Toolbox','warning')
+%
+% -----------------------------------------------
+% code by: Reza Ahmadzadeh (reza.ahmadzadeh@iit.it)
+% -----------------------------------------------
+%
+%% check the input arguments
+if nargin < 1 || nargin > 2
+    error('usage: isToolboxAvailable(toolboxName,action)');
+elseif nargin == 1
+    action = 'warning';
+end;
+
+%% Find the toolbox and give proper output
+v_=ver;
+[installedToolboxes{1:length(v_)}] = deal(v_.Name);
+result = all(ismember(toolboxName,installedToolboxes));
+switch action
+    case 'error'
+        % if you need to end the program you can use the following line of code
+        assert(result,['Error! ' toolboxName ' is not installed!']);
+    case 'warning'
+        % if you need to just give a warning you can use the following line of code
+        if ~result
+            warning([toolboxName ' is not installed!']);
+        end
+    otherwise
+end
+end
\ No newline at end of file
diff --git a/tracking/mergeResponses.m b/tracking/mergeResponses.m
new file mode 100644
index 0000000..dd343bd
--- /dev/null
+++ b/tracking/mergeResponses.m
@@ -0,0 +1,21 @@
+function [response] = mergeResponses(response_cf, response_pwp, alpha, merge_method)
+%MERGERESPONSES interpolates the two responses with the hyperparameter ALPHA
+    if strcmp(merge_method, 'const_factor')
+        response = (1 - alpha) * response_cf + alpha * response_pwp;
+    end
+
+    if strcmp(merge_method, 'fit_gaussian')
+        [~, cov2D_cf] = fitGaussian(response_cf);
+        [~, cov2D_pwp] = fitGaussian(response_pwp);
+        response_u = ones(size(response_cf));
+        [~, cov2D_u] = fitGaussian(response_u);
+        w_cf = 1 - sqrt(det(cov2D_cf)/det(cov2D_u));
+        w_pwp = 1 - sqrt(det(cov2D_pwp)/det(cov2D_u));
+        sum_cf_pwp = w_cf + w_pwp;
+        w_cf = w_cf / sum_cf_pwp;
+        w_pwp = w_pwp / sum_cf_pwp;
+
+        response = w_cf * response_cf + w_pwp * response_pwp;
+        fprintf('w_cf: %.3f    w_pwp: %.3f\n', w_cf, w_pwp);
+    end
+end
diff --git a/tracking/params.txt b/tracking/params.txt
new file mode 100644
index 0000000..888c792
--- /dev/null
+++ b/tracking/params.txt
@@ -0,0 +1,36 @@
+% configurable params (as in CVPR paper)
+params.hog_cell_size = 4;
+params.fixed_area = 150^2;           % standard area to which we resize the target
+params.n_bins = 2^5;                            % number of bins for the color histograms (bg and fg models)
+params.learning_rate_pwp = 0.04;           % bg and fg color models learning rate 
+params.inner_padding = 0.2;             % defines inner area used to sample colors from the foreground
+params.output_sigma_factor = 1/16 ;             % standard deviation for the desired translation filter output
+params.lambda = 1e-3;      
+params.lambda_scale = 1e-3;                              % regularization weight
+params.learning_rate_cf = 0.01;            % HOG model learning rate
+params.merge_factor = 0.3;              % fixed interpolation factor - how to linearly combine the two responses
+params.learning_rate_scale = 0.01  %0.025;
+params.scale_sigma_factor = 1/16;
+params.scale_sigma = 0.1;
+params.num_scales = 33;
+params.scale_step = 1.02;
+
+params.num_rotation = 1;
+params.rotation_step = 1;
+params.rotation_sigma = 0.1;
+
+
+% fixed setup
+params.hog_scale_cell_size = 4;                % Default DSST=4
+params.scale_model_factor = 1.0;
+params.scale_model_max_area = 32*16;
+params.feature_type = 'fhog';
+params.scale_adaptation = true;
+params.rotation_adaptation = 0;
+params.grayscale_sequence = false;	% suppose that sequence is colour
+params.merge_method = 'const_factor';
+params.den_per_channel = false;
+
+%% debugging stuff
+params.visualization = 0;                               % show output bbox on frame
+params.visualization_dbg = 0;                           % show also segmentation, desired response and filter output
diff --git a/tracking/precision_plot.m b/tracking/precision_plot.m
new file mode 100644
index 0000000..63ec857
--- /dev/null
+++ b/tracking/precision_plot.m
@@ -0,0 +1,46 @@
+function precisions = precision_plot(positions, ground_truth, title, show)
+%PRECISION_PLOT
+%   Calculates precision for a series of distance thresholds (percentage of
+%   frames where the distance to the ground truth is within the threshold).
+%   The results are shown in a new figure if SHOW is true.
+%
+%   Accepts positions and ground truth as Nx2 matrices (for N frames), and
+%   a title string.
+%
+%   Joao F. Henriques, 2014
+%   http://www.isr.uc.pt/~henriques/
+
+	
+	max_threshold = 50;  %used for graphs in the paper
+	
+	
+	precisions = zeros(max_threshold, 1);
+	
+	if size(positions,1) ~= size(ground_truth,1),
+% 		fprintf('%12s - Number of ground truth frames does not match number of tracked frames.\n', title)
+		
+		%just ignore any extra frames, in either results or ground truth
+		n = min(size(positions,1), size(ground_truth,1));
+		positions(n+1:end,:) = [];
+		ground_truth(n+1:end,:) = [];
+	end
+	
+	%calculate distances to ground truth over all frames
+	distances = sqrt((positions(:,1) - ground_truth(:,1)).^2 + ...
+				 	 (positions(:,2) - ground_truth(:,2)).^2);
+	distances(isnan(distances)) = [];
+
+	%compute precisions
+	for p = 1:max_threshold,
+		precisions(p) = nnz(distances <= p) / numel(distances);
+	end
+	
+	%plot the precisions
+	if show == 1,
+		figure('Number','off', 'Name',['Precisions - ' title])
+		plot(precisions, 'k-', 'LineWidth',2)
+		xlabel('Threshold'), ylabel('Precision')
+	end
+	
+end
+
diff --git a/tracking/readParams.m b/tracking/readParams.m
new file mode 100644
index 0000000..8387b36
--- /dev/null
+++ b/tracking/readParams.m
@@ -0,0 +1,26 @@
+function params = readParams(params_file_name, hp_name, hp_value)
+% read params.txt and convert into struct
+	fparams = fopen(params_file_name);
+	C = textscan(fparams, '%s', 'Delimiter', '', 'CommentStyle', '%');
+	fclose(fparams);
+	% feed the lines one by one into eval to create the variables in the MATLAB workspace:
+	cellfun(@evalc, C{1}, 'UniformOutput', false);
+	% clear all the variables not in the params file
+	clear C fparams
+	% save current workspace in a struct
+	this_workspace = evalin('caller','who');
+	% save everything in a struct
+	for i=1:size(this_workspace,1)
+	    thisvar=evalin('caller', this_workspace{i});
+	    params.(this_workspace{i})=thisvar;
+	end
+
+	if(isfield(params,'ans'))
+	    rmfield(params,'ans');
+    end
+
+    if nargin>1
+        % if hyperparam is set: override default parameters with specified one
+        params.(hp_name) = hp_value;
+    end
+end
\ No newline at end of file
diff --git a/tracking/run_KCC.m b/tracking/run_KCC.m
new file mode 100644
index 0000000..a87068f
--- /dev/null
+++ b/tracking/run_KCC.m
@@ -0,0 +1,41 @@
+function results = run_KCC(seq, res_path, bSaveImage)
+% Entry point for the Wu - CVPR2013 benchmark
+
+    %% Read params.txt
+    params = readParams('params.txt');
+    params.img_files = seq.s_frames;
+    params.img_path = '';
+
+    im = imread(params.img_files{1});
+    % grayscale sequence? --> use 1D instead of 3D histograms
+    if(size(im,3)==1)
+        params.grayscale_sequence = true;
+    end
+
+    region = seq.init_rect;
+
+    if(numel(region)==8)
+        % polygon format (VOT14, VOT15)
+        [cx, cy, w, h] = getAxisAlignedBB(region);
+    else % rectangle format (WuCVPR13)
+        x = region(1);
+        y = region(2);
+        w = region(3);
+        h = region(4);
+        cx = x+w/2;
+        cy = y+h/2;
+    end
+
+    % init_pos is the centre of the initial bounding box
+    params.init_pos = [cy cx];
+    params.target_sz = round([h w]);
+
+    [params, bg_area, fg_area, area_resize_factor] = initializeAllAreas(im, params);
+
+	% in runTracker we do not output anything because it is just for debug
+	params.fout = -1;
+
+	% start the actual tracking
+	results = trackerMain_otb_wangchen(params, im, bg_area, fg_area, area_resize_factor);
+    fclose('all');
+end
\ No newline at end of file
diff --git a/tracking/trackerMain_otb_wangchen.m b/tracking/trackerMain_otb_wangchen.m
new file mode 100644
index 0000000..ddae5fd
--- /dev/null
+++ b/tracking/trackerMain_otb_wangchen.m
@@ -0,0 +1,276 @@
+function [results] = trackerMain_otb_wangchen(p, im, bg_area, fg_area, area_resize_factor)
+%TRACKERMAIN contains the main loop of the tracker, P contains all the parameters set in runTracker
+    %% INITIALIZATION
+   num_frames = numel(p.img_files);
+    % used for OTB-13 benchmark
+    OTB_rect_positions = zeros(num_frames, 4);
+	pos = p.init_pos;
+    target_sz = p.target_sz;
+	num_frames = numel(p.img_files);
+    % patch of the target + padding
+    patch_padded = getSubwindow(im, pos, p.norm_bg_area, bg_area);
+    % initialize hist model
+    new_pwp_model = true;
+    [bg_hist, fg_hist] = updateHistModel(new_pwp_model, patch_padded, bg_area, fg_area, target_sz, p.norm_bg_area, p.n_bins, p.grayscale_sequence);
+    new_pwp_model = false;
+    % Hann (cosine) window
+    if isToolboxAvailable('Signal Processing Toolbox')
+        hann_window = single(hann(p.cf_response_size(1)) * hann(p.cf_response_size(2))');
+    else
+        hann_window = single(myHann(p.cf_response_size(1)) * myHann(p.cf_response_size(2))');
+    end
+    % gaussian-shaped desired response, centred in (1,1)
+    % bandwidth proportional to target size
+    output_sigma = sqrt(prod(p.norm_target_sz)) * p.output_sigma_factor / p.hog_cell_size;
+    y = gaussianResponse(p.cf_response_size, output_sigma);
+    yf = fft2(y);
+    %% SCALE ADAPTATION INITIALIZATION
+    if p.scale_adaptation
+        % Code from DSST
+        scale_factor = 1;
+        base_target_sz = target_sz;
+        scale_sigma = sqrt(p.num_scales) * p.scale_sigma_factor;
+        ss = (1:p.num_scales) - ceil(p.num_scales/2);
+        
+        rotation_factor = 0;
+        ys=zeros(size(ss));
+        ys(1)=1;
+        ysf = (fft(ys));
+        
+        yr = zeros(p.num_rotation,1);
+        yr(1) = 1;
+        yrf = fft(yr);
+        rotation_factors = p.rotation_step*(-floor(p.num_rotation/2):floor(p.num_rotation/2));
+
+
+        if mod(p.num_scales,2) == 0
+            scale_window = single(hann(p.num_scales+1));
+            scale_window = scale_window(2:end);
+        else
+            scale_window = single(hann(p.num_scales));
+        end;
+
+        ss = 1:p.num_scales;
+        scale_factors = p.scale_step.^(ceil(p.num_scales/2) - ss);
+
+        if p.scale_model_factor^2 * prod(p.norm_target_sz) > p.scale_model_max_area
+            p.scale_model_factor = sqrt(p.scale_model_max_area/prod(p.norm_target_sz));
+        end
+
+        scale_model_sz = floor(p.norm_target_sz * p.scale_model_factor);
+        % find maximum and minimum scales
+        min_scale_factor = p.scale_step ^ ceil(log(max(5 ./ bg_area)) / log(p.scale_step));
+        max_scale_factor = p.scale_step ^ floor(log(min([size(im,1) size(im,2)] ./ target_sz)) / log(p.scale_step));
+    end
+
+    %% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %%
+    t_imread = 0;
+    %% MAIN LOOP
+    tic;
+    for frame = 1:num_frames
+         
+        if frame>1
+            tic_imread = tic;
+            im = imread([p.img_path p.img_files{frame}]);
+            t_imread = t_imread + toc(tic_imread);
+	    %% TESTING step
+            % extract patch of size bg_area and resize to norm_bg_area
+            im_patch_cf = getSubwindow(im, pos, p.norm_bg_area, bg_area);
+            pwp_search_area = round(p.norm_pwp_search_area / area_resize_factor);
+            % extract patch of size pwp_search_area and resize to norm_pwp_search_area
+            im_patch_pwp = getSubwindow(im, pos, p.norm_pwp_search_area, pwp_search_area);
+            % compute feature map
+            xt = getFeatureMap(im_patch_cf, p.feature_type, p.cf_response_size, p.hog_cell_size);
+            % apply Hann window
+            xt_windowed = bsxfun(@times, hann_window, xt);
+            % compute FFT
+            xtf = fft2(xt_windowed);
+            % Correlation between filter and test patch gives the response
+            % Solve diagonal system per pixel.
+            if p.den_per_channel
+                hf = hf_num ./ (hf_den + p.lambda);
+            else
+                hf = bsxfun(@rdivide, hf_num, sum(hf_den, 3)+p.lambda);
+            end
+            response_cf = ensure_real(ifft2(sum(conj(hf) .* xtf, 3)));
+
+            % Crop square search region (in feature pixels).
+            response_cf = cropFilterResponse(response_cf, ...
+                floor_odd(p.norm_delta_area / p.hog_cell_size));
+            if p.hog_cell_size > 1
+                % Scale up to match center likelihood resolution.
+                response_cf = mexResize(response_cf, p.norm_delta_area,'auto');
+            end
+
+            [likelihood_map] = getColourMap(im_patch_pwp, bg_hist, fg_hist, p.n_bins, p.grayscale_sequence);
+            % (TODO) in theory it should be at 0.5 (unseen colors shoud have max entropy)
+            likelihood_map(isnan(likelihood_map)) = 0;
+
+            % each pixel of response_pwp loosely represents the likelihood that
+            % the target (of size norm_target_sz) is centred on it
+            response_pwp = getCenterLikelihood(likelihood_map, p.norm_target_sz);
+
+            %% ESTIMATION
+            response = mergeResponses(response_cf, response_pwp, p.merge_factor, p.merge_method);
+            [row, col] = find(response == max(response(:)), 1);
+            center = (1+p.norm_delta_area) / 2;
+            pos = pos + ([row, col] - center) / area_resize_factor;
+            rect_position = [pos([2,1]) - target_sz([2,1])/2, target_sz([2,1])];
+
+            %% SCALE SPACE SEARCH
+            if p.scale_adaptation
+                current_patch = getScalePatch(im, pos, base_target_sz,  scale_factor, -rotation_factor, scale_model_sz, p.hog_scale_cell_size);
+                ksf = gaussian_correlation_scale_single(im_patch_scale, current_patch, p.scale_sigma);                    
+                scale_response = abs(ifft((model_hsf.*ksf)));
+                [max_value, recovered_scale] = max(scale_response(:));
+                
+                %set the scale
+                scale_factor = scale_factor / scale_factors(recovered_scale);
+
+                if scale_factor < min_scale_factor
+                    scale_factor = min_scale_factor;
+                elseif scale_factor > max_scale_factor
+                    scale_factor = max_scale_factor;
+                end
+                % use new scale to update bboxes for target, filter, bg and fg models
+                target_sz = round(base_target_sz * scale_factor);
+                avg_dim = sum(target_sz)/2;
+                bg_area = round(target_sz + avg_dim);
+                if(bg_area(2)>size(im,2)),  bg_area(2)=size(im,2)-1;    end
+                if(bg_area(1)>size(im,1)),  bg_area(1)=size(im,1)-1;    end
+
+                bg_area = bg_area - mod(bg_area - target_sz, 2);
+                fg_area = round(target_sz - avg_dim * p.inner_padding);
+                fg_area = fg_area + mod(bg_area - fg_area, 2);
+                % Compute the rectangle with (or close to) params.fixed_area and
+                % same aspect ratio as the target bboxgetScaleSubwindow_v1
+                area_resize_factor = sqrt(p.fixed_area/prod(bg_area));
+            end
+            
+            %% ROTATION SPACE SEARCH
+            if p.rotation_adaptation
+                current_rotation_patch = getRotationPatch(im, pos, base_target_sz,  scale_factor, -rotation_factor, p.hog_scale_cell_size);
+                krf = gaussian_correlation_scale_single(im_patch_rotation, current_rotation_patch, p.rotation_sigma);                    
+                rotation_response = abs(ifft((model_hrf.*krf)));
+                [max_value, recovered_rotation] = max(rotation_response(:));
+                
+                %set the scale
+                rotation_factor = rotation_factor + rotation_factors(recovered_rotation);
+                fprintf('frame %d: current sclae factor is %.4f; current rotation %.4f:\n', frame, scale_factor, rotation_factor)
+                
+                
+                if rotation_factor <= -180
+                    rotation_factor = rotation_factor + 360;
+                elseif rotation_factor > 180
+                    rotation_factor = rotation_factor - 360;
+                end
+            end
+
+            if p.visualization_dbg==1
+                mySubplot(2,1,5,1,im_patch_cf,'FG+BG','gray');
+                mySubplot(2,1,5,2,likelihood_map,'obj.likelihood','parula');
+                mySubplot(2,1,5,3,response_cf,'CF response','parula');
+                mySubplot(2,1,5,4,response_pwp,'center likelihood','parula');
+                mySubplot(2,1,5,5,response,'merged response','parula');
+                drawnow
+            end
+        end
+
+        %% TRAINING
+        % extract patch of size bg_area and resize to norm_bg_area
+        im_patch_bg = getSubwindow(im, pos, p.norm_bg_area, bg_area);
+        % compute feature map, of cf_response_size
+        xt = getFeatureMap(im_patch_bg, p.feature_type, p.cf_response_size, p.hog_cell_size);
+        % apply Hann window
+        xt = bsxfun(@times, hann_window, xt);
+        % compute FFT
+        xtf = fft2(xt);
+        %% FILTER UPDATE
+        % Compute expectations over circular shifts,
+        % therefore divide by number of pixels.
+		new_hf_num = bsxfun(@times, conj(yf), xtf) / prod(p.cf_response_size);
+		new_hf_den = (conj(xtf) .* xtf) / prod(p.cf_response_size);
+
+        if frame == 1
+            % first frame, train with a single image
+		    hf_den = new_hf_den;
+		    hf_num = new_hf_num;
+		else
+		    % subsequent frames, update the model by linear interpolation
+        	hf_den = (1 - p.learning_rate_cf) * hf_den + p.learning_rate_cf * new_hf_den;
+	   	 	hf_num = (1 - p.learning_rate_cf) * hf_num + p.learning_rate_cf * new_hf_num;
+
+            % BG/FG MODEL UPDATE
+            % patch of the target + padding
+            [bg_hist, fg_hist] = updateHistModel(new_pwp_model, im_patch_bg, bg_area, fg_area, target_sz, p.norm_bg_area, p.n_bins, p.grayscale_sequence, bg_hist, fg_hist, p.learning_rate_pwp);
+        end
+
+        %% SCALE UPDATE
+        if p.scale_adaptation
+            im_patch_scale_frame = getScaleSubwindow_v1(im, pos, base_target_sz, scale_factors*scale_factor, -rotation_factor, scale_model_sz, p.hog_scale_cell_size);
+
+            if frame == 1
+                im_patch_scale=im_patch_scale_frame;
+            else
+                im_patch_scale=(1 - p.learning_rate_scale)*im_patch_scale + p.learning_rate_scale*im_patch_scale_frame;
+            end
+            
+            ksf = gaussian_correlation_scale_single(im_patch_scale, im_patch_scale(:,uint8(size(im_patch_scale,2)/2)), p.scale_sigma);
+            model_hsf = ysf'./(ksf+0.1);
+        end
+        
+                %% ROTATION UPDATE
+        if p.rotation_adaptation
+            im_patch_rotation_frame = getRotationSubwindow(im, pos, base_target_sz, scale_factor, rotation_factor + rotation_factors, p.hog_scale_cell_size);
+
+            if frame == 1
+                im_patch_rotation=im_patch_rotation_frame;
+            else
+                im_patch_rotation=(1 - p.learning_rate_scale)*im_patch_rotation + p.learning_rate_scale*im_patch_rotation_frame;
+            end
+            
+            krf = gaussian_correlation_scale_single(im_patch_rotation, im_patch_rotation(:,uint8(size(im_patch_rotation,2)/2)), p.rotation_sigma);
+            model_hrf = yrf./(krf+0.1);
+        end
+
+        %% update bbox position
+        if frame==1, rect_position = [pos([2,1]) - target_sz([2,1])/2, target_sz([2,1])]; end
+
+        rect_position_padded = [pos([2,1]) - bg_area([2,1])/2, bg_area([2,1])];
+
+        OTB_rect_positions(frame,:) = rect_position;
+
+        if p.fout > 0,  fprintf(p.fout,'%.2f,%.2f,%.2f,%.2f\n', rect_position(1),rect_position(2),rect_position(3),rect_position(4));   end
+
+        %% VISUALIZATION
+        if p.visualization == 1
+           
+                figure(1)
+                imshow(im)
+                rectangle('Position',rect_position, 'LineWidth',2, 'EdgeColor','g');
+                rectangle('Position',rect_position_padded, 'LineWidth',2, 'LineStyle','--', 'EdgeColor','b');
+                drawnow
+           
+        end
+    end
+    elapsed_time = toc;
+    % save data for OTB-13 benchmark
+    results.type = 'rect';
+    results.res = OTB_rect_positions;
+    results.fps = num_frames/(elapsed_time - t_imread);
+end
+
+% Reimplementation of Hann window (in case signal processing toolbox is missing)
+function H = myHann(X)
+    H = .5*(1 - cos(2*pi*(0:X-1)'/(X-1)));
+end
+
+% We want odd regions so that the central pixel can be exact
+function y = floor_odd(x)
+    y = 2*floor((x-1) / 2) + 1;
+end
+
+function y = ensure_real(x)
+    assert(norm(imag(x(:))) <= 1e-5 * norm(real(x(:))));
+    y = real(x);
+end
\ No newline at end of file
diff --git a/tracking/unfreezeColors.m b/tracking/unfreezeColors.m
new file mode 100644
index 0000000..0334f4f
--- /dev/null
+++ b/tracking/unfreezeColors.m
@@ -0,0 +1,106 @@
+function unfreezeColors(h)
+% unfreezeColors  Restore colors of a plot to original indexed color. (v2.3)
+%
+%   Useful if you want to apply a new colormap to plots whose
+%       colors were previously frozen with freezeColors.
+%
+%   Usage:
+%       unfreezeColors          unfreezes all objects in current axis, 
+%       unfreezeColors(axh)     same, but works on axis axh. axh can be vector.
+%       unfreezeColors(figh)    same, but for all objects in figure figh.
+%
+%       Has no effect on objects on which freezeColors was not already called.
+%				(Note: if colorbars were frozen using cbfreeze, use cbfreeze('off') to 
+%       unfreeze them. See freezeColors for information on cbfreeze.)
+%
+%
+%   See also freezeColors, freezeColors_pub.html, cbfreeze.
+%
+%
+%   John Iversen (iversen@nsi.edu) 3/23/05
+%
+
+%   Changes:
+%   JRI 9/1/06 now restores any object with frozen CData;
+%              can unfreeze an entire figure at once.
+%   JRI 4/7/10 Change documentation for colorbars
+
+% Free for all uses, but please retain the following:
+%
+%   Original Author:
+%   John Iversen, 2005-10
+%   john_iversen@post.harvard.edu
+
+error(nargchk(0,1,nargin,'struct'))
+
+appdatacode = 'JRI__freezeColorsData';
+
+%default: operate on gca
+if nargin < 1,
+    h = gca;
+end
+
+if ~ishandle(h),
+     error('JRI:unfreezeColors:invalidHandle',...
+            'The argument must be a valid graphics handle to a figure or axis')
+end
+
+%if h is a figure, loop on its axes
+if strcmp(get(h,'type'),'figure'),
+    h = get(h,'children');
+end
+
+for h1 = h', %loop on axes
+
+    %process all children, acting only on those with saved CData
+    %   ( in appdata JRI__freezeColorsData)
+    ch = findobj(h1);
+    
+    for hh = ch',
+        
+        %some object handles may be invalidated when their parent changes 
+        %   (e.g. restoring colors of a scattergroup unfortunately changes
+        %   the handles of all its children). So, first check to make sure
+        %   it's a valid handle
+        if ishandle(hh)
+            if isappdata(hh,appdatacode),
+                ad = getappdata(hh,appdatacode);
+                %get oroginal cdata
+                %patches have to be handled separately (see note in freezeColors)
+                if ~strcmp(get(hh,'type'),'patch'),
+                    cdata = get(hh,'CData');
+                else
+                    cdata = get(hh,'faceVertexCData');
+                    cdata = permute(cdata,[1 3 2]);
+                end
+                indexed = ad{1};
+                scalemode = ad{2};
+                
+                %size consistency check
+                if all(size(indexed) == size(cdata(:,:,1))),
+                    %ok, restore indexed cdata
+                    if ~strcmp(get(hh,'type'),'patch'),
+                        set(hh,'CData',indexed);
+                    else
+                        set(hh,'faceVertexCData',indexed);
+                    end
+                    %restore cdatamapping, if needed
+                    g = get(hh);
+                    if isfield(g,'CDataMapping'),
+                        set(hh,'CDataMapping',scalemode);
+                    end
+                    %clear appdata
+                    rmappdata(hh,appdatacode)
+                else
+                    warning('JRI:unfreezeColors:internalCdataInconsistency',...
+                        ['Could not restore indexed data: it is the wrong size. ' ...
+                        'Were the axis contents changed since the call to freezeColors?'])
+                end
+                
+            end %test if has our appdata
+        end %test ishandle
+
+    end %loop on children
+
+end %loop on axes
+
diff --git a/tracking/updateHistModel.m b/tracking/updateHistModel.m
new file mode 100644
index 0000000..dfee53c
--- /dev/null
+++ b/tracking/updateHistModel.m
@@ -0,0 +1,32 @@
+function [bg_hist_new, fg_hist_new] = updateHistModel(new_model, patch, bg_area, fg_area, target_sz, norm_area, n_bins, grayscale_sequence, bg_hist, fg_hist, learning_rate_pwp)
+%UPDATEHISTMODEL create new models for foreground and background or update the current ones
+
+	% Get BG (frame around target_sz) and FG masks (inner portion of target_sz)
+	pad_offset1 = (bg_area-target_sz)/2; % we constrained the difference to be mod2, so we do not have to round here
+% 	assert(sum(pad_offset1==round(pad_offset1))==2, 'difference between bg_area and target_sz has to be even.');
+
+	bg_mask = true(bg_area); % init bg_mask
+	pad_offset1(pad_offset1<=0)=1;
+	bg_mask(pad_offset1(1)+1:end-pad_offset1(1), pad_offset1(2)+1:end-pad_offset1(2)) = false;
+
+	pad_offset2 = (bg_area-fg_area)/2; % we constrained the difference to be mod2, so we do not have to round here
+% 	assert(sum(pad_offset2==round(pad_offset2))==2, 'difference between bg_area and fg_area has to be even.');
+	fg_mask = false(bg_area); % init fg_mask
+	pad_offset2(pad_offset2<=0)=1;
+	fg_mask(pad_offset2(1)+1:end-pad_offset2(1), pad_offset2(2)+1:end-pad_offset2(2)) = true;
+
+	fg_mask = mexResize(fg_mask, norm_area, 'auto');
+	bg_mask = mexResize(bg_mask, norm_area, 'auto');
+
+	%% (TRAIN) BUILD THE MODEL
+	if new_model
+		% from scratch (frame=1)
+		bg_hist_new = computeHistogram(patch, bg_mask, n_bins, grayscale_sequence);
+		fg_hist_new = computeHistogram(patch, fg_mask, n_bins, grayscale_sequence);
+	else
+		% update the model
+		bg_hist_new = (1 - learning_rate_pwp)*bg_hist + learning_rate_pwp*computeHistogram(patch, bg_mask, n_bins, grayscale_sequence);
+		fg_hist_new = (1 - learning_rate_pwp)*fg_hist + learning_rate_pwp*computeHistogram(patch, fg_mask, n_bins, grayscale_sequence);
+	end
+
+end
\ No newline at end of file
diff --git a/tracking/utility/affparam2geom.m b/tracking/utility/affparam2geom.m
new file mode 100644
index 0000000..bef3e6d
--- /dev/null
+++ b/tracking/utility/affparam2geom.m
@@ -0,0 +1,39 @@
+function q = affparam2geom(p)
+% function q = affparam2geom(p)
+%
+% The functions affparam2geom and affparam2mat convert a 'geometric'
+% affine parameter to/from a matrix form (2x3 matrix).
+% 
+% affparam2geom converts a 2x3 matrix to 6 affine parameters
+% (x, y, th, scale, aspect, skew), and affparam2mat does the inverse.
+%
+%    p(6) : [p(1) p(3) p(4); p(2) p(5) p(6)]
+%    q(6) : [dx dy sc th sr phi]
+%
+% Reference "Multiple View Geometry in Computer Vision" by Richard
+% Hartley and Andrew Zisserman. 
+
+% Copyright (C) Jongwoo Lim and David Ross.  All rights reserved.
+
+A = [ p(3), p(4); p(5), p(6) ];
+%%  A = USVt = (UVt)(VSVt) = R(th)R(-phi)SR(phi)
+[U,S,V] = svd(A);
+if (det(U) < 0)
+  U = U(:,2:-1:1);  V = V(:,2:-1:1);  S = S(2:-1:1,2:-1:1);
+end
+q(1) = p(1);  q(2) = p(2);
+q(4) = atan2(U(2,1)*V(1,1)+U(2,2)*V(1,2), U(1,1)*V(1,1)+U(1,2)*V(1,2));
+
+phi = atan2(V(1,2),V(1,1));
+if (phi <= -pi/2)
+  c = cos(-pi/2); s = sin(-pi/2);
+  R = [c -s; s c];  V = V * R;  S = R'*S*R;
+end
+if (phi >= pi/2)
+  c = cos(pi/2); s = sin(pi/2);
+  R = [c -s; s c];  V = V * R;  S = R'*S*R;
+end
+q(3) = S(1,1);
+q(5) = S(2,2)/S(1,1);
+q(6) = atan2(V(1,2),V(1,1));
+q = reshape(q, size(p));
diff --git a/tracking/utility/affparam2mat.m b/tracking/utility/affparam2mat.m
new file mode 100644
index 0000000..e316da6
--- /dev/null
+++ b/tracking/utility/affparam2mat.m
@@ -0,0 +1,30 @@
+function q = affparam2mat(p)
+% function q = affparam2mat(p)
+%
+% The functions affparam2geom and affparam2mat convert a 'geometric'
+% affine parameter to/from a matrix form (2x3 matrix).
+% 
+% affparam2geom converts a 2x3 matrix to 6 affine parameters
+% (x, y, th, scale, aspect, skew), and affparam2mat does the inverse.
+%
+%    p(6,n) : [dx dy sc th sr phi]'
+%    q(6,n) : [q(1) q(3) q(4); q(2) q(5) q(6)]
+%
+% Reference "Multiple View Geometry in Computer Vision" by Richard
+% Hartley and Andrew Zisserman. 
+
+% Copyright (C) Jongwoo Lim and David Ross.  All rights reserved.
+
+
+sz = size(p);
+if (length(p(:)) == 6)
+  p = p(:);
+end
+s = p(3,:);  th = p(4,:);  r = p(5,:);  phi = p(6,:);
+cth = cos(th);  sth = sin(th);  cph = cos(phi);  sph = sin(phi);
+ccc = cth.*cph.*cph;  ccs = cth.*cph.*sph;  css = cth.*sph.*sph;
+scc = sth.*cph.*cph;  scs = sth.*cph.*sph;  sss = sth.*sph.*sph;
+q(1,:) = p(1,:);  q(2,:) = p(2,:);
+q(3,:) = s.*(ccc +scs +r.*(css -scs));  q(4,:) = s.*(r.*(ccs -scc) -ccs -sss);
+q(5,:) = s.*(scc -ccs +r.*(ccs +sss));  q(6,:) = s.*(r.*(ccc +scs) -scs +css);
+q = reshape(q, sz);
diff --git a/tracking/utility/affparaminv.m b/tracking/utility/affparaminv.m
new file mode 100644
index 0000000..ef86a12
--- /dev/null
+++ b/tracking/utility/affparaminv.m
@@ -0,0 +1,17 @@
+function q = affparaminv(p,q)
+% function q = affparaminv(p[, q])
+%
+%    p(6,n) : [dx dy sc th sr phi]'
+%    q(6,n) : [q(1) q(3) q(4); q(2) q(5) q(6)]
+
+% Copyright (C) Jongwoo Lim and David Ross.  All rights reserved.
+
+if (length(p) == 6)
+  p = p(:);
+end
+if (nargin > 1)
+  q = inv([p(3) p(4); p(5) p(6)]) * [q(1)-p(1) q(3:4); q(2)-p(2) q(5:6)];
+else
+  q = inv([p(3) p(4); p(5) p(6)]) * [-p(1) 1 0; -p(2) 0 1];
+end
+q = q([1,2,3,5,4,6]);
diff --git a/tracking/utility/affwarpimg.m b/tracking/utility/affwarpimg.m
new file mode 100644
index 0000000..df6a8e3
--- /dev/null
+++ b/tracking/utility/affwarpimg.m
@@ -0,0 +1,24 @@
+function wimg = affwarpimg(img, p, sz)
+% function wimg = affwarpimg(img, p, sz)
+%
+%    img(h,w)
+%    p(6,n) : mat format affine parameter
+%    sz(th,tw)
+%
+
+%% Copyright (C) Jongwoo Lim and David Ross.
+%% All rights reserved.
+
+
+if (nargin < 3)
+    sz = size(img);
+end
+if (size(p,1) == 1)
+    p = p(:);
+end
+w = sz(2);  h = sz(1);  n = size(p,2);
+[x,y] = meshgrid(1:w, 1:h);
+pos = reshape(cat(2, ones(h*w,1),x(:),y(:)) ...
+              * [p(1,:) p(2,:); p(3:4,:) p(5:6,:)], [h,w,n,2]);
+wimg = squeeze(interp2(img, pos(:,:,:,1), pos(:,:,:,2)));
+wimg(find(isnan(wimg))) = 0;
diff --git a/tracking/utility/box_to_p.m b/tracking/utility/box_to_p.m
new file mode 100644
index 0000000..cfcee01
--- /dev/null
+++ b/tracking/utility/box_to_p.m
@@ -0,0 +1,39 @@
+function [ p ] = box_to_p( varargin )
+%BOX_TO_P Summary of this function goes here
+%   Detailed explanation goes here
+% function p_to_box(width,height, param, properties)
+%                 ([width,height], param, properties)
+%
+%
+
+%----------------------------------------------------------
+% Process the input.
+%----------------------------------------------------------
+if (length(varargin{1}) == 2)
+  w = varargin{1}(1);
+  h = varargin{1}(2);
+  varargin(1) = [];
+else
+  [w,h] = deal(varargin{1:2});
+  varargin(1:2) = [];
+end
+%%***********获得宽，高***********%%
+
+if (length(varargin) < 1 )%|| any(length(varargin{1}) ~= 2))
+  M = [0,1,0; 0,0,1];
+else
+  corners = cell2mat(varargin{1});
+  varargin(1) = [];
+%  M = [p(1) p(3) p(4); p(2) p(5) p(6)];     %%affine变换参数
+end
+
+%----------------------------------------------------------
+% Draw the box.
+%----------------------------------------------------------
+W = [ 1,-w/2,-h/2; 1,w/2,-h/2; 1,w/2,h/2; 1,-w/2,h/2; 1,-w/2,-h/2 ]';
+M = corners/W;
+
+p = [M(1,1), M(2,1),M(1,2),M(1,3),M(2,2),M(2,3) ];        %%中心
+
+end
+
diff --git a/tracking/utility/interp2.cpp b/tracking/utility/interp2.cpp
new file mode 100644
index 0000000..914d76f
--- /dev/null
+++ b/tracking/utility/interp2.cpp
@@ -0,0 +1,110 @@
+// interp2.cpp
+//
+//   Copyright (C) Jongwoo Lim.
+//   All rights reserved.
+//
+//////////////////////////////////////////////////////////////////////
+
+#include "mex.h"
+#include <math.h>
+
+//#include "cvbase/cvimage.h"
+
+#define ARR(A,r,c,nr,nc)  (((r)<0 || (r)>=(nr) || (c)<0 || (c)>=(nc))? 0 : (A)[(c)*(nr)+(r)])
+
+inline double interp(double *img, int w, int h, double x, double y)
+{
+  register int x0 = (int) x, y0 = (int) y, x1 = x0+1, y1 = y0+1;
+  register double rx0 = x-x0, rx1 = 1-rx0, ry = y-y0;
+  return ((rx1*ARR(img,y0,x0,h,w) + rx0*ARR(img,y0,x1,h,w))*(1-ry) +
+          (rx1*ARR(img,y1,x0,h,w) + rx0*ARR(img,y1,x1,h,w))*ry);
+}
+
+//-----------------------------------------------------------------------------
+
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
+{
+  // check inputs
+  if (nrhs < 3 || nrhs > 6)
+    mexErrMsgTxt("invalid number of inputs: [x,y,]z,xi,yi[,method]");
+  if (nlhs != 1)
+    mexErrMsgTxt("invalid number of outputs: zi");
+
+  const char *_argin[] = { "x","y","z","xi","yi","method" };
+  const char **argin = (nrhs < 4)? &_argin[2] : _argin;
+  char buf[1024];
+
+  // Check input data
+  const mxArray *X = NULL, *Y = NULL, *Z, *XI, *YI;
+  if (nrhs < 4)
+    Z = prhs[0], XI = prhs[1], YI = prhs[2];
+  else
+    X = prhs[0], Y = prhs[1], Z = prhs[2], XI = prhs[3], YI = prhs[4];
+
+  if (nrhs == 4 || nrhs == 6)
+    nrhs--;  // Ignore method
+
+  for (int ri = 0; ri < nrhs; ri++) {
+    if (!prhs[ri])
+    { sprintf(buf,"invalid mxArray(%s)", argin[ri]);  mexErrMsgTxt(buf); }
+    if (!mxIsDouble(prhs[ri]))
+    { sprintf(buf,"data must be double (%s)", argin[ri]);  mexErrMsgTxt(buf); }
+  }
+  if (mxGetNumberOfDimensions(Z) != 2)
+    mexErrMsgTxt("invalid Z dimension: should be 2");
+  // Need to check whether dims of X,Y,Z matche
+  if (mxGetNumberOfElements(XI) != mxGetNumberOfElements(YI))
+    mexErrMsgTxt("XI and YI does not match");
+
+  // Build output data
+  const int *dimz = mxGetDimensions(Z), *dim = mxGetDimensions(XI);
+  int ndim = mxGetNumberOfDimensions(XI);
+
+  register double dx = -1, dy = -1, sx = 1, sy = 1;
+  if (nrhs >= 4) {
+    int nx = mxGetNumberOfElements(X), ny = mxGetNumberOfElements(Y);
+    double *pX = mxGetPr(X), *pY = mxGetPr(Y);
+    if (nx == dimz[1] && ny == dimz[0])
+      dx = -pX[0], sx = (dimz[1]-1)/(pX[nx-1]-pX[0]),
+      dy = -pY[0], sy = (dimz[0]-1)/(pY[ny-1]-pY[0]);
+    else if (nx == ny && ny == mxGetNumberOfElements(Z))
+      dx = -pX[0], sx = (dimz[1]-1)/(pX[(dimz[1]-1)*dimz[0]]-pX[0]),
+      dy = -pY[0], sy = (dimz[0]-1)/(pY[dimz[0]-1]-pY[0]);
+    else
+      mexErrMsgTxt("X and Y do not match with Z");
+  }
+  if (mxIsComplex(Z)) {
+    mxArray *ZI = plhs[0] = mxCreateNumericArray(ndim, dim, mxDOUBLE_CLASS, mxCOMPLEX);
+
+    // Load input data and fill output array
+    double *pZr = mxGetPr(Z), *pZi = mxGetPi(Z), *pXI = mxGetPr(XI), *pYI = mxGetPr(YI);
+    double *pZIr = mxGetPr(ZI), *pZIi = mxGetPi(ZI);
+
+//    cvImage<double> imgr(pZr, dimz[0], dimz[1]), imgi(pZi, dimz[0], dimz[1]);
+    int h = dimz[0], w = dimz[1];
+    for (int i = 0, n = mxGetNumberOfElements(XI); i < n; i++) {
+      pZIr[i] = interp(pZr, w,h, sx*(pXI[i]+dx), sy*(pYI[i]+dy));
+      pZIi[i] = interp(pZi, w,h, sx*(pXI[i]+dx), sy*(pYI[i]+dy));
+//      cvtPoint<double> p(sy*(pYI[i]+dy),sx*(pXI[i]+dx));
+//      pZIr[i] = imgr.interp<double,double>(p);
+//      pZIi[i] = imgi.interp<double,double>(p);
+    }
+  }
+  else {
+    mxArray *ZI = plhs[0] = mxCreateNumericArray(ndim, dim, mxDOUBLE_CLASS, mxREAL);
+
+    // Load input data and fill output array
+    double *pZ = mxGetPr(Z), *pXI = mxGetPr(XI), *pYI = mxGetPr(YI);
+    double *pZI = mxGetPr(ZI);
+
+//    cvImage<double> img(pZ, dimz[0], dimz[1]);
+    int h = dimz[0], w = dimz[1];
+    for (int i = 0, n = mxGetNumberOfElements(XI); i < n; i++) {
+      pZI[i] = interp(pZ, w,h, sx*(pXI[i]+dx), sy*(pYI[i]+dy));
+//      cvtPoint<double> p(sy*(pYI[i]+dy),sx*(pXI[i]+dx));
+//      pZI[i] = img.interp<double,double>(p);
+    }
+  }
+}
+
+//-----------------------------------------------------------------------------
diff --git a/tracking/utility/interp2.dll b/tracking/utility/interp2.dll
new file mode 100644
index 0000000000000000000000000000000000000000..921b8d0b085b15a3a3cd5886356e2f4065345cc9
GIT binary patch
literal 24576
zcmeI4e{5UT700g=r*&JWc3C=QtYWX|wFH(`4oypJ)8VuZYatXbN%Nyk+Bmg?-PuX7
zpPAdG709)yc{Kur27j<gj7^$0twD%E;E&A>0haW~LMK?P##RU2_L9(wR)dC#*886K
zZ0c<M5E7dvbf0>EoO91T_q@-&=Q;Ijde`<xSSe#nkCMwVHjI`ohgyGjHHqP+-+Fl|
zdw$W0mBYHm6DwPM{;)6*>c1!C4hp;7eSQ6s(Bl<Ckv_rSCphnJ7J~gA@0#-R#bVy|
zmP`94%YR%U7s|sAT`L<Q_av^94QgwaP3rh^c?H`3#KOE>q1s)1!7ML<zRADaM?MM>
zyBZj4)RnPc{{7^}g6<4kx^jtb0b^!p9Gdpq*P<m^P-U~)hfF-^DzvN{0@>jzQELdQ
zTxj!du4n8Hj5i95c~rYksa%Xb@+m`~sTSb-GYr;9-l&A}^NUeYe1%dpM_{abO~~Vx
z+>Gs^n}I^LqqgXBfYKVxjH&x&!65mgw&-#&mQw9%7^~1>01DNP+M>%rb|G(|e>Y<z
zq<5i^KWd9E2V)zHMBe2DH~~(86W|0m0ZxDu-~>1UPJk2m5(qTQQ?aviZLQlLc=sv%
zT~`hp<R6%0mrDoBVwV>jSmLX7Ftf605n~;(D>`Y3&qwlFl3zb4-j4t(p@?19Lplg)
z1tfG2p(|+J!{~}y_YrifwC>~RR-?OrQSx!sdUhc?_e1*<3QIJK9`&1ugj4L?zf8Wf
zLaveb>g9R^6w%om<V_WFojjnIH=7fKqDPJz6Oz~|4;U~pDpBZ>(%7rHwW;{I1HnW<
z^vMeBBp40I;fh3y7?3ND#7~(HE;KOKq^w2VVR~UAekyX4jKSn#@kvCHtNoSlDFpfb
z-116+nGxX6p>0xXVICVb!_pUfRcD_Wd`_N`;a-OOR^{D=1Z%DYu!b>NSwPO)l+#cn
z_*&U14!V>RkYb|-xUkqKU^zvWQ#7+xc`<L<ru;_hb+;)`XuZ1Fs4-J!hH;@bsC}x_
z$1<Jb=JXh97fKz90~@oS;<O|;>alp$^zHcTIU?}5lsh$#0Z830WrHS3kOD5nes#>}
zQofTPqg!^07P{@^X1(mt19p{tB61HcAUA^PHTLmHha3>QA^9>5jQ%)J%w^wBpP|@d
zsa$L%m(Yu+^&&iZGAuow{ZNgSsTMBmjN!CXBrbuBa{@E77wA5w9|m|@nBoJla^D;Q
z`!&AH>6neQDhH@TlMv2XwA~bFrX*7$8xW_qdeWD}Nl@J<t6Ax}q%%He8aj$_J^4`e
zCE2s-W7&UF%;|B8D?6PY(;`dO=KxTht%cO8e3KS6y^xBhOb35~@6h(>_AH>ao3MDD
zODVyytjNvA*2GQ9bQyMEOPE&WJP9}kZOSQ4O_?J<E>co7DZk5)Zd0DWI=W4HNmI7~
zZ95il$;)KR+H(mlY1))Q497;R5zE7?$ub~Xutbzje8wqm&a<;j=|`p-aW=?7qu4p6
zKMK)-n;aC)rWc)}xttl*;40GrvGR0_D4;cik??so*fs$!h3s^Ja&{(DMhA29?c_FH
zvR)S-8>~QzwA)YaH^$Em8Yo@t=LZd$GJIiVo;5=b$d^q+gJOH>hZ+5qdG@FB_HWIz
zZ%5h%OhXq-FJ$zW;;Fp<*FMqzGSg5%Y)?w%I@}IU*ZWMJ0G<6+OrQmv2x#K!0;ZuB
zm31K5vRIeVPo$4c9#4-?o=lHr(PmF%kK(i(ogC9{6Q_JEJDxoX2;IPKQL*FEqq4?h
z3@2*E8&l3V8Ers%Y)U`z-$bt=*CNYr!#neUe)0woP>*ToeGO>whJE@5!2XMGfai;B
zfM@;&{4u`)C$RzJlkWifD|z%zv6|3h3!<476R{DU{S7=)ET+SAALVk><t0c`HZMuL
zp!zsToe(H#+7ui04vP0{c}bsJP0DI*ba0!pBHzO@E}WUo`BZG7ng758Ha03!9{v<v
zIt0jvD!8>E7pv1LOrob)m6m#1CgP!?C8k5$s`U8h<Rj(f=j5!KkH}l(WB#F<ox!P_
zN#3DLHt3LlS}vZ?AAQ7TFU~FcRJ~{#Lb{f|Pq|fWKc?nlwf*GPRH3}qj?M+mzBm`P
zqchLHmRgdCoxKUE^_H57$Ts^3;vxUoN&6qJ`Zo>HW3L7F&KqGrf>b<-y+4k9N1H_&
zqK8v3HJ(q#zhuYMRP4m0DtNOX`EVY<0v}nX!_t*pE;ISSDE-XSwR4yJ&W=}LSp30%
z2+yQ{|KPXdAke2^#gXaPFXeK%RZfAKbE_O`eXCm6sdcSdThzK*t*g{pKz)4Ph0qM~
zRAi~O%xq!i18Z{S)XL=sSjOgb@}fj}1r(9nd~_Id%Z<>-QytH&Af<k7YPL5fu%}7#
z>8ZkZ^X6wNQ7JQIuX^q1nB9Oitg2Z&C0XQhvS>}-`;q;<&^7YX_?WaX_l7zgJFlOe
zapjg*(E6*$FYDO3N7p`LWZUF-o8>M;s~j=5$@|Um4<d_V8z@>fyPAGDEQ?(!iCq7*
z01Xyd)YjUZn?T2ij<OuI5CR3fmE>k3c}16KFh-2f-_qu4RxV&q;vbNk(#UEuIzx(7
zq;zG2F}B~xq-6!)2Gutv%MHd|J9l&;xZEl;MLZ}Pbdj&cHc${OvN$J+`rLAc5$zVR
zD!(`ZPJk2O1XKbSS85>MMmvde9Ay;c2nt<~<=X>jdr{g^wxZmMBBCrt`AASr*<`*Y
z=@9hz%wrQpx}QNOKprmA4MT^|I`(XlZt3-`@ID@$`PY{juUi+4ZV5a4BRv5xV>U?N
z^GaJI!5(ku?p~)q=<N&p`}@MoK7XVk;Kh&>W|f-0$a15mafMhXdRxMq`hx=jZ<I~4
zV02T+>z2H5?+y8P*M~yxy^Pi8?fop^?+HfNcme@d$AaEyLnw4t_?{Mel}{sHZHPTr
z5*`To`=nlWS99y8rWSPq+gZ{Z^5PCmx?TQ0za)7>LH1}#&>c|UwEu+da(nhf!qTo@
zkN-aWxWk<Oup77B;zLB9Mt@Jp9ooAI)?s#ML(|rV#<e$VehcIa|AwLo`=W!gP&Pwu
zL~+r6;XVFNY8k@wh6Xej2>kX)-);#QX^3B(04Kl+Z~~kFC%_4C0-OLRzzJ{yoB$`l
z32*|O04Kl+Z~~kFC%_4C0-OLRzzJ{y|04qL-KH({-tVL6-~~0?hZ6Opj}r7szJ3qA
z;l0lt@Oy+l^;Nmh--}=FABad{o3JBVy|=n+Uv<=9z1P2^TAR>CKL2x0e?-zI?St!&
z$AeMz)p*sa@QsXq^;`(z%{BcW0v`1>y`T;h*>4v-g%|ZUA?)jq1UzJYGp*UaMR50d
zgpMtk><tU_pAiJz(r%wRrVU}_s$Te^U-AjtX_#M}04Kl+Z~~kFC%_4C0-OLRzzJ{y
zoWNI$0RG)J{eGihu~_OYcU!hwyp{pWk1RsX4K-ag57s<i^Kwn9=CzvPn|Ia5Ykyn&
z*V?bIyM5i6b$?s8VEy{_z3ZP?|L85t>gwz6t$VuePj%n03fAvhoz^?7E^Di`)7oS8
zS^KOHSs$}LX?@1}taaG>n)R4<!uqE5jP;!LlGSLdv|VR&*qUrdY~!{++s@j)T5$Q@
fZ~~kFC%_4C0-OLRzzJ{yoB$`l32*}c9|Zmd!H0Zo

literal 0
HcmV?d00001

diff --git a/tracking/utility/interp2.mexa64 b/tracking/utility/interp2.mexa64
new file mode 100644
index 0000000000000000000000000000000000000000..989d748af536adb288406acc173db745de914db4
GIT binary patch
literal 12982
zcmdT~dvH|Oc|W@=34wuKUbX@_Y$>?dmKTXv5cr`htzg#(Eb|D$U@xoPm9(+6tJ%G)
z2<iup*w!0$iP|w!8c!3_$CUh$z_j($acWl}9&+kw6*6&XXF6skalHy@47iC^V|)Ai
z&OP7i?$s(L4xP3KR^R!4-{U;*Ip?1HJ=@g0$!@a=P7d*)Anx*X1F2VnO1H{@)Qk0^
z9Dm*7PN`cyA)1t%Uu8f{g>XrL%BaAvO3ACL44mm@yJ448l<kNmyLx3;uk4t*ls?l+
z!;NTpTMw5Zuis$+R=QZ}F(q9}kHU`1yV=C$&4Iz4%4oi#p8;KnA5~$t!e%K-&2g&y
z?Vw~&by=jaPfKQgg5g}q=P4Z;fd!wSPhy9_<WTd&tnLuE73NJbK2!LpQI_NovqWo2
ze0E`+%q9->(4ZQ>V%&v!WXzvpe2!RA5}zy9m&9)umll@tbB7qW;NL<))OVbCtak!1
z<2<InYr(fG!wx(MbHueg`(gT*7X2AY@UFyf6CbL9;`#HCMSrQ%zi7eP54DnYNOgq%
zAbkfvbN$!TgNhWu&y1%fD#U7u-zJ_?j-S5G5Hwr#k0|{mpl+2nlYa$c_Wxz|0REu`
z&#3wkapN~B1p*!2@mL_K2NQZAAOe9<Uog-fjRm`+2O}czaBpCHq$8TtBZ<bYU@{p=
zitb3?rc^AXN8>TdHzg8VlN~$zbdl^yL}Plp=<f3;z426AR|MJ35q(>#yDgG<q}?0s
zj>M8MO9nnyQ&TK!>WV<9!#i1&(rigc$Bpsso~}qA@{Nf|P>;ZQBoPg1iA3<ggxi)V
zTGG+i7wCy3l5uK`ejw0`G9!JT2SQcvinew4d6IFNi4<}{U13Q*9J=>jFrIZ%rI70g
zg#t<9Y6<VrJ5i4?@ioMU3_Wy2B+;|Vg93x=)`+c5t)j8DRcr7EYCWsP=4O9GV_=nM
zm1nK!{TvLZOb)Cb4q=zSGW;!9JnbPB$aAt`XQ6**dj{>ae$ralY8rbxWG5aKg4b4h
zjwq4M9mv!CcH{ShlTX5Rg?Q2=-)54h7Dlm>m5W}J++mWtzB&zF%EXsU@(+OUDq+YS
z1v#%DybmycS*;(8*QoX9ss-P#*0&ag)84~vL{_w#aN18e?Kk1PhmkC&aH>A-J)90H
zoLXnzBc3<myta_dpb6)_i|}J6oco0E<0jnvT%R!Eyzh|yDHCpP@5?6Kyf+P*aPt;2
zV#3uCf!vEx6K>v9ef>YH%>4)oxdwDm?j!Vtum5c2aDMIfIWCc3z7u~3k@Cxl?<7}a
zm-9%=A10o<b!=4f4aC#wV{AzB8;GX~J$6d+tB9vnZS1(@?;)P1{n()7ml99IJ$6X)
z3yG)U8S9b!T;gf^jO~~F4B}}Q9&3@jop@S>$Lhg@%0%t+4V0%`F5ycL>Gu4nY#ud9
zrb7>0gGlW1W!_G0^z}cmkvM&kZ{UGfkjq8j>dGu%|4^kbQ+^6a{vy2n6*9g7moI%W
z^^@}ww!J>w-ad%`yWr!kS?~de`a<Zw=QK*5?aMxJ0Rn$!_;oG+<JX4=p__SI+pT?G
z+okQ?;d}09Rp^b(4PWM^%OB+PTINlE=AX6frjPo^9lbMQ{Y`8wxxf85pYOk6(+L<S
z;Lid0vgI>q;4)3eec1&zncZ~)C}{va1@v8+-SslicMNn0=vxLlqV<35aDHVF-RB>e
z=F6nTjeOoeup6yDS}C)eT|iyP`YxoeqI%T1$c}S1yXzEmb5!yEkIS4d41rWQUpO;;
z^nRPx|7N}^{Y&SOKY@|<Dw%43Pk(*c^l@^TNne8Dmh9qpeEl0mYNc-=O`aYen19)q
zUE4r)9C(#fUM8!Ype7qY?k_%~{Ki}t=;0x&=_0G?IjJsv<)X-I`E&mK*`~Umr=EjH
zYIW|90QZmEWe>;DAvw*kp?=6MpnPT+y3`lB3S#SqrS1snO1W}G{{5M=(3Pf-Nz-%{
zHCREv14zqUpw54GZM$=B`D^fT-8m=ih7!>gpoD*0;WY|h4>~r6hYy};G;&XWn9s{<
z(SYRQ;Izy+7WRPkXU>C_(jZuLc|Cmiv)*N8$R9&q+qiF8nKS)0^l;P0>SfOK3&^!=
z&bgmOP5yKQetiQ+C2_v;!9y5qs${ki^}J`)7aAr_@A>iT&V9qzollIy<|x$`LrDdp
z5iWvCxS=*`*(UNc;t!$HN0ov140$^PZzC|x$~ON3b@J9-NzF%N1PpVBsf^GMe067Z
z4Ufb=EC{(lqrmbWFtD7|oBI`B)SE9(Zf#n%b{wshPEa@MsR4U_L>rjDj9TtIdB%5)
z+Fo}4rkY;llbW7wkS+Ixu&{7fc}FL;oZG8o91agT4=1paeuq}yS7^!f4ZKtZU*1eU
zH%tRXQ{gDBBEF2BEx#byv=4d=ebIng45$($=QW^K19A~mYe4%6%Dk0ZYOtM(U2U*E
zV6}n8%bEtlm{Dqezw=-De0s=v<Qdda`v=rhyytvu?GLf2=SO_`A@u4fthMiwWLtLr
zM?}#JRMOV$8apP9d;+d(&ejSrM`;RNckUTQ<EVmLv(<A+$4xMFngqD4F{*noXR5E^
z_dRcRo6VcmY#1#Jk?UVDOqW3a9rVY*P9y(Q<a3$ty?^oaPnXpFqx0)S?$OKz{qD@w
z_eW3XmS}bV>HPYc!=rjt=92T|`D~fZcMLV4jzoDg&s=dfXXNuirPGu_72l{ooY(L5
zW=5okrn)x_51+|g(wrwxTU?YCxsWpiF7Oz&q7g2dSLZ;#Pc8Rm&d6GN>xNOy%xU-c
z^?R(<L?NiAmMs(JdX_B2g#z}xarkw;dQw#v8dbe1H?pqi$UUfM>lut7O%9C6oU*CM
z@PBXMXuAJ~1E*Y+?9*E|aKEaMUt{2Ym77Te*HaicfBr&q-K+F`IuE~wW8K(ZY#a1;
zJZZ&w^p{htIOD&dedhwN24i6H$oHWLKrWEnn+9}<AS^8x4EA}jlh%?KA<Lv+fQVL<
zYSc)3jcQ|@Cum)H8i&ooTJj$xqdB&Md=I3li0ci`)-o5s!a~B65RWreODfQka_bh6
z-?l!Tb0-!N+BRlhm~#18KIxg9XIehwo>OG}>n<G2tW|^i!oo47`@+n-Yf|;*o0p6!
z9ZZ^YG~DR<Ib|*_f8L)z=dF8*hFRY-wgv}1|G*tMYtSLbKd>3+QF-LV(G=&=qXKF`
ztAJL@>~1$unpwVV({Z3Vn2kM2v+B9Cbh_TF?b9CDp3nkx{Fdi$oY2S0Ttf2nLO^Ho
z9qOfF+=ioj|9Bqf`0Slmz~Fq3_lJSTabN$HN|L0m(%ZydSj+dQlkcG<dV3Ng8tcU+
zMc5sacO~w4yE_`|N$JTA?!A502dW=GSlt(`J`mko-5t?8<Kf3?x>91h$J}AQUD@DH
zcE(d(VRu`^y~@y<(3|idPw7Sh2cx1-91sUZUsN23;);X{C^+u!lUFM%mM8BP;h-LL
zcc+p%xem)49=8;<`rW};*uBROG?H}3;<~##sE0a5tAPWa0Nr1V>YeV#sAh*AGU|6#
zZuUZ4eQ2w;vsr6!Z?9eBjs^8-Zv+aqMUIX1Zb~nUmFTkjL7&G^KLyI~?zV$1+l8&q
zRq^5`hf6G)IjwS#PCH_m`klR=&&#_uZMN(8?2WT3TKcE`K}ELwNZD89B?`Yo+YU-M
zHugaN=~ssG<oi9_12ZRjrKH4|!A<at%SwTvI}u9E(;bSLpJn0cj>as%%fi=K_y;XK
z-QiFw++E;IONs9+_`X?PBNpzC=!Kios1nfWkx~^61`=L0Wcq#OA1_cc@m+%KTp*q5
zK0%p@N~O>JPA~41__}&r0gSh(4BrQxFw_gz>RcW+oj8e(sfyx_N%CAShEejPN{(@o
zJh%I*;#-v*CsHZrcBo?i#|eL*u@>Y}<MyDU+Z265(U_voD0)QEuPOR%MPF0&pA@~M
z#)t3OxLm$p6Pp_wH@H`9-nQ#*_iE1?&pP*tU2UnDo-(-mYEzP1b+<5UuJP1(Ru?oS
z$C_*ZFKK>4i?}mDoAEQ|k}(sl&kavV)vbl^J3OyVu4MLp>5p#h4S$cS%Q=4T<|#+a
zTkW&bo>1e7TG^J0c7=_<;gub9A${K}&uY&)1G-O7L?SDL-Qjg>?jzaCP%^bLl+dG|
z9pa`;%^qPZFd^L94)-cg?K)3gt*{+3DnwVCl{R`uzA3MhD-^qB6|&6{t_nxxqRQo!
z^QafC{2v3e3tWv0o9)KzS@<hI?_zzjeTupZp?XuUZ`18OKOD<!_A>Yu)9pM(996b0
z_HqNetcKY!!?q>$z{;8}c3jNcp1$=}qkdIY<EiFwq8?Nn>c-P;p+Y>#1KoPCO?Z-e
z!l0d;>kKA4g(rL<1{H(q3E}C8r98dr3yYZvgo|VWs3sy^L9$o5o-SQ@qOmCc>JfZR
z^57#>6tXy-I)o?E8E8)gyCZ?lFl?C<o{%0-BvHOWz7T>CK7-xS5DJXL4gPqfZ^KVp
z65647H+^k_YNRdI5eO!Nu?~FMQfxHV9_M_ptt}DhWv(k4i(pJ{L83c#nvJ|K@wXbL
zL#hmx^S;OX)Nesfdo9a(KW6%}641Uydp_IqK8e*=em!D)-rt$>ehxL6<eXs!v=%?w
zJK3K1OQsu?F3YX-LC9%uW}NqHrt6hG)t_{@{t6abk)geu?Rh_EdPpgfJ@uK@exJhG
zp3fIdPbhnqvpx4;yTzW*D@^%(L-uq=vfB3mqdrx~jEVRFwGePQP?3qpo%tt`G22(M
zouWgms3_ah!c3_bKeIhAKukLe2Suq~_<m<_LOf-$=XH*$`I{HnaQSTiC1ub3%W(}%
z7g{D9_aEDT*<wFcd<yGZ%RgqZ=Qs(bhuKh3Yx#es?D<^C@e~}N!~M)|toGlu*mK;)
z3N`SbWPi$H&+!@SmA$n*EB_MgY3_0T`Mk+-LL3joa<=Ds(On7+5zG1eYm545YMsAq
z&ve*g&+iAF%DxA3b7DKDZ$M}+pU<zCl)c+*j6B;hd{JRs|8aHx9ar{Nee3-B5$r8f
zRJvyS54j*kt-6<iTQRA+*f1D0NZ~K{I|J_`Q_6l&*)O#SEv&dPh5cXOX_VkLNa3&5
zemVdvM#>kJ{bGyI!pe9p?)yR^$o%}X2DiXsU}0t1ympxFD_0oYA@#nv1T0-L@ObjP
zq-C;H`OmL3>^&BtwR}UX^j`;R$NRWhgW@=Ne=f!yh4<@XyiD+ZT8x($-fxTX3c>qd
zF+NT3zE+G^3f`}ZarBPbM~d+og7=wX{Qr9#+dRZaTAwhQAM*vTYsKwdP`J-1#&13z
zj_e$QXQSB*;pRw{37MlPBE^i;X>xLM6mG|gWsX9mr4?5XUonjO&{I^L=y@nB5r;>|
zhEn)u;Pw*nj=QJOZv$?_se|`burkrUV8iOc`=x@#GZLR1r#K+-iFibQj=w7L+r?$I
z?z5k7AVc*mTx==9-zmI>FT4~!RsDw{r~11s^J-M$bA<Jt<=-Uk6fVn@xFYe%?Rp<L
z)zf7T=8?Lk`=9O;OO1<f3Vi1jco4Wv+%B%Fe&cmDO1S-6-ufIRBtH4Mcow*P6Sk}~
zX{DdD`2PlQf%(S!G;x%^CH+i}m;IY5^k0$wujTo?VR*s)rqUOSir0Z30yoEFsqIa?
zW5J!u(YOWQuJHNTtV+eL(tcp}&(rHJiBFDOt+nVcRsOveoc(ME&N`$zqV!`H{q?f_
zLS!v?TIrvZ_~f|N?^^UZZgmv6TiKe)Pb~Tzw_1+<*4z&qx4Hy)DL?o=5!LC2z!Unk
zS4?#06$*BBnF3(x@+PUL+S}nZ5NO=_$o4?9f5%Q47wawJ9&Rj%aG8LJR6P>*tXW@E
zD|(33z&soebachrf?WZ5IS~k^`h;<X5l(e?AK(I-ws|M~GqY*CwzVn192Lg-39s_z
zLwfxb_f&N0M|VLSryA<%q2O3`b%W;fB~HM*XPdUw-$)Hg#si)B&V#t#fV`;j?#3lf
zyQgHNZ6JV~79*@!P&_dLawGBvCa}5rkp`_f@W`f3JDPR|c4`gHO^7>=#KPwI-uikO
z7dsJ3TmVb~#T+!OWE5m8!-rcEtUVF#YmQWADTNxhL={&Y`CMObj?tYe*x2l-bd)n|
gm|~_?m~?Ror2?%>^=o~-DU7&8taPC<MN!@V4gT390ssI2

literal 0
HcmV?d00001

diff --git a/tracking/utility/interp2.mexglx b/tracking/utility/interp2.mexglx
new file mode 100644
index 0000000000000000000000000000000000000000..e1c5d6a6c35e56829cc3f1a42330161af3b2d1fa
GIT binary patch
literal 10719
zcmd5?eRNdEb)Wr^kb{h5EZ4N)u-L+dMl6z02r`i91A)omLk1aR6Y{d!U1>L2?W((P
zVa2gy$nsJ)k9BDiC8Wg-K8{cAlMg!w)F!Pt6^4MZ6PH6dt`tj@riqWhu8*-3n>r-w
z?>Fzw>Pe^s$Nr<UXZFtT-aB{he7-mH=HZq#D_t&^P{AWSg6O~$A?iRYW?7Ot;S)t-
zo|q#(D&<O_bk8kTEFvfuQZ9LO!DocQy@<uO+zomdc@9IySOOU%<zB?rU+SSc$W^QW
zkC8IQdmsn7DyaN7LGn>s`o}h%3CgO>K!l4jl-cx8K$!<c%5XFOZ$q4B!4Ld@<25Y&
z5oD~}t+vkeYv>-A4Jgn<H%Kk^5raoem9&(G5Imw((y9WsJYah=c<h?pUBFgLu&*q+
z)k6CUn5`qk{sOkVJW#-FF@!(0h{vu?c>&K5Pgo*xJ1BGKZM@vZZ8oNial4J%ZM+5J
zLxr2frO`)iz?xJ4zYw2=8*I<%XpbL&GF(JH$~B3EY$@VT*|-Um_=Js>y+y!nB|>z7
zL3?V<>;<27SYzwgfHMDV8<Y4DUt$8mj{!@$!+j|in%3Esh-)c5kkmB|FqjEw9np9o
z7TpmR+NZZ^8^WE@lpant#{#KTI3>EmnU(2yP>&|!Ot&PHYg3(V8C|5hlhL@|A-XcF
zQ!5hb_E;FnRbhQyx~n~$T;H)G+7*tckPSAR{FYc40^R0L)M`yigUyMq?pQd3baOHs
z(8DkqPDX=`$z-6%GS?bqWilD9JDf}<SOL99+g42<H-STUW6}1mj6ao-iEtqm&=n8V
zznIdL9twW!V+Gh53~DK&8VbAh2z&?;R}%*<JybrN?5_2rz|L4ktZmsWnm2E5Y+9|=
z_!o&)YgRWkYqiLsQsi#?1Q^!{TnG7c+tGuShp0%dQr5(Cp;!5f?AYh#RO9-x#_nhm
zH)K!VW6f_CrN%Jx+yVuTDIem$J^h6uA*Kkekk>7mY&o4Qk47<pi^b<{UTSs}3qX<3
z3;HhcF7*(7q<=&Iyke(A{s#Go43^wekaHemsEdO({kBb?wCQ(k`ixCE&#CXVX{k-0
zy%Szy!R;r-f_t143uy`DE;P<r;3=rd04WyNAyO<o=Skt*FGw+2Mo7`1uajbpdXp3r
z`&Xn`C@+(uoGYYQyzkBZGW&YT$WP(Y$ZFJiWEAV)sbgaF5T6#KkJDF%(Z{(iMwq;C
zeiY$xwn&8H3#TP!FGk3`&@VBc77_kkct&FOaKuM^P+}${^N9CJTuEF_yi;O|B1?$7
zC9WcFApV@hRE^wC+$wP$aVv2HFzU5&gb%J=$L_%Tl6~E4X0L#6^y}`?)7j(REnAP7
zvnr2K55kc*A=`EAG9G<9|CvTbdWo6+m6c(pnWY~&Z{>&q&TCV%{oY)58CaudOi9gh
z<4pR(i7a!vEaS@XH>BZj!La$&)4ziQXY^aluU@9*Gy2qgmONM<JyXbb-=>YbUY}Kj
zzStu32AQ{|JUimaUiNI8md{@PZ<%Fe09msSmAG{9w%fep;N>40$R9e3#FZZv!=bmI
z$Uao=Df`Q(K0-ScrpOtyelt31-dkaA@aB4bOs}mlSC{1Sl}Mj7D)5hd;-oRdPUX|*
zbA9Ciu;kT?@n4+wX1mH=NQtsXo`mwhFh|}PfKtvV<^Pd>{G*h7-9KeNDbkC~JS~{s
zALUk5mgKX&^ms+3*R-_JxJ%%&*(>u_RF1s&TcBb#`az_2_owU4ytF7CF};siW}dpr
z%q3<lm^lQLJ?W}{sJzUGqCL$#E1xT^9LZ5~+0E9gdIKz!^kru4V~gbK;15vlglVu<
zD4hpfe}ekY*{!{a`n}R1JRZRS%bhgyGCv9!AU5W)(gm9pv}sg<$4+^8{$REhWe@jO
zgERc-68sNl+a#w!a_+|eg<oT0nzCEy#$z;l!t7<|<WTHi$zoc8%sjh}b<0a$8+o!5
znXgCkHiMVH5B+`|iKW}UWyTZeI5W@uOUufPZx9*GIV^^I1*mEL=e<bpx&-mgoOCyr
zuLCE)62($|D0b;`>(7=|4nIoI%{~@1Tx^>#EJM=7ThfHVy5@S9QP*JM*%6o7%fwLZ
z<;(CFrt+am@9?9t37BbVAI0hEV_o=lVD!i*^4V9^_y^~@+VrvaSx+<11~XmxyGQfh
zA9MbeJw3(jlQpj@zGtlBhOGEC+D=wHFDpKUl`gZ|YAE*7Rn<1-Ol9Xi_#Wp#+3x$W
zET3FaS<1%~XbC5{*~>AOYr|a1b7XM<9DqIsNlR(2dTfAhC)4;0>rt2mk4xU?$TNJ*
zjKbDSUWmNjO7bwumP%fP2vcK$L@^s_5_MC*ESJ|wANJ;sWSA(7{EPp?!gqpa^p`!7
zfa52P2Z2xlVU790HwKTq1WwsM8YG8HJ95TPkT9j*(es|1aTw5JNH2_MkRF!9Xh;kj
z93|&YnY~h@m%itEe<(Go@jv7q$~S-xevbxU>3v?NS<2ukv(G}4Nb;2ASGsb>ZYFwM
zn|@YVIe2W~KL=kxBO}cbntS2hGe<AZu9p=r`&Pg2^d93~w1__6cqcQrzK_lD#oJ7S
z?NaviiM$u)u0U?~vL}YU&%Ne`$rpDXVJ4y1W{<Eik(qDGQtJ(AV5Tx)x#GoS8Owfr
z;Kg%C2Y)>9^K++i?^)g)k;TEA=kz<^$86h=Ta_93F~!VVCXT{`<L7>|`;6{ah0itn
z=*6VDoCTkIaqxw6$IiWQ?x>k()0=%96vMM>te!DtFOBuju)*FiWzT%zKG_3PA9|l~
zCQQ~REKK$Z6Ib<#bcKB~mi-Uz6U)T7KCuc%pUA?k5jDOiodsXrC!Ct-ldPObE^9v3
z%eh##`xQ4{nZAi<1D6#}#%nDrE$?7nO<Y!dWEy)>#@Ms3RPy$bS9pG2{x$W!1|FA?
zixTap#;v*D1EgFu4l+@Ag!W0LC&A;QaTi=YX?z>#Iu?yS>K%R6lk0mH2{jy9`jEF;
ztru!M($OJr?FU>Wz6&eAZ$Pr5lPnUfHWvvg9=AxaT_+omEKCkaCh84o0E>iWz>uzR
zK#pafY(Oe+3X6m^gP~RM1B--ZV%&hFd6B_@WTL)rl3dP$odZ&KxYZ>bknAK3NRoPk
z{j9fLXykh!-cG0IdMm(@FQ?vIeg@O&3UgaYu5T98c$aM`&B>>he3>oJ$){C*?Lbbx
zYF=1}(K9?SJn)*A|GSR!rFzR}wXMhY?ZfN&2+yc?9rfbm!`!t3Y3v9}D=pckjU(sY
z=8oW1?5Xlu`BvR}#ca5ey*#Ry<z5~=ZLYkMeZvcwetXMjx5|CQ$aiFxL&(A{#bZF$
zOFe97M32}ZGEvbJ#pYg*Btjw@--i8u$QPGK4ZcK&FOk-}(|YPI-;QWO+Le*V1@q>n
z=8I535BR#$Dc#o|_J!nWg-;4Lul5DvA>aM0!40Q;@r3T{3h2Rz*lh8^9%K5pNA-yB
z0Z}Mo%z<d!Dt}9+s;6q}j;c(wswcXo%5r!sC!KQo0bfX+jNIi*MH1;)h)&f~2XBbj
z4eb2X+QzmujZMA{HFdstK#y(<L*kn0S;|c!-(P%aAnxM0Nn~*CXx$)O50(`F?e}jI
zVrHq=yN6p};X@pGYjl+7NuqIDiQ#U#v1oO+=*!-}@PH%mE+!WQMakp&ogc9pu@>>;
zh)WPzF8A}K-1n362$yuKO?kdRTxru<oAO+ND<sbnu@x6QGvqG$EFwmh$iN%EkL`yZ
z_FBRKrOT`@9ew^zo@)rP84+3Kxhb~fB8E(;u<(Ln?8sd<PX$X5u`VZ;!8o=J3FgJr
zNa(<*Yl3BC-<tRx5_ok;SSOoczR$pOWEI|m5AG@O1MzmmTEr&A4T$ZCNyIN9{xu@c
zGQMW#_3=D(Rde%QzIm(G-8<j6$Y1AQ?3;IQdpfSCEwrE}Em7_Kahi4hYX72wrbJ5f
z<Nu#De@~52{W?VUku@rCmV-88A90+pUpa0#9(eX#39155Hr!$#{;55d?+)Zc7BR>^
zkrL}~FNj?ZMwk!_+tbll$oda0?BwCzf?EF~|6+@^Ku?Cl3!{scEM7p-!eA=BFqqV%
z{*Bj{(`gWS1^U<5xY1YZuUYJ`uMw`DmOmWjUf1{jV1MEtnm?|aM0v5NWM;|SlAAGA
z#6<iQcy7^3iR(AhZp6R*h=aCtz5k{Z!fs316-91#x~E^bmb-CO=@Lb5_P)pOnzL}u
z{he-nAR><=!804Y<?bSj!!hCch^4i>bK%13<#Xmp3wu6vXW7tBXXlLh$`R+8?owsL
zQ;t2)4KB+U;ZOB+>4A1oJ!#R1N<{*xi13Gc;!v@uo)rGhc-p_sKFD?gEi_K1K{XkU
zF}s!Mj_JZ5jYsjXhjC2r$LW3)q68m-!XJ)k9mzmfSc`;^Qz79G>WO3urCa3lK?s>a
z+D;~<TY;`<5M?J|6MueLkY%(zg{;9u7ti&f8g5T_YJp@R-ihxAHjKtQXvGq>w<p8f
z6pBUTVf55>2OdpwY;qo-wc~LcgUFSDtHB&Zt{K$hbA|CdBA+LG22qDI4=bhQQIF3d
zMm~R3P{|rF`0St_pBs!k$s}eVW&9+fdgk!?!^lHWmP<LyWmt>IXA*gQE-`W?qaJmX
z-e-`$0a50~nBuc8Vq*}g$LC!bQR(ry$jIj<^YNKVd(?{~vQ7}ob(J%U&r%gseWPF~
zy&>o@_A4b&<)iM<<6LCy=EL3&svbLll|9aEMjX=@f+}|>I84*;D(Em4D<OEw9;dg`
zbAGc0TLz`~bs(k3Z;y-(mdZH2gFw_{`^)dg&|`b4pzQrU7)p=dD);>#df#{G@f+rT
z=qVetsnUlWdi*AG0D9`Z&Z*ZAfG&{!@x9|P%Y`HKmx1=^?<qv(AHR(}i*Fc??vQ%a
zeaWH6_p9fj$8o9xbr}C0QQ4E<S}|BSZ&aWT!!RPtr$77#TMj)X1FMk!InuPJ0`;ya
z0IFoKInorDqfsJE*g^Hpd%!9WMH%S%91NA@7*c8K?Z>wq9}{*^vWf|MKZIVr5<*uc
z`QEA8o+d8gn~+Zl9eQMOEUU0=4#r=RRdzF&Y$LWO=dRP<gXoi>gW)g@&GF%^4Xs-g
z)>=mvza27gO=IA<MFy^24E(;xz%`1&BUe^yEwV629DKA0<e1zDTw54C_6*?~Vqw<8
zS}QEfmB(5WEX)<zn)?>EMkeRDg}E|YcQ7nGP2vyx&VYJv_sHD@=E$A07}N)o3W4*1
z!6U!yBZ2YEFxedgOJe<d1O|&EcmxK8BX|U#USq){zur2rb46r%_=4pKo^fCNDfbn|
zBh|dPqcDE1Um9HrKG$RhH=c-_tfc&>egG+T2g8Mtm;hFHFkIMOJp`=oV7Tx${uN-(
zZ;nS7R;q7;@|!G!3r~`RpvN8M{R8lA_GF!`{O2Iwi}j3w<^L<F+8eR_A>iFsw?JMd
zE&@CKc?Woan?3~Ecg_>p9u5Za$_aP_FxOqKxh`?Dh)$4i1Ln78&J@<~i@;}8dG4vA
z7g&uam;2J_-!Z>x5BII&JHTpwxG~Fr0L*XG(jMeL1y*-fu!@Wh0`u9(;Ks`R9~0!S
z0`pZ_%2EDrfz^0r`8Q+ASJ>s#-VA&3xg8m}m-;gsFSl`njVWWi$Hwh8?y@l;dgIdQ
zpWB!-nfea{JMX}J#g?l(Fi%X-KW5`5JO9fQ<i7xR-hp}3mhXlB?TDrJTE_P=mOtCE
zM*4sk*;w6yX`3L&K}=NVOR+zAeY?0x_kY1aEatdN!>vq8Pj_^{sHQcyt>2)nS-r7M
zUa?s*4t=`0Ai^6nF1hI8kiV{O@gmVpQj6@Ngw`2Lv<G6Ed>z*U>5Q;m+C%BCt{zoD
z%eoa~_7zyUp>b`CrVtBL=`pMFU0+&#|C=`Mf4L1(w_k$Y-9mfdftnh}J)xGmW`9i`
zZ{pY+7&MC7%4k}2@sh><&akd^2Q@vCj(@=)%m{78{p%Xnu5M=SQwc4CGdbKR(&Rp7
z#eLX?boeK|?W1Yfdsx?*1n)87Qczg#NhT|&9p$WAv%aZujkbQ}%8e~;T3ch&nigDx
z3dcjvYex+Y@}iEswlsD@slZkDoeB)cttWN6$(Hd}kTlbpRN7<rlP0^|q_mtjrml9u
saSzLxp*5<I<I>jnJ6;V9&ihFA9j-BE^`v*eCTseJ2FGou3D?U07i`^%r2qf`

literal 0
HcmV?d00001

diff --git a/tracking/utility/interp2.mexw64 b/tracking/utility/interp2.mexw64
new file mode 100644
index 0000000000000000000000000000000000000000..787c0a0274f5b62172273668b2433e21f65aeda0
GIT binary patch
literal 10752
zcmeHN4{#J!ntzkb#1IlXC?lXmnlQnzU<etMKoXoW8ST*yMhF-X9EO<=8JU01bmPQT
zH*U!4+%q1vx3%Jycc-jtx!%^Qd)_X1?wOFRN#tL^Ria!iD6ZFB*d>b=;6M6)ue&EC
zs8?QHUDeeck9zO>-uLf)?|c8go>0Gb8_6bwSWz{Nkez_6xb*fvZ&a`+&3ttdd1b=B
zSvxJ_zFCbep(Gzq#5N`T5kBCLMq@I+QQ{Lyln+Juy4zRskyubFD=5fwn5J9%zVk%Q
zb9KF$I{TBuy*Hvi`;&vc)jI6!y#es54L9oaIaNcwmAbsAcZ$yM>RqP8cTD-Wb+{rF
zXkj)pGidM-68u39+4DDFHfCreq?n(WHGz=JL78Gu<*NY;b-c~O0K>Bgv9Sgv14$b?
z2-<O%&PK;j%7Cx-4GO=6ka|c@TLgaq$p%7RVx}7i8DafQiH};yP`({3C!|Gp{OiTb
zWGN+sb+8Z(t8DI=$JoNdtFlDUFZ&6JAUp&%NgnEzsAGx?+oa48BdQGy5cDwdfMbe_
zkgl@0AvnV54XA9+#jJCoxR|{}G7$htkI4og%r^(>nBvmqq;L#^Ix7w)vsXGs2&jvY
zo5oWv^)8NpkY4A=y;EpE$)Dy5OCFj<NHcfm05?~56moN`9rQhcehQB4XgFQR&EM0^
zl@2$%xoEeibB%>j-T~!{kGYS!+?^{dr@3`1<L(CXB_jfkcsdTH=6I?GADrds7)Z^i
ztJ?G6RBxcm8yFC{7yG?|-9mOhMwoG6*a}<kjuk(`#E+Qbd&Y@d4e@Sf-2F|)g^ckj
zhV(&5JDBuLPsfqeRm`?scSW~7DxH5opz9pvr@1x1Lp_Q~d~}*?>H)qF^)Odj?Lhxc
z;QLW~xl#lRgl?g|m-qGT9ajX}E2NupSFOI&-RNHBUhTe9c=+ToY-(<9d?z<I^?V&S
z*R>su@Tfp-TgnI#oY2_n*ajJQ`?+i{{nSlIgmk0BMc;iO2c(>b@`cWV9-a{O3kbD;
zZjl+`hm0@`g3xLEKJfPQR-tp9qeV#9IhsJb2Q-Fn0KV45TY%qc;&I?hO*{pBj*f42
zv=L^K+2HDjJkB$)E2d>f<%0zxy;h{(%N1#?aCn+PZvnPgq?-ysoCYG#MynLlVTS-6
z+17qON6FIm3xPj(oOviF%=?_1=^{dV^-Q9ah%{U1wEz7~Hq>K|Z5UFsues(zPraN?
zjONc&vw>{w%4pP&)3|o5#<6T@wD;r+?SHZ;+r4zJx*L@Cb6I+{PC``uEeKw^PyMM5
zvFbzWcBWO;qsu(cWS-P{d(~D$Mx?459TjIYY*6n(A3L}v4tbo3Lfc{>m@+*Fb~!fn
z*=7${+Rr_94K8Lc-LE>KAkeU*QlR?;dQP3K8{E&v&(~R?j;IzLRo&EMtvaWychh=O
z^{%oCBVB*K2m$Jc<@4B_3Tdw6qFq9u53IU3&VjV-V-HSWr(0v3ej$LpyDo>(2H@Y_
zIF^NVzhiP9=I`D8!{4j=T#2y|J&9{n2<#E6`W2ryu!?vC^;R}uCkABZG9$7vj8Q?w
z!Q8D62ubT3Y^ctgjFCr2d)pXo)P??Gwq>Wns#?j0#CdC=4>-!157{8NcjmBdN8h3c
zwsvt_2?FQo6h2qq*CV(WkrHT~!zHArsqdL2#M-6)v)(zN?>w~&9f67tfkpEl_0Ap^
z?&mY8s#2gHs-Q-!UaIeaAA=+Yj&*c#kC!nNZn=bcr-XTDg?UG}Ap%EH-|)~k-1Kdz
z;NV||I%C;sH@BvnHAAPl$#o81edaSwW5ZWMe!n3Pj&%W~IcA2*4j8K3WYNLvy6R@8
z)!vuuZ9j?obaEY*o}Ty_g5hxvs2>0#L}RCaL?;~(=`nRs51BsSb83&?!F)$lmeY-!
zq^iEmxPhULu7~RNiJjtRZUu?m{^RGj7Y<sNp_%!W`E`G4cq+6H>z*z+r*sPxr<Wa7
z=U&yqqSyJ9gWYysdY8-1l_p2feD@emL+v1%+q%5;xZCHx+r8es&Ug3v9<8MJD(vi%
zU+eAjdixW-eN1m3(%ZP+Zq!>5t!pKFzdRIo(0d{saJL^n=RQNMB0VTy;ij(K`apyA
z?2wp!Pzdykvk!U$J-5;+g#xWDq(ZKceOvAscQz;IbH8?Hw;vl3wu{;QqVr8J9bgyw
zT_=#Sk-Zq9w4S<gV!=45rZr95iTUA%=`FqJ%E$fR+?<K0P*v~@(h_B|(A{0we70ND
zl&p_?gzkaD6W1dydkx=ru??=-@0`vyY(JJ)v#(@$+O~+j2x(gh106>cTW76hc(%|!
zo2Xafgy0O_0-&Bo_62+afcg;{XP53e95de^A>C3aq*K$x^akF1h7dBRR!@(K^hDKR
zuA>qjP+Jjvt?J_*pUwQ~ucN2Zwmgs|fm%z-*>q`B?r<TyAg|h`;(pTXv!}7)Mg;EV
zqSSSm=<(qk_itHbsBHuK?7ro8;g%-qi}3WsGhhj{U_Drxor5kmu7DA@{#dtOP^~l0
zOhuMVzYy{8Q56?>uD9xBVy&>`K_dS}vkOoYmUxK#M*(-DX77iz9<v3Ek@gVnEe!G5
zd^W@{bY8(9(a}YoLSMvfXA`SN<e!Q3bCDh&egbDov)>DrW*=q>Y1=DIA#M97fOO3h
zXzQv@aUI8*laJK=HTowG19x_vcpISoGYi+z2V^u(^wd`)^(oS4Gg8Z80Fknpfm)1i
zRiB)1#OV<%;P5g$hbnWWv)~{=pN)eI3#+;j?5aLxD+FFN1fB*^Js|hlPRuhh;i^+g
zo<RGAM{LJI5th_^0wB+VwP7IYJJ@;{aB@BEeHFJYj^{qjeh?hZcJ&;51hRqHo$O+-
z6t2$)91S83A_c&KFwQI6^P_t!uJ0{(BMu2}@Dci%NZ)5SF-bS6Y5T4cLOKp9moKru
z4ot~&9s&8pn_weuNDDs8Qk=+rzQY!;qf5@!2I&wkbH3w<e1%BgpoiI`u;@!o%aFOP
zolqWSg)kX|)vj5T>pQxXh2S|4F|~T|Mmr;4Xqrfix)H8=7@wG7xV5)BG&rJVE688W
z0ec{PcIfJX*?Z`_Y24284EXp17I6586vq(D_&S_}S*>&MeYQRaqKHC!%mB(}01>-#
z_c}dpEaw)GLlx?Yx0jw{Y4%?b4aKIuj#DMlK8h`jgaB!%K#^>*!23?n^P0xu4u5iy
zMn1|`sk2}@qN>>g_<+UatXMt~_M1Go13DX#?coJ#Qr(Dyo~#J%SW6<ig+AnF(1(FK
zILHzehPfD4itMi!(_lC+LLZCI@;?f5c!W;D!{K6qZguDmq;1>afY51n+q*C|TM;86
zpJgCz`?uMwhn$v;m$s>}-Dmq35GPnkSD$SsFk{Er6efvbTcCWeiN`hj-(7*dWGjW<
zF8vHF&&E+2(Cohi9T7NdXAyW9o%Wh62w{^Ma&sW3+2dfm%7((-s@d;hIuG|U-`V}x
zS@3gQa3``LY23vEeHXXKYanX&D3gUdi2pfFnJCabntcao!w>Dsg}(X+9G-2Eeq5J!
zjY`u3S*Mhl6ql~q#(K10)Z>h`sZ_T^Uw_PzC1~~yy0py=Rkms_&A!~ATnmaZpkG($
zF>ZJECaJezuL;<wU1x{1e$eq6UHGOGI*UHWz7RTVP67;15|CzN`(|l_dS+a8H1Q{6
zkFo95Gfr=hq1XE@^u~Vh2()`w{aWM^OLFCGXIGC95UTNk=j?iT3Hu2~$$xbsrqI0P
zj%nPFx<2mND+b)v$ZSM<NX+hbbH8y-@o>Ls%DJ^uB-8MbCSODE$QP>)CGw~{H+#2O
zHJISs8oso;OFdP45+AS;Xg51lW??ZAX}8!}Q2f58;SE>sL4Yvl;`({rUOLDqkFt<~
z#8H={Tb+4!u>F{2I7@#s5k0M>e=^p-<Q~11^>&NirqD9`?4TA3?KSsdb%>Iq*^d@x
z?#$8nPt%V!cJWVnrWACP6XErW9Y?QwlHJGbvKP({;w(G{_V5;=Q*@MT_89Ew>!I1h
zK<ns`*Lg&<w~WDqOa#)E8RlIeh`1*zMduM=CuW8M%IqOwr_&rxtU4-kbwgsh=t_*{
zJ^L20Da<YXhS!%Kv)%ziTFJ9|`<&iBueUGg?G(MuL7R!R0bQ#Q5niDGF3^KQy1rDP
z7Mh!`cc}YtSFvrne$3)PWmLY;_P59in74-1D8tz3sZqu^j%RTWJX~qRZ%f>dx+(q{
zcF3jjUx~n<^+$k^ArtO4;jc`%-GmRAu*HPSO}NyAb4+-0u`&9P342ZWnhA4Fxr~oz
zZZI&_gp*7;-h^c)yxD|JCXAZ!hbDZ=gfE!zS0?;d)9!u~zGlMbOnBae8GrujqJF~c
zH=3~8gmX+d!-SJec=~#CJ|@iA{U+Epm9()ZZG0KIv5EIbgZvsk7~|QG=n=mhXyNY*
z$t^rC*9)c1Ku2D0@)`<WXh_KmIvLSWbc;V63i50DU??I*lc892DW7bKDd8adg}x#q
z$9|>dBT7>L-98voHijj>v?Mt<!_U}>q}+*wzg1^NB)KIPB&0P&NDB3S2)Cl9P&4vT
zC9+XU@Udn-rpR$cPA=u|hv?VvLs5f$S89Ii{Pp+GPle{UhVGhg46+_7<BH2=#P661
z$1%UAW0{#tX2P+n@gsTOq8yS}krSU@a|M~V{Z=yP(c8$Zlt_wN>dEv)&ywkur^w|a
zg~uih=6B`ACv3=Vz)#Fm@#0&NujFKroQkY+OHn>4>YYf6<UCT;IDr%`$|Xe=(?_Nr
zvk#gw4b$zp#16haIDy!~x38Go<-p~KhmdT-_p29wb8ztn^rL?}`f>DGv2Zzr3tmPz
z_{vpGj9c=vNWNt<nQV^F?6Cb)i=R?{j+*9~7HE0a(fLiVT8Q-~{8nzt9>&j8dB9)d
z4QMg>hEA>?lc&e*-!|-Rej(;PhV~QC3PB@dXfJ`*R*t%ed0>=W_Vo<fH<e5G&-lpv
zPu)KsWJglBNOD7nYz9`Fko>Z=T!~1DP(VM(zQyHim8_AJH(48t#KTgG95;0m<Yfb^
z%Um5vUHm<R3Ntv`w>Nvj5^febNruMLFVeE=W7O$a3nAgq#z?9x7!DIBiAX6=B2k~*
z)X4rs#W<3|B6eFco(M(dW>UX$b?u5J<+{dVi%*IM@CQX-JO-weAj`9S5q}te81e-Y
zGTCqOg`y!Dj0ib!0ryGd7uM#4BoR`Xu@Lk{75H$S#TS&uS^R+b{E_4)pOgy0r!~aq
ztMw)0QXteE3iw)ZM~5*iS<C40cvB!Q-y8}vPuf9&4=}e;*|bSY_*!DI&E&g`8kYP~
zpWGUk;9Yaf7xhP^Bv5EY{b53$&Gso#Okp$WfHWc#i|T{dWchUJhgm+8Qes^hmZUgY
zXl=w4Ls2YgLn0QCl1V~duwb=ol|%w7k|Ew?T@^K44ti1nDK6s?Fx@57tt%z@8`*oZ
zyvaIqCJhNG38UnR>=kl&btoY#{_v{keW55`fNhK>Vx~&`W+fV619+lx0*^tOaccc>
zSxHE>F-#UqyqDC4l774}jV;)ZK{2#3;ZL;IVgZw6dG;NOlxS^`63wwh#2*bnOVQW3
zISVsu3<Wmp1VVl}W~PQAFFg5&F;sK@Am7e(SGFc)Dbk3&>`r3(62r)Co)ybIqO(HZ
zSmTf6=--T2;4Wi#ZKew7i>M8#8zvLt2V^fARx(~uhRu854U5>O&V=G3v;Ql)%67wz
z%%{n|_57jc-?L8$qp$>>l}!H<{67|fLPAzx?~Wh?WPAAWEAa?83F!sA7P(F(aF+9U
zP_1~^k*|_JnDk?SOOUU)AX5!^FY=ZkaQ1&FUPgTy_$z?tXCT)Ao{#+GHPoHJ8GH+s
z+5A0VK5~*FCIeV!;)8&bk*9nJdO6@qR3_gDco`2_;4E(npss^V9FXNd3}<kMiN6GR
z(8S55^0g5_N{eD4)6k-}u^uZKI?Gj&Aha<WD;YX|RwDxXO&d@l^r68Z)&68sifjzG
z@~KETnye{S649l}K#LUdC$Ecy0*P2M)+}Eah((tAlaaD53yblM6Ad*>NqM!DV4o*7
z#S6>Ii<cGT@qD$Mz!!-ZujRM2<}~!6nUn%bA|$u!c6gNPgmkZhQKVo)BD4k9|0XFp
z${fq~q@c=PFj3kfh50aRYl{6zZ*)s+vy>?2m5@8Y9^0DYW`8&-6)&q^Flu|e@~HhW
z5?>j;dcinfs~3z;+K53=ic8TT5}wRlsxuVXZN~b-IGN~sp&(=WJj*v~it$Q2D~frR
za?}*pXI7zXbuBWg+7*q({FY4Ymz0&4RV*y8s3@x}X8)s$O*EEpCju=YSqiXM9I=Q5
zZ&(B_J`zKL+oaY;#Dy(ebD(mgbK?^KLTO>8<lM*tbFn|Br_~F-ehwGt`fM-gFKcku
k)p?iS($r)km-_#`c(6l_CxEm))3u|22arqU|5XJ312|k8kN^Mx

literal 0
HcmV?d00001

diff --git a/tracking/utility/p_to_box.m b/tracking/utility/p_to_box.m
new file mode 100644
index 0000000..20989c3
--- /dev/null
+++ b/tracking/utility/p_to_box.m
@@ -0,0 +1,38 @@
+function [ center corners ] = p_to_box(varargin)
+% function p_to_box(width,height, param, properties)
+%                 ([width,height], param, properties)
+%
+%
+
+%----------------------------------------------------------
+% Process the input.
+%----------------------------------------------------------
+if (length(varargin{1}) == 2)
+  w = varargin{1}(1);
+  h = varargin{1}(2);
+  varargin(1) = [];
+else
+  [w,h] = deal(varargin{1:2});
+  varargin(1:2) = [];
+end
+%%***********获得宽，高***********%%
+
+if (length(varargin) < 1 || any(length(varargin{1}) ~= 6))
+  M = [0,1,0; 0,0,1];
+else
+  p = varargin{1};
+  if (length(varargin) > 1 && strcmp(varargin{2},'geom'))
+    p = affparam2mat(p);
+    varargin(1:2) = [];
+  else
+    varargin(1) = [];
+  end
+  M = [p(1) p(3) p(4); p(2) p(5) p(6)];     %%affine变换参数
+end
+
+%----------------------------------------------------------
+% Draw the box.
+%----------------------------------------------------------
+corners = [ 1,-w/2,-h/2; 1,w/2,-h/2; 1,w/2,h/2; 1,-w/2,h/2; 1,-w/2,-h/2 ]';
+corners = M * corners;                  %%顶点
+center = mean(corners(:,1:4),2);        %%中心
diff --git a/tracking/utility/warpimg.m b/tracking/utility/warpimg.m
new file mode 100644
index 0000000..109f02d
--- /dev/null
+++ b/tracking/utility/warpimg.m
@@ -0,0 +1,29 @@
+function wimg = warpimg(img, p, sz)
+% function wimg = warpimg(img, p, sz)
+%
+%    img(h,w)
+%    p(6,n) : mat format
+%    sz(th,tw)
+%
+
+%% Copyright (C) 2005 Jongwoo Lim and David Ross.
+%% All rights reserved.
+
+
+if (nargin < 3)
+    sz = size(img);
+end
+if (size(p,1) == 1)
+    p = p(:);
+end
+w = sz(2);  h = sz(1);  n = size(p,2);
+%[x,y] = meshgrid(1:w, 1:h);
+[x,y] = meshgrid([1:w]-w/2+0.5, [1:h]-h/2);
+pos = reshape(cat(2, ones(h*w,1),x(:),y(:)) ...
+              * [p(1,:) p(2,:); p(3:4,:) p(5:6,:)], [h,w,n,2]);
+cn=size(img,3);
+wimg=zeros([sz cn]);
+for i=1:cn
+    wimg(:,:,i) = squeeze(interp2(img(:,:,i), pos(:,:,:,1), pos(:,:,:,2)));
+end
+wimg(find(isnan(wimg))) = 0;
diff --git a/tracking/wrappers.hpp b/tracking/wrappers.hpp
new file mode 100644
index 0000000..3009281
--- /dev/null
+++ b/tracking/wrappers.hpp
@@ -0,0 +1,42 @@
+/*******************************************************************************
+* Piotr's Computer Vision Matlab Toolbox      Version 3.00
+* Copyright 2014 Piotr Dollar.  [pdollar-at-gmail.com]
+* Licensed under the Simplified BSD License [see external/bsd.txt]
+*******************************************************************************/
+#ifndef _WRAPPERS_HPP_
+#define _WRAPPERS_HPP_
+#ifdef MATLAB_MEX_FILE
+
+// wrapper functions if compiling from Matlab
+#include "mex.h"
+inline void wrError(const char *errormsg) { mexErrMsgTxt(errormsg); }
+inline void* wrCalloc( size_t num, size_t size ) { return mxCalloc(num,size); }
+inline void* wrMalloc( size_t size ) { return mxMalloc(size); }
+inline void wrFree( void * ptr ) { mxFree(ptr); }
+
+#else
+
+// wrapper functions if compiling from C/C++
+inline void wrError(const char *errormsg) { throw errormsg; }
+inline void* wrCalloc( size_t num, size_t size ) { return calloc(num,size); }
+inline void* wrMalloc( size_t size ) { return malloc(size); }
+inline void wrFree( void * ptr ) { free(ptr); }
+
+#endif
+
+// platform independent aligned memory allocation (see also alFree)
+void* alMalloc( size_t size, int alignment ) {
+  const size_t pSize = sizeof(void*), a = alignment-1;
+  void *raw = wrMalloc(size + a + pSize);
+  void *aligned = (void*) (((size_t) raw + pSize + a) & ~a);
+  *(void**) ((size_t) aligned-pSize) = raw;
+  return aligned;
+}
+
+// platform independent alignned memory de-allocation (see also alMalloc)
+void alFree(void* aligned) {
+  void* raw = *(void**)((char*)aligned-sizeof(void*));
+  wrFree(raw);
+}
+
+#endif