lorenz.c

学习c++必备

/*
 * Lorenz Attractor Demo
 *
 * Adapted from code originally written for the 4D60GT by
 * Aaron T. Ferrucci ([e-mail]), 7/3/92.
 *
 * Description:
 *
 * This program shows some particles stuck in a Lorenz attractor (the parameters
 * used are r=28, b=8/3, sigma=10). The eye is attracted to the red particle,
 * with a force directly proportionate to distance. A command line
 * puts the whole mess inside a box made of hexagons. I think this helps to
 * maintain the illusion of 3 dimensions, but it can slow things down.
 * Other options allow you to play with the redraw rate and the number of new
 * lines per redraw. So you can customize it to the speed of your machine.
 * 
 * For general info on Lorenz attractors I recommend "An Introduction to
 * the Lorenz Equations", IEEE Transactions on Circuits and Systems, August '83.
 *
 * Bugs: hidden surface removal doesn't apply to hexagons, and
 * works poorly on lines when they are too close together.
 *
 * Notes on OpenGL port:
 * 
 * The timer functions do not exist in OpenGL, so the drawing occurs in a
 * continuous loop, controlled by step, stop and go input from the keyboard.
 * Perhaps system function could be called to control timing.
 *
 */

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include <GL/glut.h>

float
random_float(void)
{
  return (float)rand()/RAND_MAX;
}

void
seed_random_float(long seed)
{
  srand(seed);
}

static GLuint asphere;

#define POINTMASK 511
#define G (0.002)	/* eyept to red sphere gravity */
#define LG (0.3)
#define CUBESIDE (120.)
#define CUBESCALE (23.)
#define CUBEOFFX (-4.)
#define CUBEOFFY (0.)
#define CUBEOFFZ (57.)
#define FALSE 0
#define TRUE 1

/* globals */
float sigma = 10., r = 28., b = 8./3., dt = 0.003;
int rp = 0, bp = 0, gp = 0, yp = 0, mp = 0;
long xmax, ymax, zmax, zmin;
float rv[POINTMASK+1][3],			/* red points */
	bv[POINTMASK+1][3],		/* blue points */
	gv[POINTMASK+1][3],		/* green points */
	yv[POINTMASK+1][3],		/* yellow points */
	mv[POINTMASK+1][3];		/* magenta points */

int lpf;				/* number of new lines per frame */

float eyex[3],	/* eye location */
	 eyev[3],	/* eye velocity */
	 eyel[3];	/* lookat point location */
GLint fovy = 600;
float dx, dy, dz;
GLUquadricObj *quadObj;

float cubeoffx = CUBEOFFX;
float cubeoffy = CUBEOFFY;
float cubeoffz = CUBEOFFZ;
float farplane = 80.;

int animate = 1;

/* option flags */
GLboolean hexflag,		/* hexagons? */
	sflag, 			
	fflag,			
	wflag,
	gflag,
	debug;

/* option values */
short hexbright;	/* brightness for hexagon color */
int speed,		/* speed (number of new line segs per redraw) */
    frame;		/* frame rate (actually noise value for TIMER0) */
float a = 0,
    da;			/* hexagon rotational velocity (.1 degree/redraw) */
float gravity;

/* function declarations */
void init_3d(void);
void init_graphics(void);
void draw_hexcube(void);
void draw_hexplane(void);
void draw_hexagon(void);
void move_eye(void);
void redraw(void);
void next_line(float v[][3], int *p);
void parse_args(int argc, char **argv);
void print_usage(char*);
void print_info(void);
void sphdraw(float args[4]);
void setPerspective(int angle, float aspect, float zNear, float zFar);


static void Reshape(int width, int height)
{

      glViewport(0,0,width,height);
      glClear(GL_COLOR_BUFFER_BIT);
      xmax = width;
      ymax = height;
}

/* ARGSUSED1 */
static void Key(unsigned char key, int x, int y)
{

    switch (key) {
      case 'g':
	animate = 1;
	glutPostRedisplay();
	break;
      case 's':
	 animate = 0;
	 glutPostRedisplay();
	 break;
      case 27:
	gluDeleteQuadric(quadObj);
	exit(0);
    }
}

static void Draw(void)
{
    int i;

    if (animate) {
	i = speed;
	while (i--) {
	    next_line(rv, &rp);
	    next_line(bv, &bp);
	    next_line(gv, &gp);
	    next_line(yv, &yp);
	    next_line(mv, &mp);
	}
	glPushMatrix();
	move_eye();
	redraw();
	glPopMatrix();
    }
}

