/************************************************************************
* *
* julia_cgi.c - By Mark McClure - mcmcclur@bulldog.unca.edu *
* Originally written Feburary, 1996 *
* Updated April, 1996. *
* *
************************************************************************/
/************************************************************************
* *
* This program is meant to be executed as a cgi. It generates a *
* Julia Set (or sets) for functions of the form f(z) = z^2 + c and *
* returns the result as a gif image using the gd gif image *
* generating library available from *
* *
* The exact behavior is dependent upon the brower and the query *
* string. The Julia Set is generated using the inverse *
* iteration algorithm. *
* *
* The URL: *
* *
* *
* The cgi may be called three ways. When called with a query string *
* of the form c=a+bi or c=a-bi where a and are doubles, a single gif *
* image of the corresponding Julia Set is returned. Examples: *
* *
* *
* *
* *
* If the cgi is called with the query string is c=random, then a *
* single randomly generated image is returned. *
* *
* If the cgi is called without a query string, then the behavior *
* depends upon the browser. If the HTTP_USER_AGENT environment *
* variable starts with "Mozilla" (ie: a Netscape compatible *
* browser), then a randomly generated sequence of gifs is returned *
* using Netscape's server push scheme. This results in an *
* animation. The sequence is chosen by letting the parameter c *
* go for a constrained random walk in the complex plane. *
* *
* If some other browser calls the cgi without a query string, then *
* a single randomly chosen Julia Set is returned. *
* *
* To compile with gcc (assuming libgd.a has been built): *
* gcc julia_cgi.c -o julia.cgi -L./ -lgd -lm *
* *
************************************************************************/
/* some standard includes */
#include
#include
#include
#include
/* header file for the gd gif library */
#include "gd.h"
/* HALF_RAND_MAX isn't defined in stdlib.h */
#define HALF_RAND_MAX 16384
/* Pi */
#define PI 3.1415926
/* number of frames for server push animation */
#define NUM_FRAMES 200
/* step size for c's random walk */
#define STEP_SIZE 0.05
/* size in pixels of the square Julia Set image */
#define IMG_SIZE 200
/* size of the Julia Set in the complex plane */
#define COMP_SIZE 4.0
/* An initial point */
#define ZINITIAL_RE 0.78
#define ZINITIAL_IM 0.52
/* How many times to plot a pixel before */
/* giving up on that region */
#define BAILOUT 3
#define BRUSHSIZE 10
/*----------Structures------------*/
/*---Complex Numbers---*/
struct complex {
double re;
double im;
};
typedef struct complex complex;
/*---The Julia set is a tree with these nodes---*/
struct julia_tree_node {
complex z;
int x, y;
struct julia_tree_node *left, *right;
};
typedef struct julia_tree_node julia_tree_node;
/*---A queue to facilitate level traversal---*/
/*---of the julia_tree---*/
struct queue_node {
julia_tree_node *val;
struct queue_node *next;
};
typedef struct queue_node queue_node;
/*-----------end Structures-----------*/
/*-----------Function Prototypes------*/
/* Gif construction functions */
/* General */
void makeGif( void );
void mandelInit( void );
void singleGif(void);
void plot(int, int);
/* Server Push stuff */
void server_pushGifs(void);
double xDisplace(complex);
double yDisplace(complex);
/* Functions related to data structures */
/* The Julia set is stored in a binary tree */
void construct_julia_tree( void );
void inv(julia_tree_node *);
/* The queue */
void queueinit( void );
void put(julia_tree_node *);
julia_tree_node *get( void );
int queueempty( void );
/* Complex number functions */
/* general complex */
complex cMinus(complex);
complex cSubtract(complex, complex);
complex cSqrt(complex);
complex random_complex(void);
/* Julia set specific */
complex finv_1(complex);
complex finv_2(complex);
/*---------end Function Prototypes------*/
/*----------Global Variables------------*/
complex c; /*The Julia constant*/
julia_tree_node *top; /*The top of the Julia Tree*/
queue_node *head, *tail; /*Pointers to the ends of*/
/*the queue*/
/*This matrix keeps track of which points have been plotted*/
int plot_record[IMG_SIZE + 1][IMG_SIZE + 1];
/*gd Variables*/
gdImagePtr im; /*The gdImage*/
gdImagePtr mandel_im; /*The Mandelbrot Set*/
gdImagePtr brush; /*The brush*/
int black, white, blue; /*The colors*/
FILE *mandelIN;
/*-------------end Global Variables-----------*/
/*-----------------The Functions----------------*/
/* main() simply checks the environment variables and */
/* calls the appropriate functions */
int main( void ) {
char query_string[50];
if(getenv("QUERY_STRING") != NULL) {
strcpy(query_string, getenv("QUERY_STRING"));
}
else {
printf("This program is meant to be executed as a cgi.\n");
exit(1);
}
mandelInit();
if(strcmp(query_string, "c=random") == 0) {
c = random_complex();
singleGif();
}
else if(sscanf(query_string, "c=%lf+%lfi", &c.re, &c.im) == 2) {
singleGif();
}
else if(strlen(query_string) > 0) {
exit(1);
}
else if(strncmp(getenv("HTTP_USER_AGENT"), "Mozilla", 7) == 0) {
server_pushGifs();
}
else {
c = random_complex();
singleGif();
}
gdImageDestroy(mandel_im);
gdImageDestroy(brush);
exit(0);
}
/* end main() */
void mandelInit() {
mandelIN = fopen("mandel.gif", "rb");
mandel_im = gdImageCreateFromGif(mandelIN);
fclose(mandelIN);
/* Create the brush */
brush = gdImageCreate(BRUSHSIZE,BRUSHSIZE);
white = gdImageColorAllocate(brush,255,255,255);
gdImageColorTransparent(brush, white);
blue = gdImageColorAllocate(brush,0,0,255);
gdImageArc(brush,BRUSHSIZE/2,BRUSHSIZE/2,BRUSHSIZE/2,BRUSHSIZE/2,0,360,blue);
gdImageFillToBorder(brush,BRUSHSIZE/2,BRUSHSIZE/2,blue,blue);
}
/*------------------Gif construction functions--------------*/
/* singleGif let's the browser know that a gif is coming and */
/* calls makeGif() */
/* called by main() */
void singleGif( void ) {
printf("Content-type: image/gif\n\n");
makeGif();
}
/* end single_gif() */
/* server_pushGifs() sets coordinates server pushing */
/* called by main() */
void server_pushGifs() {
int j;
complex c_old;
c = random_complex();
printf("Content-type: multipart/x-mixed-replace;boundary=ThisRandomString\n");
for(j=0; j= 0.5 )
return -STEP_SIZE;
if( z.re <= -2 )
return STEP_SIZE;
if( (z.re <= -1) && (fabs(z.im) >= .3) )
return STEP_SIZE;
if( rand()/HALF_RAND_MAX == 0 )
return STEP_SIZE;
else
return -STEP_SIZE;
}
double yDisplace(complex z) {
if( z.im >= 1.2 )
return -STEP_SIZE;
if( z.im <= -1.2 )
return STEP_SIZE;
if( (z.im >= 0.3) && (z.re <= -1) )
return -STEP_SIZE;
if( (z.im <= 0.3) && (z.re <= -1) )
return STEP_SIZE;
if( rand()/HALF_RAND_MAX == 0 )
return STEP_SIZE;
else
return -STEP_SIZE;
}
/* end Displacements */
/* This is the main point of entry to generate a Julia Set and */
/* make the gif */
/* called by singleGif() and server_pushGifs() */
void makeGif( void ) {
int cx,cy;
/* Bounds for mandel.gif */
float xmin = -2.0;
float xmax = 0.6;
float ymin = -1.3;
float ymax = 1.3;
cx = ( IMG_SIZE*(c.re-xmin)/(xmax-xmin) );
cy = ( IMG_SIZE*(c.im-ymax)/(ymin-ymax) );
/* Set up the gdImage */
im = gdImageCreate(2*IMG_SIZE, IMG_SIZE);
white = gdImageColorAllocate(im, 255, 255, 255);
black = gdImageColorAllocate(im, 0, 0, 0);
gdImageColorTransparent(im, white);
/* Draw the Mandelbrot Set */
gdImageCopy(im, mandel_im,0,0,0,0,IMG_SIZE,IMG_SIZE);
gdImageCopy(im, brush,cx - BRUSHSIZE/2, cy - BRUSHSIZE/2,0,0,
BRUSHSIZE,BRUSHSIZE);
/* Construct the Julia Set */
/* The points will be plotted in subsequent functions */
construct_julia_tree();
/* Write the gif to stdout and clean up */
gdImageGif(im, stdout);
gdImageDestroy(im);
}
/* end makeGif() */
/* plot() plots the points in the gif */
/* called by construct_julia_tree() */
void plot(int x, int y) {
gdImageLine(im, IMG_SIZE + x, y, IMG_SIZE + x, y, black);
}
/* end plot() */
/*-----------end Gif construction functions--------------*/
/*--------Julia Set Functions---------*/
/* construct_julia_tree() is the main function responsible for */
/* setting up and constructing the binary tree which holds the */
/* the Julia Set */
/* called by makeGif() */
void construct_julia_tree() {
julia_tree_node *tree_item;
int x, y; /* points to plot */
complex z0; /* The starting point */
int i, j; /* Dummies */
/* various initializations */
for(i=0; i <= IMG_SIZE; i++) {
for(j=0; j <= IMG_SIZE; j++) {
plot_record[i][j] = 0;
}
}
z0.re = ZINITIAL_RE;
z0.im = ZINITIAL_IM;
for(i=0; i<10; i++){
z0 = finv_1(z0);
z0 = finv_2(z0);
}
queueinit();
top = (julia_tree_node *) malloc(sizeof( *top ));
top->z = z0;
put(top);
/* level construction of the tree using the queue to */
/* keep track of which nodes have been constructed and */
/* plot record keep track of which points have been plotted */
while(!queueempty()) {
tree_item = get();
inv(tree_item);
x = (tree_item->left)->x;
y = (tree_item->left)->y;
if(plot_record[x][y] < BAILOUT) {
plot(x,y);
plot_record[x][y] = plot_record[x][y] + 1;
put(tree_item->left);
}
x = (tree_item->right)->x;
y = (tree_item->right)->y;
if(plot_record[x][y] < BAILOUT) {
plot(x,y);
plot_record[x][y] = plot_record[x][y] + 1;
put(tree_item->right);
}
}
}
/* end construct_julia_tree() */
/* inv() constructs the two children of tree_item which contain */
/* the two points in the inverse image of the point in tree_item */
/* called by construct_julia_tree() */
void inv(julia_tree_node *tree_item) {
julia_tree_node *left_child, *right_child;
left_child = (julia_tree_node *) malloc(sizeof(*left_child));
right_child = (julia_tree_node *) malloc(sizeof(*right_child));
tree_item->left = left_child;
tree_item->right = right_child;
left_child->z = finv_1(tree_item->z);
right_child->z = finv_2(tree_item->z);
left_child->x = (int)
( (IMG_SIZE/COMP_SIZE)*((left_child->z).re + COMP_SIZE/2) );
left_child->y = (int)
( (IMG_SIZE/COMP_SIZE)*(COMP_SIZE/2 - (left_child->z).im) );
right_child->x = (int)
( (IMG_SIZE/COMP_SIZE)*((right_child->z).re + COMP_SIZE/2) );
right_child->y = (int)
( (IMG_SIZE/COMP_SIZE)*(COMP_SIZE/2 - (right_child->z).im) );
}
/* end inv() */
/*---------end Julia Set Functions------*/
/*---------queue Functions-------------*/
/* These are standard manipulation functions for a linked list */
/* representation for a queue */
/* all are called by construct_julia_tree() */
/* queueinit() sets up an empty queue */
void queueinit( void ) {
head = (queue_node *) malloc(sizeof *head);
tail = (queue_node *) malloc(sizeof *tail);
head->next = tail;
tail->next=head;
}
/* end queueinit() */
/* put() puts a julia_tree_node * at the end of the queue */
void put(julia_tree_node *v) {
queue_node *queue_item;
queue_item = (queue_node *) malloc(sizeof *queue_item);
queue_item->val = v;
queue_item->next = tail;
(tail->next)->next = queue_item;
tail->next = queue_item;
}
/* end put() */
/* get() gets a julia_tree_node * from the front of the queue */
julia_tree_node *get() {
queue_node *queue_item;
julia_tree_node *x;
queue_item = head->next;
head->next = queue_item->next;
if(queueempty()) {
tail->next = head;
}
x = queue_item->val;
free(queue_item);
return x;
}
/* end get() */
/* queueempty() tests to see if anything is in the queue */
int queueempty( void ) {
return head->next == tail;
}
/* end queueempty() */
/*--------end queue Functions---------*/
/*-------Complex Number Functions------*/
/* random_complex() returns a random complex number from */
/* a region which gives reasonable Julia Sets */
/* called by main() and server_pushGifs() */
complex random_complex() {
complex z;
srand(time(NULL));
while(1) {
z.re = (double) ( 2.5*rand()/RAND_MAX - 2 );
z.im = (double) ( 2.4*rand()/RAND_MAX - 1.2 );
if( (z.re > -1) || (fabs(z.im) < .3) )
break;
}
return z;
}
/* end random_complex */
/* finv_1() and finv_2() are the two inverse images of a */
/* complex number */
/* called by inv() */
complex finv_1(complex z) {
return cSqrt(cSubtract(z,c));
}
complex finv_2(complex z) {
return cMinus(cSqrt(cSubtract(z,c)));
}
/* end inverses */
/* Standard complex number operations */
/* called by finv_1() and finv_2() */
complex cMinus(complex z) {
complex z_out;
z_out.re = -z.re;
z_out.im = -z.im;
return z_out;
}
complex cSubtract(complex z1, complex z2) {
complex z_out;
z_out.re = z1.re - z2.re;
z_out.im = z1.im - z2.im;
return z_out;
}
complex cSqrt(complex z) {
double theta;
double r = sqrt(sqrt(z.re*z.re + z.im*z.im));
complex z_out;
if(z.re >= 0) {
theta = atan(z.im/z.re)/2;
}
else {
if(z.im >= 0) {
theta = (atan(z.im/z.re) + PI)/2;
}
else {
theta = (atan(z.im/z.re) - PI)/2;
}
}
z_out.re = r*(cos(theta));
z_out.im = r*(sin(theta));
return z_out;
}
/*--------end Complex Number Functions--------*/
/*------------end File---------------------*/