SPC/distance.c at main · eytandomany/SPC · GitHub

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
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
#include "SW.h"

/* -------------------------------------------------------------------- */
/*  Distance: returns the euclidean distance between data poins i and j */
/* -------------------------------------------------------------------- */
double Distance(int D, double *X, double *Y) {
  int d;
  double dist = 0.0;
  double diff;
  for (d = 0; d <  D; d++){
    diff = X[d] - Y[d];
    dist += diff * diff;
  }
  return( sqrt(dist) );
}

/* -------------------------------------------------------------------- */
/*  Distance_Linf: returns the L-infinity distance between data poins i */
/*                 and j                                                */
/* -------------------------------------------------------------------- */
double Distance_Linf(int D, double *X, double *Y) {
  int d;
  double dist;
  double diff;
  dist = 0.0;
  for (d = 0; d <  D; d++){
    diff = fabs(X[d] - Y[d]);
    if (dist < diff) dist = diff;
  }
  return( dist );
}

RaggedArray EdgeDistance( int D, UIRaggedArray NK, double **X ){
   RaggedArray M;
   int i,k,metric,nedges;
   float chd;

   M = InitRaggedArray( NK );
   metric = ( GetParam( "InfMetric" ) == NULL );

   chd = 0.0;
   nedges = 0;
   for(i=0; i < M.n; i++){
     for(k = 0; k<M.c[i]; k++ ) {
         M.p[i][k] = D? ( metric? Distance(D,X[i],X[NK.p[i][k]])
			  : Distance_Linf(D,X[i],X[NK.p[i][k]]) )
	                : X[i][NK.p[i][k]];
         if( M.p[i][k] < INF ) chd+=M.p[i][k], nedges++;
      }
   }
   nedges /= 2;

   ISetParam( "NumberOfEdges", nedges );
   FSetParam( "NearestNeighbrs", 2.0 * (float)nedges / (float)M.n );
   FSetParam( "CharDist", chd / (2.0 * (float)nedges) );
   return M;
}

void DistanceToInteraction( RaggedArray M, UIRaggedArray NK, UIRaggedArray KN )
{
   int i, k;
   float dd, lambdaq, *z, nn, chd;

   if( (nn=FGetParam( "ForceNN" ))==0. ) nn=FGetParam( "NearestNeighbrs" );
   if( (chd=FGetParam( "ForceChD" ))==0. ) chd=FGetParam( "CharDist" );

   /* the interactions are calculated */
   for(i = 0; i < M.n; i++)
      for( k = M.c[i]-1; NK.p[i][k]>i && k>=0; k-- ) {
	 dd = (M.p[i][k]*M.p[ NK.p[i][k] ][ KN.p[i][k] ])/(chd*chd);
	 M.p[i][k] = M.p[ NK.p[i][k] ][ KN.p[i][k] ] = exp(-dd/2.0) / nn;
      }

   /* calc z_i */
   if ( GetParam( "UseZ" ) ) {
   /* allocate z */
     z=InitVector(M.n);
     for(i=0; i < M.n; i++)
       for(k= 0; k<M.c[i]; k++) {
	 z[i] += M.p[i][k];
       }
     lambdaq = FGetParam( "Lambda" ) * (float) IGetParam( "PottsSpins" );
     for(i=0; i < M.n; i++)
       for(k = 0; k<M.c[i]; k++) {
	 M.p[i][k] *= ( 1.0/(z[i]+lambdaq) + 1.0/(z[NK.p[i][k]]+lambdaq) );
       }
   free(z);
   }
}