git.hungrycats.org Git - xscreensaver/blob - hacks/xlyap.c

   1 /* Lyap - calculate and display Lyapunov exponents */
   2
   3 /* Written by Ron Record (rr@sco) 03 Sep 1991 */
   4
   5 /* The idea here is to calculate the Lyapunov exponent for a periodically
   6  * forced logistic map (later i added several other nonlinear maps of the unit
   7  * interval). In order to turn the 1-dimensional parameter space of the
   8  * logistic map into a 2-dimensional parameter space, select two parameter
   9  * values ('a' and 'b') then alternate the iterations of the logistic map using
  10  * first 'a' then 'b' as the parameter. This program accepts an argument to
  11  * specify a forcing function, so instead of just alternating 'a' and 'b', you
  12  * can use 'a' as the parameter for say 6 iterations, then 'b' for 6 iterations
  13  * and so on. An interesting forcing function to look at is abbabaab (the
  14  * Morse-Thue sequence, an aperiodic self-similar, self-generating sequence).
  15  * Anyway, step through all the values of 'a' and 'b' in the ranges you want,
  16  * calculating the Lyapunov exponent for each pair of values. The exponent
  17  * is calculated by iterating out a ways (specified by the variable "settle")
  18  * then on subsequent iterations calculating an average of the logarithm of
  19  * the absolute value of the derivative at that point. Points in parameter
  20  * space with a negative Lyapunov exponent are colored one way (using the
  21  * value of the exponent to index into a color map) while points with a
  22  * non-negative exponent are colored differently.
  23  *
  24  * The algorithm was taken from the September 1991 Scientific American article
  25  * by A. K. Dewdney who gives credit to Mario Markus of the Max Planck
  26  * Institute for its creation. Additional information and ideas were gleaned
  27  * from the discussion on alt.fractals involving Stephen Hall, Ed Kubaitis,
  28  * Dave Platt and Baback Moghaddam. Assistance with colormaps and spinning
  29  * color wheels and X was gleaned from Hiram Clawson. Rubber banding code was
  30  * adapted from an existing Mandelbrot program written by Stacey Campbell.
  31  */
  32
  33 #define LYAP_PATCHLEVEL 4
  34 #define LYAP_VERSION "#(@) lyap 2.3 2/20/92"
  35
  36 #include <assert.h>
  37 #include <math.h>
  38
  39 #include "screenhack.h"
  40 #include "yarandom.h"
  41 #include "hsv.h"
  42
  43 #ifndef HAVE_JWXYZ
  44 # include <X11/cursorfont.h>
  45 #endif
  46
  47 static const char *xlyap_defaults [] = {
  48   ".background:         black",
  49   ".foreground:         white",
  50 /*  ".lowrez:           true", */
  51   "*fpsSolid:           true",
  52   "*randomize:          true",
  53   "*builtin:            -1",
  54   "*minColor:           1",
  55   "*maxColor:           256",
  56   "*dwell:              50",
  57   "*useLog:             false",
  58   "*colorExponent:      1.0",
  59   "*colorOffset:        0",
  60   "*randomForce:        ",              /* 0.5 */
  61   "*settle:             50",
  62   "*minA:               2.0",
  63   "*minB:               2.0",
  64   "*wheels:             7",
  65   "*function:           10101010",
  66   "*forcingFunction:    abbabaab",
  67   "*bRange:             ",              /* 2.0 */
  68   "*startX:             0.65",
  69   "*mapIndex:           ",              /* 0 */
  70   "*outputFile:         ",
  71   "*beNegative:         false",
  72   "*rgbMax:             65000",
  73   "*spinLength:         256",
  74   "*show:               false",
  75   "*aRange:             ",              /* 2.0 */
  76   "*delay:              10000",
  77   "*linger:             5",
  78   "*colors:             200",
  79 #ifdef HAVE_MOBILE
  80   "*ignoreRotation:     True",
  81 #endif
  82   0
  83 };
  84
  85 static XrmOptionDescRec xlyap_options [] = {
  86   { "-randomize", ".randomize", XrmoptionNoArg, "true" },
  87   { "-builtin",   ".builtin",   XrmoptionSepArg, 0 },
  88   { "-C", ".minColor",          XrmoptionSepArg, 0 },   /* n */
  89   { "-D", ".dwell",             XrmoptionSepArg, 0 },   /* n */
  90   { "-L", ".useLog",            XrmoptionNoArg, "true" },
  91   { "-M", ".colorExponent",     XrmoptionSepArg, 0 },   /* r */
  92   { "-O", ".colorOffset",       XrmoptionSepArg, 0 },   /* n */
  93   { "-R", ".randomForce",       XrmoptionSepArg, 0 },   /* p */
  94   { "-S", ".settle",            XrmoptionSepArg, 0 },   /* n */
  95   { "-a", ".minA",              XrmoptionSepArg, 0 },   /* r */
  96   { "-b", ".minB",              XrmoptionSepArg, 0 },   /* n */
  97   { "-c", ".wheels",            XrmoptionSepArg, 0 },   /* n */
  98   { "-F", ".function",          XrmoptionSepArg, 0 },   /* 10101010 */
  99   { "-f", ".forcingFunction",   XrmoptionSepArg, 0 },   /* abbabaab */
 100   { "-h", ".bRange",            XrmoptionSepArg, 0 },   /* r */
 101   { "-i", ".startX",            XrmoptionSepArg, 0 },   /* r */
 102   { "-m", ".mapIndex",          XrmoptionSepArg, 0 },   /* n */
 103   { "-o", ".outputFile",        XrmoptionSepArg, 0 },   /* filename */
 104   { "-p", ".beNegative",        XrmoptionNoArg, "true" },
 105   { "-r", ".rgbMax",            XrmoptionSepArg, 0 },   /* n */
 106   { "-s", ".spinLength",        XrmoptionSepArg, 0 },   /* n */
 107   { "-v", ".show",              XrmoptionNoArg, "true" },
 108   { "-w", ".aRange",            XrmoptionSepArg, 0 },   /* r */
 109   { "-delay", ".delay",         XrmoptionSepArg, 0 },   /* delay */
 110   { "-linger", ".linger",       XrmoptionSepArg, 0 },   /* linger */
 111   { 0, 0, 0, 0 }
 112 };
 113
 114
 115 #define ABS(a)  (((a)<0) ? (0-(a)) : (a) )
 116 #define Min(x,y) ((x < y)?x:y)
 117 #define Max(x,y) ((x > y)?x:y)
 118
 119 #ifdef SIXTEEN_COLORS
 120 # define MAXPOINTS  128
 121 # ifdef BIGMEM
 122 #  define MAXFRAMES 4
 123 # else  /* !BIGMEM */
 124 #  define MAXFRAMES 2
 125 # endif /* !BIGMEM */
 126 # define MAXCOLOR 16
 127 #else  /* !SIXTEEN_COLORS */
 128 # define MAXPOINTS  256
 129 # ifdef BIGMEM
 130 #  define MAXFRAMES 8
 131 # else  /* !BIGMEM */
 132 #  define MAXFRAMES 2
 133 # endif /* !BIGMEM */
 134 # define MAXCOLOR 256
 135 #endif /* !SIXTEEN_COLORS */
 136
 137
 138 #define MAXINDEX 64
 139 #define FUNCMAXINDEX 16
 140 #define MAXWHEELS 7
 141 #define NUMMAPS 5
 142 #define NBUILTINS 22
 143
 144 #ifndef TRUE
 145 # define TRUE 1
 146 # define FALSE 0
 147 #endif
 148
 149
 150 typedef struct {
 151   int x, y;
 152 } xy_t;
 153
 154 #if 0
 155 typedef struct {
 156   int start_x, start_y;
 157   int last_x, last_y;
 158 } rubber_band_data_t;
 159 #endif
 160
 161 typedef struct {
 162 # ifndef HAVE_JWXYZ
 163   Cursor band_cursor;
 164 # endif
 165   double p_min, p_max, q_min, q_max;
 166 /*  rubber_band_data_t rubber_band;*/
 167 } image_data_t;
 168
 169 typedef struct points_t {
 170   XPoint data[MAXCOLOR][MAXPOINTS];
 171   int npoints[MAXCOLOR];
 172 } points_t;
 173
 174
 175 typedef double (*PFD)(double,double);
 176
 177 /* #### What was this for?  Everything was drawn twice, to the window and
 178    to this, and this was never displayed! */
 179 /*#define BACKING_PIXMAP*/
 180
 181 struct state {
 182   Display *dpy;
 183   Screen *screen;
 184   Visual *visual;
 185   Colormap cmap;
 186
 187   unsigned long foreground, background;
 188
 189   Window canvas;
 190   int delay, linger;
 191
 192   unsigned int maxcolor, startcolor, mincolindex;
 193   int color_offset;
 194   int dwell, settle;
 195   int width, height, xposition, yposition;
 196
 197   points_t Points;
 198 /*  image_data_t rubber_data;*/
 199
 200   GC Data_GC[MAXCOLOR]/*, RubberGC*/;
 201   PFD map, deriv;
 202
 203   int aflag, bflag, wflag, hflag, Rflag;
 204
 205   int   maxindex;
 206   int   funcmaxindex;
 207   double  min_a, min_b, a_range, b_range, minlyap;
 208   double  max_a, max_b;
 209   double  start_x, lyapunov, a_inc, b_inc, a, b;
 210   int   numcolors, numfreecols, lowrange;
 211   xy_t  point;
 212 #ifdef BACKING_PIXMAP
 213   Pixmap  pixmap;
 214 #endif
 215 /*  XColor  Colors[MAXCOLOR];*/
 216   double  *exponents[MAXFRAMES];
 217   double  a_minimums[MAXFRAMES], b_minimums[MAXFRAMES];
 218   double  a_maximums[MAXFRAMES], b_maximums[MAXFRAMES];
 219   double  minexp, maxexp, prob;
 220   int     expind[MAXFRAMES], resized[MAXFRAMES];
 221   int     numwheels, force, Force, negative;
 222   int     rgb_max, nostart, stripe_interval;
 223   int     save, show, useprod, spinlength;
 224   int     maxframe, frame, dorecalc, mapindex, run;
 225   char    *outname;
 226
 227   int sendpoint_index;
 228
 229   int   forcing[MAXINDEX];
 230   int   Forcing[FUNCMAXINDEX];
 231
 232   int reset_countdown;
 233
 234   int ncolors;
 235   XColor colors[MAXCOLOR];
 236 };
 237
 238
 239 static const double pmins[NUMMAPS] = { 2.0, 0.0, 0.0, 0.0, 0.0 };
 240 static const double pmaxs[NUMMAPS] = { 4.0, 1.0, 6.75, 6.75, 16.0 };
 241 static const double amins[NUMMAPS] = { 2.0, 0.0, 0.0, 0.0, 0.0 };
 242 static const double aranges[NUMMAPS] = { 2.0, 1.0, 6.75, 6.75, 16.0 };
 243 static const double bmins[NUMMAPS] = { 2.0, 0.0, 0.0, 0.0, 0.0 };
 244 static const double branges[NUMMAPS] = { 2.0, 1.0, 6.75, 6.75, 16.0 };
 245
 246 /****************************************************************************/
 247
 248 /* callback function declarations
 249  */
 250
 251 static double logistic(double,double);
 252 static double circle(double,double);
 253 static double leftlog(double,double);
 254 static double rightlog(double,double);
 255 static double doublelog(double,double);
 256 static double dlogistic(double,double);
 257 static double dcircle(double,double);
 258 static double dleftlog(double,double);
 259 static double drightlog(double,double);
 260 static double ddoublelog(double,double);
 261
 262 static const PFD Maps[NUMMAPS] = { logistic, circle, leftlog, rightlog,
 263                                    doublelog };
 264 static const PFD Derivs[NUMMAPS] = { dlogistic, dcircle, dleftlog,
 265                                      drightlog, ddoublelog };
 266
 267
 268 /****************************************************************************/
 269
 270 /* other function declarations
 271  */
 272
 273 static void resize(struct state *);
 274 /*static void Spin(struct state *);*/
 275 static void show_defaults(struct state *);
 276 /*static void StartRubberBand(struct state *, image_data_t *, XEvent *);
 277 static void TrackRubberBand(struct state *, image_data_t *, XEvent *);
 278 static void EndRubberBand(struct state *, image_data_t *, XEvent *);*/
 279 /*static void CreateXorGC(struct state *);*/
 280 static void InitBuffer(struct state *);
 281 static void BufferPoint(struct state *, int color, int x, int y);
 282 static void FlushBuffer(struct state *);
 283 static void init_data(struct state *);
 284 static void init_color(struct state *);
 285 static void parseargs(struct state *);
 286 static void Clear(struct state *);
 287 static void setupmem(struct state *);
 288 static int complyap(struct state *);
 289 static Bool Getkey(struct state *, XKeyEvent *);
 290 static int sendpoint(struct state *, double expo);
 291 /*static void save_to_file(struct state *);*/
 292 static void setforcing(struct state *);
 293 static void check_params(struct state *, int mapnum, int parnum);
 294 static void usage(struct state *);
 295 static void Destroy_frame(struct state *);
 296 static void freemem(struct state *);
 297 static void Redraw(struct state *);
 298 static void redraw(struct state *, double *exparray, int index, int cont);
 299 static void recalc(struct state *);
 300 /*static void SetupCorners(XPoint *, image_data_t *);
 301 static void set_new_params(struct state *, image_data_t *);*/
 302 static void go_down(struct state *);
 303 static void go_back(struct state *);
 304 static void go_init(struct state *);
 305 static void jumpwin(struct state *);
 306 static void print_help(struct state *);
 307 static void print_values(struct state *);
 308
 309
 310 /****************************************************************************/
 311
 312
 313 /* complyap() is the guts of the program. This is where the Lyapunov exponent
 314  * is calculated. For each iteration (past some large number of iterations)
 315  * calculate the logarithm of the absolute value of the derivative at that
 316  * point. Then average them over some large number of iterations. Some small
 317  * speed up is achieved by utilizing the fact that log(a*b) = log(a) + log(b).
 318  */
 319 static int
 320 complyap(struct state *st)
 321 {
 322   int i, bindex;
 323   double total, prod, x, dx, r;
 324
 325   if (st->maxcolor > MAXCOLOR)
 326     abort();
 327
 328   if (!st->run)
 329     return TRUE;
 330   st->a += st->a_inc;
 331   if (st->a >= st->max_a) {
 332     if (sendpoint(st, st->lyapunov) == TRUE)
 333       return FALSE;
 334     else {
 335       FlushBuffer(st);
 336       /*      if (savefile)
 337               save_to_file(); */
 338       return TRUE;
 339     }
 340   }
 341   if (st->b >= st->max_b) {
 342     FlushBuffer(st);
 343     /*    if (savefile)
 344           save_to_file();*/
 345     return TRUE;
 346   }
 347   prod = 1.0;
 348   total = 0.0;
 349   bindex = 0;
 350   x = st->start_x;
 351   r = (st->forcing[bindex]) ? st->b : st->a;
 352 #ifdef MAPS
 353   findex = 0;
 354   map = Maps[st->Forcing[findex]];
 355 #endif
 356   for (i=0;i<st->settle;i++) {     /* Here's where we let the thing */
 357     x = st->map (x, r);  /* "settle down". There is usually */
 358     if (++bindex >= st->maxindex) { /* some initial "noise" in the */
 359       bindex = 0;    /* iterations. How can we optimize */
 360       if (st->Rflag)      /* the value of settle ??? */
 361         setforcing(st);
 362     }
 363     r = (st->forcing[bindex]) ? st->b : st->a;
 364 #ifdef MAPS
 365     if (++findex >= funcmaxindex)
 366       findex = 0;
 367     map = Maps[st->Forcing[findex]];
 368 #endif
 369   }
 370 #ifdef MAPS
 371   deriv = Derivs[st->Forcing[findex]];
 372 #endif
 373   if (st->useprod) {      /* using log(a*b) */
 374     for (i=0;i<st->dwell;i++) {
 375       x = st->map (x, r);
 376       dx = st->deriv (x, r); /* ABS is a macro, so don't be fancy */
 377       dx = ABS(dx);
 378       if (dx == 0.0) /* log(0) is nasty so break out. */
 379         {
 380           i++;
 381           break;
 382         }
 383       prod *= dx;
 384       /* we need to prevent overflow and underflow */
 385       if ((prod > 1.0e12) || (prod < 1.0e-12)) {
 386         total += log(prod);
 387         prod = 1.0;
 388       }
 389       if (++bindex >= st->maxindex) {
 390         bindex = 0;
 391         if (st->Rflag)
 392           setforcing(st);
 393       }
 394       r = (st->forcing[bindex]) ? st->b : st->a;
 395 #ifdef MAPS
 396       if (++findex >= funcmaxindex)
 397         findex = 0;
 398       map = Maps[st->Forcing[findex]];
 399       deriv = Derivs[st->Forcing[findex]];
 400 #endif
 401     }
 402     total += log(prod);
 403     st->lyapunov = (total * M_LOG2E) / (double)i;
 404   }
 405   else {        /* use log(a) + log(b) */
 406     for (i=0;i<st->dwell;i++) {
 407       x = st->map (x, r);
 408       dx = st->deriv (x, r); /* ABS is a macro, so don't be fancy */
 409       dx = ABS(dx);
 410       if (x == 0.0)  /* log(0) check */
 411         {
 412           i++;
 413           break;
 414         }
 415       total += log(dx);
 416       if (++bindex >= st->maxindex) {
 417         bindex = 0;
 418         if (st->Rflag)
 419           setforcing(st);
 420       }
 421       r = (st->forcing[bindex]) ? st->b : st->a;
 422 #ifdef MAPS
 423       if (++findex >= funcmaxindex)
 424         findex = 0;
 425       map = Maps[st->Forcing[findex]];
 426       deriv = Derivs[st->Forcing[findex]];
 427 #endif
 428     }
 429     st->lyapunov = (total * M_LOG2E) / (double)i;
 430   }
 431
 432   if (sendpoint(st, st->lyapunov) == TRUE)
 433     return FALSE;
 434   else {
 435     FlushBuffer(st);
 436     /*    if (savefile)
 437           save_to_file();*/
 438     return TRUE;
 439   }
 440 }
 441
 442 static double
 443 logistic(double x, double r)        /* the familiar logistic map */
 444 {
 445   return(r * x * (1.0 - x));
 446 }
 447
 448 static double
 449 dlogistic(double x, double r)       /* the derivative of logistic map */
 450 {
 451   return(r - (2.0 * r * x));
 452 }
 453
 454 static double
 455 circle(double x, double r)        /* sin() hump or sorta like the circle map */
 456 {
 457   return(r * sin(M_PI * x));
 458 }
 459
 460 static double
 461 dcircle(double x, double r)       /* derivative of the "sin() hump" */
 462 {
 463   return(r * M_PI * cos(M_PI * x));
 464 }
 465
 466 static double
 467 leftlog(double x, double r)       /* left skewed logistic */
 468 {
 469   double d;
 470
 471   d = 1.0 - x;
 472   return(r * x * d * d);
 473 }
 474
 475 static double
 476 dleftlog(double x, double r)    /* derivative of the left skewed logistic */
 477 {
 478   return(r * (1.0 - (4.0 * x) + (3.0 * x * x)));
 479 }
 480
 481 static double
 482 rightlog(double x, double r)    /* right skewed logistic */
 483 {
 484   return(r * x * x * (1.0 - x));
 485 }
 486
 487 static double
 488 drightlog(double x, double r)    /* derivative of the right skewed logistic */
 489 {
 490   return(r * ((2.0 * x) - (3.0 * x * x)));
 491 }
 492
 493 static double
 494 doublelog(double x, double r)    /* double logistic */
 495 {
 496   double d;
 497
 498   d = 1.0 - x;
 499   return(r * x * x * d * d);
 500 }
 501
 502 static double
 503 ddoublelog(double x, double r)   /* derivative of the double logistic */
 504 {
 505   double d;
 506
 507   d = x * x;
 508   return(r * ((2.0 * x) - (6.0 * d) + (4.0 * x * d)));
 509 }
 510
 511 static void
 512 init_data(struct state *st)
 513 {
 514   st->numcolors = get_integer_resource (st->dpy, "colors", "Integer");
 515   if (st->numcolors < 2)
 516     st->numcolors = 2;
 517   if (st->numcolors > st->maxcolor)
 518     st->numcolors = st->maxcolor;
 519   st->numfreecols = st->numcolors - st->mincolindex;
 520   st->lowrange = st->mincolindex - st->startcolor;
 521   st->a_inc = st->a_range / (double)st->width;
 522   st->b_inc = st->b_range / (double)st->height;
 523   st->point.x = -1;
 524   st->point.y = 0;
 525   st->a = /*st->rubber_data.p_min = */st->min_a;
 526   st->b = /*st->rubber_data.q_min = */st->min_b;
 527 /*  st->rubber_data.p_max = st->max_a;
 528   st->rubber_data.q_max = st->max_b;*/
 529   if (st->show)
 530     show_defaults(st);
 531   InitBuffer(st);
 532 }
 533
 534 #if 0
 535 static void
 536 hls2rgb(int hue_light_sat[3],
 537         int rgb[3])             /*      Each in range [0..65535]        */
 538 {
 539   unsigned short r, g, b;
 540   hsv_to_rgb((int) (hue_light_sat[0] / 10),             /* 0-3600 -> 0-360 */
 541              (int) ((hue_light_sat[2]/1000.0) * 64435), /* 0-1000 -> 0-65535 */
 542              (int) ((hue_light_sat[1]/1000.0) * 64435), /* 0-1000 -> 0-65535 */
 543              &r, &g, &b);
 544   rgb[0] = r;
 545   rgb[1] = g;
 546   rgb[2] = b;
 547 }
 548 #endif /* 0 */
 549
 550
 551 static void
 552 init_color(struct state *st)
 553 {
 554   int i;
 555   if (st->ncolors)
 556     free_colors (st->screen, st->cmap, st->colors, st->ncolors);
 557   st->ncolors = st->maxcolor;
 558   make_smooth_colormap(st->screen, st->visual, st->cmap,
 559                        st->colors, &st->ncolors, True, NULL, True);
 560
 561   for (i = 0; i < st->maxcolor; i++) {
 562     if (! st->Data_GC[i]) {
 563       XGCValues gcv;
 564       gcv.background = BlackPixelOfScreen(st->screen);
 565       st->Data_GC[i] = XCreateGC(st->dpy, st->canvas, GCBackground, &gcv);
 566     }
 567     XSetForeground(st->dpy, st->Data_GC[i],
 568                    st->colors[((int) ((i / ((float)st->maxcolor)) *
 569                                       st->ncolors))].pixel);
 570   }
 571 }
 572
 573
 574 static void
 575 parseargs(struct state *st)
 576 {
 577   int i;
 578   int bindex=0, findex;
 579   char *s, *ch;
 580
 581   st->map = Maps[0];
 582   st->deriv = Derivs[0];
 583   st->maxexp=st->minlyap; st->minexp= -1.0 * st->minlyap;
 584
 585   st->mincolindex = get_integer_resource(st->dpy, "minColor", "Integer");
 586   st->dwell = get_integer_resource(st->dpy, "dwell", "Integer");
 587 #ifdef MAPS
 588   {
 589     char *optarg = get_string_resource(st->dpy, "function", "String");
 590     funcmaxindex = strlen(optarg);
 591     if (funcmaxindex > FUNCMAXINDEX)
 592       usage();
 593     ch = optarg;
 594     st->Force++;
 595     for (findex=0;findex<funcmaxindex;findex++) {
 596       st->Forcing[findex] = (int)(*ch++ - '0');;
 597       if (st->Forcing[findex] >= NUMMAPS)
 598         usage();
 599     }
 600     if (optarg) free (optarg);
 601   }
 602 #endif
 603   if (get_boolean_resource(st->dpy, "useLog", "Boolean"))
 604     st->useprod=0;
 605
 606   st->minlyap=ABS(get_float_resource(st->dpy, "colorExponent", "Float"));
 607   st->maxexp=st->minlyap;
 608   st->minexp= -1.0 * st->minlyap;
 609
 610   st->color_offset = get_integer_resource(st->dpy, "colorOffset", "Integer");
 611
 612   st->maxcolor=ABS(get_integer_resource(st->dpy, "maxColor", "Integer"));
 613   if ((st->maxcolor - st->startcolor) <= 0)
 614     st->startcolor = get_pixel_resource(st->dpy, st->cmap,
 615                                         "background", "Background");
 616   if ((st->maxcolor - st->mincolindex) <= 0) {
 617     st->mincolindex = 1;
 618     st->color_offset = 0;
 619   }
 620
 621   s = get_string_resource(st->dpy, "randomForce", "Float");
 622   if (s && *s) {
 623     st->prob=atof(s); st->Rflag++; setforcing(st);
 624   }
 625   if (s) free (s);
 626
 627   st->settle = get_integer_resource(st->dpy, "settle", "Integer");
 628
 629 #if 0
 630   s = get_string_resource(st->dpy, "minA", "Float");
 631   if (s && *s) {
 632     st->min_a = atof(s);
 633     st->aflag++;
 634   }
 635   if (s) free (s);
 636
 637   s = get_string_resource(st->dpy, "minB", "Float");
 638   if (s && *s) {
 639     st->min_b=atof(s); st->bflag++;
 640   }
 641   if (s) free (s);
 642 #else
 643   st->min_a = get_float_resource (st->dpy, "minA", "Float");
 644   st->aflag++;
 645   st->min_b = get_float_resource (st->dpy, "minB", "Float");
 646   st->bflag++;
 647 #endif
 648
 649
 650   st->numwheels = get_integer_resource(st->dpy, "wheels", "Integer");
 651
 652   s = get_string_resource(st->dpy, "forcingFunction", "String");
 653   if (s && *s) {
 654     st->maxindex = strlen(s);
 655     if (st->maxindex > MAXINDEX)
 656       usage(st);
 657     ch = s;
 658     st->force++;
 659     while (bindex < st->maxindex) {
 660       if (*ch == 'a')
 661         st->forcing[bindex++] = 0;
 662       else if (*ch == 'b')
 663         st->forcing[bindex++] = 1;
 664       else
 665         usage(st);
 666       ch++;
 667     }
 668   }
 669   if (s) free (s);
 670
 671   s = get_string_resource(st->dpy, "bRange", "Float");
 672   if (s && *s) {
 673     st->b_range = atof(s);
 674     st->hflag++;
 675   }
 676   if (s) free (s);
 677
 678   st->start_x = get_float_resource(st->dpy, "startX", "Float");
 679
 680   s = get_string_resource(st->dpy, "mapIndex", "Integer");
 681   if (s && *s) {
 682     st->mapindex=atoi(s);
 683     if ((st->mapindex >= NUMMAPS) || (st->mapindex < 0))
 684       usage(st);
 685     st->map = Maps[st->mapindex];
 686     st->deriv = Derivs[st->mapindex];
 687     if (!st->aflag)
 688       st->min_a = amins[st->mapindex];
 689     if (!st->wflag)
 690       st->a_range = aranges[st->mapindex];
 691     if (!st->bflag)
 692       st->min_b = bmins[st->mapindex];
 693     if (!st->hflag)
 694       st->b_range = branges[st->mapindex];
 695     if (!st->Force)
 696       for (i=0;i<FUNCMAXINDEX;i++)
 697         st->Forcing[i] = st->mapindex;
 698   }
 699   if (s) free (s);
 700
 701   st->outname = get_string_resource(st->dpy, "outputFile", "Integer");
 702
 703   if (get_boolean_resource(st->dpy, "beNegative", "Boolean"))
 704     st->negative--;
 705
 706   st->rgb_max = get_integer_resource(st->dpy, "rgbMax", "Integer");
 707   st->spinlength = get_integer_resource(st->dpy, "spinLength", "Integer");
 708   st->show = get_boolean_resource(st->dpy, "show", "Boolean");
 709
 710   s = get_string_resource(st->dpy, "aRange", "Float");
 711   if (s && *s) {
 712     st->a_range = atof(s); st->wflag++;
 713   }
 714   if (s) free (s);
 715
 716   st->max_a = st->min_a + st->a_range;
 717   st->max_b = st->min_b + st->b_range;
 718
 719   st->a_minimums[0] = st->min_a; st->b_minimums[0] = st->min_b;
 720   st->a_maximums[0] = st->max_a; st->b_maximums[0] = st->max_b;
 721
 722   if (st->Force)
 723     if (st->maxindex == st->funcmaxindex)
 724       for (findex=0;findex<st->funcmaxindex;findex++)
 725         check_params(st, st->Forcing[findex],st->forcing[findex]);
 726     else
 727       fprintf(stderr, "Warning! Unable to check parameters\n");
 728   else
 729     check_params(st, st->mapindex,2);
 730 }
 731
 732 static void
 733 check_params(struct state *st, int mapnum, int parnum)
 734 {
 735
 736   if (parnum != 1) {
 737     if ((st->max_a > pmaxs[mapnum]) || (st->min_a < pmins[mapnum])) {
 738       fprintf(stderr, "Warning! Parameter 'a' out of range.\n");
 739       fprintf(stderr, "You have requested a range of (%f,%f).\n",
 740               st->min_a,st->max_a);
 741       fprintf(stderr, "Valid range is (%f,%f).\n",
 742               pmins[mapnum],pmaxs[mapnum]);
 743     }
 744   }
 745   if (parnum != 0) {
 746     if ((st->max_b > pmaxs[mapnum]) || (st->min_b < pmins[mapnum])) {
 747       fprintf(stderr, "Warning! Parameter 'b' out of range.\n");
 748       fprintf(stderr, "You have requested a range of (%f,%f).\n",
 749               st->min_b,st->max_b);
 750       fprintf(stderr, "Valid range is (%f,%f).\n",
 751               pmins[mapnum],pmaxs[mapnum]);
 752     }
 753   }
 754 }
 755
 756 static void
 757 usage(struct state *st)
 758 {
 759   fprintf(stderr,"lyap [-BLs][-W#][-H#][-a#][-b#][-w#][-h#][-x xstart]\n");
 760   fprintf(stderr,"\t[-M#][-S#][-D#][-f string][-r#][-O#][-C#][-c#][-m#]\n");
 761 #ifdef MAPS
 762   fprintf(stderr,"\t[-F string]\n");
 763 #endif
 764   fprintf(stderr,"\tWhere: -C# specifies the minimum color index\n");
 765   fprintf(stderr,"\t       -r# specifies the maxzimum rgb value\n");
 766   fprintf(stderr,"\t       -u displays this message\n");
 767   fprintf(stderr,"\t       -a# specifies the minimum horizontal parameter\n");
 768   fprintf(stderr,"\t       -b# specifies the minimum vertical parameter\n");
 769   fprintf(stderr,"\t       -w# specifies the horizontal parameter range\n");
 770   fprintf(stderr,"\t       -h# specifies the vertical parameter range\n");
 771   fprintf(stderr,"\t       -D# specifies the dwell\n");
 772   fprintf(stderr,"\t       -S# specifies the settle\n");
 773   fprintf(stderr,"\t       -H# specifies the initial window height\n");
 774   fprintf(stderr,"\t       -W# specifies the initial window width\n");
 775   fprintf(stderr,"\t       -O# specifies the color offset\n");
 776   fprintf(stderr,"\t       -c# specifies the desired color wheel\n");
 777   fprintf(stderr,"\t       -m# specifies the desired map (0-4)\n");
 778   fprintf(stderr,"\t       -f aabbb specifies a forcing function of 00111\n");
 779 #ifdef MAPS
 780   fprintf(stderr,"\t       -F 00111 specifies the function forcing function\n");
 781 #endif
 782   fprintf(stderr,"\t       -L indicates use log(x)+log(y) rather than log(xy)\n");
 783   fprintf(stderr,"\tDuring display :\n");
 784   fprintf(stderr,"\t       Use the mouse to zoom in on an area\n");
 785   fprintf(stderr,"\t       e or E recalculates color indices\n");
 786   fprintf(stderr,"\t       f or F saves exponents to a file\n");
 787   fprintf(stderr,"\t       KJmn increase/decrease minimum negative exponent\n");
 788   fprintf(stderr,"\t       r or R redraws\n");
 789   fprintf(stderr,"\t       s or S spins the colorwheel\n");
 790   fprintf(stderr,"\t       w or W changes the color wheel\n");
 791   fprintf(stderr,"\t       x or X clears the window\n");
 792   fprintf(stderr,"\t       q or Q exits\n");
 793   exit(1);
 794 }
 795
 796 static void
 797 Cycle_frames(struct state *st)
 798 {
 799   int i;
 800   for (i=0;i<=st->maxframe;i++)
 801     redraw(st, st->exponents[i], st->expind[i], 1);
 802 }
 803
 804 #if 0
 805 static void
 806 Spin(struct state *st)
 807 {
 808   int i, j;
 809   long tmpxcolor;
 810
 811   if (!mono_p) {
 812     for (j=0;j<st->spinlength;j++) {
 813       tmpxcolor = st->Colors[st->mincolindex].pixel;
 814       for (i=st->mincolindex;i<st->numcolors-1;i++)
 815         st->Colors[i].pixel = st->Colors[i+1].pixel;
 816       st->Colors[st->numcolors-1].pixel = tmpxcolor;
 817       XStoreColors(st->dpy, st->cmap, st->Colors, st->numcolors);
 818     }
 819     for (j=0;j<st->spinlength;j++) {
 820       tmpxcolor = st->Colors[st->numcolors-1].pixel;
 821       for (i=st->numcolors-1;i>st->mincolindex;i--)
 822         st->Colors[i].pixel = st->Colors[i-1].pixel;
 823       st->Colors[st->mincolindex].pixel = tmpxcolor;
 824       XStoreColors(st->dpy, st->cmap, st->Colors, st->numcolors);
 825     }
 826   }
 827 }
 828 #endif
 829
 830 static Bool
 831 Getkey(struct state *st, XKeyEvent *event)
 832 {
 833   unsigned char key;
 834   int i;
 835   if (XLookupString(event, (char *)&key, sizeof(key), (KeySym *)0,
 836                     (XComposeStatus *) 0) > 0) {
 837
 838     if (st->reset_countdown)
 839       st->reset_countdown = st->linger;
 840
 841     switch (key) {
 842     case '<': st->dwell /= 2; if (st->dwell < 1) st->dwell = 1; return True;
 843     case '>': st->dwell *= 2; return True;
 844     case '[': st->settle /= 2; if (st->settle < 1) st->settle = 1; return True;
 845     case ']': st->settle *= 2; return True;
 846     case 'd': go_down(st); return True;
 847     case 'D': FlushBuffer(st); return True;
 848     case 'e':
 849     case 'E': FlushBuffer(st);
 850       st->dorecalc = (!st->dorecalc);
 851       if (st->dorecalc)
 852         recalc(st);
 853       else {
 854         st->maxexp = st->minlyap; st->minexp = -1.0 * st->minlyap;
 855       }
 856       redraw(st, st->exponents[st->frame], st->expind[st->frame], 1);
 857       return True;
 858     case 'f':
 859       /*  case 'F': save_to_file(); return True;*/
 860     case 'i': if (st->stripe_interval > 0) {
 861         st->stripe_interval--;
 862         if (!mono_p) {
 863           init_color(st);
 864         }
 865       }
 866       return True;
 867     case 'I': st->stripe_interval++;
 868       if (!mono_p) {
 869         init_color(st);
 870       }
 871       return True;
 872     case 'K': if (st->minlyap > 0.05)
 873         st->minlyap -= 0.05;
 874       return True;
 875     case 'J': st->minlyap += 0.05;
 876       return True;
 877     case 'm': st->mapindex++;
 878       if (st->mapindex >= NUMMAPS)
 879         st->mapindex=0;
 880       st->map = Maps[st->mapindex];
 881       st->deriv = Derivs[st->mapindex];
 882       if (!st->aflag)
 883         st->min_a = amins[st->mapindex];
 884       if (!st->wflag)
 885         st->a_range = aranges[st->mapindex];
 886       if (!st->bflag)
 887         st->min_b = bmins[st->mapindex];
 888       if (!st->hflag)
 889         st->b_range = branges[st->mapindex];
 890       if (!st->Force)
 891         for (i=0;i<FUNCMAXINDEX;i++)
 892           st->Forcing[i] = st->mapindex;
 893       st->max_a = st->min_a + st->a_range;
 894       st->max_b = st->min_b + st->b_range;
 895       st->a_minimums[0] = st->min_a; st->b_minimums[0] = st->min_b;
 896       st->a_maximums[0] = st->max_a; st->b_maximums[0] = st->max_b;
 897       st->a_inc = st->a_range / (double)st->width;
 898       st->b_inc = st->b_range / (double)st->height;
 899       st->point.x = -1;
 900       st->point.y = 0;
 901       st->a = /*st->rubber_data.p_min = */st->min_a;
 902       st->b = /*st->rubber_data.q_min = */st->min_b;
 903 /*      st->rubber_data.p_max = st->max_a;
 904       st->rubber_data.q_max = st->max_b;*/
 905       Clear(st);
 906       return True;
 907     case 'M': if (st->minlyap > 0.005)
 908         st->minlyap -= 0.005;
 909       return True;
 910     case 'N': st->minlyap += 0.005;
 911       return True;
 912     case 'p':
 913     case 'P': st->negative = (!st->negative);
 914       FlushBuffer(st); redraw(st, st->exponents[st->frame],
 915                               st->expind[st->frame], 1);
 916       return True;
 917     case 'r': FlushBuffer(st); redraw(st, st->exponents[st->frame],
 918                                       st->expind[st->frame], 1);
 919       return True;
 920     case 'R': FlushBuffer(st); Redraw(st); return True;
 921     case 's':
 922       st->spinlength=st->spinlength/2;
 923 #if 0
 924     case 'S': if (!mono_p)
 925         Spin(st);
 926       st->spinlength=st->spinlength*2; return True;
 927 #endif
 928     case 'u': go_back(st); return True;
 929     case 'U': go_init(st); return True;
 930     case 'v':
 931     case 'V': print_values(st); return True;
 932     case 'W': if (st->numwheels < MAXWHEELS)
 933         st->numwheels++;
 934       else
 935         st->numwheels = 0;
 936       if (!mono_p) {
 937         init_color(st);
 938       }
 939       return True;
 940     case 'w': if (st->numwheels > 0)
 941         st->numwheels--;
 942       else
 943         st->numwheels = MAXWHEELS;
 944       if (!mono_p) {
 945         init_color(st);
 946       }
 947       return True;
 948     case 'x': Clear(st); return True;
 949     case 'X': Destroy_frame(st); return True;
 950     case 'z': Cycle_frames(st); redraw(st, st->exponents[st->frame],
 951                                        st->expind[st->frame], 1);
 952       return True;
 953 #if 0
 954     case 'Z': while (!XPending(st->dpy)) Cycle_frames(st);
 955       redraw(st, st->exponents[st->frame], st->expind[st->frame], 1);
 956       return True;
 957 #endif
 958     case 'q':
 959     case 'Q': exit(0); return True;
 960     case '?':
 961     case 'h':
 962     case 'H': print_help(st); return True;
 963     default:  return False;
 964     }
 965   }
 966
 967   return False;
 968 }
 969
 970 /* Here's where we index into a color map. After the Lyapunov exponent is
 971  * calculated, it is used to determine what color to use for that point.  I
 972  * suppose there are a lot of ways to do this. I used the following : if it's
 973  * non-negative then there's a reserved area at the lower range of the color
 974  * map that i index into. The ratio of some "minimum exponent value" and the
 975  * calculated value is used as a ratio of how high to index into this reserved
 976  * range. Usually these colors are dark red (see init_color).  If the exponent
 977  * is negative, the same ratio (expo/minlyap) is used to index into the
 978  * remaining portion of the colormap (which is usually some light shades of
 979  * color or a rainbow wheel). The coloring scheme can actually make a great
 980  * deal of difference in the quality of the picture. Different colormaps bring
 981  * out different details of the dynamics while different indexing algorithms
 982  * also greatly effect what details are seen. Play around with this.
 983  */
 984 static int
 985 sendpoint(struct state *st, double expo)
 986 {
 987   double tmpexpo;
 988
 989   if (st->maxcolor > MAXCOLOR)
 990     abort();
 991
 992 #if 0
 993   /* The relationship st->minexp <= expo <= maxexp should always be true. This
 994      test enforces that. But maybe not enforcing it makes better pictures. */
 995   if (expo < st->minexp)
 996     expo = st->minexp;
 997   else if (expo > maxexp)
 998     expo = maxexp;
 999 #endif
1000
1001   st->point.x++;
1002   tmpexpo = (st->negative) ? expo : -1.0 * expo;
1003   if (tmpexpo > 0) {
1004     if (!mono_p) {
1005       st->sendpoint_index = (int)(tmpexpo*st->lowrange/st->maxexp);
1006       st->sendpoint_index = ((st->sendpoint_index % st->lowrange) +
1007                              st->startcolor);
1008     }
1009     else
1010       st->sendpoint_index = 0;
1011   }
1012   else {
1013     if (!mono_p) {
1014       st->sendpoint_index = (int)(tmpexpo*st->numfreecols/st->minexp);
1015       st->sendpoint_index = ((st->sendpoint_index % st->numfreecols)
1016                              + st->mincolindex);
1017     }
1018     else
1019       st->sendpoint_index = 1;
1020   }
1021   BufferPoint(st, st->sendpoint_index, st->point.x, st->point.y);
1022   if (st->save) {
1023     if (st->frame > MAXFRAMES)
1024       abort();
1025     st->exponents[st->frame][st->expind[st->frame]++] = expo;
1026   }
1027   if (st->point.x >= st->width) {
1028     st->point.y++;
1029     st->point.x = 0;
1030     if (st->save) {
1031       st->b += st->b_inc;
1032       st->a = st->min_a;
1033     }
1034     if (st->point.y >= st->height)
1035       return FALSE;
1036     else
1037       return TRUE;
1038   }
1039   return TRUE;
1040 }
1041
1042
1043 static void
1044 resize(struct state *st)
1045 {
1046   Window r;
1047   int n, x, y;
1048   unsigned int bw, d, new_w, new_h;
1049
1050   XGetGeometry(st->dpy,st->canvas,&r,&x,&y,&new_w,&new_h,&bw,&d);
1051   if ((new_w == st->width) && (new_h == st->height))
1052     return;
1053   st->width = new_w; st->height = new_h;
1054   XClearWindow(st->dpy, st->canvas);
1055 #ifdef BACKING_PIXMAP
1056   if (st->pixmap)
1057     XFreePixmap(st->dpy, st->pixmap);
1058   st->pixmap = XCreatePixmap(st->dpy, st->canvas, st->width, st->height, d);
1059 #endif
1060   st->a_inc = st->a_range / (double)st->width;
1061   st->b_inc = st->b_range / (double)st->height;
1062   st->point.x = -1;
1063   st->point.y = 0;
1064   st->run = 1;
1065   st->a = /*st->rubber_data.p_min = */st->min_a;
1066   st->b = /*st->rubber_data.q_min = */st->min_b;
1067 /*  st->rubber_data.p_max = st->max_a;
1068   st->rubber_data.q_max = st->max_b;*/
1069   freemem(st);
1070   setupmem(st);
1071   for (n=0;n<MAXFRAMES;n++)
1072     if ((n <= st->maxframe) && (n != st->frame))
1073       st->resized[n] = 1;
1074   InitBuffer(st);
1075   Clear(st);
1076   Redraw(st);
1077 }
1078
1079 static void
1080 redraw(struct state *st, double *exparray, int index, int cont)
1081 {
1082   int i, x_sav, y_sav;
1083
1084   x_sav = st->point.x;
1085   y_sav = st->point.y;
1086
1087   st->point.x = -1;
1088   st->point.y = 0;
1089
1090   st->save=0;
1091   for (i=0;i<index;i++)
1092     sendpoint(st, exparray[i]);
1093   st->save=1;
1094
1095   if (cont) {
1096     st->point.x = x_sav;
1097     st->point.y = y_sav;
1098   }
1099   else {
1100     st->a = st->point.x * st->a_inc + st->min_a;
1101     st->b = st->point.y * st->b_inc + st->min_b;
1102   }
1103   FlushBuffer(st);
1104 }
1105
1106 static void
1107 Redraw(struct state *st)
1108 {
1109   FlushBuffer(st);
1110   st->point.x = -1;
1111   st->point.y = 0;
1112   st->run = 1;
1113   st->a = st->min_a;
1114   st->b = st->min_b;
1115   st->expind[st->frame] = 0;
1116   st->resized[st->frame] = 0;
1117 }
1118
1119 static void
1120 recalc(struct state *st)
1121 {
1122   int i;
1123
1124   st->minexp = st->maxexp = 0.0;
1125   for (i=0;i<st->expind[st->frame];i++) {
1126     if (st->exponents[st->frame][i] < st->minexp)
1127       st->minexp = st->exponents[st->frame][i];
1128     if (st->exponents[st->frame][i] > st->maxexp)
1129       st->maxexp = st->exponents[st->frame][i];
1130   }
1131 }
1132
1133 static void
1134 Clear(struct state *st)
1135 {
1136   XClearWindow(st->dpy, st->canvas);
1137 #ifdef BACKING_PIXMAP
1138   XCopyArea(st->dpy, st->canvas, st->pixmap, st->Data_GC[0],
1139             0, 0, st->width, st->height, 0, 0);
1140 #endif
1141   InitBuffer(st);
1142 }
1143
1144 static void
1145 show_defaults(struct state *st)
1146 {
1147
1148   printf("Width=%d  Height=%d  numcolors=%d  settle=%d  dwell=%d\n",
1149          st->width,st->height,st->numcolors,st->settle,st->dwell);
1150   printf("min_a=%f  a_range=%f  max_a=%f\n", st->min_a,st->a_range,st->max_a);
1151   printf("min_b=%f  b_range=%f  max_b=%f\n", st->min_b,st->b_range,st->max_b);
1152   printf("minlyap=%f  minexp=%f  maxexp=%f\n", st->minlyap,st->minexp,
1153          st->maxexp);
1154   exit(0);
1155 }
1156
1157 #if 0
1158 static void
1159 CreateXorGC(struct state *st)
1160 {
1161   XGCValues values;
1162
1163   values.foreground = st->foreground;
1164   values.function = GXxor;
1165   st->RubberGC = XCreateGC(st->dpy, st->canvas,
1166                            GCForeground | GCFunction, &values);
1167 }
1168
1169 static void
1170 StartRubberBand(struct state *st, image_data_t *data, XEvent *event)
1171 {
1172   XPoint corners[5];
1173
1174   st->nostart = 0;
1175   data->rubber_band.last_x = data->rubber_band.start_x = event->xbutton.x;
1176   data->rubber_band.last_y = data->rubber_band.start_y = event->xbutton.y;
1177   SetupCorners(corners, data);
1178   XDrawLines(st->dpy, st->canvas, st->RubberGC,
1179              corners, sizeof(corners) / sizeof(corners[0]), CoordModeOrigin);
1180 }
1181
1182 static void
1183 SetupCorners(XPoint *corners, image_data_t *data)
1184 {
1185   corners[0].x = data->rubber_band.start_x;
1186   corners[0].y = data->rubber_band.start_y;
1187   corners[1].x = data->rubber_band.start_x;
1188   corners[1].y = data->rubber_band.last_y;
1189   corners[2].x = data->rubber_band.last_x;
1190   corners[2].y = data->rubber_band.last_y;
1191   corners[3].x = data->rubber_band.last_x;
1192   corners[3].y = data->rubber_band.start_y;
1193   corners[4] = corners[0];
1194 }
1195
1196 static void
1197 TrackRubberBand(struct state *st, image_data_t *data, XEvent *event)
1198 {
1199   XPoint corners[5];
1200   int xdiff, ydiff;
1201
1202   if (st->nostart)
1203     return;
1204   SetupCorners(corners, data);
1205   XDrawLines(st->dpy, st->canvas, st->RubberGC,
1206              corners, sizeof(corners) / sizeof(corners[0]), CoordModeOrigin);
1207   ydiff = event->xbutton.y - data->rubber_band.start_y;
1208   xdiff = event->xbutton.x - data->rubber_band.start_x;
1209   data->rubber_band.last_x = data->rubber_band.start_x + xdiff;
1210   data->rubber_band.last_y = data->rubber_band.start_y + ydiff;
1211   if (data->rubber_band.last_y < data->rubber_band.start_y ||
1212       data->rubber_band.last_x < data->rubber_band.start_x)
1213     {
1214       data->rubber_band.last_y = data->rubber_band.start_y;
1215       data->rubber_band.last_x = data->rubber_band.start_x;
1216     }
1217   SetupCorners(corners, data);
1218   XDrawLines(st->dpy, st->canvas, st->RubberGC,
1219              corners, sizeof(corners) / sizeof(corners[0]), CoordModeOrigin);
1220 }
1221
1222 static void
1223 EndRubberBand(struct state *st, image_data_t *data, XEvent *event)
1224 {
1225   XPoint corners[5];
1226   XPoint top, bot;
1227   double delta, diff;
1228
1229   st->nostart = 1;
1230   SetupCorners(corners, data);
1231   XDrawLines(st->dpy, st->canvas, st->RubberGC,
1232              corners, sizeof(corners) / sizeof(corners[0]), CoordModeOrigin);
1233   if (data->rubber_band.start_x >= data->rubber_band.last_x ||
1234       data->rubber_band.start_y >= data->rubber_band.last_y)
1235     return;
1236   top.x = data->rubber_band.start_x;
1237   bot.x = data->rubber_band.last_x;
1238   top.y = data->rubber_band.start_y;
1239   bot.y = data->rubber_band.last_y;
1240   diff = data->q_max - data->q_min;
1241   delta = (double)top.y / (double)st->height;
1242   data->q_min += diff * delta;
1243   delta = (double)(st->height - bot.y) / (double)st->height;
1244   data->q_max -= diff * delta;
1245   diff = data->p_max - data->p_min;
1246   delta = (double)top.x / (double)st->width;
1247   data->p_min += diff * delta;
1248   delta = (double)(st->width - bot.x) / (double)st->width;
1249   data->p_max -= diff * delta;
1250   set_new_params(st, data);
1251 }
1252
1253 static void
1254 set_new_params(struct state *st, image_data_t *data)
1255 {
1256   st->frame = (st->maxframe + 1) % MAXFRAMES;
1257   if (st->frame > st->maxframe)
1258     st->maxframe = st->frame;
1259   st->a_range = data->p_max - data->p_min;
1260   st->b_range = data->q_max - data->q_min;
1261   st->a_minimums[st->frame] = st->min_a = data->p_min;
1262   st->b_minimums[st->frame] = st->min_b = data->q_min;
1263   st->a_inc = st->a_range / (double)st->width;
1264   st->b_inc = st->b_range / (double)st->height;
1265   st->point.x = -1;
1266   st->point.y = 0;
1267   st->run = 1;
1268   st->a = st->min_a;
1269   st->b = st->min_b;
1270   st->a_maximums[st->frame] = st->max_a = data->p_max;
1271   st->b_maximums[st->frame] = st->max_b = data->q_max;
1272   st->expind[st->frame] = 0;
1273   Clear(st);
1274 }
1275 #endif
1276
1277 static void
1278 go_down(struct state *st)
1279 {
1280   st->frame++;
1281   if (st->frame > st->maxframe)
1282     st->frame = 0;
1283   jumpwin(st);
1284 }
1285
1286 static void
1287 go_back(struct state *st)
1288 {
1289   st->frame--;
1290   if (st->frame < 0)
1291     st->frame = st->maxframe;
1292   jumpwin(st);
1293 }
1294
1295 static void
1296 jumpwin(struct state *st)
1297 {
1298   /*st->rubber_data.p_min =*/ st->min_a = st->a_minimums[st->frame];
1299   /*st->rubber_data.q_min =*/ st->min_b = st->b_minimums[st->frame];
1300   /*st->rubber_data.p_max =*/ st->max_a = st->a_maximums[st->frame];
1301   /*st->rubber_data.q_max =*/ st->max_b = st->b_maximums[st->frame];
1302   st->a_range = st->max_a - st->min_a;
1303   st->b_range = st->max_b - st->min_b;
1304   st->a_inc = st->a_range / (double)st->width;
1305   st->b_inc = st->b_range / (double)st->height;
1306   st->point.x = -1;
1307   st->point.y = 0;
1308   st->a = st->min_a;
1309   st->b = st->min_b;
1310   Clear(st);
1311   if (st->resized[st->frame])
1312     Redraw(st);
1313   else
1314     redraw(st, st->exponents[st->frame], st->expind[st->frame], 0);
1315 }
1316
1317 static void
1318 go_init(struct state *st)
1319 {
1320   st->frame = 0;
1321   jumpwin(st);
1322 }
1323
1324 static void
1325 Destroy_frame(struct state *st)
1326 {
1327   int i;
1328
1329   for (i=st->frame; i<st->maxframe; i++) {
1330     st->exponents[st->frame] = st->exponents[st->frame+1];
1331     st->expind[st->frame] = st->expind[st->frame+1];
1332     st->a_minimums[st->frame] = st->a_minimums[st->frame+1];
1333     st->b_minimums[st->frame] = st->b_minimums[st->frame+1];
1334     st->a_maximums[st->frame] = st->a_maximums[st->frame+1];
1335     st->b_maximums[st->frame] = st->b_maximums[st->frame+1];
1336   }
1337   st->maxframe--;
1338   go_back(st);
1339 }
1340
1341 static void
1342 InitBuffer(struct state *st)
1343 {
1344   int i;
1345
1346   for (i = 0 ; i < st->maxcolor; ++i)
1347     st->Points.npoints[i] = 0;
1348 }
1349
1350 static void
1351 BufferPoint(struct state *st, int color, int x, int y)
1352 {
1353   if (st->maxcolor > MAXCOLOR)
1354     abort();
1355
1356   /* Guard against bogus color values. Shouldn't be necessary but paranoia
1357      is good. */
1358   if (color < 0)
1359     color = 0;
1360   else if (color >= st->maxcolor)
1361     color = st->maxcolor - 1;
1362
1363   if (st->Points.npoints[color] == MAXPOINTS)
1364     {
1365       XDrawPoints(st->dpy, st->canvas, st->Data_GC[color],
1366                   st->Points.data[color], st->Points.npoints[color],
1367                   CoordModeOrigin);
1368 #ifdef BACKING_PIXMAP
1369       XDrawPoints(st->dpy, st->pixmap, st->Data_GC[color],
1370                   st->Points.data[color], st->Points.npoints[color],
1371                   CoordModeOrigin);
1372 #endif
1373       st->Points.npoints[color] = 0;
1374     }
1375   st->Points.data[color][st->Points.npoints[color]].x = x;
1376   st->Points.data[color][st->Points.npoints[color]].y = y;
1377   ++st->Points.npoints[color];
1378 }
1379
1380 static void
1381 FlushBuffer(struct state *st)
1382 {
1383   int color;
1384
1385   for (color = 0; color < st->maxcolor; ++color)
1386     if (st->Points.npoints[color])
1387       {
1388         XDrawPoints(st->dpy, st->canvas, st->Data_GC[color],
1389                     st->Points.data[color], st->Points.npoints[color],
1390                     CoordModeOrigin);
1391 #ifdef BACKING_PIXMAP
1392         XDrawPoints(st->dpy, st->pixmap, st->Data_GC[color],
1393                     st->Points.data[color], st->Points.npoints[color],
1394                     CoordModeOrigin);
1395 #endif
1396         st->Points.npoints[color] = 0;
1397       }
1398 }
1399
1400 static void
1401 print_help(struct state *st)
1402 {
1403   printf("During run-time, interactive control can be exerted via : \n");
1404   printf("Mouse buttons allow rubber-banding of a zoom box\n");
1405   printf("< halves the 'dwell', > doubles the 'dwell'\n");
1406   printf("[ halves the 'settle', ] doubles the 'settle'\n");
1407   printf("D flushes the drawing buffer\n");
1408   printf("e or E recalculates color indices\n");
1409   printf("f or F saves exponents to a file\n");
1410   printf("h or H or ? displays this message\n");
1411   printf("i decrements, I increments the stripe interval\n");
1412   printf("KJMN increase/decrease minimum negative exponent\n");
1413   printf("m increments the map index, changing maps\n");
1414   printf("p or P reverses the colormap for negative/positive exponents\n");
1415   printf("r redraws without recalculating\n");
1416   printf("R redraws, recalculating with new dwell and settle values\n");
1417   printf("s or S spins the colorwheel\n");
1418   printf("u pops back up to the last zoom\n");
1419   printf("U pops back up to the first picture\n");
1420   printf("v or V displays the values of various settings\n");
1421   printf("w decrements, W increments the color wheel index\n");
1422   printf("x or X clears the window\n");
1423   printf("q or Q exits\n");
1424 }
1425
1426 static void
1427 print_values(struct state *st)
1428 {
1429   int i;
1430   printf("\nminlyap=%f minexp=%f maxexp=%f\n",
1431          st->minlyap,st->minexp, st->maxexp);
1432   printf("width=%d height=%d\n",st->width,st->height);
1433   printf("settle=%d  dwell=%d st->start_x=%f\n",
1434          st->settle,st->dwell, st->start_x);
1435   printf("min_a=%f  a_rng=%f        max_a=%f\n",
1436          st->min_a,st->a_range,st->max_a);
1437   printf("min_b=%f  b_rng=%f        max_b=%f\n",
1438          st->min_b,st->b_range,st->max_b);
1439   if (st->Rflag)
1440     printf("pseudo-random forcing\n");
1441   else if (st->force) {
1442     printf("periodic forcing=");
1443     for (i=0;i<st->maxindex;i++)
1444       printf("%d",st->forcing[i]);
1445     printf("\n");
1446   }
1447   else
1448     printf("periodic forcing=01\n");
1449   if (st->Force) {
1450     printf("function forcing=");
1451     for (i=0;i<st->funcmaxindex;i++) {
1452       printf("%d",st->Forcing[i]);
1453     }
1454     printf("\n");
1455   }
1456   printf("numcolors=%d\n",st->numcolors-1);
1457 }
1458
1459 static void
1460 freemem(struct state *st)
1461 {
1462   int i;
1463   for (i=0;i<MAXFRAMES;i++)
1464     free(st->exponents[i]);
1465 }
1466
1467 static void
1468 setupmem(struct state *st)
1469 {
1470   int i;
1471   for (i=0;i<MAXFRAMES;i++) {
1472     if((st->exponents[i]=
1473         (double *)malloc(sizeof(double)*st->width*(st->height+1)))==NULL){
1474       fprintf(stderr,"Error malloc'ing exponent array.\n");
1475       exit(-1);
1476     }
1477   }
1478 }
1479
1480 static void
1481 setforcing(struct state *st)
1482 {
1483   int i;
1484   for (i=0;i<MAXINDEX;i++)
1485     st->forcing[i] = (random() > st->prob) ? 0 : 1;
1486 }
1487
1488 /****************************************************************************/
1489
1490 static void
1491 do_defaults (struct state *st)
1492 {
1493   int i;
1494
1495   memset (st->expind,  0, sizeof(st->expind));
1496   memset (st->resized, 0, sizeof(st->resized));
1497
1498   st->aflag = 0;
1499   st->bflag = 0;
1500   st->hflag = 0;
1501   st->wflag = 0;
1502   st->minexp = 0;
1503   st->mapindex = 0;
1504
1505 # ifdef SIXTEEN_COLORS
1506   st->maxcolor=16;
1507   st->startcolor=0;
1508   st->color_offset=0;
1509   st->mincolindex=1;
1510   st->dwell=50;
1511   st->settle=25;
1512   st->xposition=128;
1513   st->yposition=128;
1514 # else  /* !SIXTEEN_COLORS */
1515   st->maxcolor=256;
1516   st->startcolor=17;
1517   st->color_offset=96;
1518   st->mincolindex=33;
1519   st->dwell=100;
1520   st->settle=50;
1521 # endif /* !SIXTEEN_COLORS */
1522
1523   st->maxindex = MAXINDEX;
1524   st->funcmaxindex = FUNCMAXINDEX;
1525   st->min_a=2.0;
1526   st->min_b=2.0;
1527   st->a_range=2.0;
1528   st->b_range=2.0;
1529   st->minlyap=1.0;
1530   st->max_a=4.0;
1531   st->max_b=4.0;
1532   st->numcolors=16;
1533   st->prob=0.5;
1534   st->numwheels=MAXWHEELS;
1535   st->negative=1;
1536   st->rgb_max=65000;
1537   st->nostart=1;
1538   st->stripe_interval=7;
1539   st->save=1;
1540   st->useprod=1;
1541   st->spinlength=256;
1542   st->run=1;
1543
1544   for (i = 0; i < countof(st->forcing); i++)
1545     st->forcing[i] = (i & 1) ? 1 : 0;
1546 }
1547
1548 static void
1549 do_preset (struct state *st, int builtin)
1550 {
1551   char *ff = 0;
1552   switch (builtin) {
1553   case 0:
1554     st->min_a = 3.75; st->aflag++;
1555     st->min_b = 3.299999; st->bflag++;
1556     st->a_range = 0.05; st->wflag++;
1557     st->b_range = 0.05; st->hflag++;
1558     st->dwell = 200;
1559     st->settle = 100;
1560     ff = "abaabbaaabbb";
1561     break;
1562
1563   case 1:
1564     st->min_a = 3.8; st->aflag++;
1565     st->min_b = 3.2; st->bflag++;
1566     st->b_range = .05; st->hflag++;
1567     st->a_range = .05; st->wflag++;
1568     ff = "bbbbbaaaaa";
1569     break;
1570
1571   case 2:
1572     st->min_a =  3.4; st->aflag++;
1573     st->min_b =  3.04; st->bflag++;
1574     st->a_range =  .5; st->wflag++;
1575     st->b_range =  .5; st->hflag++;
1576     ff = "abbbbbbbbb";
1577     st->settle = 500;
1578     st->dwell = 1000;
1579     break;
1580
1581   case 3:
1582     st->min_a = 3.5; st->aflag++;
1583     st->min_b = 3.0; st->bflag++;
1584     st->a_range = 0.2; st->wflag++;
1585     st->b_range = 0.2; st->hflag++;
1586     st->dwell = 600;
1587     st->settle = 300;
1588     ff = "aaabbbab";
1589     break;
1590
1591   case 4:
1592     st->min_a = 3.55667; st->aflag++;
1593     st->min_b = 3.2; st->bflag++;
1594     st->b_range = .05; st->hflag++;
1595     st->a_range = .05; st->wflag++;
1596     ff = "bbbbbaaaaa";
1597     break;
1598
1599   case 5:
1600     st->min_a = 3.79; st->aflag++;
1601     st->min_b = 3.22; st->bflag++;
1602     st->b_range = .02999; st->hflag++;
1603     st->a_range = .02999; st->wflag++;
1604     ff = "bbbbbaaaaa";
1605     break;
1606
1607   case 6:
1608     st->min_a = 3.7999; st->aflag++;
1609     st->min_b = 3.299999; st->bflag++;
1610     st->a_range = 0.2; st->wflag++;
1611     st->b_range = 0.2; st->hflag++;
1612     st->dwell = 300;
1613     st->settle = 150;
1614     ff = "abaabbaaabbb";
1615     break;
1616
1617   case 7:
1618     st->min_a = 3.89; st->aflag++;
1619     st->min_b = 3.22; st->bflag++;
1620     st->b_range = .028; st->hflag++;
1621     st->a_range = .02999; st->wflag++;
1622     ff = "bbbbbaaaaa";
1623     st->settle = 600;
1624     st->dwell = 1000;
1625     break;
1626
1627   case 8:
1628     st->min_a = 3.2; st->aflag++;
1629     st->min_b = 3.7; st->bflag++;
1630     st->a_range = 0.05; st->wflag++;
1631     st->b_range = .005; st->hflag++;
1632     ff = "abbbbaa";
1633     break;
1634
1635   case 9:
1636     ff = "aaaaaabbbbbb";
1637     st->mapindex = 1;
1638     st->dwell =  400;
1639     st->settle =  200;
1640     st->minlyap = st->maxexp = ABS(-0.85);
1641     st->minexp = -1.0 * st->minlyap;
1642     break;
1643
1644   case 10:
1645     ff = "aaaaaabbbbbb";
1646     st->mapindex = 1;
1647     st->dwell =  400;
1648     st->settle = 200;
1649     st->minlyap = st->maxexp = ABS(-0.85);
1650     st->minexp = -1.0 * st->minlyap;
1651     break;
1652
1653   case 11:
1654     st->mapindex = 1;
1655     st->dwell =  400;
1656     st->settle = 200;
1657     st->minlyap = st->maxexp = ABS(-0.85);
1658     st->minexp = -1.0 * st->minlyap;
1659     break;
1660
1661   case 12:
1662     ff = "abbb";
1663     st->mapindex = 1;
1664     st->dwell =  400;
1665     st->settle = 200;
1666     st->minlyap = st->maxexp = ABS(-0.85);
1667     st->minexp = -1.0 * st->minlyap;
1668     break;
1669
1670   case 13:
1671     ff = "abbabaab";
1672     st->mapindex = 1;
1673     st->dwell =  400;
1674     st->settle = 200;
1675     st->minlyap = st->maxexp = ABS(-0.85);
1676     st->minexp = -1.0 * st->minlyap;
1677     break;
1678
1679   case 14:
1680     ff = "abbabaab";
1681     st->dwell =  800;
1682     st->settle = 200;
1683     st->minlyap = st->maxexp = ABS(-0.85);
1684     st->minexp = -1.0 * st->minlyap;
1685     /* ####  -x 0.05 */
1686     st->min_a = 3.91; st->aflag++;
1687     st->a_range =  0.0899999999; st->wflag++;
1688     st->min_b =  3.28; st->bflag++;
1689     st->b_range =  0.35; st->hflag++;
1690     break;
1691
1692   case 15:
1693     ff = "aaaaaabbbbbb";
1694     st->dwell =  400;
1695     st->settle = 200;
1696     st->minlyap = st->maxexp = ABS(-0.85);
1697     st->minexp = -1.0 * st->minlyap;
1698     break;
1699
1700   case 16:
1701     st->dwell =  400;
1702     st->settle = 200;
1703     st->minlyap = st->maxexp = ABS(-0.85);
1704     st->minexp = -1.0 * st->minlyap;
1705     break;
1706
1707   case 17:
1708     ff = "abbb";
1709     st->dwell =  400;
1710     st->settle = 200;
1711     st->minlyap = st->maxexp = ABS(-0.85);
1712     st->minexp = -1.0 * st->minlyap;
1713     break;
1714
1715   case 18:
1716     ff = "abbabaab";
1717     st->dwell =  400;
1718     st->settle = 200;
1719     st->minlyap = st->maxexp = ABS(-0.85);
1720     st->minexp = -1.0 * st->minlyap;
1721     break;
1722
1723   case 19:
1724     st->mapindex = 2;
1725     ff = "aaaaaabbbbbb";
1726     st->dwell =  400;
1727     st->settle = 200;
1728     st->minlyap = st->maxexp = ABS(-0.85);
1729     st->minexp = -1.0 * st->minlyap;
1730     break;
1731
1732   case 20:
1733     st->mapindex = 2;
1734     st->dwell =  400;
1735     st->settle = 200;
1736     st->minlyap = st->maxexp = ABS(-0.85);
1737     st->minexp = -1.0 * st->minlyap;
1738     break;
1739
1740   case 21:
1741     st->mapindex = 2;
1742     ff = "abbb";
1743     st->dwell =  400;
1744     st->settle = 200;
1745     st->minlyap = st->maxexp = ABS(-0.85);
1746     st->minexp = -1.0 * st->minlyap;
1747     break;
1748
1749   case 22:
1750     st->mapindex = 2;
1751     ff = "abbabaab";
1752     st->dwell =  400;
1753     st->settle = 200;
1754     st->minlyap = st->maxexp = ABS(-0.85);
1755     st->minexp = -1.0 * st->minlyap;
1756     break;
1757
1758   default:
1759     abort();
1760     break;
1761   }
1762
1763   if (ff) {
1764     char *ch;
1765     int bindex = 0;
1766     st->maxindex = strlen(ff);
1767     if (st->maxindex > MAXINDEX)
1768       usage(st);
1769     ch = ff;
1770     st->force++;
1771     while (bindex < st->maxindex) {
1772       if (*ch == 'a')
1773         st->forcing[bindex++] = 0;
1774       else if (*ch == 'b')
1775         st->forcing[bindex++] = 1;
1776       else
1777         usage(st);
1778       ch++;
1779     }
1780   }
1781 }
1782
1783
1784 static void *
1785 xlyap_init (Display *d, Window window)
1786 {
1787   struct state *st = (struct state *) calloc (1, sizeof(*st));
1788   XWindowAttributes xgwa;
1789   int builtin = -1;
1790   XGetWindowAttributes (d, window, &xgwa);
1791   st->dpy = d;
1792   st->width = xgwa.width;
1793   st->height = xgwa.height;
1794   st->visual = xgwa.visual;
1795   st->screen = xgwa.screen;
1796   st->cmap = xgwa.colormap;
1797
1798   do_defaults(st);
1799   parseargs(st);
1800
1801   if (get_boolean_resource(st->dpy, "randomize", "Boolean"))
1802     builtin = random() % NBUILTINS;
1803   else {
1804     char *s = get_string_resource(st->dpy, "builtin", "Integer");
1805     if (s && *s)
1806       builtin = atoi(s);
1807     if (s) free (s);
1808   }
1809
1810   if (builtin >= 0)
1811     do_preset (st, builtin);
1812
1813   st->background = BlackPixelOfScreen(st->screen);
1814   setupmem(st);
1815   init_data(st);
1816   if (!mono_p)
1817     st->foreground = st->startcolor;
1818   else
1819     st->foreground = WhitePixelOfScreen(st->screen);
1820
1821   /*
1822    * Create the window to display the Lyapunov exponents
1823    */
1824   st->canvas = window;
1825   init_color(st);
1826
1827 #ifdef BACKING_PIXMAP
1828   st->pixmap = XCreatePixmap(st->dpy, window, st->width, st->height,
1829                              xgwa.depth);
1830 #endif
1831 /*  st->rubber_data.band_cursor = XCreateFontCursor(st->dpy, XC_hand2);*/
1832 /*  CreateXorGC(st);*/
1833   Clear(st);
1834
1835   st->delay  = get_integer_resource(st->dpy, "delay", "Delay");
1836   st->linger = get_integer_resource(st->dpy, "linger", "Linger");
1837   if (st->linger < 1) st->linger = 1;
1838
1839   return st;
1840 }
1841
1842
1843 static unsigned long
1844 xlyap_draw (Display *dpy, Window window, void *closure)
1845 {
1846   struct state *st = (struct state *) closure;
1847   int i;
1848
1849   if (!st->run && st->reset_countdown) {
1850     st->reset_countdown--;
1851     if (st->reset_countdown)
1852       return 1000000;
1853     else {
1854       do_defaults (st);
1855       do_preset (st, (random() % NBUILTINS));
1856       Clear (st);
1857       init_data(st);
1858       init_color(st);
1859       resize (st);
1860       st->frame = 0;
1861       st->run = 1;
1862     }
1863   }
1864
1865   for (i = 0; i < 2000; i++)
1866     if (complyap(st) == TRUE)
1867       {
1868         st->run = 0;
1869         st->reset_countdown = st->linger;
1870         break;
1871       }
1872   return st->delay;
1873 }
1874
1875 static void
1876 xlyap_reshape (Display *dpy, Window window, void *closure,
1877                  unsigned int w, unsigned int h)
1878 {
1879   struct state *st = (struct state *) closure;
1880   resize(st);
1881 }
1882
1883 static Bool
1884 xlyap_event (Display *dpy, Window window, void *closure, XEvent *event)
1885 {
1886   struct state *st = (struct state *) closure;
1887
1888   switch(event->type)
1889     {
1890     case KeyPress:
1891       if (Getkey(st, &event->xkey))
1892         return True;
1893       break;
1894 #if 0
1895     case ButtonPress:
1896       StartRubberBand(st, &st->rubber_data, event);
1897       return True;
1898     case MotionNotify:
1899       TrackRubberBand(st, &st->rubber_data, event);
1900       return True;
1901     case ButtonRelease:
1902       EndRubberBand(st, &st->rubber_data, event);
1903       return True;
1904 #endif
1905     default:
1906       break;
1907     }
1908
1909   if (screenhack_event_helper (dpy, window, event))
1910     {
1911       Clear(st);
1912       return True;
1913     }
1914
1915   return False;
1916 }
1917
1918 static void
1919 xlyap_free (Display *dpy, Window window, void *closure)
1920 {
1921   int i;
1922   struct state *st = (struct state *) closure;
1923
1924   freemem (st);
1925
1926 #ifdef BACKING_PIXMAP
1927   XFreePixmap (st->dpy, st->pixmap);
1928 #endif
1929 /*  XFreeGC (st->dpy, st->RubberGC);*/
1930   for (i = 0; i < st->maxcolor; i++)
1931     XFreeGC (st->dpy, st->Data_GC[i]);
1932   if (st->outname) free (st->outname);
1933
1934   free (st);
1935 }
1936
1937
1938 XSCREENSAVER_MODULE ("XLyap", xlyap)