void main(int argc, char **argv)
{
    glutInitWindowSize(600, 600);

    glutInit(&argc, argv);

    parse_args(argc, argv);

    glutInitDisplayMode(GLUT_RGB | GLUT_DOUBLE | GLUT_DEPTH);

    glutCreateWindow("Lorenz Attractors");

    init_3d();
    init_graphics();

    /* draw the first POINTMASK points in each color */
    while(rp < POINTMASK) {
	next_line(rv, &rp);
	next_line(bv, &bp);
	next_line(gv, &gp);
	next_line(yv, &yp);
	next_line(mv, &mp);
    }

    eyex[0] = eyex[1] = eyex[2] = 0.;
    eyel[0] = rv[rp][0];
    eyel[1] = rv[rp][1];
    eyel[2] = rv[rp][2];
	
    glPushMatrix();
    move_eye();
    redraw();
    glPopMatrix();

    glutReshapeFunc(Reshape);
    glutKeyboardFunc(Key);
    glutIdleFunc(Draw);
    glutDisplayFunc(Draw);
    glutMainLoop();
}

/* compute the next point on the path according to Lorenz' equations. */
void next_line(float v[][3], int *p)
{

    dx = sigma * (v[*p][1] - v[*p][0]) * dt;
    dy = (r*v[*p][0] - v[*p][1] + v[*p][0]*v[*p][2]) * dt;
    dz = (v[*p][0] *v[*p][1] + b*v[*p][2]) * dt;	
	
    v[(*p + 1) & POINTMASK][0] = v[*p][0] + dx;
    v[(*p + 1) & POINTMASK][1] = v[*p][1] + dy;
    v[(*p + 1) & POINTMASK][2] = v[*p][2] - dz;
    *p = (*p + 1) & POINTMASK;
}

void drawLines(int index, float array[POINTMASK][3])
{
    int p;
    int i;

#define LINE_STEP 4

    p = (index+1)&POINTMASK;
    i = LINE_STEP-(p % LINE_STEP);
    if (i == LINE_STEP) i=0;
    glBegin(GL_LINE_STRIP);
	/* draw points in order from oldest to newest */
	while(p != index) {
	    if (i == 0) {
		glVertex3fv(array[p]);
		i = LINE_STEP;
	    } 
	    i--;
	    p = (p+1) & POINTMASK;
	}
	glVertex3fv(array[index]);
    glEnd();
}

void redraw(void)
{
    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);

    if(hexflag)
	draw_hexcube();

    glColor3f(1.0, 0.0, 0.0);
    drawLines(rp, rv);
    sphdraw(rv[rp]);

    glColor3f(0.0, 0.0, 1.0);
    drawLines(bp, bv);
    sphdraw(bv[bp]);
	
    glColor3f(0.0, 1.0, 0.0);
    drawLines(gp, gv);
    sphdraw(gv[gp]);

    glColor3f(1.0, 0.0, 1.0);
    drawLines(yp, yv);
    sphdraw(yv[yp]);

    glColor3f(0.0, 1.0, 1.0);
    drawLines(mp, mv);
    sphdraw(mv[mp]);

    glutSwapBuffers();
}

void move_eye(void)
{
    /* first move the eye */
    eyev[0] += gravity * (rv[rp][0] - eyex[0]);
    eyev[1] += gravity * (rv[rp][1] - eyex[1]);
    eyev[2] += gravity * (rv[rp][2] - eyex[2]);

    /* adjust position using new velocity */
    eyex[0] += eyev[0] * dt;
    eyex[1] += eyev[1] * dt;
    eyex[2] += eyev[2] * dt;

    /* move the lookat point */
    /* it catches up to the red point if it's moving slowly enough */
    eyel[0] += LG * (rv[rp][0] - eyel[0]);
    eyel[1] += LG * (rv[rp][1] - eyel[1]);
    eyel[2] += LG * (rv[rp][2] - eyel[2]);

    /* change view */
    gluLookAt(eyex[0], eyex[1], eyex[2], eyel[0], eyel[1], eyel[2],
	      0, 1, 0);
}

void draw_hexcube(void)
{

    a += da;
    if(a >= 720.)		/* depends on slowest rotation factor */
	a = 0.;

    /* draw hexplanes, without changing z-values */
    glDepthMask(GL_FALSE); 
    glDisable(GL_DEPTH_TEST);

    /* x-y plane */
    glColor3f(0.2, 0.2, 0.6);
    glPushMatrix();
    glTranslatef(cubeoffx, cubeoffy, cubeoffz);
    glScalef(CUBESCALE, CUBESCALE, CUBESCALE);
    draw_hexplane();
    glPopMatrix();

    /* x-y plane, translated */
    glPushMatrix();
    glTranslatef(cubeoffx, cubeoffy, cubeoffz - 2*CUBESIDE);