NNP.zip_最近点对

#include <iostream> #include <fstream> #include <string> #include <stdio.h> #include <stdlib.h> #include <math.h> #include <cmath> #include <cassert> #include "RTree.h" using namespace std; #define MAXINT 32767 #define sqr(a) ((a)*(a)) #define DIM 10 int NUMBER=0; int Statistics_Count; double* Statistics_SumX; double* Statistics_SumXsquared; double* Statistics_SumProduct; void InitStatistics(int Dimensions) // ============== // initializes variables for the statistics { Statistics_SumX = new double[Dimensions]; Statistics_SumXsquared = new double[Dimensions]; Statistics_SumProduct = new double[Dimensions*Dimensions]; Statistics_Count = 0; for (int d=0; d<Dimensions; d++) { Statistics_SumX[d]=0.0; Statistics_SumXsquared[d]=0.0; for (int dd=0; dd<Dimensions; dd++) Statistics_SumProduct[d*Dimensions+dd] = 0.0; } } void EnterStatistics(int Dimensions,double* x) // =============== // counts the vector "x" in the statistics { Statistics_Count++; for (int d=0; d<Dimensions; d++) { Statistics_SumX[d] += x[d]; Statistics_SumXsquared[d] += x[d]*x[d]; for (int dd=0; dd<Dimensions; dd++) Statistics_SumProduct[d*Dimensions+dd] += x[d]*x[dd]; } } void OutputStatistics(int Dimensions) // ================ // outputs the statistics { for (int d=0; d<Dimensions; d++) { double E = Statistics_SumX[d] / Statistics_Count; double V = Statistics_SumXsquared[d]/Statistics_Count - E*E; double s = sqrt(V); } for ( d=0; d<Dimensions; d++) { for (int dd=0; dd<Dimensions; dd++) { double Kov = (Statistics_SumProduct[d*Dimensions+dd]/Statistics_Count) - (Statistics_SumX[d]/Statistics_Count) * (Statistics_SumX[dd]/Statistics_Count); double Cor = Kov / sqrt(Statistics_SumXsquared[d]/Statistics_Count - sqr(Statistics_SumX[d] / Statistics_Count)) / sqrt(Statistics_SumXsquared[dd]/Statistics_Count - sqr(Statistics_SumX[dd] / Statistics_Count)); } } } double RandomEqual(double min,double max) // =========== // generates a random number equally distributed in the interval [min,max[ { // double x = (double)rand()/MAXINT; <-- patch double x = (double)rand()/RAND_MAX; return x*(max-min)+min; } double RandomPeak(double min,double max,int dim) // ========== // generates a random number in the interval [min,max[ // as a sum of "dim" equally distributes random numbers { double sum = 0.0; for (int d=0; d<dim; d++) sum += RandomEqual(0,1); sum /= dim; return sum*(max-min)+min; } double RandomNormal(double med,double var) // ============ // generates a normally distributed random number with mean "med" // in the interval ]med-var,med+var[ { return RandomPeak(med-var,med+var,12); } void GenerateDataEqually(FILE* f,int Count,int Dimensions) // =================== // generates "Count" independently equally distributed vectors // in the file "f" { InitStatistics(Dimensions); for (int i=0; i<Count; i++) { double x[DIM]; for (int d=0; d<Dimensions; d++) { x[d] = RandomEqual(0,1); fprintf(f,"%8.6f ",x[d]); } EnterStatistics(Dimensions,x); fprintf(f,"\n"); } OutputStatistics(Dimensions); } void GenerateDataCorrelated(FILE* f,int Count,int Dimensions) // ====================== // generates "Count" correlated vectors in the file "f" { InitStatistics(Dimensions); double x[DIM]; for (int i=0; i<Count; i++) { again: double v = RandomPeak(0,1,Dimensions); for (int d=0; d<Dimensions; d++) x[d] = v; double l = v<=0.5 ? v:1.0-v; for ( d=0; d<Dimensions; d++) { double h = RandomNormal(0,l); x[d] += h; x[(d+1)%Dimensions] -= h; } for ( d=0; d<Dimensions; d++) if (x[d]<0 || x[d]>=1) goto again; for ( d=0; d<Dimensions; d++) fprintf(f,"%8.6f ",x[d]); EnterStatistics(Dimensions,x); fprintf(f,"\n"); } OutputStatistics(Dimensions); } void GenerateDataAnticorrelated(FILE* f,int Count,int Dimensions) // ========================== // generates "Count" anti-correlated vectors in the file "f" { InitStatistics(Dimensions); double x[DIM]; for (int i=0; i<Count; i++) { again: double v = RandomNormal(0.5,0.25); for (int d=0; d<Dimensions; d++) x[d] = v; double l = v<=0.5 ? v:1.0-v; for ( d=0; d<Dimensions; d++) { double h = RandomEqual(-l,l); x[d] += h; x[(d+1)%Dimensions] -= h; } for (d=0; d<Dimensions; d++) if (x[d]<0 || x[d]>=1) goto again; for ( d=0; d<Dimensions; d++) fprintf(f,"%8.6f ",x[d]); EnterStatistics(Dimensions,x); fprintf(f,"\n"); } OutputStatistics(Dimensions); } void GenerateData(int Dimensions,char Distribution,int Count,char* FileName) // ============ // generates a file with random data { if (Count <= 0) { printf("Ungtige Anzahl von Punkten.\n"); return; } if (Dimensions < 2) { printf("Ungtige Anzahl von Dimensionen.\n"); return; } switch (Distribution) { case 'E': case 'e': Distribution = 'E'; break; case 'C': case 'c': Distribution = 'C'; break; case 'A': case 'a': Distribution = 'A'; break; default: printf("Invalid distribution.\n"); return; } FILE* f = fopen(FileName,"wt"); if (f == NULL) { printf("Cannot create file \"%s\".\n",FileName); return; } fprintf(f,"%d %d\n",Count,Dimensions); switch (Distribution) { case 'E': GenerateDataEqually(f,Count,Dimensions); break; case 'C': GenerateDataCorrelated(f,Count,Dimensions); break; case 'A': GenerateDataAnticorrelated(f,Count,Dimensions); break; } fclose(f); printf("%d vectors generated, filename \"%s\".\n",Count,FileName); } class NNTree:public RTree<long,int,2>{ public: int MINDIST(Branch *branch,int[]); int MINMAXDIST(Branch *branch,int[]); void Push(Node* node,int[]); void NNSEARCH(int[]); private: struct Heap { Branch* branch; int level; }; int E[10000][2]; Heap heap[1000]; }; struct info { int MinDist; int no; }; struct info nn[100]; int NNTree::MINDIST(Branch* branch,int lookup[]) { int r,distant=0; for(int i=0;i<2;i++) { if(lookup[i]<branch->m_rect.m_min[i]) r=branch->m_rect.m_min[i]; else if(lookup[i]>branch->m_rect.m_max[i]) r=branch->m_rect.m_max[i]; else r=lookup[i]; distant=distant+(lookup[i]-r)*(lookup[i]-r); } return distant; } int NNTree::MINMAXDIST(Branch* branch,int lookup[]) { int r,distant=0,dis[2]; for(int i=0;i<2;i++) { distant=0; for(int j=0;j<2;j++) { if(j==i) { if(lookup[j]<=(branch->m_rect.m_min[j]+branch->m_rect.m_max[j])/2) r=branch->m_rect.m_min[j]; else r=branch->m_rect.m_max[j]; } else if(j!=i) { if(lookup[j]>=(branch->m_rect.m_min[j]+branch->m_rect.m_max[j])/2) r=branch->m_rect.m_min[j]; else r=branch->m_rect.m_max[j]; } distant=distant+(lookup[j]-r)*(lookup[j]-r); } dis[i]=distant; } if(dis[0]<=dis[1]) distant=dis[0]; else distant=dis[1]; return distant; } void NNTree::Push(Node* node,int lookup[]) { int i,j=0,min_dist,node_no,k,num,change; for(i=0;i<MAXNODES;i++) { if((node->m_branch[i].m_rect.m_min[0]<0)&&(node->m_branch[i].m_rect.m_max[0]<0)) continue; nn[j].MinDist=MINDIST(node->m_branch+i,lookup); nn[j].no=i;//结点标号 j++; } num=j;//数组nn中结点数 for(i=num-1,change=1;i>=1&&change;--i) { change=0; for(j=0;j<i;++j) { if(nn[j].MinDist>nn[j+1].MinDist) { min_dist=nn[j].MinDist; node_no=nn[j].no; nn[j].MinDist=nn[j+1].MinDist; nn[j].