Применение теста глубины OpenGL к орбитальным сферам

Я пишу гравитационное моделирование n-тела на C ++, которое анимируется с использованием OpenGL и GLUT (это хобби-проект). По большей части анимация работает нормально, однако у меня возникли две основные проблемы, которые я не могу решить:

• Хотя глубинное тестирование включено, оно не работает должным образом, и
• Поверхности … ну, грязные.

Мой вопрос, как я могу исправить эти проблемы?

Обе проблемы можно увидеть на следующих изображениях (извинения за ссылки, но у меня недостаточно репутации для публикации изображений). Это снимки с симуляции простой орбиты, видимой на грани.

Вот желтая сфера рисуется перед фиолетовой, как и должно быть.

После половины орбиты желтая сфера все еще нарисована впереди, хотя она еще дальше.

Код, использованный для создания анимации, приведен ниже.

#include <GL/glut.h>
#include "Cluster.h" // My own class.

// Scale for animation.  Each unit in the animation = 1/SCALE m.
const double SCALE = 1e-10;
// Size of spheres for animation.
const double SPHERE_SIZE = 2e10*SCALE;

// Cluster object contains bodies and updates their positions.
Cluster cluster();

// Array of rgb colors for spheres.
GLfloat colorArr[4][4] =
{
{0.7, 0.7, 0.0, 1.0},
{0.73, 0.24, 0.95, 1.0},
}

void display()
{
glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);
glMatrixMode(GL_MODELVIEW);
glLoadIdentity();
gluLookAt( 0.0*SCALE, 5e11*SCALE, 0.0*SCALE  // eye is on y-axis outside orbit.
, 0.0, 0.0, 0.0
, 0.0, 0.0, 1.0 );

for (int i=0; i<N; i++) // N is the number of bodies in cluster.
{
glPushMatrix();
glTranslated( SCALE*cluster.getX(i)     // Get coordinate of ith body.
, SCALE*cluster.getY(i)
, SCALE*cluster.getZ(i) );
glMaterialfv(GL_FRONT, GL_AMBIENT_AND_DIFFUSE, colorArr[i]);
glutSolidSphere(SPHERE_SIZE, 50, 50);
glCullFace(GL_BACK);
glPopMatrix();
}

glutSwapBuffers();
}

void reshape(GLint w, GLint h)
{
glViewport(0, 0, w, h);
glMatrixMode(GL_PROJECTION);
glLoadIdentity();
gluPerspective(30.0, (GLfloat)w/(GLfloat)h, SCALE, 5e11*SCALE);
glMatrixMode(GL_MODELVIEW);
}

void animate()
{
// Update positions and redraw.
cluster.update();
display();
}

void init()
{
GLfloat black[] = {0.0, 0.0, 0.0, 1.0};
GLfloat white[] = {1.0, 1.0, 1.0, 0.5};
GLfloat direction[] = {1.0, 1.0, 1.0, 0.0};

glMaterialfv(GL_FRONT, GL_SPECULAR, white);
glMaterialf(GL_FRONT, GL_SHININESS, 10);

glLightfv(GL_LIGHT0, GL_AMBIENT, black);
glLightfv(GL_LIGHT0, GL_DIFFUSE, white);
glLightfv(GL_LIGHT0, GL_SPECULAR, white);
glLightfv(GL_LIGHT0, GL_POSITION, direction);

glEnable(GL_LIGHTING);
glEnable(GL_LIGHT0);
glEnable(GL_DEPTH_TEST);
}

int main(int argc, char** argv)
{
glutInit(&argc, argv);
glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGB | GLUT_DEPTH);
glutInitWindowSize(800, 600);
glutCreateWindow("Test Orbit");
glutDisplayFunc(display);
glutReshapeFunc(reshape);
glutIdleFunc(animate);
init();
glutMainLoop();
}