Graphviz  2.41.20170921.2350
layout_similarity.c
Go to the documentation of this file.
1 /* $Id$ $Revision$ */
2 #include <stdio.h>
3 #include <stdlib.h>
4 #include <getopt.h>
5 #include "general.h"
6 #include "SparseMatrix.h"
7 #include "call_tri.h"
8 
9 static double boundingbox_area(int n, double *x){
10  double xmax, xmin, ymax, ymin;
11  int i;
12  xmax = xmin = x[0];
13  ymax = ymin = x[1];
14  for (i = 0; i < n; i++){
15  xmax = MAX(xmax, x[2*i]);
16  xmin = MIN(xmin, x[2*i]);
17  ymax = MAX(ymax, x[2*i+1]);
18  ymin = MIN(ymin, x[2*i+1]);
19  }
20  return (xmax-xmin)*(ymax-ymin);
21 
22 }
23 
24 
25 #define MINDIST 0.000000001
26 
27 real distance_cropped(real *x, int dim, int i, int j){
28  int k;
29  real dist = 0.;
30  for (k = 0; k < dim; k++) dist += (x[i*dim+k] - x[j*dim + k])*(x[i*dim+k] - x[j*dim + k]);
31  dist = sqrt(dist);
32  return MAX(dist, MINDIST);
33 }
34 
35 real distance(real *x, int dim, int i, int j){
36  int k;
37  real dist = 0.;
38  for (k = 0; k < dim; k++) dist += (x[i*dim+k] - x[j*dim + k])*(x[i*dim+k] - x[j*dim + k]);
39  dist = sqrt(dist);
40  return dist;
41 }
42 
43 static real get_r(int n, real *x, real *y, real theta, real *r, real *p, real *q){
44  /* find the best scaling factor r* such that
45  f[theta, r] = Sum[(x[i,1] + p - r(Sin[theta] y[i,1] + Cos[theta] y[i,2]))^2+(x[i,2] + q - r(Cos[theta] y[i,1] -Sin[theta] y[i,2]))^2, i = 1, n] is minimized. Return the value f[theta, r*],
46  also return r.
47  */
48  real top = 0, bot = 0, c, s, f = 0, yy, yx, xxavg = 0, xyavg = 0, yxavg = 0, yyavg = 0;
49  int i;
50 
51  for (i = 0; i < n; i++){
52  xxavg += x[2*i];
53  xyavg += x[2*i+1];
54  }
55  xxavg /= n;
56  xyavg /= n;
57 
58  for (i = 0; i < n; i++){
59  s = sin(theta);
60  c = cos(theta);
61  yx = c*y[2*i] - s*y[2*i+1];
62  yy = s*y[2*i] + c*y[2*i+1];
63  top += (x[2*i] - xxavg)*yx + (x[2*i+1] - xyavg)*yy;
64  bot += yx*yx+yy*yy;
65  yxavg += yx;
66  yyavg += yy;
67  }
68  *r = top/(bot - yxavg*yxavg/n - yyavg*yyavg/n);
69 
70  *p = yxavg/n*(*r)-xxavg;
71  *q = yyavg/n*(*r)-xyavg;
72 
73  for (i = 0; i < n; i++){
74  s = sin(theta);
75  c = cos(theta);
76  yx = c*y[2*i] - s*y[2*i+1];
77  yy = s*y[2*i] + c*y[2*i+1];
78  f += (x[2*i] + *p - (*r)*yx)*(x[2*i] + *p - (*r)*yx)+(x[2*i+1] + *q - (*r)*yy)*(x[2*i+1] + *q - (*r)*yy);
79  }
80  return f;
81 }
82 
83 static real dispacement(int n, real *x, real *y, real *rmin, real *thetamin, real *pmin, real *qmin){
84  /* find the distance between two sets of cordinates by finding a suitable
85  rotation and scaling to minimize the distance */
86  real dist = 0, theta, normx = 0, r, d, p, q;
87  int i;
88 
89  *thetamin = 0.;
90  dist = get_r(n, x, y, 0., &r, &p, &q);
91  *rmin = r; *pmin = p; *qmin = q;
92 
93  for (i = 0; i < 180; i++){
94  theta = 3.1415626*i/180;
95  d = get_r(n, x, y, theta, &r, &p, &q);
96  if (d < dist){
97  dist = d;
98  *rmin = r; *pmin = p; *qmin = q;
99  *thetamin = theta;
100  }
101  }
102 
103 
104  for (i = 0; i < n; i++) normx += x[2*i]*x[2*i]+x[2*i+1]*x[2*i+1];
105  return sqrt(dist/normx);
106 }
107 
108 int main(int argc, char *argv[])
109 {
110 
111 
112  char *infile;
113 
114  SparseMatrix B = NULL;
115  int dim = 2, n = 0, i, *ia, *ja, j, k, nz = 0;
116  real *x, *y, mean, dev, ratio, *z, r, theta, p, q, xx;
117 
118  FILE *fp;
119 
120 
121  if (argc < 2){
122  fprintf(stderr,"Usage: similarity graph layout1 layout2\n");
123  }
124 
125  infile = argv[1];
126  fp = fopen(infile,"r");
127  while (fscanf(fp,"%lf %lf\n",&xx, &xx) == 2) n++;
128  fclose(fp);
129 
130  infile = argv[1];
131  fp = fopen(infile,"r");
132  x = N_GNEW(n*2,real);
133  for (i = 0; i < n; i++){
134  fscanf(fp,"%lf %lf\n",&(x[2*i]), &(x[2*i+1]));
135  }
136  fclose(fp);
137 
138 
139  infile = argv[2];
140  fp = fopen(infile,"r");
141  y = N_GNEW(n*2,real);
142  nz = 0;
143  for (i = 0; i < n; i++){
144  if (fscanf(fp,"%lf %lf\n",&(y[2*i]), &(y[2*i+1])) != 2) goto ERROR;
145  }
146  fclose(fp);
147 
148  B = call_tri(n, 2, x);
149 
150  z = N_GNEW(B->nz,real);
151  ia = B->ia; ja = B->ja;
152  nz = 0;
153  for (i = 0; i < n; i++){
154  for (j = ia[i]; j < ia[i+1]; j++){
155  k = ja[j];
156  if (k == i) continue;
157  z[nz++] = distance(y,dim,i,k)/distance_cropped(x,dim,i,k);
158  }
159  }
160 
161  /* mean */
162  mean = 0;
163  for (i = 0; i < nz; i++) mean += z[i];
164  mean /= nz;
165 
166  /* deviation*/
167  dev = 0;
168  for (i = 0; i < nz; i++) {
169  dev += (z[i]-mean)*(z[i]-mean);
170  }
171  dev /= nz;
172  dev = sqrt(dev);
173 
174  /* bounding box area */
175  ratio = boundingbox_area(n, y);
176  /*/MAX(boundingbox_area(n, x), 0.001);*/
177 
178 
179  fprintf(stderr, "mean = %f std = %f disimilarity = %f area = %f badness = %f displacement = %f\n",
180  mean, dev, dev/mean, ratio, (dev/mean+1)*ratio, dispacement(n, x, y, &r, &theta, &p, &q));
181  fprintf(stderr, "theta = %f scaling = %f, shift = {%f, %f}\n",theta, 1/r, p, q);
182 
183 
184  printf("%.2f %.2f %.2f\n",dev/mean, dispacement(n, x, y, &r, &theta, &p, &q),ratio/1000000.);
185 
187  FREE(x);
188  FREE(y);
189  FREE(z);
190 
191 
192  return 0;
193 
194  ERROR:
195  printf("- - -\n");
196  return 0;
197 }
198 
#define MAX(a, b)
Definition: agerror.c:17
#define MINDIST
double xmax
Definition: geometry.c:20
#define MIN(a, b)
Definition: arith.h:35
#define FREE
Definition: PriorityQueue.c:23
double xmin
Definition: geometry.c:20
double ymax
Definition: geometry.c:20
real distance(real *x, int dim, int i, int j)
double ymin
Definition: geometry.c:20
void SparseMatrix_delete(SparseMatrix A)
Definition: SparseMatrix.c:411
Definition: grammar.c:79
SparseMatrix call_tri(int n, int dim, real *x)
Definition: call_tri.c:21
#define NULL
Definition: logic.h:39
#define top(sp)
Definition: stack.h:35
double dist(Site *s, Site *t)
Definition: site.c:41
real distance_cropped(real *x, int dim, int i, int j)
#define N_GNEW(n, t)
Definition: agxbuf.c:20
int main(int argc, char **argv)
Definition: dot.c:95
#define real
Definition: general.h:34