source: trunk/gsdl3/extensions/vishnu/src/vishnu/util/Sammon.java@ 8417

Last change on this file since 8417 was 8282, checked in by kjdon, 20 years ago

moved from sammon package to util package
  • Property svn:keywords set to Author Date Id Revision
File size: 6.3 KB
Line 
1package vishnu.util;
2
3import java.io.*;
4import java.util.*;
5import vishnu.datablock.Point2D;
6
7public class Sammon
8{
9
10 public static final double MF = 0.3; /* Magic Factor */
11 public static final double _MAX_Err = 0.00005; /* Max Mapping Error */
12 public static final int _MAX_M = 1000; /* Max Iteration Number */
13 public static final int _MAX_N = 1000; /* Max Number of Entries */
14 public static final int _MAX_L = 128; /* Max Dimension for X */
15
16
17 double[][] Xs;
18 double[][] Ys;
19 int number;
20 int dimension;
21 double total_Xs_dists;
22
23/* SS298: Added function for scaling values between 0 and 1
24* n : number of entries
25*/
26 void scale_values(int n)
27 {
28
29 // get the maximum and the minimum of Ys array
30
31 // get the maximum and the minimum of Ys array
32
33 double maxX = Ys[0][0];
34 double minX = Ys[0][0];
35 double maxY = Ys[0][1];
36 double minY = Ys[0][1];
37
38 int i;
39 for (i = 1 ; i < n ; i++)
40 {
41
42 if (maxX < Ys[i][0]) maxX = Ys[i][0];
43 if (minX > Ys[i][0]) minX = Ys[i][0];
44 if (maxY < Ys[i][1]) maxY = Ys[i][1];
45 if (minY > Ys[i][1]) minY = Ys[i][1];
46 }
47 maxY -= minY;
48 maxX -= minX;
49 for (i = 0 ; i < n ; i++)
50 {
51 if (maxX > 0)
52 {
53 Ys[i][0] -= minX;
54 Ys[i][0] /= maxX;
55 Ys[i][0] *=0.8;
56 Ys[i][0] +=0.1;
57 }
58 else Ys[i][0]=0.5;
59 if (maxY > 0)
60 {
61 Ys[i][1] -= minY;
62 Ys[i][1] /= maxY;
63 Ys[i][1] *= 0.8;
64 Ys[i][1] += 0.1;
65 }
66 else Ys[i][1]=0.5;
67 }
68
69 }
70
71
72
73/* return a Euclidean distance between two 'd'-dimensional vectors 'xs' and 'ys'
74*/
75 double eucl_dist(int d, double[] xs, double[] ys)
76 {
77 double ans=0.0;
78 double tmp;
79 int i;
80 for(i=0;i<d;i++)
81 {
82 tmp = xs[i]-ys[i];
83 ans += tmp * tmp;
84 }
85 return Math.sqrt(ans);
86 }
87
88
89/* n - number of entries
90* d - dimension
91* post:
92* return an array of eucidean distances for X
93* and assign its total distance (will be used in the mapping error) to 'total'
94*/
95 double[] compute_Xs_eucl_dists()
96 {
97 double tmp[];
98 int i,j,k;
99 tmp = new double[number*(number-1)/2];
100
101 total_Xs_dists=0.0;
102 k=0;
103 for(i=1;i<number;i++)
104 for(j=0;j<i;j++)
105 {
106 tmp[k] = eucl_dist(dimension, Xs[i], Xs[j]);
107 total_Xs_dists+=tmp[k];
108 k++;
109 }
110 return tmp;
111 }
112
113
114 public Point2D [] getMapping()
115 {
116 Point2D [] sam = new Point2D[Ys.length];
117 for (int c=0;c<Ys.length;c++)
118 {
119 float fx = (float) Ys[c][0];
120 float fy = (float) Ys[c][1];
121 sam[c]=new Point2D(fx,fy);
122 }
123 return sam;
124 }
125
126
127
128/*----------------------------------------------------------------*/
129
130 public Sammon(double [][] centroids)
131 {
132
133 int i, j, m, p;
134 double[] Xs_dists;
135
136 double dist_star,dist, alpha, beta, gamma, zeta, err;
137 double uppDer0, uppDer1, lowDer0, lowDer1;
138 double maxX=0.0;
139 double minX=0.0;
140 double maxY=0.0;
141 double minY=0.0;
142
143 number = centroids.length;
144 dimension = centroids[0].length;
145
146 Xs = centroids;
147 Ys = new double[number][2];
148
149 Xs_dists = compute_Xs_eucl_dists();
150
151 p=0;
152 for (i=0;i<number;i++)
153 {
154 Ys[i][0] = Xs[i][0];
155 double value = Xs_dists.length>0?Xs_dists[(i==(number-1))?i:p]:
156 0.0;
157 Ys[i][1] = Xs[i][0]+value;
158 p+=(number-(i+1));
159 }
160
161 m=0; /* iteration counter */
162 do{
163
164 for(p=0;p<number;p++)
165 {
166 uppDer0=uppDer1=lowDer0=lowDer1=0.0;
167
168
169 /* calcuate partial derivatives */
170
171 for(j=0;j<number;j++)
172 {
173 if (p==j) continue;
174 dist_star = Xs_dists[(j>p) ? (j*(j-1)/2+p) : (p*(p-1)/2+j)];
175
176 dist = eucl_dist(2,Ys[p],Ys[j]);
177 if (dist == 0) dist = 0.00001;
178 alpha = dist_star-dist;
179 beta = dist_star*dist;
180 if (beta == 0) beta = 0.00000000000001;
181 gamma = alpha/beta;
182
183 /* 1st dimension */
184 zeta = Ys[p][0]-Ys[j][0];
185 uppDer0 += gamma * zeta;
186 lowDer0 += (alpha-(zeta*zeta/dist) * (1+alpha/dist))/beta;
187
188 /* 2nd dimension */
189 zeta = Ys[p][1]-Ys[j][1];
190 uppDer1 += gamma * zeta;
191 lowDer1 += (alpha-(zeta*zeta/dist) * (1+alpha/dist))/beta;
192 } /* for */
193
194 Ys[p][0] = Ys[p][0]+MF*uppDer0/Math.abs(lowDer0);
195 Ys[p][1] = Ys[p][1]+MF*uppDer1/Math.abs(lowDer1);
196 if(p==0)
197 {
198 maxX = Ys[0][0];
199 minX = Ys[0][0];
200 maxY = Ys[0][1];
201 minY = Ys[0][1];
202 }
203 else
204 {
205 if (maxX < Ys[p][0]) maxX = Ys[p][0];
206 if (minX > Ys[p][0]) minX = Ys[p][0];
207 if (maxY < Ys[p][1]) maxY = Ys[p][1];
208 if (minY > Ys[p][1]) minY = Ys[p][1];
209 }
210
211
212 } /* for */
213
214
215 /* calcuate the mapping error */
216
217 err=0.0;
218 p=0;
219 for(j=1;j<number;j++)
220 for(i=0;i<j;i++)
221 {
222 dist_star = Xs_dists[p++];
223 alpha = dist_star - eucl_dist(2,Ys[j],Ys[i]);
224 err += (alpha*alpha/dist_star);
225 }
226 err=err / total_Xs_dists;
227
228 m++;
229 } while ((err>_MAX_Err) && (m!=_MAX_M));
230
231
232 if(maxY==minY)
233 {
234 System.out.println("Warning: simulated 'y' values in mapping");
235 for(p=0;p<number;p++)
236 Ys[p][1]=Math.random();
237 }
238 if(maxX==minX)
239 {
240 System.out.println("Warning: simulated 'x' values in mapping");
241 for(p=0;p<number;p++)
242 Ys[p][0]=Math.random();
243 }
244
245 scale_values(number);
246 }
247}
Note: See TracBrowser for help on using the repository browser.