bzr branch
http://suren.me/webbzr/normxcorr/trunk
1
by Suren A. Chilingaryan
Initial import |
1 |
function [xpeak, ypeak, max_f] = findpeak(f,subpixel) |
2 |
%FINDPEAK Find extremum of matrix.
|
|
3 |
% [XPEAK,YPEAK,MAX_F] = FINDPEAK(F,SUBPIXEL) finds the extremum of F,
|
|
4 |
% MAX_F, and its location (XPEAK, YPEAK). F is a matrix. MAX_F is the maximum
|
|
5 |
% absolute value of F, or an estimate of the extremum if a subpixel
|
|
6 |
% extremum is requested.
|
|
7 |
%
|
|
8 |
% SUBPIXEL is a boolean that controls if FINDPEAK attempts to estimate the
|
|
9 |
% extremum location to subpixel precision. If SUBPIXEL is false, FINDPEAK
|
|
10 |
% returns the coordinates of the maximum absolute value of F and MAX_F is
|
|
11 |
% max(abs(F(:))). If SUBPIXEL is true, FINDPEAK fits a 2nd order
|
|
12 |
% polynomial to the 9 points surrounding the maximum absolute value of
|
|
13 |
% F. In this case, MAX_F is the absolute value of the polynomial evaluated
|
|
14 |
% at its extremum.
|
|
15 |
%
|
|
16 |
% Note: Even if SUBPIXEL is true, there are some cases that result
|
|
17 |
% in FINDPEAK returning the coordinates of the maximum absolute value
|
|
18 |
% of F:
|
|
19 |
% * When the maximum absolute value of F is on the edge of matrix F.
|
|
20 |
% * When the coordinates of the estimated polynomial extremum would fall
|
|
21 |
% outside the coordinates of the points used to constrain the estimate.
|
|
22 |
||
23 |
% Copyright 1993-2004 The MathWorks, Inc.
|
|
24 |
% $Revision $ $Date: 2004/10/20 17:54:47 $
|
|
25 |
||
26 |
% get absolute peak pixel
|
|
27 |
[max_f, imax] = max(abs(f(:))); |
|
28 |
[ypeak, xpeak] = ind2sub(size(f),imax(1)); |
|
29 |
||
30 |
if ~subpixel || ... |
|
31 |
xpeak==1 || xpeak==size(f,2) || ypeak==1 || ypeak==size(f,1) % on edge |
|
32 |
return % return absolute peak |
|
33 |
||
34 |
else
|
|
35 |
% fit a 2nd order polynomial to 9 points |
|
36 |
% using 9 pixels centered on irow,jcol
|
|
37 |
u = f(ypeak-1:ypeak+1, xpeak-1:xpeak+1); |
|
38 |
u = u(:); |
|
39 |
x = [-1 -1 -1 0 0 0 1 1 1]'; |
|
40 |
y = [-1 0 1 -1 0 1 -1 0 1]'; |
|
41 |
||
42 |
% u(x,y) = A(1) + A(2)*x + A(3)*y + A(4)*x*y + A(5)*x^2 + A(6)*y^2
|
|
43 |
X = [ones(9,1), x, y, x.*y, x.^2, y.^2]; |
|
44 |
||
45 |
% u = X*A
|
|
46 |
A = X\u; |
|
47 |
||
48 |
% get absolute maximum, where du/dx = du/dy = 0
|
|
49 |
x_offset = (-A(3)*A(4)+2*A(6)*A(2)) / (A(4)^2-4*A(5)*A(6)); |
|
50 |
y_offset = -1 / ( A(4)^2-4*A(5)*A(6))*(A(4)*A(2)-2*A(5)*A(3)); |
|
51 |
||
52 |
if abs(x_offset)>1 || abs(y_offset)>1 |
|
53 |
% adjusted peak falls outside set of 9 points fit,
|
|
54 |
return % return absolute peak |
|
55 |
end
|
|
56 |
||
57 |
% return only one-tenth of a pixel precision |
|
58 |
x_offset = round(10*x_offset)/10; |
|
59 |
y_offset = round(10*y_offset)/10; |
|
60 |
||
61 |
xpeak = xpeak + x_offset; |
|
62 |
ypeak = ypeak + y_offset; |
|
63 |
||
64 |
% Calculate extremum of fitted function
|
|
65 |
max_f = [1 x_offset y_offset x_offset*y_offset x_offset^2 y_offset^2] * A; |
|
66 |
max_f = abs(max_f); |
|
67 |
||
68 |
end
|