开通博客赚积分
发布资源赚积分
分类

源码开发语言/平台
当前位置：首页 > 源码/资料 > 图形图像 > 3D图形编程 > 查看源码
Mesh.cpp
资源名称：3dSimplifier.zip　[点击查看]
上传用户：lin85885
上传日期：2013-04-27
资源大小：796k
文件大小：33k
源码类别：

3D图形编程

开发平台：

Visual C++

Mesh.cpp：源码内容
  1. // Mesh.cpp: implementation of the CMesh class.
  2. //
  3. //////////////////////////////////////////////////////////////////////
  5. #include "stdafx.h"
  6. #include <gl/gl.h>
  7. #include "Mesh.h"
  8. #include "Heap.h"
  10. #ifdef _DEBUG
  11. #undef THIS_FILE
  12. static char THIS_FILE[]=__FILE__;
  13. #define new DEBUG_NEW
  14. #endif
  16. //////////////////////////////////////////////////////////////////////
  17. // Construction/Destruction
  18. //////////////////////////////////////////////////////////////////////
  20. CMesh::CMesh()
  21. {
  22. m_Triangle = NULL;
  23. m_sTriangle = NULL;
  24. m_Vertex = NULL;
  25. m_sVertex = NULL;
  26. m_nTriangle = 0 ;
  27. m_nVertex = 0;
  28. m_VNormal = NULL;
  29. m_TNormal = NULL;
  30. m_sVNormal = NULL;
  31. m_sTNormal = NULL;
  33. m_pProgress = NULL;
  34. }
  36. CMesh::~CMesh()
  37. {
  38. if (m_Vertex)
  39. delete []m_Vertex;
  40. if (m_sVertex)
  41. delete []m_sVertex;
  42. if (m_Triangle)
  43. delete []m_Triangle;
  44. if (m_sTriangle)
  45. delete []m_sTriangle;
  46. if (m_VNormal)
  47. delete []m_VNormal;
  48. if (m_TNormal)
  49. delete []m_TNormal;
  50. if (m_sVNormal)
  51. delete []m_sVNormal;
  52. if (m_sTNormal)
  53. delete []m_sTNormal;
  54. }
  56. int CMesh::ReadData(FILE *fp)
  57. {
  58. int i,j;
  59. float box[2][3];
  60. float center[3];
  62. fscanf(fp,"%d%d",&m_nVertex,&m_nTriangle);
  64. if (m_pProgress)
  65. m_pProgress->SetRange(0, m_nVertex+m_nTriangle);
  66. m_Vertex = new Coord3[m_nVertex];
  67. m_Triangle = new Triangle[m_nTriangle];
  69. for (i=0;i<m_nVertex;i++)
  70. {
  71. fscanf(fp,"%f%f%f", &(m_Vertex[i][0]), &(m_Vertex[i][1]), &(m_Vertex[i][2]));
  72. m_pProgress->SetPos(i);
  73. }
  74. for (i=0;i<m_nTriangle;i++)
  75. {
  76. fscanf(fp,"%d%d%d", &(m_Triangle[i][0]), &(m_Triangle[i][1]), &(m_Triangle[i][2]));
  77. m_pProgress->SetPos(i+m_nVertex);
  78. }
  80. box[0][0] = box[0][1] = box[0][2] = MAXFLOAT;
  81. box[1][0] = box[1][1] = box[1][2] = -MAXFLOAT;
  82. for (i=0;i<m_nVertex;i++)
  83. for (j=0;j<3;j++)
  84. {
  85. if (box[0][j]>m_Vertex[i][j])
  86. box[0][j] = m_Vertex[i][j];
  87. if (box[1][j]<m_Vertex[i][j])
  88. box[1][j] = m_Vertex[i][j];
  89. }
  90. center[0] = (box[0][0]+box[1][0])/2.0;
  91. center[1] = (box[0][1]+box[1][1])/2.0;
  92. center[2] = (box[0][2]+box[1][2])/2.0;
  93. float scale;
  94. scale = FMAX(box[1][0]-box[0][0],FMAX(box[1][1]-box[0][1],box[1][2]-box[0][2]));
  95. scale /=2.0;
  96. for (i=0;i<m_nVertex;i++)
  97. {
  98. VEC3_VOPV_OP_S(m_Vertex[i],m_Vertex[i],-,center,/,scale);
  99. }
  100. VEC3_VOPV_OP_S(box[0],box[0],-,center,/,scale);
  101. VEC3_VOPV_OP_S(box[1],box[1],-,center,/,scale);
  102. VEC3_ZERO(center);
  104. m_VNormal = new Coord3[m_nVertex];
  105. m_TNormal = new Coord3[m_nTriangle];
  106. m_nsTriangle = m_nTriangle;
  107. m_nsVertex = m_nVertex;
  108. m_sTriangle = new Triangle[m_nTriangle];
  109. m_sVertex = new Coord3[m_nVertex];
  111. for (i=0;i<m_nVertex;i++)
  112. {
  113. VEC3_ASN_OP(m_sVertex[i],=,m_Vertex[i]);
  114. }
  115. for (i=0;i<m_nTriangle;i++)
  116. {
  117. VEC3_ASN_OP(m_sTriangle[i],=,m_Triangle[i]);
  118. }
  120. for (i=0;i<m_nVertex;i++)
  121. {
  122. VEC3_ZERO(m_VNormal[i]);
  123. }
  124. for (i=0;i<m_nTriangle;i++)
  125. {
  126. Coord3 vect[3];
  127. VEC3_V_OP_V(vect[0],m_Vertex[m_Triangle[i][1]],-,m_Vertex[m_Triangle[i][0]]);
  128. VEC3_V_OP_V(vect[1],m_Vertex[m_Triangle[i][2]],-,m_Vertex[m_Triangle[i][0]]);
  129. CROSSPROD3(vect[2],vect[0],vect[1]);
  130. V3Standarize(vect[2]);
  131. VEC3_ASN_OP(m_TNormal[i],=,vect[2]);
  132. for (j=0;j<3;j++)
  133. {
  134. VEC3_V_OP_V(m_VNormal[m_Triangle[i][j]],m_VNormal[m_Triangle[i][j]],+,vect[2]);
  135. }
  136. }
  137. for (i=0;i<m_nVertex;i++)
  138. V3Standarize(m_VNormal[i]);
  140. m_sVNormal = new Coord3[m_nVertex];
  141. m_sTNormal = new Coord3[m_nTriangle];
  142. for (i=0;i<m_nVertex;i++)
  143. {
  144. VEC3_ASN_OP(m_sVNormal[i],=,m_VNormal[i]);
  145. }
  147. for (i=0;i<m_nTriangle;i++)
  148. {
  149. VEC3_ASN_OP(m_sTNormal[i],=,m_TNormal[i]);
  150. }
  152. return m_nTriangle;
  153. }
  155. void CMesh::SaveData(FILE *fp)
  156. {
  157. int i;
  159. fprintf(fp,"%d %dn",m_nsVertex,m_nsTriangle);
  161. if (m_pProgress)
  162. m_pProgress->SetRange(0, m_nsVertex+m_nsTriangle);
  164. for (i=0;i<m_nsVertex;i++)
  165. {
  166. fprintf(fp,"%f %f %fn", m_sVertex[i][0], m_sVertex[i][1], m_sVertex[i][2]);
  167. m_pProgress->SetPos(i);
  168. }
  169. for (i=0;i<m_nsTriangle;i++)
  170. {
  171. fprintf(fp,"%d %d %dn", m_sTriangle[i][0], m_sTriangle[i][1], m_sTriangle[i][2]);
  172. m_pProgress->SetPos(i+m_nsVertex);
  173. }
  174. }
  176. int CMesh::DrawModel(int RenderMode, int RenderModel)
  177. {
  178. int i;
  180. if (RenderMode == 0)
  181. {
  182. glDisable(GL_LIGHTING);
  183. if (RenderModel == 0)
  184. {
  185. for (i=0; i<m_nTriangle; i++)
  186. {
  187. glBegin(GL_LINE_LOOP);
  188. glVertex3fv(m_Vertex[m_Triangle[i][0]]);
  189. glVertex3fv(m_Vertex[m_Triangle[i][1]]);
  190. glVertex3fv(m_Vertex[m_Triangle[i][2]]);
  191. glEnd();
  192. }
  193. }
  194. else if (RenderModel == 1)
  195. {
  196. for (i=0; i<m_nsTriangle; i++)
  197. {
  198. glBegin(GL_LINE_LOOP);
  199. glVertex3fv(m_sVertex[m_sTriangle[i][0]]);
  200. glVertex3fv(m_sVertex[m_sTriangle[i][1]]);
  201. glVertex3fv(m_sVertex[m_sTriangle[i][2]]);
  202. glEnd();
  203. }
  204. }
  205. glEnable(GL_LIGHTING);
  206. }
  207. else if (RenderMode == 1)
  208. {
  209. glBegin(GL_TRIANGLES);
  210. if (RenderModel == 0)
  211. {
  212. for (i=0; i<m_nTriangle; i++)
  213. {
  214. glNormal3fv(m_TNormal[i]);
  215. glVertex3fv(m_Vertex[m_Triangle[i][0]]);
  216. glVertex3fv(m_Vertex[m_Triangle[i][1]]);
  217. glVertex3fv(m_Vertex[m_Triangle[i][2]]);
  218. }
  219. }
  220. else if (RenderModel == 1)
  221. {
  222. for (i=0; i<m_nsTriangle; i++)
  223. {
  224. glNormal3fv(m_sTNormal[i]);
  225. glVertex3fv(m_sVertex[m_sTriangle[i][0]]);
  226. glVertex3fv(m_sVertex[m_sTriangle[i][1]]);
  227. glVertex3fv(m_sVertex[m_sTriangle[i][2]]);
  228. }
  229. }
  230. glEnd();
  231. }
  233. glEnable(GL_LIGHTING);
  234. if (RenderModel == 0)
  235. return m_nTriangle;
  236. else if (RenderModel == 1)
  237. return m_nsTriangle;
  238. else
  239. return 0;
  240. }
  242. int CMesh::Simplify(float ratio)
  243. {
  244.   Matrix4x4 *tri_Q;
  245.   Matrix4x4 *vex_Q;
  246.   Coord3 *Vertex;
  247.   Triangle *Mesh;
  248.   
  249.   int **vex_tri;
  250.   int **vRing;
  251.   int *Vex_map;
  252.   int i,j,k;
  253.   int tmp,tmp1,tmp2;
  255.   int **cad_pair;
  256.   float **pair_cost;
  257.   Coord3 **pair_pos;
  259.   BOOL *bMark;
  260.   int  *nMap;
  261.   int  *mMap;
  263.   if (m_nTriangle == 0)
  264.   return 0;
  266.   delete []m_sVertex;
  267.   delete []m_sTriangle;
  268.   delete []m_sVNormal;
  269.   delete []m_sTNormal;
  271.   Vertex= new Coord3[m_nVertex];
  272.   for (i=0;i<m_nVertex;i++)
  273.   {
  274.   VEC3_ASN_OP(Vertex[i],=,m_Vertex[i]);
  275.   }
  276.   Mesh = new Triangle[m_nTriangle];
  277.   for (i=0; i<m_nTriangle; i++)
  278.   {
  279.   VEC3_ASN_OP(Mesh[i],=,m_Triangle[i]);
  280.   }
  281.   
  282.   vex_tri=(int **)malloc(m_nVertex*sizeof(int *));
  283.   vRing=(int **)malloc(m_nVertex*sizeof(int *));
  284.   for (i=0;i<m_nVertex;i++) {
  285. vex_tri[i]=(int *)malloc(6*Size_int);
  286. vex_tri[i][0]=0;
  287. vRing[i]=(int *)malloc(6*Size_int);
  288. vRing[i][0]=0;
  289.   }
  290.   cad_pair=(int **)malloc(m_nVertex*sizeof(int *));
  291.   for (i=0;i<m_nVertex;i++) {
  292. cad_pair[i]=(int *)malloc(6*Size_int);
  293. cad_pair[i][0]=0;
  294.   }
  295.   pair_cost=(float **)malloc(m_nVertex*sizeof(float *));
  296.   pair_pos=(Coord3 **)malloc(m_nVertex*sizeof(Coord3 *));
  297.   Vex_map= new int[m_nVertex];
  298.   for (i=0;i<m_nVertex;i++)
  299. Vex_map[i]=i;
  301.   for (k=0;k<m_nTriangle;k++)
  302.   {
  303. for (i=0;i<3;i++) {
  304. tmp=tmp1=Mesh[k][i];
  305. tmp2=Mesh[k][(i+2)%3];
  306. for (j=1;j<vRing[tmp][0]+1;j++)
  307. if (vRing[tmp][j] == tmp2)
  308. break;
  309. if (j==vRing[tmp][0]+1)
  310. {
  311. vRing[tmp][j]=tmp2;
  312. vRing[tmp][0] += 1;
  313. if (!(vRing[tmp][0]%5))
  314. vRing[tmp] = (int *)realloc(vRing[tmp],(vRing[tmp][0]+6)*sizeof(int));
  315. }
  316. tmp2=Mesh[k][(i+1)%3];
  317. for (j=1;j<vRing[tmp][0]+1;j++)
  318. if (vRing[tmp][j] == tmp2)
  319. break;
  320. if (j==vRing[tmp][0]+1)
  321. {
  322. vRing[tmp][j]=tmp2;
  323. vRing[tmp][0] += 1;
  324. if (!(vRing[tmp][0]%5))
  325. vRing[tmp] = (int *)realloc(vRing[tmp],(vRing[tmp][0]+6)*sizeof(int));
  326. }
  327. if (vex_tri[tmp]==NULL) {
  328. vex_tri[tmp]=(int *)malloc(6*sizeof(int));
  329. vex_tri[tmp][0]=1;
  330. vex_tri[tmp][1]=k;
  331. }
  332. else {
  333. if (!(vex_tri[tmp][0]%5))
  334. vex_tri[tmp]=(int *)realloc(vex_tri[tmp],(vex_tri[tmp][0]+6)*sizeof(int));
  335. vex_tri[tmp][0]=vex_tri[tmp][0]+1;
  336. vex_tri[tmp][vex_tri[tmp][0]]=k;
  337. }
  338. if (tmp1>tmp2) {
  339. tmp=tmp1;
  340. tmp1=tmp2;
  341. tmp2=tmp;
  342. }
  343. for (j=1;j<=cad_pair[tmp1][0];j++)
  344. if (cad_pair[tmp1][j]==tmp2)
  345. break;
  346. if (j>cad_pair[tmp1][0]) {
  347. if (!(cad_pair[tmp1][0]%5))
  348. cad_pair[tmp1]=(int *)realloc(cad_pair[tmp1],(cad_pair[tmp1][0]+6)*sizeof(int));
  349. cad_pair[tmp1][0]=j;
  350. cad_pair[tmp1][j]=tmp2;
  351. }
  352. }
  354.   }
  356.   BOOL *border;
  357.   border = new BOOL[m_nVertex];
  358.   for (i=0;i<m_nVertex;i++)
  359.   {
  360.   if (vex_tri[i][0]==vRing[i][0])
  361.   border[i]=FALSE;
  362.   else
  363.   border[i]=TRUE;
  364.   free(vRing[i]);
  365.   }
  366.   free(vRing);
  368.   tri_Q=(Matrix4x4 *)malloc(m_nTriangle*Size_Matrix4x4);
  369.   cal_tri_Q(m_nTriangle,Mesh,Vertex,m_TNormal,tri_Q);
  370.   vex_Q=(Matrix4x4 *)malloc(m_nVertex*Size_Matrix4x4);
  371.   cal_vex_Q(m_nVertex,border,vex_tri,Mesh,Vertex,m_TNormal,tri_Q,vex_Q);
  372.   free(tri_Q);
  375.   bMark = new BOOL[m_nVertex];
  376.   for (i=0;i<m_nVertex;i++)
  377.   bMark[i]=FALSE;
  378.   int nvex = 0;
  379.   for (i=0;i<m_nTriangle;i++)
  380.   for (j=0;j<3;j++)
  381.   if (!bMark[Mesh[i][j]])
  382.   {
  383.   nvex++;
  384.   bMark[Mesh[i][j]] = TRUE;
  385.   }
  387.   for (i=0;i<m_nVertex;i++)
  388.   if (vex_tri[i])
  389. free(vex_tri[i]);
  390.   free(vex_tri);
  392.   tmp=(int)(m_nTriangle*ratio);
  394.   if (m_pProgress)
  395.   m_pProgress->SetRange(0,tmp);
  396.   
  397.   pair_contract(m_nVertex,border,cad_pair,pair_cost,pair_pos,
  398.         vex_Q,Vertex,Vex_map,tmp);
  400.   free(vex_Q);
  401.   delete []border;
  403.   nMap = new int[m_nVertex];
  404.   mMap = new int[m_nVertex];
  405.   for (i=0;i<m_nVertex;i++)
  406.   bMark[i]=FALSE;
  407.   for (i=0;i<m_nVertex;i++)
  408.   {
  409. tmp=i;
  410. while (Vex_map[tmp]!=tmp)
  411. tmp=Vex_map[tmp];
  412. Vex_map[i]=tmp;
  413.   }
  414.   m_nsTriangle=0;
  415.   for (i=0;i<m_nTriangle;i++)
  416.   if ((Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][1]])&&
  417.   (Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][2]])&&
  418.   (Vex_map[Mesh[i][1]]!=Vex_map[Mesh[i][2]]))
  419.   {
  420.   bMark[Vex_map[Mesh[i][0]]] = TRUE;
  421.   bMark[Vex_map[Mesh[i][1]]] = TRUE;
  422.   bMark[Vex_map[Mesh[i][2]]] = TRUE;
  423.   m_nsTriangle++;
  424.   }
  425.   m_nsVertex=0;
  426.   for (i=0;i<m_nVertex;i++)
  427.   if (bMark[i])
  428.   {
  429.   nMap[i] = m_nsVertex;
  430.   mMap[m_nsVertex] = i;
  431.   m_nsVertex++;
  432.   }
  433.   m_sVertex = new Coord3[m_nsVertex];
  434.   for (i=0;i<m_nsVertex;i++)
  435.   {
  436.   VEC3_ASN_OP(m_sVertex[i],=,Vertex[mMap[i]]);
  437.   }
  438.   m_sTriangle = new Triangle[m_nsTriangle];
  439.   m_nsTriangle = 0;
  440.   for (i=0;i<m_nTriangle;i++)
  441.   if ((Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][1]])&&
  442.   (Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][2]])&&
  443.   (Vex_map[Mesh[i][1]]!=Vex_map[Mesh[i][2]]))
  444.   {
  445.   m_sTriangle[m_nsTriangle][0] = nMap[Vex_map[Mesh[i][0]]];
  446.   m_sTriangle[m_nsTriangle][1] = nMap[Vex_map[Mesh[i][1]]];
  447.   m_sTriangle[m_nsTriangle][2] = nMap[Vex_map[Mesh[i][2]]];
  448.   m_nsTriangle++;
  449.   }
  451. m_sVNormal = new Coord3[m_nsVertex];
  452. m_sTNormal = new Coord3[m_nsTriangle];
  453. for (i=0;i<m_nsVertex;i++)
  454. {
  455. VEC3_ZERO(m_sVNormal[i]);
  456. }
  457. for (i=0;i<m_nsTriangle;i++)
  458. {
  459. Coord3 vect[3];
  460. VEC3_V_OP_V(vect[0],m_sVertex[m_sTriangle[i][1]],-,m_sVertex[m_sTriangle[i][0]]);
  461. VEC3_V_OP_V(vect[1],m_sVertex[m_sTriangle[i][2]],-,m_sVertex[m_sTriangle[i][0]]);
  462. CROSSPROD3(vect[2],vect[0],vect[1]);
  463. V3Standarize(vect[2]);
  464. VEC3_ASN_OP(m_sTNormal[i],=,vect[2]);
  465. for (j=0;j<3;j++)
  466. {
  467. VEC3_V_OP_V(m_sVNormal[m_sTriangle[i][j]],m_sVNormal[m_sTriangle[i][j]],+,vect[2]);
  468. }
  469. }
  470. for (i=0;i<m_nsVertex;i++)
  471. V3Standarize(m_sVNormal[i]);
  473.   delete []Vex_map;
  474.   delete []Mesh;
  475.   delete []Vertex;
  476.   delete []mMap;
  477.   delete []nMap;
  478.   delete []bMark;
  480.   return m_nsTriangle;
  481. }
  483. int CMesh::Simplify(int nFaceNo)
  484. {
  485.   Matrix4x4 *tri_Q;
  486.   Matrix4x4 *vex_Q;
  487.   Coord3 *Vertex;
  488.   Triangle *Mesh;
  489.   
  490.   int **vex_tri;
  491.   int **vRing;
  492.   int *Vex_map;
  493.   int i,j,k;
  494.   int tmp,tmp1,tmp2;
  496.   int **cad_pair;
  497.   float **pair_cost;
  498.   Coord3 **pair_pos;
  500.   BOOL *bMark;
  501.   int  *nMap;
  502.   int  *mMap;
  504.   if (m_nTriangle == 0)
  505.   return 0;
  507.   delete []m_sVertex;
  508.   delete []m_sTriangle;
  509.   delete []m_sVNormal;
  510.   delete []m_sTNormal;
  512.   Vertex= new Coord3[m_nVertex];
  513.   for (i=0;i<m_nVertex;i++)
  514.   {
  515.   VEC3_ASN_OP(Vertex[i],=,m_Vertex[i]);
  516.   }
  517.   Mesh = new Triangle[m_nTriangle];
  518.   for (i=0; i<m_nTriangle; i++)
  519.   {
  520.   VEC3_ASN_OP(Mesh[i],=,m_Triangle[i]);
  521.   }
  522.   
  523.   vex_tri=(int **)malloc(m_nVertex*sizeof(int *));
  524.   vRing=(int **)malloc(m_nVertex*sizeof(int *));
  525.   for (i=0;i<m_nVertex;i++) {
  526. vex_tri[i]=(int *)malloc(6*Size_int);
  527. vex_tri[i][0]=0;
  528. vRing[i]=(int *)malloc(6*Size_int);
  529. vRing[i][0]=0;
  530.   }
  531.   cad_pair=(int **)malloc(m_nVertex*sizeof(int *));
  532.   for (i=0;i<m_nVertex;i++) {
  533. cad_pair[i]=(int *)malloc(6*Size_int);
  534. cad_pair[i][0]=0;
  535.   }
  536.   pair_cost=(float **)malloc(m_nVertex*sizeof(float *));
  537.   pair_pos=(Coord3 **)malloc(m_nVertex*sizeof(Coord3 *));
  538.   Vex_map= new int[m_nVertex];
  539.   for (i=0;i<m_nVertex;i++)
  540. Vex_map[i]=i;
  542.   for (k=0;k<m_nTriangle;k++)
  543.   {
  544. for (i=0;i<3;i++) {
  545. tmp=tmp1=Mesh[k][i];
  546. tmp2=Mesh[k][(i+2)%3];
  547. for (j=1;j<vRing[tmp][0]+1;j++)
  548. if (vRing[tmp][j] == tmp2)
  549. break;
  550. if (j==vRing[tmp][0]+1)
  551. {
  552. vRing[tmp][j]=tmp2;
  553. vRing[tmp][0] += 1;
  554. if (!(vRing[tmp][0]%5))
  555. vRing[tmp] = (int *)realloc(vRing[tmp],(vRing[tmp][0]+6)*sizeof(int));
  556. }
  557. tmp2=Mesh[k][(i+1)%3];
  558. for (j=1;j<vRing[tmp][0]+1;j++)
  559. if (vRing[tmp][j] == tmp2)
  560. break;
  561. if (j==vRing[tmp][0]+1)
  562. {
  563. vRing[tmp][j]=tmp2;
  564. vRing[tmp][0] += 1;
  565. if (!(vRing[tmp][0]%5))
  566. vRing[tmp] = (int *)realloc(vRing[tmp],(vRing[tmp][0]+6)*sizeof(int));
  567. }
  568. if (vex_tri[tmp]==NULL) {
  569. vex_tri[tmp]=(int *)malloc(6*sizeof(int));
  570. vex_tri[tmp][0]=1;
  571. vex_tri[tmp][1]=k;
  572. }
  573. else {
  574. if (!(vex_tri[tmp][0]%5))
  575. vex_tri[tmp]=(int *)realloc(vex_tri[tmp],(vex_tri[tmp][0]+6)*sizeof(int));
  576. vex_tri[tmp][0]=vex_tri[tmp][0]+1;
  577. vex_tri[tmp][vex_tri[tmp][0]]=k;
  578. }
  579. if (tmp1>tmp2) {
  580. tmp=tmp1;
  581. tmp1=tmp2;
  582. tmp2=tmp;
  583. }
  584. for (j=1;j<=cad_pair[tmp1][0];j++)
  585. if (cad_pair[tmp1][j]==tmp2)
  586. break;
  587. if (j>cad_pair[tmp1][0]) {
  588. if (!(cad_pair[tmp1][0]%5))
  589. cad_pair[tmp1]=(int *)realloc(cad_pair[tmp1],(cad_pair[tmp1][0]+6)*sizeof(int));
  590. cad_pair[tmp1][0]=j;
  591. cad_pair[tmp1][j]=tmp2;
  592. }
  593. }
  595.   }
  597.   BOOL *border;
  598.   border = new BOOL[m_nVertex];
  599.   for (i=0;i<m_nVertex;i++)
  600.   {
  601.   if (vex_tri[i][0]==vRing[i][0])
  602.   border[i]=FALSE;
  603.   else
  604.   border[i]=TRUE;
  605.   free(vRing[i]);
  606.   }
  607.   free(vRing);
  609.   tri_Q=(Matrix4x4 *)malloc(m_nTriangle*Size_Matrix4x4);
  610.   cal_tri_Q(m_nTriangle,Mesh,Vertex,m_TNormal,tri_Q);
  611.   vex_Q=(Matrix4x4 *)malloc(m_nVertex*Size_Matrix4x4);
  612.   cal_vex_Q(m_nVertex,border,vex_tri,Mesh,Vertex,m_TNormal,tri_Q,vex_Q);
  613.   free(tri_Q);
  616.   bMark = new BOOL[m_nVertex];
  617.   for (i=0;i<m_nVertex;i++)
  618.   bMark[i]=FALSE;
  619.   int nvex = 0;
  620.   for (i=0;i<m_nTriangle;i++)
  621.   for (j=0;j<3;j++)
  622.   if (!bMark[Mesh[i][j]])
  623.   {
  624.   nvex++;
  625.   bMark[Mesh[i][j]] = TRUE;
  626.   }
  628.   for (i=0;i<m_nVertex;i++)
  629.   if (vex_tri[i])
  630. free(vex_tri[i]);
  631.   free(vex_tri);
  633.   tmp = m_nTriangle - nFaceNo;
  635.   if (m_pProgress)
  636.   m_pProgress->SetRange(0,tmp);
  638.   pair_contract(m_nVertex,border,cad_pair,pair_cost,pair_pos,
  639.         vex_Q,Vertex,Vex_map,tmp);
  641.   free(vex_Q);
  642.   delete []border;
  644.   nMap = new int[m_nVertex];
  645.   mMap = new int[m_nVertex];
  646.   for (i=0;i<m_nVertex;i++)
  647.   bMark[i]=FALSE;
  648.   for (i=0;i<m_nVertex;i++)
  649.   {
  650. tmp=i;
  651. while (Vex_map[tmp]!=tmp)
  652. tmp=Vex_map[tmp];
  653. Vex_map[i]=tmp;
  654.   }
  655.   m_nsTriangle=0;
  656.   for (i=0;i<m_nTriangle;i++)
  657.   if ((Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][1]])&&
  658.   (Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][2]])&&
  659.   (Vex_map[Mesh[i][1]]!=Vex_map[Mesh[i][2]]))
  660.   {
  661.   bMark[Vex_map[Mesh[i][0]]] = TRUE;
  662.   bMark[Vex_map[Mesh[i][1]]] = TRUE;
  663.   bMark[Vex_map[Mesh[i][2]]] = TRUE;
  664.   m_nsTriangle++;
  665.   }
  666.   m_nsVertex=0;
  667.   for (i=0;i<m_nVertex;i++)
  668.   if (bMark[i])
  669.   {
  670.   nMap[i] = m_nsVertex;
  671.   mMap[m_nsVertex] = i;
  672.   m_nsVertex++;
  673.   }
  674.   m_sVertex = new Coord3[m_nsVertex];
  675.   for (i=0;i<m_nsVertex;i++)
  676.   {
  677.   VEC3_ASN_OP(m_sVertex[i],=,Vertex[mMap[i]]);
  678.   }
  679.   m_sTriangle = new Triangle[m_nsTriangle];
  680.   m_nsTriangle = 0;
  681.   for (i=0;i<m_nTriangle;i++)
  682.   if ((Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][1]])&&
  683.   (Vex_map[Mesh[i][0]]!=Vex_map[Mesh[i][2]])&&
  684.   (Vex_map[Mesh[i][1]]!=Vex_map[Mesh[i][2]]))
  685.   {
  686.   m_sTriangle[m_nsTriangle][0] = nMap[Vex_map[Mesh[i][0]]];
  687.   m_sTriangle[m_nsTriangle][1] = nMap[Vex_map[Mesh[i][1]]];
  688.   m_sTriangle[m_nsTriangle][2] = nMap[Vex_map[Mesh[i][2]]];
  689.   m_nsTriangle++;
  690.   }
  692. m_sVNormal = new Coord3[m_nsVertex];
  693. m_sTNormal = new Coord3[m_nsTriangle];
  694. for (i=0;i<m_nsVertex;i++)
  695. {
  696. VEC3_ZERO(m_sVNormal[i]);
  697. }
  698. for (i=0;i<m_nsTriangle;i++)
  699. {
  700. Coord3 vect[3];
  701. VEC3_V_OP_V(vect[0],m_sVertex[m_sTriangle[i][1]],-,m_sVertex[m_sTriangle[i][0]]);
  702. VEC3_V_OP_V(vect[1],m_sVertex[m_sTriangle[i][2]],-,m_sVertex[m_sTriangle[i][0]]);
  703. CROSSPROD3(vect[2],vect[0],vect[1]);
  704. V3Standarize(vect[2]);
  705. VEC3_ASN_OP(m_sTNormal[i],=,vect[2]);
  706. for (j=0;j<3;j++)
  707. {
  708. VEC3_V_OP_V(m_sVNormal[m_sTriangle[i][j]],m_sVNormal[m_sTriangle[i][j]],+,vect[2]);
  709. }
  710. }
  711. for (i=0;i<m_nsVertex;i++)
  712. V3Standarize(m_sVNormal[i]);
  714.   delete []Vex_map;
  715.   delete []Mesh;
  716.   delete []Vertex;
  717.   delete []mMap;
  718.   delete []nMap;
  719.   delete []bMark;
  721.    return m_nsTriangle;
  722. }
  724. void CMesh::cal_tri_Q(int Num_triangle,Triangle *mesh,Coord3 *vertex,Vector3 *tri_normal,Matrix4x4 *tri_Q)
  725. {
  726.   int i;
  727.   float tmp;
  729.   for (i=0;i<Num_triangle;i++) {
  730. tmp=-DOTPROD3(tri_normal[i],vertex[mesh[i][0]]);
  731. tri_Q[i][0][0]=tri_normal[i][0]*tri_normal[i][0];
  732. tri_Q[i][1][1]=tri_normal[i][1]*tri_normal[i][1];
  733. tri_Q[i][2][2]=tri_normal[i][2]*tri_normal[i][2];
  734. tri_Q[i][3][3]=tmp*tmp;
  735. tri_Q[i][0][1]=tri_Q[i][1][0]=tri_normal[i][0]*tri_normal[i][1];
  736. tri_Q[i][0][2]=tri_Q[i][2][0]=tri_normal[i][0]*tri_normal[i][2];
  737. tri_Q[i][0][3]=tri_Q[i][3][0]=tri_normal[i][0]*tmp;
  738. tri_Q[i][1][2]=tri_Q[i][2][1]=tri_normal[i][1]*tri_normal[i][2];
  739. tri_Q[i][1][3]=tri_Q[i][3][1]=tri_normal[i][1]*tmp;
  740. tri_Q[i][2][3]=tri_Q[i][3][2]=tri_normal[i][2]*tmp;
  741.   }
  742. }
  744. void CMesh::cal_vex_Q(int Num_vertex,BOOL *border,int **vex_tri,Triangle *mesh,Coord3 *vertex,Vector3 *tri_normal,Matrix4x4 *tri_Q,Matrix4x4 *vex_Q)
  745. {
  746.   int i,j,k,m,p;
  747.   Matrix4x4 tmp_Q;
  748.   float tmp;
  749.   Vector3 vect[2];
  751.   for (i=0;i<Num_vertex;i++) {
  752. for (j=0;j<4;j++)
  753. VEC4_ZERO(vex_Q[i][j]);
  754. for (j=1;j<=vex_tri[i][0];j++) {
  755. p=vex_tri[i][j];
  756. for (k=0;k<4;k++)
  757. VEC4_V_OP_V(vex_Q[i][k],vex_Q[i][k],+,tri_Q[p][k]);
  758. for (k=0;k<3;k++)
  759. if (mesh[p][k]==i)
  760. break;
  761. int a,b;
  762. if (mesh[p][k]>mesh[p][(k+2)%3])
  763. {
  764. a=mesh[p][k];
  765. b=mesh[p][(k+2)%3];
  766. }
  767. else
  768. {
  769. b=mesh[p][k];
  770. a=mesh[p][(k+2)%3];
  771. }
  772. if (border[a]&&border[b]) {
  773. VEC3_V_OP_V(vect[0],vertex[mesh[p][(k+2)%3]],-,vertex[i]);
  774. V3Standarize(vect[0]);
  775. CROSSPROD3(vect[1],vect[0],tri_normal[p]);
  776. V3Standarize(vect[1]);
  777. tmp=-DOTPROD3(vect[1],vertex[i]);
  778. tmp_Q[0][0]=vect[1][0]*vect[1][0];
  779. tmp_Q[1][1]=vect[1][1]*vect[1][1];
  780. tmp_Q[2][2]=vect[1][2]*vect[1][2];
  781. tmp_Q[3][3]=tmp*tmp;
  782. tmp_Q[0][1]=tmp_Q[1][0]=vect[1][0]*vect[1][1];
  783. tmp_Q[0][2]=tmp_Q[2][0]=vect[1][0]*vect[1][2];
  784. tmp_Q[0][3]=tmp_Q[3][0]=vect[1][0]*tmp;
  785. tmp_Q[1][2]=tmp_Q[2][1]=vect[1][1]*vect[1][2];
  786. tmp_Q[1][3]=tmp_Q[3][1]=vect[1][1]*tmp;
  787. tmp_Q[2][3]=tmp_Q[3][2]=vect[1][2]*tmp;
  788. for (m=0;m<4;m++)
  789. VEC4_V_OP_V(vex_Q[i][m],vex_Q[i][m],+,tmp_Q[m]);
  790. }
  791. if (mesh[p][k]>mesh[p][(k+1)%3])
  792. {
  793. a=mesh[p][k];
  794. b=mesh[p][(k+1)%3];
  795. }
  796. else
  797. {
  798. b=mesh[p][k];
  799. a=mesh[p][(k+1)%3];
  800. }
  801. if (border[a]&&border[b]) {
  802. VEC3_V_OP_V(vect[0],vertex[mesh[p][(k+1)%3]],-,vertex[i]);
  803. V3Standarize(vect[0]);
  804. CROSSPROD3(vect[1],vect[0],tri_normal[p]);
  805. V3Standarize(vect[1]);
  806. tmp=-DOTPROD3(vect[1],vertex[i]);
  807. tmp_Q[0][0]=vect[1][0]*vect[1][0];
  808. tmp_Q[1][1]=vect[1][1]*vect[1][1];
  809. tmp_Q[2][2]=vect[1][2]*vect[1][2];
  810. tmp_Q[3][3]=tmp*tmp;
  811. tmp_Q[0][1]=tmp_Q[1][0]=vect[1][0]*vect[1][1];
  812. tmp_Q[0][2]=tmp_Q[2][0]=vect[1][0]*vect[1][2];
  813. tmp_Q[0][3]=tmp_Q[3][0]=vect[1][0]*tmp;
  814. tmp_Q[1][2]=tmp_Q[2][1]=vect[1][1]*vect[1][2];
  815. tmp_Q[1][3]=tmp_Q[3][1]=vect[1][1]*tmp;
  816. tmp_Q[2][3]=tmp_Q[3][2]=vect[1][2]*tmp;
  817. for (m=0;m<4;m++)
  818. VEC4_V_OP_V(vex_Q[i][m],vex_Q[i][m],+,tmp_Q[m]);
  819. }
  820. }
  821.   }
  822. }
  824. void CMesh::pair_contract(int Num_vertex,BOOL *border,int **cad_pair,float **pair_cost,Coord3 **pair_pos,Matrix4x4 *vex_Q,Coord3 *vertex,int *vex_map,int require)
  825. {
  826.   int i,j,k,m,n,min_i,min_j;
  827.   float **tmp_Q;
  828.   float M3x3[3][3];
  829.   float vec_b[4];
  830.   int remain,deletenum;
  831.   float min_cost;
  833.   tmp_Q=(float **)malloc(4*sizeof(float *));
  834.   for (i=0;i<4;i++)
  835. tmp_Q[i]=(float *)malloc(4*Size_float);
  837.   for (i=0;i<Num_vertex;i++) {
  838. pair_cost[i]=(float *)malloc(cad_pair[i][0]*Size_float);
  839. pair_pos[i]=(Coord3 *)malloc(cad_pair[i][0]*Size_Coord3);
  840. for (k=1;k<=cad_pair[i][0];k++) {
  841. j=cad_pair[i][k];
  842. for (m=0;m<3;m++)
  843. VEC4_V_OP_V(tmp_Q[m],vex_Q[i][m],+,vex_Q[j][m]);
  844. tmp_Q[3][0]=tmp_Q[3][1]=tmp_Q[3][2]=0.0;
  845. tmp_Q[3][3]=1.0;
  846. vec_b[0]=vec_b[1]=vec_b[2]=0.0;
  847. vec_b[3]=1.0;
  848. for (m=0;m<3;m++)
  849. VEC3_ASN_OP(M3x3[m],=,tmp_Q[m]);
  850. if ((fabs(deta3(M3x3))>=0.01)&&(Gauss(4,tmp_Q,vec_b))) {
  851. VEC3_ASN_OP(pair_pos[i][k-1],=,vec_b);
  852. for (m=0;m<4;m++)
  853. VEC4_V_OP_V(tmp_Q[m],vex_Q[i][m],+,vex_Q[j][m]);
  854. pair_cost[i][k-1]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  855.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  856.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  857.                   2*tmp_Q[0][3]*vec_b[0]+
  858.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  859.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  860.                   2*tmp_Q[1][3]*vec_b[1]+
  861.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  862.                   2*tmp_Q[2][3]*vec_b[2]+
  863.   tmp_Q[3][3];
  864. }
  865. else {
  867. for (m=0;m<4;m++)
  868. VEC4_V_OP_V(tmp_Q[m],vex_Q[i][m],+,vex_Q[j][m]);
  869. VEC3_ASN_OP(vec_b,=,vertex[i]);
  870. VEC3_ASN_OP(pair_pos[i][k-1],=,vec_b);
  871. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  872.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  873.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  874.                   2*tmp_Q[0][3]*vec_b[0]+
  875.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  876.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  877.                   2*tmp_Q[1][3]*vec_b[1]+
  878.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  879.                   2*tmp_Q[2][3]*vec_b[2]+
  880.   tmp_Q[3][3];
  881. min_cost=vec_b[3];
  882. VEC3_ASN_OP(vec_b,=,vertex[j]);
  883. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  884.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  885.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  886.                   2*tmp_Q[0][3]*vec_b[0]+
  887.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  888.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  889.                   2*tmp_Q[1][3]*vec_b[1]+
  890.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  891.                   2*tmp_Q[2][3]*vec_b[2]+
  892.   tmp_Q[3][3];
  893. if (min_cost>vec_b[3]) {
  894. min_cost=vec_b[3];
  895. VEC3_ASN_OP(pair_pos[i][k-1],=,vec_b);
  896. }
  897. VEC3_V_OP_V(vec_b,vertex[i],+,vertex[j]);
  898. VEC3_V_OP_S(vec_b,vec_b,/,2.0f);
  899. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  900.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  901.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  902.                   2*tmp_Q[0][3]*vec_b[0]+
  903.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  904.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  905.                   2*tmp_Q[1][3]*vec_b[1]+
  906.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  907.                   2*tmp_Q[2][3]*vec_b[2]+
  908.   tmp_Q[3][3];
  909. if (min_cost>vec_b[3]) {
  910. min_cost=vec_b[3];
  911. VEC3_ASN_OP(pair_pos[i][k-1],=,vec_b);
  912. }
  913. pair_cost[i][k-1]=min_cost;
  914. }
  915. }
  916.   }
  919.   // 99,6,1
  920.   int numEdge = 0;
  921.   Heap *heap;
  922.   Heapable *pair;
  923.   int **VexEdge;
  924.   //VexEdge = new int51[Num_vertex];
  925.   VexEdge = new Pint[Num_vertex];
  926.   for (i=0;i<Num_vertex;i++)
  927.   {
  928.   VexEdge[i] = (int *)malloc(51*sizeof(int));
  929.   VexEdge[i][0] = 0;
  930.   numEdge += cad_pair[i][0];
  931.   }
  932.   heap = new Heap(numEdge);
  933.   pair = new Heapable[numEdge];
  935.   k = 0;
  936.   for (i=0;i<Num_vertex;i++)
  937. for (j=1; j<=cad_pair[i][0]; j++)
  938. {
  940. pair[k].vert[0] = i;
  941. m = pair[k].vert[1] = cad_pair[i][j];
  942. pair[k].cost = -pair_cost[i][j-1];
  943. VEC3_ASN_OP(pair[k].newposition,=,pair_pos[i][j-1])
  944. heap->insert(&(pair[k]), pair[k].cost);
  945. if (!(pair[k].isInHeap()))
  946. break;
  947. n = VexEdge[i][0]+1;
  948. if (n>50)
  949. {
  950. VexEdge[i] = (int *)realloc(VexEdge[i],(n+50)*sizeof(int));
  951. }
  952. VexEdge[i][n] = k;
  953. VexEdge[i][0] = n;
  954. n = VexEdge[m][0]+1;
  955. if (n>50)
  956. {
  957. VexEdge[m] = (int *)realloc(VexEdge[m],(n+50)*sizeof(int));
  958. }
  959. VexEdge[m][n] = k;
  960. VexEdge[m][0] = n;
  961. k++;
  962. }
  964.   for (i=0;i<Num_vertex;i++)
  965.   {
  966.   if (pair_cost[i])
  967.   free(pair_cost[i]);
  968.   if (pair_pos[i])
  969.   free(pair_pos[i]);
  970.   if (cad_pair[i])
  971.   free(cad_pair[i]);
  972.   }
  973.   free(pair_cost);
  974.   free(pair_pos);
  975.   free(cad_pair);
  977.   // 99,6,1
  979.   deletenum=0;
  980.   remain=Num_vertex;  
  981.   int delEdge = 0;
  982.   
  983.   while (1) {
  985. heap_node *top;
  986. Heapable *dpair;
  988. top = heap->extract();
  989. if (!top)
  990. break;
  991. dpair = top->obj;
  993. min_i = dpair->vert[0];
  994. min_j = dpair->vert[1];
  995. if (border[min_i]&&border[min_j])
  996. deletenum ++;
  997. else
  998. deletenum +=2;
  999. border[min_i] = border[min_i]||border[min_j];
  1001. VEC3_ASN_OP(vertex[min_i],=,dpair->newposition);
  1002. VEC3_ASN_OP(vertex[min_j],=,dpair->newposition);
  1003. vex_map[min_j]=min_i;
  1005. delEdge++;
  1007. for (i=1;i<=VexEdge[min_i][0];i++)
  1008. {
  1009. m = VexEdge[min_i][i];
  1010. if ((pair[m].vert[0] == min_j)||(pair[m].vert[1] == min_j))
  1011. {
  1012. for (j=i; j<VexEdge[min_i][0]; j++)
  1013. VexEdge[min_i][j] = VexEdge[min_i][j+1];
  1014. VexEdge[min_i][0] = VexEdge[min_i][0]-1;
  1015. break;
  1016. }
  1017. }
  1019. for (i=1;i<=VexEdge[min_j][0];i++)
  1020. {
  1021. int p,q;
  1023. m = VexEdge[min_j][i];
  1024. if ((pair[m].vert[0] == min_i)||(pair[m].vert[1] == min_i))
  1025. continue;
  1026. if (pair[m].vert[0] == min_j)
  1027. {
  1028. pair[m].vert[0] = min_i;
  1029. for (j=1;j<=VexEdge[min_i][0]; j++)
  1030. {
  1031. n = VexEdge[min_i][j];
  1032. if (pair[n].vert[0] == pair[m].vert[1])
  1033. {
  1034. k = pair[n].vert[0];
  1035. heap->kill(pair[m].getHeapPos());
  1036. for (p=1; p<=VexEdge[k][0]; p++)
  1037. {
  1038. if (VexEdge[k][p] == m)
  1039. {
  1040. for (q=p; q<VexEdge[k][0]; q++)
  1041. VexEdge[k][q]=VexEdge[k][q+1];
  1042. VexEdge[k][0] = VexEdge[k][0]-1;
  1043. break;
  1044. }
  1045. }
  1046. break;
  1047. }
  1048. else if (pair[n].vert[1] == pair[m].vert[1])
  1049. {
  1050. k = pair[n].vert[1];
  1051. heap->kill(pair[m].getHeapPos());
  1052. for (p=1; p<=VexEdge[k][0]; p++)
  1053. {
  1054. if (VexEdge[k][p] == m)
  1055. {
  1056. for (q=p; q<VexEdge[k][0]; q++)
  1057. VexEdge[k][q]=VexEdge[k][q+1];
  1058. VexEdge[k][0] = VexEdge[k][0]-1;
  1059. break;
  1060. }
  1061. }
  1062. break;
  1063. }
  1064. }
  1066. if (j>VexEdge[min_i][0])
  1067. {
  1068. if (j>50)
  1069. {
  1070. VexEdge[min_i] = (int *)realloc(VexEdge[min_i],(j+50)*sizeof(int));
  1071. }
  1072. VexEdge[min_i][0] = j;
  1073. VexEdge[min_i][j] = m;
  1074. }
  1075. }
  1076. else if (pair[m].vert[1] == min_j)
  1077. {
  1078. pair[m].vert[1] = min_i;
  1079. for (j=1;j<=VexEdge[min_i][0]; j++)
  1080. {
  1081. n = VexEdge[min_i][j];
  1082. if (pair[n].vert[0] == pair[m].vert[0])
  1083. {
  1084. k = pair[n].vert[0];
  1085. heap->kill(pair[m].getHeapPos());
  1086. for (p=1; p<=VexEdge[k][0]; p++)
  1087. {
  1088. if (VexEdge[k][p] == m)
  1089. {
  1090. for (q=p; q<VexEdge[k][0]; q++)
  1091. VexEdge[k][q]=VexEdge[k][q+1];
  1092. VexEdge[k][0] = VexEdge[k][0]-1;
  1093. break;
  1094. }
  1095. }
  1096. break;
  1097. }
  1098. else if (pair[n].vert[1] == pair[m].vert[0])
  1099. {
  1100. k = pair[n].vert[1];
  1101. heap->kill(pair[m].getHeapPos());
  1102. for (p=1; p<=VexEdge[k][0]; p++)
  1103. {
  1104. if (VexEdge[k][p] == m)
  1105. {
  1106. for (q=p; q<VexEdge[k][0]; q++)
  1107. VexEdge[k][q]=VexEdge[k][q+1];
  1108. VexEdge[k][0] = VexEdge[k][0]-1;
  1109. break;
  1110. }
  1111. }
  1112. break;
  1113. }
  1114. }
  1116. if (j>VexEdge[min_i][0])
  1117. {
  1118. if (j>50)
  1119. {
  1120. VexEdge[min_i] = (int *)realloc(VexEdge[min_i],(j+50)*sizeof(int));
  1121. }
  1122. VexEdge[min_i][0] = j;
  1123. VexEdge[min_i][j] = m;
  1124. }
  1125. }
  1126. }
  1128. for (i=0;i<4;i++)
  1129. VEC4_V_OP_V(vex_Q[min_i][i],vex_Q[min_i][i],+,vex_Q[min_j][i]);
  1131. for (i=1; i<=VexEdge[min_i][0];i++)
  1132. {
  1133. n = VexEdge[min_i][i];
  1134. j = pair[n].vert[0];
  1135. if (j == min_i )
  1136. j = pair[n].vert[1];
  1137. for (m=0;m<3;m++)
  1138. VEC4_V_OP_V(tmp_Q[m],vex_Q[min_i][m],+,vex_Q[j][m]);
  1139. tmp_Q[3][0]=tmp_Q[3][1]=tmp_Q[3][2]=0.0;
  1140. tmp_Q[3][3]=1.0;
  1141. vec_b[0]=vec_b[1]=vec_b[2]=0.0;
  1142. vec_b[3]=1.0;
  1143. for (m=0;m<3;m++)
  1144. VEC3_ASN_OP(M3x3[m],=,tmp_Q[m]);
  1145. if ((fabs(deta3(M3x3))>=0.01)&&(Gauss(4,tmp_Q,vec_b))) {
  1146. VEC3_ASN_OP(pair[n].newposition,=,vec_b);
  1147. for (m=0;m<4;m++)
  1148. VEC4_V_OP_V(tmp_Q[m],vex_Q[min_i][m],+,vex_Q[j][m]);
  1149. pair[n].cost  =   tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  1150.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  1151.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  1152.                   2*tmp_Q[0][3]*vec_b[0]+
  1153.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  1154.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  1155.                   2*tmp_Q[1][3]*vec_b[1]+
  1156.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  1157.                   2*tmp_Q[2][3]*vec_b[2]+
  1158.   tmp_Q[3][3];
  1159. }
  1160. else {
  1162. for (m=0;m<4;m++)
  1163. VEC4_V_OP_V(tmp_Q[m],vex_Q[min_i][m],+,vex_Q[j][m]);
  1164. VEC3_ASN_OP(vec_b,=,vertex[min_i]);
  1165. VEC3_ASN_OP(pair[n].newposition,=,vec_b);
  1166. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  1167.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  1168.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  1169.                   2*tmp_Q[0][3]*vec_b[0]+
  1170.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  1171.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  1172.                   2*tmp_Q[1][3]*vec_b[1]+
  1173.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  1174.                   2*tmp_Q[2][3]*vec_b[2]+
  1175.   tmp_Q[3][3];
  1176. min_cost=vec_b[3];
  1177. VEC3_ASN_OP(vec_b,=,vertex[j]);
  1178. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  1179.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  1180.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  1181.                   2*tmp_Q[0][3]*vec_b[0]+
  1182.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  1183.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  1184.                   2*tmp_Q[1][3]*vec_b[1]+
  1185.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  1186.                   2*tmp_Q[2][3]*vec_b[2]+
  1187.   tmp_Q[3][3];
  1188. if (min_cost>vec_b[3]) {
  1189. min_cost=vec_b[3];
  1190. VEC3_ASN_OP(pair[n].newposition,=,vec_b);
  1191. }
  1192. VEC3_V_OP_V(vec_b,vertex[min_i],+,vertex[j]);
  1193. VEC3_V_OP_S(vec_b,vec_b,/,2.0);
  1194. vec_b[3]=tmp_Q[0][0]*vec_b[0]*vec_b[0]+
  1195.                   2*tmp_Q[0][1]*vec_b[0]*vec_b[1]+
  1196.                   2*tmp_Q[0][2]*vec_b[0]*vec_b[2]+
  1197.                   2*tmp_Q[0][3]*vec_b[0]+
  1198.                   tmp_Q[1][1]*vec_b[1]*vec_b[1]+
  1199.                   2*tmp_Q[1][2]*vec_b[1]*vec_b[2]+
  1200.                   2*tmp_Q[1][3]*vec_b[1]+
  1201.                   tmp_Q[2][2]*vec_b[2]*vec_b[2]+
  1202.                   2*tmp_Q[2][3]*vec_b[2]+
  1203.   tmp_Q[3][3];
  1204. if (min_cost>vec_b[3]) {
  1205. min_cost=vec_b[3];
  1206. VEC3_ASN_OP(pair[n].newposition,=,vec_b);
  1207. }
  1208. pair[n].cost=min_cost;
  1209. }
  1210. pair[n].cost = -pair[n].cost;
  1211. heap->update(&(pair[n]),pair[n].cost);
  1212. }
  1214. if (m_pProgress)
  1215. m_pProgress->SetPos(deletenum);
  1217. if (deletenum >= require)
  1218. break;
  1220.   } /* while (1) */
  1222.   delete heap;
  1223.   delete []pair;
  1224.   for (i=0;i<Num_vertex;i++)
  1225.   {
  1226.   free(VexEdge[i]);
  1227.   }
  1228.   delete []VexEdge;
  1230.   for (i=0;i<4;i++)
  1231. free(tmp_Q[i]);
  1232.   free(tmp_Q);
  1233. }
  1235. void CMesh::V3Standarize(Vector3 v)
  1236. {
  1237. float len;
  1238. len=sqrt(DOTPROD3(v,v));
  1239. if (len!=0.0)
  1240. {
  1241. v[X]=v[X]/len;
  1242. v[Y]=v[Y]/len;
  1243. v[Z]=v[Z]/len;
  1244. }
  1245. }
  1247. float CMesh::deta3(float matrix[3][3])
  1248. {
  1249.   float deta;
  1250.   deta=matrix[0][0]*(matrix[1][1]*matrix[2][2]-matrix[1][2]*matrix[2][1])-
  1251.        matrix[0][1]*(matrix[1][0]*matrix[2][2]-matrix[1][2]*matrix[2][0])+
  1252.        matrix[0][2]*(matrix[1][0]*matrix[2][1]-matrix[1][1]*matrix[2][0]);
  1253.   return deta;
  1254. }
  1256. int CMesh::Gauss(int n,float **A,float *b)
  1257. {
  1258.   int *js,i,j,k,is;
  1259.   float d,t;
  1261.   js=(int *)malloc(n*sizeof(int));
  1262.   for (k=0;k<n;k++) {
  1263. d=0.0;
  1264. for (i=k;i<n;i++)
  1265. for (j=k;j<n;j++)
  1266. if (A[i][j]!=0.0) {
  1267. t=fabs(A[i][j]);
  1268. if (t>d) {
  1269. d=t;
  1270. js[k]=j;
  1271. is=i;
  1272. }
  1273. }
  1275. if (d+1.0==1.0) {
  1276. free(js);
  1277. return 0;
  1278. }
  1280. if (is!=k) {
  1281. for (j=k;j<n;j++) {
  1282. t=A[k][j];
  1283. A[k][j]=A[is][j];
  1284. A[is][j]=t;
  1286. }
  1287. t=b[k];
  1288. b[k]=b[is];
  1289. b[is]=t;
  1290. }
  1292. if (js[k]!=k)
  1293. for (i=0;i<n;i++) {
  1294. t=A[i][k];
  1295. A[i][k]=A[i][js[k]];
  1296. A[i][js[k]]=t;
  1297. }
  1299. t=A[k][k];
  1300. for (j=k+1;j<n;j++)
  1301. if (A[k][j]!=0.0)
  1302. A[k][j]/=t;
  1304. b[k]=b[k]/t;
  1305. for (j=k+1;j<n;j++)
  1306. if (A[k][j]!=0.0)
  1307. for (i=0;i<n;i++)
  1308. if ((i!=k)&&(A[i][k]!=0.0))
  1309. A[i][j]=A[i][j]-A[i][k]*A[k][j];
  1311. for (i=0;i<n;i++)
  1312. if ((i!=k)&&(A[i][k]!=0.0))
  1313. b[i]=b[i]-A[i][k]*b[k];
  1315.   }
  1317.   for (k=n-1;k>=0;k--)
  1318. if (k!=js[k]) {
  1319. t=b[k];
  1320. b[k]=b[js[k]];
  1321. b[js[k]]=t;
  1322. }
  1324.   free(js);
  1326.   return 1;
  1327. }