1 /* anemotaxis, Copyright (c) 2004 Eugene Balkovski
3 * Permission to use, copy, modify, distribute, and sell this software
4 * and its documentation for any purpose is hereby granted without
5 * fee, provided that the above copyright notice appear in all copies
6 * and that both that copyright notice and this permission notice
7 * appear in supporting documentation. No representations are made
8 * about the suitability of this software for any purpose. It is
9 * provided "as is" without express or implied warranty.
12 /*------------------------------------------------------------------------
16 | DESCRIPTION Anemotaxis
18 | This code illustrates an optimal algorithm designed
19 | for searching a source of particles on a plane.
20 | The particles drift in one direction and walk randomly
21 | in the other. The only information available to the
22 | searcher is the presence of a particle at its location
23 | and the local direction from where particle arrived.
24 | The algorithm "explains" the behavior
25 | of some animals and insects
26 | who use olfactory and directional cues to find
27 | sources of odor (mates, food, home etc) in
28 | turbulent atmosphere (odor-modulated anemotaxis),
29 | e.g. male moths locating females who release
30 | pheromones to attract males. The search trajectories
31 | resemble the trajectories that the animals follow.
34 | WRITTEN BY Eugene Balkovski
36 | MODIFICATIONS june 2004 started
38 +----------------------------------------------------------------------*/
43 -distance <arg> size of the lattice
44 -sources <arg> number of sources
45 -searhers <arg> number of searcher */
48 #include "screenhack.h"
50 #ifdef HAVE_DOUBLE_BUFFER_EXTENSION
55 /*-----------------------------------------------------------------------+
57 +-----------------------------------------------------------------------*/
62 #define MAX_INV_RATE 5
64 #define RND(x) (random() % (x))
65 #define X(x) ((int) (st->ax * (x) + st->bx))
66 #define Y(x) ((int) (st->ay * (x) + st->by))
75 short y; /* y-coordinate of the particle (relative to the source) */
77 short v; /* velocity of the particle. Allowed values are -1, 0, 1.
78 2 means the particle is not valid */
83 int n; /* size of array xv */
85 YV *yv; /* yv[i] keeps velocity and y-coordinate of the
86 particle at (i + 1, yv[i].y) relative to the
89 int inv_rate; /* Inverse rate of particle emission (if 0 then
90 source doesn't emit new particles (although
91 old ones can still exist )*/
93 Point r; /* Position of the source */
100 typedef struct PList {
105 typedef enum {UP_LEFT, UP_RIGHT, LEFT, RIGHT, DONE} State_t;
109 Point r; /* Current position */
111 Point vertex; /* Position of the vertex of the most recent
112 cone, which is the region where the source
113 is located. We do exhaustive search in the
114 cone until we encounter a new particle,
115 which gives us a new cone. */
117 State_t state; /* Internal state variable */
119 unsigned char c; /* Concentration at r */
121 short v; /* Velocity at r (good only when c != 0) */
123 PList *hist; /* Trajectory */
125 int rs; /* small shift of x-coordinate to avoid
126 painting over the same x */
136 int max_dist, max_src, max_searcher;
138 double ax, ay, bx, by;
146 #ifdef HAVE_DOUBLE_BUFFER_EXTENSION
147 XdbeBackBuffer backb;
150 long delay; /* usecs to wait between updates. */
152 int scrWidth, scrHeight;
163 /*-----------------------------------------------------------------------+
165 +-----------------------------------------------------------------------*/
169 /*-----------------------------------------------------------------------+
170 | PRIVATE FUNCTIONS |
171 +-----------------------------------------------------------------------*/
173 static void *emalloc(size_t size)
175 void *ret = malloc(size);
178 fprintf(stderr, "out of memory\n");
184 static Searcher *new_searcher(struct state *st)
186 Searcher *m = (Searcher *) emalloc(sizeof(Searcher));
188 m->r.x = m->vertex.x = st->max_dist;
191 m->r.y = RND(2 * st->max_dist);
192 } while(m->r.y < MIN_DIST || m->r.y > 2 * st->max_dist - MIN_DIST);
194 m->vertex.y = m->r.y;
196 m->state = (RND(2) == 0 ? UP_RIGHT : UP_LEFT);
198 m->color = st->colors[random() % st->ncolors].pixel;
205 static void destroy_searcher(Searcher *m)
207 PList *p = m->hist, *q;
218 static void write_hist(Searcher *m)
223 m->hist = (PList *) emalloc(sizeof(PList));
229 static void move_searcher(Searcher *m)
237 m->state = (RND(2) == 0 ? UP_LEFT : UP_RIGHT);
255 if(m->vertex.x - m->r.x == m->vertex.y - m->r.y) {
274 if(m->vertex.x - m->r.x == m->r.y - m->vertex.y) {
286 static void evolve_source(Source *s)
291 /* propagate existing particles */
293 for(i = s->n - 1; i > 0; i--) {
295 if(s->yv[i - 1].v == 2)
298 s->yv[i].v = RND(3) - 1;
299 s->yv[i].y = s->yv[i - 1].y + s->yv[i].v;
305 if(s->inv_rate > 0 && (RND(s->inv_rate) == 0)) /* emit a particle */
306 s->yv[0].y = s->yv[0].v = RND(3) - 1;
312 static Source *new_source(struct state *st)
316 Source *s = (Source *) emalloc(sizeof(Source));
317 if(st->max_searcher == 0) {
319 s->r.y = RND(2 * st->max_dist);
322 s->r.x = RND(st->max_dist / 3);
324 s->r.y = RND(2 * st->max_dist);
325 } while(s->r.y < MIN_DIST || s->r.y > 2 * st->max_dist - MIN_DIST);
328 s->n = st->max_dist - s->r.x;
329 s->yv = emalloc(sizeof(YV) * s->n);
331 for(i = 0; i < s->n; i++)
334 s->inv_rate = RND(MAX_INV_RATE);
339 s->color = st->colors[random() % st->ncolors].pixel;
344 static void destroy_source(Source *s)
350 static void get_v(const Source *s, Searcher *m)
352 int x = m->r.x - s->r.x - 1;
356 if(x < 0 || x >= s->n)
359 if((s->yv[x].v == 2) || (s->yv[x].y != m->r.y - s->r.y))
369 anemotaxis_init (Display *disp, Window win)
371 struct state *st = (struct state *) calloc (1, sizeof(*st));
372 XWindowAttributes wa;
377 XGetWindowAttributes(st->dpy, st->window, &wa);
379 st->ncolors = get_integer_resource (st->dpy, "colors", "Colors");
381 st->colors = (XColor *) malloc(sizeof(*st->colors) * (st->ncolors+1));
382 make_random_colormap (wa.screen, wa.visual, wa.colormap,
383 st->colors, &st->ncolors,
384 True, True, 0, True);
386 st->delay = get_integer_resource(st->dpy, "delay", "Integer");
387 st->max_dist = get_integer_resource(st->dpy, "distance", "Integer");
389 if(st->max_dist < MIN_DIST)
390 st->max_dist = MIN_DIST + 1;
391 if(st->max_dist > MAX_DIST)
392 st->max_dist = MAX_DIST;
394 st->max_src = get_integer_resource(st->dpy, "sources", "Integer");
399 st->max_searcher = get_integer_resource(st->dpy, "searchers", "Integer");
401 if(st->max_searcher < 0)
402 st->max_searcher = 0;
406 # ifdef HAVE_COCOA /* Don't second-guess Quartz's double-buffering */
410 st->source = (Source **) emalloc(sizeof(Source *) * st->max_src);
411 memset(st->source, 0, st->max_src * sizeof(Source *));
413 st->source[0] = new_source(st);
415 st->searcher = (Searcher **) emalloc(st->max_searcher * sizeof(Searcher *));
417 memset(st->searcher, 0, st->max_searcher * sizeof(Searcher *));
419 st->b = st->ba = st->bb = 0; /* double-buffer to reduce flicker */
421 #ifdef HAVE_DOUBLE_BUFFER_EXTENSION
422 st->b = st->backb = xdbe_get_backbuffer (st->dpy, st->window, XdbeUndefined);
425 st->scrWidth = wa.width;
426 st->scrHeight = wa.height;
427 st->cmap = wa.colormap;
428 st->gcDraw = XCreateGC(st->dpy, st->window, 0, &st->gcv);
429 st->gcv.foreground = get_pixel_resource(st->dpy, st->cmap,
430 "background", "Background");
431 st->gcClear = XCreateGC(st->dpy, st->window, GCForeground, &st->gcv);
435 st->ba = XCreatePixmap (st->dpy, st->window, st->scrWidth, st->scrHeight, wa.depth);
436 st->bb = XCreatePixmap (st->dpy, st->window, st->scrWidth, st->scrHeight, wa.depth);
444 if (st->ba) XFillRectangle (st->dpy, st->ba, st->gcClear, 0, 0, st->scrWidth, st->scrHeight);
445 if (st->bb) XFillRectangle (st->dpy, st->bb, st->gcClear, 0, 0, st->scrWidth, st->scrHeight);
447 st->ax = st->scrWidth / (double) st->max_dist;
448 st->ay = st->scrHeight / (2. * st->max_dist);
452 if((st->dx = st->scrWidth / (2 * st->max_dist)) == 0)
454 if((st->dy = st->scrHeight / (4 * st->max_dist)) == 0)
456 XSetLineAttributes(st->dpy, st->gcDraw, st->dx / 3 + 1, LineSolid, CapRound, JoinRound);
457 XClearWindow(st->dpy, st->window);
462 static void draw_searcher(struct state *st, Drawable curr_window, int i)
467 if(st->searcher[i] == NULL)
470 XSetForeground(st->dpy, st->gcDraw, st->searcher[i]->color);
472 r1.x = X(st->searcher[i]->r.x) + st->searcher[i]->rs;
473 r1.y = Y(st->searcher[i]->r.y);
475 XFillRectangle(st->dpy, curr_window, st->gcDraw, r1.x - 2, r1.y - 2, 4, 4);
477 for(l = st->searcher[i]->hist; l != NULL; l = l->next) {
479 r2.x = X(l->r.x) + st->searcher[i]->rs;
482 XDrawLine(st->dpy, curr_window, st->gcDraw, r1.x, r1.y, r2.x, r2.y);
489 static void draw_image(struct state *st, Drawable curr_window)
494 for(i = 0; i < st->max_src; i++) {
496 if(st->source[i] == NULL)
499 XSetForeground(st->dpy, st->gcDraw, st->source[i]->color);
501 if(st->source[i]->inv_rate > 0) {
503 if(st->max_searcher > 0) {
504 x = (int) X(st->source[i]->r.x);
505 y = (int) Y(st->source[i]->r.y);
506 j = st->dx * (MAX_INV_RATE + 1 - st->source[i]->inv_rate) / (2 * MAX_INV_RATE);
509 XFillArc(st->dpy, curr_window, st->gcDraw, x - j, y - j, 2 * j, 2 * j, 0, 360 * 64);
512 for(j = 0; j < st->source[i]->n; j++) {
514 if(st->source[i]->yv[j].v == 2)
517 /* Move the particles slightly off lattice */
518 x = X(st->source[i]->r.x + 1 + j) + RND(st->dx);
519 y = Y(st->source[i]->r.y + st->source[i]->yv[j].y) + RND(st->dy);
520 XFillArc(st->dpy, curr_window, st->gcDraw, x - 2, y - 2, 4, 4, 0, 360 * 64);
525 for(i = 0; i < st->max_searcher; i++)
526 draw_searcher(st, curr_window, i);
530 static void animate_anemotaxis(struct state *st, Drawable curr_window)
535 for(i = 0; i < st->max_src; i++) {
537 if(st->source[i] == NULL)
540 evolve_source(st->source[i]);
542 /* reap dead sources for which all particles are gone */
543 if(st->source[i]->inv_rate == 0) {
547 for(j = 0; j < st->source[i]->n; j++) {
548 if(st->source[i]->yv[j].v != 2) {
555 destroy_source(st->source[i]);
556 st->source[i] = NULL;
561 /* Decide if we want to add new sources */
563 for(i = 0; i < st->max_src; i++) {
564 if(st->source[i] == NULL && RND(st->max_dist * st->max_src) == 0)
565 st->source[i] = new_source(st);
568 if(st->max_searcher == 0) { /* kill some sources when searchers don't do that */
569 for(i = 0; i < st->max_src; i++) {
570 if(st->source[i] != NULL && RND(st->max_dist * st->max_src) == 0) {
571 destroy_source(st->source[i]);
577 /* Now deal with searchers */
579 for(i = 0; i < st->max_searcher; i++) {
581 if((st->searcher[i] != NULL) && (st->searcher[i]->state == DONE)) {
582 destroy_searcher(st->searcher[i]);
583 st->searcher[i] = NULL;
586 if(st->searcher[i] == NULL) {
588 if(RND(st->max_dist * st->max_searcher) == 0) {
590 st->searcher[i] = new_searcher(st);
595 if(st->searcher[i] == NULL)
598 /* Check if searcher found a source or missed all of them */
599 for(j = 0; j < st->max_src; j++) {
601 if(st->source[j] == NULL || st->source[j]->inv_rate == 0)
604 if(st->searcher[i]->r.x < 0) {
605 st->searcher[i]->state = DONE;
609 if((st->source[j]->r.y == st->searcher[i]->r.y) &&
610 (st->source[j]->r.x == st->searcher[i]->r.x)) {
611 st->searcher[i]->state = DONE;
612 st->source[j]->inv_rate = 0; /* source disappears */
615 st->searcher[i]->color = WhitePixel(st->dpy, DefaultScreen(st->dpy));
621 st->searcher[i]->c = 0; /* set it here in case we don't get to get_v() */
623 /* Find the concentration at searcher's location */
625 if(st->searcher[i]->state != DONE) {
626 for(j = 0; j < st->max_src; j++) {
628 if(st->source[j] == NULL)
631 get_v(st->source[j], st->searcher[i]);
633 if(st->searcher[i]->c == 1)
638 move_searcher(st->searcher[i]);
642 draw_image(st, curr_window);
646 anemotaxis_draw (Display *disp, Window window, void *closure)
648 struct state *st = (struct state *) closure;
649 XWindowAttributes wa;
653 XGetWindowAttributes(st->dpy, st->window, &wa);
658 if(w != st->scrWidth || h != st->scrHeight) {
661 st->ax = st->scrWidth / (double) st->max_dist;
662 st->ay = st->scrHeight / (2. * st->max_dist);
666 if((st->dx = st->scrWidth / (2 * st->max_dist)) == 0)
668 if((st->dy = st->scrHeight / (4 * st->max_dist)) == 0)
670 XSetLineAttributes(st->dpy, st->gcDraw, st->dx / 3 + 1, LineSolid, CapRound, JoinRound);
673 XFillRectangle (st->dpy, st->b, st->gcClear, 0, 0, st->scrWidth, st->scrHeight);
674 animate_anemotaxis(st, st->b);
676 #ifdef HAVE_DOUBLE_BUFFER_EXTENSION
678 XdbeSwapInfo info[1];
679 info[0].swap_window = st->window;
680 info[0].swap_action = XdbeUndefined;
681 XdbeSwapBuffers (st->dpy, info, 1);
686 XCopyArea (st->dpy, st->b, st->window, st->gcClear, 0, 0,
687 st->scrWidth, st->scrHeight, 0, 0);
688 st->b = (st->b == st->ba ? st->bb : st->ba);
697 anemotaxis_reshape (Display *dpy, Window window, void *closure,
698 unsigned int w, unsigned int h)
703 anemotaxis_event (Display *dpy, Window window, void *closure, XEvent *event)
709 anemotaxis_free (Display *dpy, Window window, void *closure)
711 struct state *st = (struct state *) closure;
714 for (i = 0; i < st->max_src; i++)
715 if (st->source[i]) destroy_source (st->source[i]);
719 for (i = 0; i < st->max_searcher; i++)
720 if (st->searcher[i]) destroy_searcher (st->searcher[i]);
729 static const char *anemotaxis_defaults [] = {
730 ".background: black",
736 #ifdef HAVE_DOUBLE_BUFFER_EXTENSION
743 static XrmOptionDescRec anemotaxis_options [] = {
744 { "-delay", ".delay", XrmoptionSepArg, 0 },
745 { "-distance", ".distance", XrmoptionSepArg, 0 },
746 { "-sources", ".sources", XrmoptionSepArg, 0 },
747 { "-searchers", ".searchers", XrmoptionSepArg, 0 },
748 { "-colors", ".colors", XrmoptionSepArg, 0 },
753 XSCREENSAVER_MODULE ("Anemotaxis", anemotaxis)