Cilindro conectando dois pontos em OpenGL

Como parte do meu TCC, eu tive que quebrar um pouco a cabeça para conseguir dois pontos com um cilindro. Parece fácil, mas a única implementação que encontrei tinha um erro que demorei um pouco a encontrar.

Segue o minha classe C++ para isso:

#include<GL/glut.h>
#include<GL/gl.h>
#include<cmath>

 class Cylinder{

public:

GLdouble _x1;

GLdouble _y1;

GLdouble _z1;

GLdouble _x2;

GLdouble _y2;

GLdouble _z2;

Cylinder(double x1, double y1, double z1, double x2, double y2, double z2);

void render(GLfloat red, GLfloat green, GLfloat blue, GLdouble angle_x, GLdouble angle_y);

};

#define CYLINDER_RADIUS 0.75
#define CYLINDER_SUBDIVISIONS 36

Cylinder::Cylinder(double x1, double y1, double z1, double x2, double y2, double z2){

_x1 = x1;

_y1 = y1;

_z1 = z1;

_x2 = x2;

_y2 = y2;

_z2 = z2;

}

void Cylinder::render(GLfloat red, GLfloat green, GLfloat blue, GLdouble angle_x, GLdouble angle_y){

GLdouble vx, vy, vz, v, ax;

GLUquadricObj *_quadric;

_quadric = gluNewQuadric();

gluQuadricNormals(_quadric, GLU_SMOOTH);

vx = _x2 - _x1;

vy = _y2 - _y1;

vz = _z2 - _z1;

//handle the degenerate case of z1 == z2 with an approximation

if (vz > -0.000001 && vz < 0.000001){

if (vz >= 0.0){

vz = 0.000001;

}else{

vz = -0.000001;

}

}

if (vy > -0.000001 && vy < 0.000001){

if (vy >= 0.0){

vy = -0.000001;

}else{

vy = 0.000001;

}

}

if (vx > -0.000001 && vx < 0.000001){

if (vx >= 0.0){

vx = 0.000001;

}else{

vx = -0.000001;

}

}

v = sqrt( vx*vx + vy*vy + vz*vz );

ax = (180.0/M_PI)*acos( vz/v );

if( vz < 0.0 ) ax = -ax;

glPushMatrix();

glColor3f(red, green, blue);

glRotated(angle_x, 1.0, 0.0, 0.0);

glRotated(angle_y, 0.0, 1.0, 0.0);

//draw the cylinder body

glTranslated( _x1, _y1, _z1 );

glRotated(ax, -vy, vx, 0.0);

gluQuadricOrientation(_quadric,GLU_OUTSIDE);

gluCylinder(_quadric, CYLINDER_RADIUS, CYLINDER_RADIUS, v, CYLINDER_SUBDIVISIONS, 1);

//draw the first cap

gluQuadricOrientation(_quadric,GLU_INSIDE);

gluDisk(_quadric, 0.0, CYLINDER_RADIUS, CYLINDER_SUBDIVISIONS, 1);

glTranslated( 0,0,v );

//draw the second cap

gluQuadricOrientation(_quadric,GLU_OUTSIDE);

gluDisk(_quadric, 0.0, CYLINDER_RADIUS, CYLINDER_SUBDIVISIONS, 1);

glPopMatrix();

gluDeleteQuadric(_quadric);

}

A matemática por trás de todo este código é a seguinte:

Com os dois pontos que desejo conectar obtenho um vetor direção para o meu cilindro. Descubro o angulo entre a direção atual do cilindro e a que quero que ele sigo. Depois calculo o vetor ortogonal à direção original e à final para ter o eixo de rotação. Enfim rotaciono sobre o eixo que descobri pelo ângulo que descobri.

O código original foi escrito por Curran Kelleher e pode ser encontrado em http://lifeofaprogrammergeek.blogspot.com.br/2008/08/opengl-example-rendering-cylinders.html.
A descrição em pseudocódigo que me levou a entender a correção necessária foi escrita por Toby H. J. Smith e pode ser encontrada em http://www.thjsmith.com/40/cylinder-between-two-points-opengl-c