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
| /**
* @param image, width, height : binary image
* @return double[] { mean(x), mean(y), eigenvalue 1, eigenvalue 2,
* vector1.x, vector1.y, vector2.x, vector2.y }
*/
double[] binaryPCA(boolean[][] image, int width, int height) {
// compute sum of x,x2,y,y2,x*y of marked pixels
double sum_x=0, sum_x2=0, sum_y=0, sum_y2=0, sum_xy=0, count=0;
for(int y=0;y<height;y++) {
for(int x=0;x<width;x++) {
boolean value = image[x][y];
if (value==false) continue;
sum_x += x; sum_x2 += x*x;
sum_y += y; sum_y2 += y*y;
sum_xy += x*y;
count++;
}
}
// cov(X,Y) = E(X*Y) - E(X)*E(Y)
double Ex = sum_x/count, Ex2 = sum_x2/count;
double Ey = sum_y/count, Ey2 = sum_y2/count;
double Exy = sum_xy/count;
double COVxx = Ex2-Ex*Ex;
double COVyy = Ey2-Ey*Ey;
double COVxy = Exy-Ex*Ey;
// covariance matrix is
// | a b | | COVxx COVxy |
// | c d | = | COVxy COVyy |
double a=COVxx, b=COVxy, c=COVxy, d=COVyy;
// characteristic polynomial : P(L) = L^2 - Trace.L + Determinant
double T = a+d;
double Det = a*d-b*c;
// eingenvalues L1,L1 = roots of characteristic polynomial
double sqrt = Math.sqrt(T*T - 4*Det);
double L1 = 0.5 * ( T + sqrt );
double L2 = 0.5 * ( T - sqrt );
// eigenvectors: A.v = Li.v
// | a.vx + b.vy = Li.vx
// | c.vx + d.vy = Li.vy
double V1x,V1y,V1n,V2x,V2y,V2n, EPSILON=1E-8;
V1x=-(b+d-L1); V1y=(a+c-L1);
V2x=-(b+d-L2); V2y=(a+c-L2);
// normalize eigenvectors
V1n=Math.sqrt(V1x*V1x+V1y*V1y);
if (V1n>EPSILON) {V1x/=V1n;V1y/=V1n;}
V2n=Math.sqrt(V2x*V2x+V2y*V2y);
if (V2n>EPSILON) {V2x/=V2n;V2y/=V2n;}
// return data
return new double[] {Ex, Ey, L1, L2, V1x, V1y, V2x, V2y};
} |
