Follow the paper of Sweeping stage.

README.md

@@ -36,7 +36,7 @@ Two caps and many strips composites a cylinder.
 Almost all 3D editor have a infinite grid ground, I just made a finite one, in the future, I should expand the grid outside of the screen to make it infinite.  
 Now, for just beginning, I think it's a not bad start.
 <img src="screenshot/dust3d_sphere_cylinder.png">
-- [ ] Implement B-Mesh algorithm (Dec 18, 2016 ~ Dec 20, 2016)   
+- [ ] Implement B-Mesh algorithm (Dec 18, 2016 ~ Dec 23, 2016)   
 *Drawing Skeletal Shape Balls*  
 Draw shape ball is easy, no need to rotate, I just need scale it along the ball's radius.
 Draw the cylinder which connects two shape balls is more difficult, I need do some math to rotate it. [Here](http://www.thjsmith.com/40/cylinder-between-two-points-opengl-c) described it.  

src/bmesh.c

@@ -155,38 +155,21 @@ static int bmeshAddInbetweenNodeBetween(bmesh *bm,
   return newNode.index;
 }
 
-static void floatsToQuad(float *floats, quad *q) {
-  int i;
-  int offset = 0;
-  for (i = 0; i < 4; ++i) {
-    q->pt[i].x = floats[offset++];
-    q->pt[i].y = floats[offset++];
-    q->pt[i].z = floats[offset++];
-  }
-}
-
 static int bmeshGenerateNodeQuad(bmesh *bm, bmeshNode *node,
-    matrix *matRotate, int connectWithQuad) {
-  quad q;
-  matrix matTranslate;
-  matrix matFinal;
+    vec3 *localYaxis, vec3 *localZaxis, int connectWithQuad) {
   int i;
-  float floats[4][3] = {
-    {-node->radius, +node->radius, 0},
-    {-node->radius, -node->radius, 0},
-    {+node->radius, -node->radius, 0},
-    {+node->radius, +node->radius, 0},
-  };
-  matrixTranslate(&matTranslate, node->position.x, node->position.y,
-    node->position.z);
-  matrixLoadIdentity(&matFinal);
-  matrixAppend(&matFinal, &matTranslate);
-  matrixAppend(&matFinal, matRotate);
-  matrixTransformVector(&matFinal, floats[0]);
-  matrixTransformVector(&matFinal, floats[1]);
-  matrixTransformVector(&matFinal, floats[2]);
-  matrixTransformVector(&matFinal, floats[3]);
-  floatsToQuad(&floats[0][0], &q);
+  quad q;
+  vec3 z, y;
+  vec3Scale(localYaxis, node->radius, &y);
+  vec3Scale(localZaxis, node->radius, &z);
+  vec3Sub(&node->position, &y, &q.pt[0]);
+  vec3Add(&q.pt[0], &z, &q.pt[0]);
+  vec3Sub(&node->position, &y, &q.pt[1]);
+  vec3Sub(&q.pt[1], &z, &q.pt[1]);
+  vec3Add(&node->position, &y, &q.pt[2]);
+  vec3Sub(&q.pt[2], &z, &q.pt[2]);
+  vec3Add(&node->position, &y, &q.pt[3]);
+  vec3Add(&q.pt[3], &z, &q.pt[3]);
   if (-1 == bmeshAddQuad(bm, &q)) {
     fprintf(stderr, "%s:meshAddQuad failed.\n", __FUNCTION__);
     return -1;
@@ -221,11 +204,10 @@ static int bmeshGenerateInbetweenNodesBetween(bmesh *bm,
   float step = 0.5;
   float distance;
   int parentNodeIndex = firstNodeIndex;
-  float rotateAngle = 0;
-  vec3 rotateAround = {0, 0, 0};
-  vec3 p;
-  vec3 zAxis = {0, 0, 1};
-  matrix matRotate;
+  vec3 localZaxis;
+  vec3 localYaxis;
+  vec3 edgeDirection;
+  vec3 worldZaxis = {0, 0, 1};
   int lastQuadIndex = -1;
   
   bmeshNode *firstNode = bmeshGetNode(bm, firstNodeIndex);
@@ -234,16 +216,15 @@ static int bmeshGenerateInbetweenNodesBetween(bmesh *bm,
   if (secondNode->roundColor == bm->roundColor) {
     return 0;
   }
-  vec3Sub(&firstNode->position, &secondNode->position, &p);
-  vec3CrossProduct(&zAxis, &p, &rotateAround);
-  vec3Normalize(&rotateAround);
+  vec3Sub(&firstNode->position, &secondNode->position, &edgeDirection);
+  vec3CrossProduct(&worldZaxis, &edgeDirection, &localZaxis);
+  vec3Normalize(&localZaxis);
+  vec3CrossProduct(&localZaxis, &edgeDirection, &localYaxis);
+  vec3Normalize(&localYaxis);
   
-  distance = vec3Length(&p);
+  distance = vec3Length(&edgeDirection);
   if (distance > 0) {
     float offset = step;
-    rotateAngle = 180 / M_PI * acos(vec3DotProduct(&zAxis, &p) / distance);
-    matrixRotate(&matRotate,
-      rotateAngle, rotateAround.x, rotateAround.y, rotateAround.z);
     if (offset + step <= distance) {
       while (offset + step <= distance) {
         float frac = offset / distance;
@@ -253,7 +234,7 @@ static int bmeshGenerateInbetweenNodesBetween(bmesh *bm,
           return -1;
         }
         newNode = bmeshGetNode(bm, parentNodeIndex);
-        bmeshGenerateNodeQuad(bm, newNode, &matRotate,
+        bmeshGenerateNodeQuad(bm, newNode, &localYaxis, &localZaxis,
           lastQuadIndex);
         lastQuadIndex = -1 == lastQuadIndex ? bmeshGetQuadNum(bm) - 1 :
           bmeshGetQuadNum(bm) - 5;
@@ -266,7 +247,7 @@ static int bmeshGenerateInbetweenNodesBetween(bmesh *bm,
         return -1;
       }
       newNode = bmeshGetNode(bm, parentNodeIndex);
-      bmeshGenerateNodeQuad(bm, newNode, &matRotate, lastQuadIndex);
+      bmeshGenerateNodeQuad(bm, newNode, &localYaxis, &localZaxis, lastQuadIndex);
       lastQuadIndex = -1 == lastQuadIndex ? bmeshGetQuadNum(bm) - 1 :
           bmeshGetQuadNum(bm) - 5;
     }

src/vector3d.c

@@ -18,6 +18,12 @@ void vec3Normalize(vec3 *p) {
   }
 }
 
+void vec3Scale(vec3 *a, float scale, vec3 *result) {
+  result->x = a->x * scale;
+  result->y = a->y * scale;
+  result->z = a->z * scale;
+}
+
 void vec3Midpoint(vec3 *a, vec3 *b, vec3 *mid) {
   vec3Lerp(a, b, 0.5, mid);
 }
@@ -34,6 +40,12 @@ void vec3Sub(vec3 *a, vec3 *b, vec3 *result) {
   result->z = a->z - b->z;
 }
 
+void vec3Add(vec3 *a, vec3 *b, vec3 *result) {
+  result->x = a->x + b->x;
+  result->y = a->y + b->y;
+  result->z = a->z + b->z;
+}
+
 void vec3CrossProduct(vec3 *a, vec3 *b, vec3 *result) {
   result->x = a->y * b->z - a->z * b->y;
   result->y = a->z * b->x - a->x * b->z;

src/vector3d.h

@@ -16,6 +16,8 @@ void vec3Midpoint(vec3 *a, vec3 *b, vec3 *mid);
 void vec3Lerp(vec3 *a, vec3 *b, float frac, vec3 *result);
 void vec3CrossProduct(vec3 *a, vec3 *b, vec3 *result);
 void vec3Sub(vec3 *a, vec3 *b, vec3 *result);
+void vec3Add(vec3 *a, vec3 *b, vec3 *result);
+void vec3Scale(vec3 *a, float scale, vec3 *result);
 float vec3DotProduct(vec3 *a, vec3 *b);
 float vec3Length(vec3 *p);
 float vec3Distance(vec3 *a, vec3 *b);