1488 lines
36 KiB
C
1488 lines
36 KiB
C
/*
|
|
flame - cosmic recursive fractal flames
|
|
Copyright (C) 1992 Scott Draves <spot@cs.cmu.edu>
|
|
|
|
This program is free software: you can redistribute it and/or modify
|
|
it under the terms of the GNU General Public License as published by
|
|
the Free Software Foundation; either version 3 of the License, or
|
|
(at your option) any later version.
|
|
|
|
This program is distributed in the hope that it will be useful,
|
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
GNU General Public License for more details.
|
|
|
|
You should have received a copy of the GNU General Public License
|
|
along with this program. If not, see <https://www.gnu.org/licenses/>.
|
|
*/
|
|
|
|
#include "config.h"
|
|
|
|
#include <stdlib.h>
|
|
#include <string.h> /* strcmp */
|
|
|
|
#include "libpika/pika.h"
|
|
|
|
#include "libifs.h"
|
|
|
|
#define CHOOSE_XFORM_GRAIN 100
|
|
|
|
static int flam3_random_bit (void);
|
|
static double flam3_random01 (void);
|
|
|
|
/*
|
|
* run the function system described by CP forward N generations.
|
|
* store the n resulting 3 vectors in POINTS. the initial point is passed
|
|
* in POINTS[0]. ignore the first FUSE iterations.
|
|
*/
|
|
|
|
void
|
|
iterate (control_point *cp,
|
|
int n,
|
|
int fuse,
|
|
point *points)
|
|
{
|
|
int i, j, count_large = 0, count_nan = 0;
|
|
int xform_distrib[CHOOSE_XFORM_GRAIN];
|
|
double p[3], t, r, dr;
|
|
p[0] = points[0][0];
|
|
p[1] = points[0][1];
|
|
p[2] = points[0][2];
|
|
|
|
/*
|
|
* first, set up xform, which is an array that converts a uniform random
|
|
* variable into one with the distribution dictated by the density
|
|
* fields
|
|
*/
|
|
dr = 0.0;
|
|
for (i = 0; i < NXFORMS; i++)
|
|
dr += cp->xform[i].density;
|
|
dr = dr / CHOOSE_XFORM_GRAIN;
|
|
|
|
j = 0;
|
|
t = cp->xform[0].density;
|
|
r = 0.0;
|
|
for (i = 0; i < CHOOSE_XFORM_GRAIN; i++)
|
|
{
|
|
while (r >= t)
|
|
{
|
|
j++;
|
|
t += cp->xform[j].density;
|
|
}
|
|
xform_distrib[i] = j;
|
|
r += dr;
|
|
}
|
|
|
|
for (i = -fuse; i < n; i++)
|
|
{
|
|
/* FIXME: the following is supported only by gcc and c99 */
|
|
int fn = xform_distrib[g_random_int_range (0, CHOOSE_XFORM_GRAIN)];
|
|
double tx, ty, v;
|
|
|
|
if (p[0] > 100.0 || p[0] < -100.0 ||
|
|
p[1] > 100.0 || p[1] < -100.0)
|
|
count_large++;
|
|
if (p[0] != p[0])
|
|
count_nan++;
|
|
|
|
#define coef cp->xform[fn].c
|
|
#define vari cp->xform[fn].var
|
|
|
|
/* first compute the color coord */
|
|
p[2] = (p[2] + cp->xform[fn].color) / 2.0;
|
|
|
|
/* then apply the affine part of the function */
|
|
tx = coef[0][0] * p[0] + coef[1][0] * p[1] + coef[2][0];
|
|
ty = coef[0][1] * p[0] + coef[1][1] * p[1] + coef[2][1];
|
|
|
|
p[0] = p[1] = 0.0;
|
|
/* then add in proportional amounts of each of the variations */
|
|
v = vari[0];
|
|
if (v > 0.0)
|
|
{
|
|
/* linear */
|
|
double nx, ny;
|
|
nx = tx;
|
|
ny = ty;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[1];
|
|
if (v > 0.0)
|
|
{
|
|
/* sinusoidal */
|
|
double nx, ny;
|
|
nx = sin (tx);
|
|
ny = sin (ty);
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[2];
|
|
if (v > 0.0)
|
|
{
|
|
/* spherical */
|
|
double nx, ny;
|
|
double r2 = tx * tx + ty * ty + 1e-6;
|
|
nx = tx / r2;
|
|
ny = ty / r2;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[3];
|
|
if (v > 0.0)
|
|
{
|
|
/* swirl */
|
|
double r2 = tx * tx + ty * ty; /* /k here is fun */
|
|
double c1 = sin (r2);
|
|
double c2 = cos (r2);
|
|
double nx = c1 * tx - c2 * ty;
|
|
double ny = c2 * tx + c1 * ty;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[4];
|
|
if (v > 0.0)
|
|
{
|
|
/* horseshoe */
|
|
double a, c1, c2, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
a = atan2(tx, ty); /* times k here is fun */
|
|
else
|
|
a = 0.0;
|
|
c1 = sin (a);
|
|
c2 = cos (a);
|
|
nx = c1 * tx - c2 * ty;
|
|
ny = c2 * tx + c1 * ty;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[5];
|
|
if (v > 0.0)
|
|
{
|
|
/* polar */
|
|
double nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
nx = atan2 (tx, ty) / G_PI;
|
|
else
|
|
nx = 0.0;
|
|
|
|
ny = sqrt (tx * tx + ty * ty) - 1.0;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[6];
|
|
if (v > 0.0)
|
|
{
|
|
/* bent */
|
|
double nx, ny;
|
|
nx = tx;
|
|
ny = ty;
|
|
if (nx < 0.0) nx = nx * 2.0;
|
|
if (ny < 0.0) ny = ny / 2.0;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[7];
|
|
if (v > 0.0)
|
|
{
|
|
/* folded handkerchief */
|
|
double theta, r2, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2( tx, ty );
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt (tx * tx + ty * ty);
|
|
nx = sin (theta + r2) * r2;
|
|
ny = cos (theta - r2) * r2;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[8];
|
|
if (v > 0.0)
|
|
{
|
|
/* heart */
|
|
double theta, r2, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2( tx, ty );
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt (tx * tx + ty * ty);
|
|
theta *= r2;
|
|
nx = sin (theta) * r2;
|
|
ny = cos (theta) * -r2;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[9];
|
|
if (v > 0.0)
|
|
{
|
|
/* disc */
|
|
double theta, r2, nx, ny;
|
|
if ( tx < -EPS || tx > EPS ||
|
|
ty < - EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
nx = tx * G_PI;
|
|
ny = ty * G_PI;
|
|
r2 = sqrt (nx * nx * ny * ny);
|
|
p[0] += v * sin(r2) * theta / G_PI;
|
|
p[1] += v * cos(r2) * theta / G_PI;
|
|
}
|
|
|
|
v = vari[10];
|
|
if (v > 0.0)
|
|
{
|
|
/* spiral */
|
|
double theta, r2;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2( tx, ty );
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt (tx * tx + ty * ty) + 1e-6;
|
|
p[0] += v * (cos (theta) + sin (r2)) / r2;
|
|
p[1] += v * (cos (theta) + cos (r2)) / r2;
|
|
}
|
|
|
|
v = vari[11];
|
|
if (v > 0.0)
|
|
{
|
|
/* hyperbolic */
|
|
double theta, r2;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt (tx * tx + ty * ty) + 1e-6;
|
|
p[0] += v * sin (theta) / r2;
|
|
p[1] += v * cos (theta) * r2;
|
|
}
|
|
|
|
v = vari[12];
|
|
if (v > 0.0 )
|
|
{
|
|
double theta, r2;
|
|
/* diamond */
|
|
if ( tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt( tx * tx + ty * ty );
|
|
p[0] += v * sin (theta) * cos (r2);
|
|
p[1] += v * cos (theta) * sin (r2);
|
|
}
|
|
|
|
v = vari[13];
|
|
if (v > 0.0)
|
|
{
|
|
/* ex */
|
|
double theta, r2, n0, n1, m0, m1;
|
|
if ( tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt( tx * tx + ty * ty );
|
|
n0 = sin(theta + r2);
|
|
n1 = cos(theta - r2);
|
|
m0 = n0 * n0 * n0 * r2;
|
|
m1 = n1 * n1 * n1 * r2;
|
|
p[0] += v * (m0 + m1);
|
|
p[1] += v * (m0 - m1);
|
|
}
|
|
|
|
v = vari[14];
|
|
if ( v > 0.0)
|
|
{
|
|
double theta, r2, nx, ny;
|
|
/* julia */
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
if (flam3_random_bit ())
|
|
theta += G_PI;
|
|
r2 = pow (tx * tx + ty * ty, 0.25);
|
|
nx = r2 * cos (theta);
|
|
ny = r2 * sin (theta);
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[15];
|
|
if (v > 0.0)
|
|
{
|
|
/* waves */
|
|
double dx, dy, nx, ny;
|
|
dx = coef[2][0];
|
|
dy = coef[2][1];
|
|
nx = tx + coef[1][0] * sin (ty / ((dx * dx) + EPS));
|
|
ny = ty + coef[1][1] * sin (tx / ((dy * dy) + EPS));
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[16];
|
|
if (v > 0.0)
|
|
{
|
|
/* fisheye */
|
|
double theta, r2, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
r2 = sqrt (tx * tx + ty * ty);
|
|
r2 = 2 * r2 / (r2 + 1);
|
|
nx = r2 * cos (theta);
|
|
ny = r2 * sin (theta);
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[17];
|
|
if (v > 0.0)
|
|
{
|
|
/* popcorn */
|
|
double dx, dy, nx, ny;
|
|
dx = tan (3 * ty);
|
|
dy = tan (3 * tx);
|
|
nx = tx + coef[2][0] * sin (dx);
|
|
ny = ty + coef[2][1] * sin (dy);
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[18];
|
|
if (v > 0.0)
|
|
{
|
|
/* exponential */
|
|
double dx, dy, nx, ny;
|
|
dx = exp (tx - 1.0);
|
|
dy = G_PI * ty;
|
|
nx = cos (dy) * dx;
|
|
ny = sin (dy) * dx;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[19];
|
|
if (v > 0.0)
|
|
{
|
|
/* power */
|
|
double theta, r2, tsin, tcos, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
tsin = sin (theta);
|
|
tcos = cos (theta);
|
|
r2 = sqrt (tx * tx + ty * ty);
|
|
r2 = pow (r2, tsin);
|
|
nx = r2 * tcos;
|
|
ny = r2 * tsin;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[20];
|
|
if (v > 0.0)
|
|
{
|
|
/* cosine */
|
|
double nx, ny;
|
|
nx = cos (tx * G_PI) * cosh (ty);
|
|
ny = -sin (tx * G_PI) * sinh (ty);
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[21];
|
|
if (v > 0.0)
|
|
{
|
|
/* rings */
|
|
double theta, r2, dx, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0;
|
|
dx = coef[2][0];
|
|
dx = dx * dx + EPS;
|
|
r2 = sqrt (tx * tx + ty * ty);
|
|
r2 = fmod (r2 + dx, 2 * dx) - dx + r2 * (1 - dx);
|
|
nx = cos (theta) * r2;
|
|
ny = sin (theta) * r2;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[22];
|
|
if (v > 0.0)
|
|
{
|
|
/* fan */
|
|
double theta, r2, dx, dy, dx2, nx, ny;
|
|
if (tx < -EPS || tx > EPS ||
|
|
ty < -EPS || ty > EPS)
|
|
theta = atan2 (tx, ty);
|
|
else
|
|
theta = 0.0;
|
|
dx = coef[2][0];
|
|
dy = coef[2][1];
|
|
dx = G_PI * (dx * dx + EPS);
|
|
dx2 = dx / 2;
|
|
r2 = sqrt (tx * tx + ty * ty );
|
|
theta += (fmod (theta + dy, dx) > dx2) ? -dx2: dx2;
|
|
nx = cos (theta) * r2;
|
|
ny = sin (theta) * r2;
|
|
p[0] += v * nx;
|
|
p[1] += v * ny;
|
|
}
|
|
|
|
v = vari[23];
|
|
if (v > 0.0)
|
|
{
|
|
/* eyefish */
|
|
double r2;
|
|
r2 = 2.0 * v / (sqrt(tx * tx + ty * ty) + 1.0);
|
|
p[0] += r2 * tx;
|
|
p[1] += r2 * ty;
|
|
}
|
|
|
|
v = vari[24];
|
|
if (v > 0.0)
|
|
{
|
|
/* bubble */
|
|
double r2;
|
|
r2 = v / ((tx * tx + ty * ty) / 4 + 1);
|
|
p[0] += r2 * tx;
|
|
p[1] += r2 * ty;
|
|
}
|
|
|
|
v = vari[25];
|
|
if (v > 0.0)
|
|
{
|
|
/* cylinder */
|
|
double nx;
|
|
nx = sin (tx);
|
|
p[0] += v * nx;
|
|
p[1] += v * ty;
|
|
}
|
|
|
|
v = vari[26];
|
|
if (v > 0.0)
|
|
{
|
|
/* noise */
|
|
double rx, sinr, cosr, nois;
|
|
rx = flam3_random01 () * 2 * G_PI;
|
|
sinr = sin (rx);
|
|
cosr = cos (rx);
|
|
nois = flam3_random01 ();
|
|
p[0] += v * nois * tx * cosr;
|
|
p[1] += v * nois * ty * sinr;
|
|
}
|
|
|
|
v = vari[27];
|
|
if (v > 0.0)
|
|
{
|
|
/* blur */
|
|
double rx, sinr, cosr, nois;
|
|
rx = flam3_random01 () * 2 * G_PI;
|
|
sinr = sin (rx);
|
|
cosr = cos (rx);
|
|
nois = flam3_random01 ();
|
|
p[0] += v * nois * cosr;
|
|
p[1] += v * nois * sinr;
|
|
}
|
|
|
|
v = vari[28];
|
|
if (v > 0.0)
|
|
{
|
|
/* gaussian */
|
|
double ang, sina, cosa, r2;
|
|
ang = flam3_random01 () * 2 * G_PI;
|
|
sina = sin (ang);
|
|
cosa = cos (ang);
|
|
r2 = v * (flam3_random01 () + flam3_random01 () + flam3_random01 () +
|
|
flam3_random01 () - 2.0);
|
|
p[0] += r2 * cosa;
|
|
p[1] += r2 * sina;
|
|
}
|
|
|
|
/* if fuse over, store it */
|
|
if (i >= 0)
|
|
{
|
|
points[i][0] = p[0];
|
|
points[i][1] = p[1];
|
|
points[i][2] = p[2];
|
|
}
|
|
}
|
|
}
|
|
|
|
/* args must be non-overlapping */
|
|
void
|
|
mult_matrix (double s1[2][2],
|
|
double s2[2][2],
|
|
double d[2][2])
|
|
{
|
|
d[0][0] = s1[0][0] * s2[0][0] + s1[1][0] * s2[0][1];
|
|
d[1][0] = s1[0][0] * s2[1][0] + s1[1][0] * s2[1][1];
|
|
d[0][1] = s1[0][1] * s2[0][0] + s1[1][1] * s2[0][1];
|
|
d[1][1] = s1[0][1] * s2[1][0] + s1[1][1] * s2[1][1];
|
|
}
|
|
|
|
static
|
|
double det_matrix (double s[2][2])
|
|
{
|
|
return s[0][0] * s[1][1] - s[0][1] * s[1][0];
|
|
}
|
|
|
|
static void
|
|
interpolate_angle (double t,
|
|
double s,
|
|
double *v1,
|
|
double *v2,
|
|
double *v3,
|
|
int tie,
|
|
int cross)
|
|
{
|
|
double x = *v1;
|
|
double y = *v2;
|
|
double d;
|
|
static double lastx, lasty;
|
|
|
|
/* take the shorter way around the circle... */
|
|
if (x > y)
|
|
{
|
|
d = x - y;
|
|
if (d > G_PI + EPS ||
|
|
(d > G_PI - EPS && tie))
|
|
y += 2 * G_PI;
|
|
}
|
|
else
|
|
{
|
|
d = y - x;
|
|
if (d > G_PI + EPS ||
|
|
(d > G_PI - EPS && tie))
|
|
x += 2 * G_PI;
|
|
}
|
|
/* unless we are supposed to avoid crossing */
|
|
if (cross)
|
|
{
|
|
if (lastx > x)
|
|
{
|
|
if (lasty < y)
|
|
y -= 2 * G_PI;
|
|
}
|
|
else
|
|
{
|
|
if (lasty > y)
|
|
y += 2 * G_PI;
|
|
}
|
|
}
|
|
else
|
|
{
|
|
lastx = x;
|
|
lasty = y;
|
|
}
|
|
|
|
*v3 = s * x + t * y;
|
|
}
|
|
|
|
static void
|
|
interpolate_complex (double t,
|
|
double s,
|
|
double *r1,
|
|
double *r2,
|
|
double *r3,
|
|
int flip,
|
|
int tie,
|
|
int cross)
|
|
{
|
|
double c1[2], c2[2], c3[2];
|
|
double a1, a2, a3, d1, d2, d3;
|
|
|
|
c1[0] = r1[0];
|
|
c1[1] = r1[1];
|
|
c2[0] = r2[0];
|
|
c2[1] = r2[1];
|
|
if (flip)
|
|
{
|
|
double t = c1[0];
|
|
c1[0] = c1[1];
|
|
c1[1] = t;
|
|
t = c2[0];
|
|
c2[0] = c2[1];
|
|
c2[1] = t;
|
|
}
|
|
|
|
/* convert to log space */
|
|
a1 = atan2 (c1[1], c1[0]);
|
|
a2 = atan2 (c2[1], c2[0]);
|
|
d1 = 0.5 * log (c1[0] * c1[0] + c1[1] * c1[1]);
|
|
d2 = 0.5 * log (c2[0] * c2[0] + c2[1] * c2[1]);
|
|
|
|
/* interpolate linearly */
|
|
interpolate_angle (t, s, &a1, &a2, &a3, tie, cross);
|
|
d3 = s * d1 + t * d2;
|
|
|
|
/* convert back */
|
|
d3 = exp (d3);
|
|
c3[0] = cos (a3) * d3;
|
|
c3[1] = sin (a3) * d3;
|
|
|
|
if (flip)
|
|
{
|
|
r3[1] = c3[0];
|
|
r3[0] = c3[1];
|
|
}
|
|
else
|
|
{
|
|
r3[0] = c3[0];
|
|
r3[1] = c3[1];
|
|
}
|
|
}
|
|
|
|
static void
|
|
interpolate_matrix (double t,
|
|
double m1[3][2],
|
|
double m2[3][2],
|
|
double m3[3][2])
|
|
{
|
|
double s = 1.0 - t;
|
|
|
|
interpolate_complex (t, s, &m1[0][0], &m2[0][0], &m3[0][0], 0, 0, 0);
|
|
interpolate_complex (t, s, &m1[1][0], &m2[1][0], &m3[1][0], 1, 1, 0);
|
|
|
|
/* handle the translation part of the xform linearly */
|
|
m3[2][0] = s * m1[2][0] + t * m2[2][0];
|
|
m3[2][1] = s * m1[2][1] + t * m2[2][1];
|
|
}
|
|
|
|
#define INTERP(x) result->x = c0 * cps[i1].x + c1 * cps[i2].x
|
|
|
|
/*
|
|
* create a control point that interpolates between the control points
|
|
* passed in CPS. for now just do linear. in the future, add control
|
|
* point types and other things to the cps. CPS must be sorted by time.
|
|
*/
|
|
void
|
|
interpolate (control_point cps[],
|
|
int ncps,
|
|
double time,
|
|
control_point *result)
|
|
{
|
|
int i, j, i1, i2;
|
|
double c0, c1, t;
|
|
|
|
g_return_if_fail (ncps > 0);
|
|
|
|
if (ncps == 1)
|
|
{
|
|
*result = cps[0];
|
|
return;
|
|
}
|
|
if (cps[0].time >= time)
|
|
{
|
|
i1 = 0;
|
|
i2 = 1;
|
|
}
|
|
else if (cps[ncps - 1].time <= time)
|
|
{
|
|
i1 = ncps - 2;
|
|
i2 = ncps - 1;
|
|
}
|
|
else
|
|
{
|
|
i1 = 0;
|
|
while (i1 < ncps && cps[i1].time < time)
|
|
i1++;
|
|
i1--;
|
|
i2 = i1 + 1;
|
|
|
|
if (i2 == ncps ||
|
|
(time - cps[i1].time > -1e-7 &&
|
|
time - cps[i1].time < 1e-7))
|
|
{
|
|
*result = cps[i1];
|
|
return;
|
|
}
|
|
}
|
|
|
|
c0 = (cps[i2].time - time) / (cps[i2].time - cps[i1].time);
|
|
c1 = 1.0 - c0;
|
|
|
|
result->time = time;
|
|
|
|
if (cps[i1].cmap_inter)
|
|
{
|
|
for (i = 0; i < 256; i++)
|
|
{
|
|
double spread = 0.15;
|
|
double d0, d1, e0, e1, c = 2 * G_PI * i / 256.0;
|
|
c = cos(c * cps[i1].cmap_inter) + 4.0 * c1 - 2.0;
|
|
if (c > spread) c = spread;
|
|
if (c < -spread) c = -spread;
|
|
d1 = (c + spread) * 0.5 / spread;
|
|
d0 = 1.0 - d1;
|
|
e0 = (d0 < 0.5) ? (d0 * 2) : (d1 * 2);
|
|
e1 = 1.0 - e0;
|
|
for (j = 0; j < 3; j++)
|
|
{
|
|
result->cmap[i][j] = (d0 * cps[i1].cmap[i][j] +
|
|
d1 * cps[i2].cmap[i][j]);
|
|
#define bright_peak 2.0
|
|
result->cmap[i][j] = (e1 * result->cmap[i][j] +
|
|
e0 * 1.0);
|
|
}
|
|
}
|
|
}
|
|
else
|
|
{
|
|
for (i = 0; i < 256; i++)
|
|
{
|
|
double t[3], s[3];
|
|
rgb2hsv (cps[i1].cmap[i], s);
|
|
rgb2hsv (cps[i2].cmap[i], t);
|
|
for (j = 0; j < 3; j++)
|
|
t[j] = c0 * s[j] + c1 * t[j];
|
|
hsv2rgb (t, result->cmap[i]);
|
|
}
|
|
}
|
|
|
|
result->cmap_index = -1;
|
|
INTERP(brightness);
|
|
INTERP(contrast);
|
|
INTERP(gamma);
|
|
INTERP(width);
|
|
INTERP(height);
|
|
INTERP(spatial_oversample);
|
|
INTERP(center[0]);
|
|
INTERP(center[1]);
|
|
INTERP(pixels_per_unit);
|
|
INTERP(spatial_filter_radius);
|
|
INTERP(sample_density);
|
|
INTERP(zoom);
|
|
INTERP(nbatches);
|
|
INTERP(white_level);
|
|
for (i = 0; i < 2; i++)
|
|
for (j = 0; j < 2; j++)
|
|
{
|
|
INTERP(pulse[i][j]);
|
|
INTERP(wiggle[i][j]);
|
|
}
|
|
|
|
for (i = 0; i < NXFORMS; i++)
|
|
{
|
|
double r;
|
|
double rh_time;
|
|
|
|
INTERP(xform[i].density);
|
|
if (result->xform[i].density > 0)
|
|
result->xform[i].density = 1.0;
|
|
INTERP(xform[i].color);
|
|
for (j = 0; j < NVARS; j++)
|
|
INTERP(xform[i].var[j]);
|
|
t = 0.0;
|
|
for (j = 0; j < NVARS; j++)
|
|
t += result->xform[i].var[j];
|
|
t = 1.0 / t;
|
|
for (j = 0; j < NVARS; j++)
|
|
result->xform[i].var[j] *= t;
|
|
|
|
interpolate_matrix(c1, cps[i1].xform[i].c, cps[i2].xform[i].c,
|
|
result->xform[i].c);
|
|
|
|
rh_time = time * 2 * G_PI / (60.0 * 30.0);
|
|
|
|
/* apply pulse factor. */
|
|
r = 1.0;
|
|
for (j = 0; j < 2; j++)
|
|
r += result->pulse[j][0] * sin(result->pulse[j][1] * rh_time);
|
|
for (j = 0; j < 3; j++)
|
|
{
|
|
result->xform[i].c[j][0] *= r;
|
|
result->xform[i].c[j][1] *= r;
|
|
}
|
|
|
|
/* apply wiggle factor */
|
|
for (j = 0; j < 2; j++)
|
|
{
|
|
double tt = result->wiggle[j][1] * rh_time;
|
|
double m = result->wiggle[j][0];
|
|
result->xform[i].c[0][0] += m * cos(tt);
|
|
result->xform[i].c[1][0] += m * -sin(tt);
|
|
result->xform[i].c[0][1] += m * sin(tt);
|
|
result->xform[i].c[1][1] += m * cos(tt);
|
|
}
|
|
} /* for i */
|
|
}
|
|
|
|
|
|
|
|
/*
|
|
* split a string passed in ss into tokens on whitespace.
|
|
* # comments to end of line. ; terminates the record
|
|
*/
|
|
void
|
|
tokenize (char **ss,
|
|
char *argv[],
|
|
int *argc)
|
|
{
|
|
char *s = *ss;
|
|
int i = 0, state = 0;
|
|
gint len = 0;
|
|
|
|
len = strlen (s);
|
|
while (*s != ';' && len > 0)
|
|
{
|
|
char c = *s;
|
|
switch (state)
|
|
{
|
|
case 0:
|
|
if (c == '#')
|
|
state = 2;
|
|
else if (!g_ascii_isspace (c))
|
|
{
|
|
argv[i] = s;
|
|
i++;
|
|
state = 1;
|
|
}
|
|
break;
|
|
case 1:
|
|
if (g_ascii_isspace (c))
|
|
{
|
|
*s = 0;
|
|
state = 0;
|
|
}
|
|
break;
|
|
case 2:
|
|
if (c == '\n')
|
|
state = 0;
|
|
break;
|
|
}
|
|
s++;
|
|
len--;
|
|
}
|
|
*s = 0;
|
|
*ss = s + 1;
|
|
*argc = i;
|
|
}
|
|
|
|
static int
|
|
compare_xforms (const void *va,
|
|
const void *vb)
|
|
{
|
|
double aa[2][2];
|
|
double bb[2][2];
|
|
double ad, bd;
|
|
const xform *a = va;
|
|
const xform *b = vb;
|
|
|
|
aa[0][0] = a->c[0][0];
|
|
aa[0][1] = a->c[0][1];
|
|
aa[1][0] = a->c[1][0];
|
|
aa[1][1] = a->c[1][1];
|
|
bb[0][0] = b->c[0][0];
|
|
bb[0][1] = b->c[0][1];
|
|
bb[1][0] = b->c[1][0];
|
|
bb[1][1] = b->c[1][1];
|
|
ad = det_matrix (aa);
|
|
bd = det_matrix (bb);
|
|
if (ad < bd)
|
|
return -1;
|
|
if (ad > bd)
|
|
return 1;
|
|
return 0;
|
|
}
|
|
|
|
#define MAXARGS 1000
|
|
#define streql(x,y) (!strcmp(x,y))
|
|
|
|
/*
|
|
* given a pointer to a string SS, fill fields of a control point CP.
|
|
*/
|
|
|
|
void
|
|
parse_control_point (char **ss,
|
|
control_point *cp)
|
|
{
|
|
char *argv[MAXARGS];
|
|
int argc, i, j;
|
|
gint64 xf_index = 0;
|
|
gint parse_errors = 0;
|
|
double *slot = NULL, t;
|
|
|
|
for (i = 0; i < NXFORMS; i++)
|
|
{
|
|
cp->xform[i].density = 0.0;
|
|
cp->xform[i].color = (i == 0);
|
|
cp->xform[i].var[0] = 1.0;
|
|
for (j = 1; j < NVARS; j++)
|
|
cp->xform[i].var[j] = 0.0;
|
|
cp->xform[i].c[0][0] = 1.0;
|
|
cp->xform[i].c[0][1] = 0.0;
|
|
cp->xform[i].c[1][0] = 0.0;
|
|
cp->xform[i].c[1][1] = 1.0;
|
|
cp->xform[i].c[2][0] = 0.0;
|
|
cp->xform[i].c[2][1] = 0.0;
|
|
}
|
|
for (j = 0; j < 2; j++)
|
|
{
|
|
cp->pulse[j][0] = 0.0;
|
|
cp->pulse[j][1] = 60.0;
|
|
cp->wiggle[j][0] = 0.0;
|
|
cp->wiggle[j][1] = 60.0;
|
|
}
|
|
|
|
tokenize (ss, argv, &argc);
|
|
|
|
i = 0;
|
|
while (i < argc)
|
|
{
|
|
gint itoken;
|
|
|
|
itoken = i;
|
|
if (i < argc)
|
|
{
|
|
/* First value belonging to token. */
|
|
i++;
|
|
}
|
|
else
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
|
|
if (streql ("xform", argv[itoken]))
|
|
{
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 0, NXFORMS-1, &xf_index, NULL))
|
|
{
|
|
g_printerr ("Invalid xform index '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
xf_index = 0;
|
|
}
|
|
}
|
|
else if (streql ("density", argv[itoken]))
|
|
{
|
|
cp->xform[xf_index].density = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("color", argv[itoken]))
|
|
{
|
|
cp->xform[xf_index].color = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("coefs", argv[itoken]))
|
|
{
|
|
/* We need 6 coef values and we know are at the first */
|
|
if (i + 5 >= argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
slot = cp->xform[xf_index].c[0];
|
|
for (j = 0; j < 6; j++)
|
|
{
|
|
*slot++ = g_strtod (argv[i++], NULL);
|
|
}
|
|
cp->xform[xf_index].density = 1.0;
|
|
}
|
|
else if (streql ("var", argv[itoken]))
|
|
{
|
|
/* We need NVARS var values and we know are at the first */
|
|
if (i + NVARS > argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
slot = cp->xform[xf_index].var;
|
|
for (j = 0; j < NVARS; j++)
|
|
{
|
|
*slot++ = g_strtod (argv[i++], NULL);
|
|
}
|
|
}
|
|
else if (streql ("time", argv[itoken]))
|
|
{
|
|
cp->time = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("brightness", argv[itoken]))
|
|
{
|
|
cp->brightness = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("contrast", argv[itoken]))
|
|
{
|
|
cp->contrast = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("gamma", argv[itoken]))
|
|
{
|
|
cp->gamma = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("zoom", argv[itoken]))
|
|
{
|
|
cp->zoom = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("image_size", argv[itoken]))
|
|
{
|
|
gint64 w, h;
|
|
|
|
/* We need 2 values and we know are at the first */
|
|
if (i + 1 >= argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 1, PIKA_MAX_IMAGE_SIZE, &w, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid image width '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else if (! g_ascii_string_to_signed (argv[i++], 10, 1, PIKA_MAX_IMAGE_SIZE, &h, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid image_size heigth '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->width = w;
|
|
cp->height = h;
|
|
}
|
|
}
|
|
else if (streql ("center", argv[itoken]))
|
|
{
|
|
/* We need 2 values and we know are at the first */
|
|
if (i + 1 >= argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
cp->center[0] = g_strtod (argv[i++], NULL);
|
|
cp->center[1] = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("pixels_per_unit", argv[itoken]))
|
|
{
|
|
cp->pixels_per_unit = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("pulse", argv[itoken]))
|
|
{
|
|
/* We need 4 values and we know are at the first */
|
|
if (i + 3 >= argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
slot = &cp->pulse[0][0];
|
|
for (j = 0; j < 4; j++)
|
|
{
|
|
*slot++ = g_strtod (argv[i++], NULL);
|
|
}
|
|
}
|
|
else if (streql ("wiggle", argv[itoken]))
|
|
{
|
|
/* We need 4 values and we know are at the first */
|
|
if (i + 3 >= argc)
|
|
{
|
|
g_printerr ("Not enough parameters. File may be corrupt!\n");
|
|
parse_errors++;
|
|
break;
|
|
}
|
|
slot = &cp->wiggle[0][0];
|
|
for (j = 0; j < 4; j++)
|
|
{
|
|
*slot++ = g_strtod (argv[i++], NULL);
|
|
}
|
|
}
|
|
else if (streql ("spatial_oversample", argv[itoken]))
|
|
{
|
|
gint64 oversample;
|
|
|
|
/* Values in the gui seem to be between 1 and 4 */
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 1, 4, &oversample, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid spatial oversample value '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->spatial_oversample = oversample;
|
|
}
|
|
}
|
|
else if (streql ("spatial_filter_radius", argv[itoken]))
|
|
{
|
|
cp->spatial_filter_radius = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("sample_density", argv[itoken]))
|
|
{
|
|
cp->sample_density = g_strtod (argv[i++], NULL);
|
|
}
|
|
else if (streql ("nbatches", argv[itoken]))
|
|
{
|
|
gint64 nbatches;
|
|
|
|
/* Not sure what the maximum should be. It always seems to be set to 1. */
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 0, 2, &nbatches, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid nbatches value '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->nbatches = nbatches;
|
|
}
|
|
}
|
|
else if (streql ("white_level", argv[itoken]))
|
|
{
|
|
gint64 wl;
|
|
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 0, 255, &wl, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid white level value '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->white_level = wl;
|
|
}
|
|
}
|
|
else if (streql ("cmap", argv[itoken]))
|
|
{
|
|
gint64 cmi;
|
|
|
|
/* -1 = random */
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, -1, 255, &cmi, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid color map value '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->cmap_index = cmi;
|
|
}
|
|
}
|
|
else if (streql ("cmap_inter", argv[itoken]))
|
|
{
|
|
gint64 cmi;
|
|
|
|
/* 0 or 1 */
|
|
if (! g_ascii_string_to_signed (argv[i++], 10, 0, 1, &cmi, NULL))
|
|
{
|
|
g_printerr ("Ignoring invalid color interpolate value '%s'\n", argv[i-1]);
|
|
parse_errors++;
|
|
}
|
|
else
|
|
{
|
|
cp->cmap_inter = cmi;
|
|
}
|
|
}
|
|
else
|
|
{
|
|
g_printerr ("Invalid token '%s'. File may be corrupt!\n", argv[itoken]);
|
|
parse_errors++;
|
|
}
|
|
}
|
|
|
|
if (parse_errors > 0)
|
|
g_warning ("Input file contains %d errors. File may be corrupt!", parse_errors);
|
|
|
|
for (i = 0; i < NXFORMS; i++)
|
|
{
|
|
t = 0.0;
|
|
for (j = 0; j < NVARS; j++)
|
|
t += cp->xform[i].var[j];
|
|
t = 1.0 / t;
|
|
for (j = 0; j < NVARS; j++)
|
|
cp->xform[i].var[j] *= t;
|
|
}
|
|
qsort ((char *) cp->xform, NXFORMS, sizeof(xform), compare_xforms);
|
|
}
|
|
|
|
void
|
|
print_control_point (FILE *f,
|
|
control_point *cp,
|
|
int quote)
|
|
{
|
|
int i, j;
|
|
char *q = quote ? "# " : "";
|
|
fprintf (f, "%stime %g\n", q, cp->time);
|
|
if (cp->cmap_index != -1)
|
|
fprintf (f, "%scmap %d\n", q, cp->cmap_index);
|
|
fprintf (f, "%simage_size %d %d center %g %g pixels_per_unit %g\n",
|
|
q, cp->width, cp->height, cp->center[0], cp->center[1],
|
|
cp->pixels_per_unit);
|
|
fprintf (f, "%sspatial_oversample %d spatial_filter_radius %g",
|
|
q, cp->spatial_oversample, cp->spatial_filter_radius);
|
|
fprintf (f, " sample_density %g\n", cp->sample_density);
|
|
fprintf (f, "%snbatches %d white_level %d\n",
|
|
q, cp->nbatches, cp->white_level);
|
|
fprintf (f, "%sbrightness %g gamma %g cmap_inter %d\n",
|
|
q, cp->brightness, cp->gamma, cp->cmap_inter);
|
|
|
|
for (i = 0; i < NXFORMS; i++)
|
|
if (cp->xform[i].density > 0.0)
|
|
{
|
|
fprintf (f, "%sxform %d density %g color %g\n",
|
|
q, i, cp->xform[i].density, cp->xform[i].color);
|
|
fprintf (f, "%svar", q);
|
|
for (j = 0; j < NVARS; j++)
|
|
fprintf (f, " %g", cp->xform[i].var[j]);
|
|
fprintf (f, "\n%scoefs", q);
|
|
for (j = 0; j < 3; j++)
|
|
fprintf (f, " %g %g", cp->xform[i].c[j][0], cp->xform[i].c[j][1]);
|
|
fprintf (f, "\n");
|
|
}
|
|
fprintf (f, "%s;\n", q);
|
|
}
|
|
|
|
/* returns a uniform variable from 0 to 1 */
|
|
double
|
|
random_uniform01 (void)
|
|
{
|
|
return g_random_double ();
|
|
}
|
|
|
|
double random_uniform11 (void)
|
|
{
|
|
return g_random_double_range (-1, 1);
|
|
}
|
|
|
|
/* returns a mean 0 variance 1 random variable
|
|
see numerical recipes p 217 */
|
|
double random_gaussian(void)
|
|
{
|
|
static int iset = 0;
|
|
static double gset;
|
|
double fac, r, v1, v2;
|
|
|
|
if (iset == 0)
|
|
{
|
|
do
|
|
{
|
|
v1 = random_uniform11 ();
|
|
v2 = random_uniform11 ();
|
|
r = v1 * v1 + v2 * v2;
|
|
}
|
|
while (r >= 1.0 || r == 0.0);
|
|
fac = sqrt (-2.0 * log (r) / r);
|
|
gset = v1 * fac;
|
|
iset = 1;
|
|
return v2 * fac;
|
|
}
|
|
iset = 0;
|
|
return gset;
|
|
}
|
|
|
|
void
|
|
copy_variation (control_point *cp0,
|
|
control_point *cp1)
|
|
{
|
|
int i, j;
|
|
for (i = 0; i < NXFORMS; i++)
|
|
{
|
|
for (j = 0; j < NVARS; j++)
|
|
cp0->xform[i].var[j] = cp1->xform[i].var[j];
|
|
}
|
|
}
|
|
|
|
#define random_distrib(v) ((v)[g_random_int_range (0, vlen(v))])
|
|
|
|
void
|
|
random_control_point (control_point *cp,
|
|
int ivar)
|
|
{
|
|
int i, nxforms, var;
|
|
static int xform_distrib[] =
|
|
{
|
|
2, 2, 2,
|
|
3, 3, 3,
|
|
4, 4,
|
|
5
|
|
};
|
|
static int var_distrib[] =
|
|
{
|
|
-1, -1, -1,
|
|
0, 0, 0, 0,
|
|
1, 1, 1,
|
|
2, 2, 2,
|
|
3, 3,
|
|
4, 4,
|
|
5
|
|
};
|
|
|
|
static int mixed_var_distrib[] =
|
|
{
|
|
0, 0, 0,
|
|
1, 1, 1,
|
|
2, 2, 2,
|
|
3, 3,
|
|
4, 4,
|
|
5, 5
|
|
};
|
|
|
|
get_cmap (cmap_random, cp->cmap, 256);
|
|
cp->time = 0.0;
|
|
nxforms = random_distrib (xform_distrib);
|
|
var = (0 > ivar) ?
|
|
random_distrib(var_distrib) :
|
|
ivar;
|
|
for (i = 0; i < nxforms; i++)
|
|
{
|
|
int j, k;
|
|
cp->xform[i].density = 1.0 / nxforms;
|
|
cp->xform[i].color = i == 0;
|
|
for (j = 0; j < 3; j++)
|
|
for (k = 0; k < 2; k++)
|
|
cp->xform[i].c[j][k] = random_uniform11();
|
|
for (j = 0; j < NVARS; j++)
|
|
cp->xform[i].var[j] = 0.0;
|
|
if (var >= 0)
|
|
cp->xform[i].var[var] = 1.0;
|
|
else
|
|
cp->xform[i].var[random_distrib(mixed_var_distrib)] = 1.0;
|
|
|
|
}
|
|
for (; i < NXFORMS; i++)
|
|
cp->xform[i].density = 0.0;
|
|
}
|
|
|
|
/*
|
|
* find a 2d bounding box that does not enclose eps of the fractal density
|
|
* in each compass direction. works by binary search.
|
|
* this is stupid, it should just use the find nth smallest algorithm.
|
|
*/
|
|
void
|
|
estimate_bounding_box (control_point *cp,
|
|
double eps,
|
|
double *bmin,
|
|
double *bmax)
|
|
{
|
|
int i, j, batch = (eps == 0.0) ? 10000 : 10.0/eps;
|
|
int low_target = batch * eps;
|
|
int high_target = batch - low_target;
|
|
point min, max, delta;
|
|
point *points = g_malloc0 (sizeof (point) * batch);
|
|
|
|
iterate (cp, batch, 20, points);
|
|
|
|
min[0] = min[1] = 1e10;
|
|
max[0] = max[1] = -1e10;
|
|
|
|
for (i = 0; i < batch; i++)
|
|
{
|
|
if (points[i][0] < min[0]) min[0] = points[i][0];
|
|
if (points[i][1] < min[1]) min[1] = points[i][1];
|
|
if (points[i][0] > max[0]) max[0] = points[i][0];
|
|
if (points[i][1] > max[1]) max[1] = points[i][1];
|
|
}
|
|
|
|
if (low_target == 0)
|
|
{
|
|
bmin[0] = min[0];
|
|
bmin[1] = min[1];
|
|
bmax[0] = max[0];
|
|
bmax[1] = max[1];
|
|
return;
|
|
}
|
|
|
|
delta[0] = (max[0] - min[0]) * 0.25;
|
|
delta[1] = (max[1] - min[1]) * 0.25;
|
|
|
|
bmax[0] = bmin[0] = min[0] + 2.0 * delta[0];
|
|
bmax[1] = bmin[1] = min[1] + 2.0 * delta[1];
|
|
|
|
for (i = 0; i < 14; i++)
|
|
{
|
|
int n, s, e, w;
|
|
n = s = e = w = 0;
|
|
for (j = 0; j < batch; j++)
|
|
{
|
|
if (points[j][0] < bmin[0]) n++;
|
|
if (points[j][0] > bmax[0]) s++;
|
|
if (points[j][1] < bmin[1]) w++;
|
|
if (points[j][1] > bmax[1]) e++;
|
|
}
|
|
bmin[0] += (n < low_target) ? delta[0] : -delta[0];
|
|
bmax[0] += (s < high_target) ? delta[0] : -delta[0];
|
|
bmin[1] += (w < low_target) ? delta[1] : -delta[1];
|
|
bmax[1] += (e < high_target) ? delta[1] : -delta[1];
|
|
delta[0] = delta[0] / 2.0;
|
|
delta[1] = delta[1] / 2.0;
|
|
}
|
|
g_free (points);
|
|
}
|
|
|
|
/* this has serious flaws in it */
|
|
|
|
double
|
|
standard_metric (control_point *cp1,
|
|
control_point *cp2)
|
|
{
|
|
int i, j, k;
|
|
double t;
|
|
double dist = 0.0;
|
|
|
|
for (i = 0; i < NXFORMS; i++)
|
|
{
|
|
double var_dist = 0.0;
|
|
double coef_dist = 0.0;
|
|
for (j = 0; j < NVARS; j++)
|
|
{
|
|
t = cp1->xform[i].var[j] - cp2->xform[i].var[j];
|
|
var_dist += t * t;
|
|
}
|
|
for (j = 0; j < 3; j++)
|
|
for (k = 0; k < 2; k++)
|
|
{
|
|
t = cp1->xform[i].c[j][k] - cp2->xform[i].c[j][k];
|
|
coef_dist += t *t;
|
|
}
|
|
/* weight them equally for now. */
|
|
dist += var_dist + coef_dist;
|
|
}
|
|
return dist;
|
|
}
|
|
|
|
static int
|
|
flam3_random_bit (void)
|
|
{
|
|
static int n = 0;
|
|
static int l;
|
|
if (n == 0)
|
|
{
|
|
l = g_random_int ();
|
|
n = 20;
|
|
}
|
|
else
|
|
{
|
|
l = l >> 1;
|
|
n--;
|
|
}
|
|
return l & 1;
|
|
}
|
|
|
|
static double
|
|
flam3_random01 (void)
|
|
{
|
|
return (g_random_int () & 0xfffffff) / (double) 0xfffffff;
|
|
}
|
