1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
|
/*Calling Sequence
c = xcorr2 (a)
c = xcorr2 (a, b)
c = xcorr2 (a, b, scale)
Description:
Compute the 2D cross-correlation of matrices a and b.
If b is not specified, computes autocorrelation of a, i.e., same as xcorr (a, a).
The optional argument scale, defines the type of scaling applied to the cross-correlation matrix. Possible values are:
"none" (default)
No scaling.
"biased"
Scales the raw cross-correlation by the maximum number of elements of a and b involved in the generation of any element of c.
"unbiased"
Scales the raw correlation by dividing each element in the cross-correlation matrix by the number of products a and b used to generate that element.
"coeff"
Scales the normalized cross-correlation on the range of [0 1] so that a value of 1 corresponds to a correlation coefficient of 1.
Examples
xcorr2(5,0.8,'coeff')
ans = 1 */
function c = xcorr2 (a, b, scale)
funcprot(0);
nargin=argn(2);
if nargin < 3 then
scale = "none"
end
if (nargin < 1 || nargin > 3)
error("Wrong number of inputs")
elseif (nargin == 2 && type (b) == 10 )
scale = b;
b = a;
elseif (nargin == 1)
// we have to set this case here instead of the function line, because if
// someone calls the function with zero argument, since a is never set, we
// will fail with "`a' undefined" error rather that print_usage
b = a;
end
if (ndims (a) ~= 2 || ndims (b) ~= 2)
error ("xcorr2: input matrices must must have only 2 dimensions");
end
// compute correlation
[ma,na] = size(a);
[mb,nb] = size(b);
c = conv2 (a, conj (b (mb:-1:1, nb:-1:1)));
// bias routines by Dave Cogdell (cogdelld@asme.org)
// optimized by Paul Kienzle (pkienzle@users.sf.net)
// coeff routine by Carnë Draug (carandraug+dev@gmail.com)
switch (scale)
case {"none"}
// do nothing, it's all done
case {"biased"}
c = c / ( min ([ma, mb]) * min ([na, nb]) );
case {"unbiased"}
lo = min ([na,nb]);
hi = max ([na, nb]);
row = [ 1:(lo-1), lo*ones(1,hi-lo+1), (lo-1):-1:1 ];
lo = min ([ma,mb]);
hi = max ([ma, mb]);
col = [ 1:(lo-1), lo*ones(1,hi-lo+1), (lo-1):-1:1 ]';
bias = col*row;
c = c./bias;
case {"coeff"}
a = double (a);
b = double (b);
a = conv2 (a.^2, ones (size (b,1) , size( b ,2)));
b = sum(b(:).* conj(b(:)));
c(:,:) = c(:,:) ./ sqrt (a(:,:) * b);
else
error ("xcorr2: invalid type of scale %s", scale);
end
endfunction
|