00001 #include "vuUDSphere.h"
00002 #include <stdio.h>
00003 #include <stdlib.h>
00004 #include <assert.h>
00005
00006 #define TWO_PI 6.283185307179586
00007
00008 vuUDSphere::vuUDSphere(int _size)
00009 {
00010 m_Tolerance = 0.001;
00011 m_Size = _size + 1;
00012 m_P = new Point3d[m_Size];
00013 m_dt = TWO_PI / ((double)m_Size);
00014
00015 initPoints();
00016 }
00017
00018 vuUDSphere::~vuUDSphere()
00019 {
00020 delete [] m_P;
00021 }
00022
00023
00024
00025
00026
00027 void vuUDSphere::makeRandomPoint(Point3d& p)
00028 {
00029 double x, y, z;
00030 double dist;
00031
00032
00033 x = ((double)(rand()%RAND_MAX))/((double)(RAND_MAX-1))*2.0 - 1.0;
00034 y = ((double)(rand()%RAND_MAX))/((double)(RAND_MAX-1))*2.0 - 1.0;
00035 z = ((double)(rand()%RAND_MAX))/((double)(RAND_MAX-1))*2.0 - 1.0;
00036 dist = sqrt(x*x + y*y + z*z);
00037 x /= dist;
00038 y /= dist;
00039 z /= dist;
00040
00041
00042 p.set(x, y, z);
00043 }
00044
00045 void vuUDSphere::initPoints(void)
00046 {
00047 int i;
00048
00049 m_P[0].x = 0.0;
00050 m_P[0].y = 0.0;
00051 m_P[0].z = 0.0;
00052 for(i=1;i<m_Size;++i)
00053 makeRandomPoint(m_P[i]);
00054 }
00055
00056 void vuUDSphere::computeForces(Point3d& p0, Point3d& p1)
00057 {
00058 double x, y, z;
00059 double d;
00060
00061 x = p1.x - p0.x;
00062 y = p1.y - p0.y;
00063 z = p1.z - p0.z;
00064
00065 d = x*x + y*y + z*z;
00066
00067 if (d == 0.0)
00068 d = 0.0001;
00069 d = d*sqrt(d);
00070
00071 x /= d;
00072 y /= d;
00073 z /= d;
00074
00075 p0.fx -= x;
00076 p0.fy -= y;
00077 p0.fz -= z;
00078
00079 p1.fx += x;
00080 p1.fy += y;
00081 p1.fz += z;
00082 }
00083
00084
00085
00086
00087
00088 bool vuUDSphere::advancePoints(void)
00089 {
00090 double sx, sy, sz;
00091 double min;
00092 double t;
00093 double dt2;
00094 int i;
00095 int j;
00096 bool result;
00097
00098 result = false;
00099
00100
00101 for(i=1;i<m_Size;++i)
00102 {
00103 m_P[i].fx = 0.0;
00104 m_P[i].fy = 0.0;
00105 m_P[i].fz = 0.0;
00106 }
00107
00108
00109 for(i=1;i<m_Size;++i)
00110 for(j=i+1;j<m_Size;++j)
00111 computeForces(m_P[i], m_P[j]);
00112
00113 dt2 = m_dt*m_dt;
00114
00115
00116 for(i=1;i<m_Size;++i)
00117 {
00118 sx = m_P[i].x;
00119 sy = m_P[i].y;
00120 sz = m_P[i].z;
00121
00122 m_P[i].x += m_P[i].fx * dt2;
00123 m_P[i].y += m_P[i].fy * dt2;
00124 m_P[i].z += m_P[i].fz * dt2;
00125
00126
00127 t = sqrt(m_P[i].x*m_P[i].x + m_P[i].y*m_P[i].y + m_P[i].z*m_P[i].z);
00128 m_P[i].x /= t;
00129 m_P[i].y /= t;
00130 m_P[i].z /= t;
00131
00132
00133 if (!result)
00134 {
00135 sx -= m_P[i].x;
00136 sy -= m_P[i].y;
00137 sz -= m_P[i].z;
00138 t = sx*sx + sy*sy + sz*sz;
00139 if (t > m_Tolerance*m_Tolerance)
00140 result = true;
00141 }
00142 }
00143
00144
00145 min = m_P[1].x*m_P[1].x + (m_P[1].y-1.0)*(m_P[1].y-1.0) + m_P[1].z*m_P[1].z;
00146 j = 1;
00147 for(i=2;i<m_Size;++i)
00148 {
00149 t = m_P[i].x*m_P[i].x + (m_P[i].y-1.0)*(m_P[i].y-1.0) + m_P[i].z*m_P[i].z;
00150 if (t < min)
00151 {
00152 min = t;
00153 j = i;
00154 }
00155 }
00156 m_P[j].x = 0.0;
00157 m_P[j].y = 1.0;
00158 m_P[j].z = 0.0;
00159
00160 return result;
00161 }
00162
00163 bool vuUDSphere::doesFileExist() {
00164 FILE *fp;
00165 char name[32];
00166
00167 sprintf(name, "%i.nrm", m_Size-1);
00168 fp = fopen(name, "rb");
00169 if (fp) {
00170 fclose(fp);
00171 return true;
00172 }
00173 else
00174 return false;
00175 }
00176
00177 void vuUDSphere::writeToFile()
00178 {
00179 FILE* fp;
00180 float n[3];
00181 char name[32];
00182 int i;
00183 sprintf(name, "%i.nrm", m_Size-1);
00184 fp = fopen(name, "wb");
00185 if (fp) {
00186 for(i=1;i<m_Size;++i) {
00187 n[0] = (float)(m_P[i].x);
00188 n[1] = (float)(m_P[i].y);
00189 n[2] = (float)(m_P[i].z);
00190 fwrite(n, sizeof(float)*3, 1, fp);
00191 }
00192 fclose(fp);
00193 }
00194 }
00195
00196 void vuUDSphere::readFromFile()
00197 {
00198 FILE* fp;
00199 float n[3];
00200 char name[32];
00201 int i;
00202 sprintf(name, "%i.nrm", m_Size-1);
00203 fp = fopen(name, "rb");
00204 if (fp) {
00205 for(i=1;i<m_Size;++i) {
00206 fread(n, sizeof(float)*3, 1, fp);
00207 m_P[i].set(n[0], n[1], n[2]);
00208 }
00209 fclose(fp);
00210 }
00211 }
00212
00213 int vuUDSphere::getSize(void)
00214 {
00215 return m_Size-1;
00216 }
00217
00218 Point3d vuUDSphere::getPointAt(int idx)
00219 {
00220 assert(idx <= m_Size-1);
00221 return m_P[idx+1];
00222 }
00223
00224 void vuUDSphere::calculate(void) {
00225 while(advancePoints());
00226 }
00227
00228 void vuUDSphere::lazyCalculate(void) {
00229 if (doesFileExist())
00230 readFromFile();
00231 else {
00232 calculate();
00233 writeToFile();
00234 }
00235 }