Je reporte ici un fragment de code Java posté sur le forum qui effectue une analyse en composante principale (ACP) sur une image binaire.

Dans cet exemple, on utilise l'ACP pour trouver l'axe principal (rouge) et l'axe secondaire (bleu) d'une forme binaire :


Code java : Sélectionner tout - Visualiser dans une fenêtre à part
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};
}