partibrains.c

static char local_id[] = "$Id$";
/*
 * $Log$
 * Revision 1.52  2001/07/09 23:46:56  slevy
 * A bit more ripening of selection system...
 * And, report match counts from "sel", "thresh", "emph".
 *
 * Revision 1.51  2001/07/07 15:35:49  slevy
 * Implement more (most?) of selection framework:
 *     selsrc(), seldest(), selname(), seltoken().
 * "thresh" now accepts "-s 'destexpr'".
 * Add new selthresh SelOp.
 *
 * Revision 1.50  2001/07/04 03:44:12  slevy
 * Add framework of set-selection stuff; not really usable yet.
 * Next: implement "see" in drawspecks().
 * Add "picked" callback which gets notified of all pick results.
 *
 * Revision 1.49  2001/06/30 18:14:22  slevy
 * Provide for storing another bit in rgba alpha field: EMPHBIT, for
 * emphasis tag.
 * Use new sfStrDrawTJ (with transform & justification) for text.
 * Yeow, handle polymax properly -- default == infinity!
 * Add "vcmap" for per-datafield colormaps.
 * Enable GL_SMOOTH mode to get texture colors to work right (huh? why do we need to?).
 * "tfm c"/"tfm w" choose cam vs. world (default) location for object tfms.
 *
 * Revision 1.48  2001/05/30 14:32:55  slevy
 * Add subcameras ("subcam" command).
 * popen() needs "r" not "rb" -- no such thing as binary mode.
 * kira_open() gets original filename if we can't find it, for error messaging.
 *
 * Revision 1.47  2001/05/15 12:18:57  slevy
 * New interface to dynamic-data routines.
 * Now, the only #ifdef KIRA/WARP needed in partibrains.c are
 * the data-command initialization routines.  All others,
 * including control-command parsing and specialized drawing,
 * is now via a function table.
 *
 * Revision 1.46  2001/05/14 15:51:34  slevy
 * Initialize new speckseq value whenever we make a new specklist.
 * Put colorseq/sizeseq/threshseq invalidation in kira_parti, not here,
 * since we needn't do it when warping.
 * By default, clock range is whatever it had been set to before.
 *
 * Revision 1.45  2001/05/12 07:22:24  slevy
 * Add get-time-range func for dynamic data.
 * Initialize all dynamic-data stuff explicitly.
 *
 * Revision 1.44  2001/05/11 10:05:38  slevy
 * Add "warp" command if -DUSE_WARP.
 * Ellipsoids allow 3-component (Rx Ry Rz) orientation numbers too.
 * st->dyndata doesn't mean we should trash the anima[][] specklist.
 * Need to keep it intact for warping.
 *
 * Revision 1.43  2001/05/02 09:51:01  slevy
 * Add "ellipsoids" and "meshes" commands to toggle their display.
 * Finally parse snapshot arguments properly...?
 *
 * Revision 1.42  2001/04/26 08:52:46  slevy
 * Add hook for off-screen rendering: "-w width[xheight]" option for
 * snapshot/snapset.  Not actually implemented yet.
 * Add "clipbox hide".  This (a) doesn't display the big yellow box
 * and (b) doesn't enable OpenGL clipping -- just whatever object culling
 * the draw routine does anyway.
 * specks_draw_mesh() now pays attention to drawing style (solid/line/plane/point).
 * When drawing ellipsoids, disable blending altogether if bgcolor != 0.
 * Small, distant labels now disappear by default.  To get the old
 * draw-tiny-line behavior, use a negative "labelmin" value.
 *
 * Revision 1.41  2001/04/20 13:47:39  slevy
 * Allow *cment commands to alter multiple colormap entries: new editcmap().
 *
 * Revision 1.40  2001/04/10 19:18:33  slevy
 * Open sdb files (and while we're at it, other files too) in binary mode.
 *
 * Revision 1.39  2001/04/04 20:33:02  slevy
 * Use findfile() for "read" ctrl-command.
 *
 * Revision 1.38  2001/03/30 16:49:00  slevy
 * Change enum SurfStyle to avoid mentioning POINT, which Windows uses too.
 *
 * Revision 1.37  2001/03/30 14:00:42  slevy
 * Allow cment/boxcment/textcment as data commands too.
 * Add "ghosts" ctl command (not really implemented).
 * Add "winsize" ctl command.
 * Ellipsoids now accept -l levelno option, in which case
 * "hide levelno"/"show levelno" applies to them.  Default is -1,
 * which means "show if any level is shown", but hidden by "hide all".
 * Scrap ghost maintenance.  Do that elsewhere instead.
 * Make "maxcomment" work properly.
 *
 * Revision 1.36  2001/03/19 11:55:08  slevy
 * New version.c (derived from ../VERSION) contains current version string.
 * New "version" command in partibrains reports that and the partibrains.c CVS ver no.
 *
 * Revision 1.35  2001/03/19 10:39:07  slevy
 * Handle speck comments properly.
 *
 * Revision 1.34  2001/03/15 18:19:06  slevy
 * Don't include comments in argc/argv -- it complicates things.
 * Have a separate "comment" pointer.  Pass it to specks_read_ellipsoid too.
 *
 * Revision 1.33  2001/03/15 15:37:06  slevy
 * Yeow -- handle VIRDIR prefix properly!
 * Complain of unrecognized data commands.
 *
 * Revision 1.32  2001/03/14 17:27:05  slevy
 * Make rejoinargs() work -- don't omit last arg.
 * "object" command alone reports our currently selected object.
 *
 * Revision 1.31  2001/03/13 22:45:08  slevy
 * Add "verbose" flag to parti_allobjs, so "gall -v ..."
 * reports object name before invoking each command.
 *
 * Revision 1.30  2001/03/13 08:23:39  slevy
 * Use parti_object's new "create-if-not-present" flag.
 * Data-language references can create, command-language ones can't.
 * Switch specks_read() to argc/argv style, using new "tokenize()" function.
 * Some commands still want the original string, so we also have rejoinargs().
 * Allow for adjustable comment length with "maxcomment" data command.
 * New "ellipsoid" and "mesh" objects.  Only quadmeshes implemented right now.
 * Data-language "tfm" command is now quiet.  Command-language tfm still verbose.
 * Check at run time for endian-config errors and refuse to run if wrong!
 *
 * Revision 1.29  2001/03/08 22:02:29  slevy
 * Disable the parti menu unless PARTIMENU envar set.
 * Allow alpha to adjust brightness of mesh objects.
 *
 * Revision 1.28  2001/03/05 03:04:41  slevy
 * Add mesh objects.  Or, quad meshes, anyway.
 *
 * Revision 1.27  2001/03/04 16:47:35  slevy
 * Make "add" and "eval" behave consistently.
 * polyorivar and texturevar now accept field names as well as numbers.
 *
 * Revision 1.26  2001/02/22 20:04:34  slevy
 * CONSTify.
 *
 * Revision 1.25  2001/02/19 22:01:16  slevy
 * Satisfy windows C compiler: pull enum FadeModel outside struct speck, etc.
 *
 * Revision 1.24  2001/02/19 20:50:33  slevy
 * Oops, "object" should always invoke parti_object() whether there's an
 * alias or not!
 *
 * Revision 1.23  2001/02/17 22:02:54  slevy
 * Enlarge polygons so that unit circle is inscribed, not circumscribed.
 * Then "txscale .5" always shows entire texture regardless of polysides.
 * Add new "ptsize" command -- makes more sense than "fast".  "fast" still works.
 *
 * Revision 1.22  2001/02/17 17:44:05  slevy
 * For polygons, rotate circle of vertices by 1/2 step.
 * Then, for "polysides 4" and "txscale .707",
 * the vertices coincide with the corners of the 0..1 texture.
 *
 * Revision 1.21  2001/02/17 05:39:45  slevy
 * Allow (in data language) "object gN=NAME".
 *
 * Revision 1.20  2001/02/15 05:41:12  slevy
 * new textcmap, textcment commands.  Regularize color-gamma-mapping.
 * Object aliases: "object gN=ALIAS", or in command mode, "gN=ALIAS".
 * Accept "ellipsoid" data tag; not yet implemented.
 *
 * Revision 1.19  2001/02/05 00:41:52  slevy
 * Accept "time" as synonym for "step".
 * Add "pickrange".
 * Mention jump, center in help msg.
 *
 * Revision 1.18  2001/02/03 23:43:10  slevy
 * Don't let "kira tree" hide "kira track" -- demand 3 chars!
 *
 * Revision 1.17  2001/02/03 16:49:42  slevy
 * Update cookedcmap when "cment" changes cmap.
 *
 * Revision 1.16  2001/02/03 15:29:03  slevy
 * Toss unused variable.
 *
 * Revision 1.15  2001/02/03 14:50:44  slevy
 * Add "setgamma" (abbr. "setgam" or "cgam") command to adjust colors.
 * Add "kira tree {off|on|cross|tick} [tickscale]" subcommand
 * for showing tree structure of interacting groups.
 *
 * Revision 1.14  2001/01/31 17:11:54  slevy
 * Ensure that, for starlab, clock is always in "continuous" mode.
 *
 * Revision 1.13  2001/01/31 17:07:12  slevy
 * Add RCS Id and Log strings.
 *
 */

#define __USE_MISC	/* makes <math.h> define sqrtf() on GNU libc */

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <ctype.h>

#if !sgi && !defined(sqrtf)
# define sqrtf(x)  sqrt(x)	/* if no sqrtf() */
#endif

#undef isspace		/* hack for irix 6.5 back-compat */
#undef isdigit
#undef isalnum

#if unix
# include <unistd.h>
# include <sys/types.h>
# include <netinet/in.h>  /* for htonl */
# include <time.h>
# include <alloca.h>
# ifndef WORDS_BIGENDIAN
#  include "config.h"	/* for WORDS_BIGENDIAN */
# endif
#else /*WIN32*/
# include "winjunk.h"
# define WORDS_BIGENDIAN 0
#endif
#include <string.h>
#include <errno.h>


#if WORDS_BIGENDIAN
#define RGBALPHA(rgb, alpha)	((rgb) | (alpha))
#define RGBWHITE		0xFFFFFF00
#define	PACKRGBA(r,g,b,a)	((r)<<24 | (g)<<16 | (b)<<8 | (a))
#else
#define RGBALPHA(rgb, alpha)	((rgb) | ((alpha)<<24))
#define RGBWHITE		0x00FFFFFF
#define	PACKRGBA(r,g,b,a)	((a)<<24 | (b)<<16 | (g)<<8 | (r))
#endif

#define	THRESHBIT		PACKRGBA(0,0,0,1)
#define	EMPHBIT			PACKRGBA(0,0,0,2)
#define	EXTRABITS		PACKRGBA(0,0,0,0xff)
#define	RGBBITS			PACKRGBA(0xff,0xff,0xff,0)

#if CAVE
# define  CAVEMENU  1
#endif

#include "geometry.h"
#include <GL/gl.h>

#if CAVEMENU
#include <cave_ogl.h>
#include "cavemenu.h"
#include "vd_util.h"
#define IFMENU(x)  (x)
#else
#define IFMENU(x)
#endif

#include "shmem.h"	/* NewN(), etc. */
#include "futil.h"

#include "specks.h"

#include "textures.h"
#include "findfile.h"
#include "partiviewc.h"
#include "sfont.h"

#include <sys/types.h>
#include <signal.h>

#ifdef sgi
#include <malloc.h>	/* for mallinfo(), amallinfo() */
#endif

#ifdef USE_KIRA
#include "kira_parti.h"
#endif

  /* Star Renderer (.sdb) structure -- from stardef.h */
typedef enum {ST_POINT, ST_BRIGHT_CLOUD ,ST_DARK_CLOUD, ST_BOTH_CLOUD, ST_OFF} stype;

typedef struct {
        float  x, y, z;
        float  dx, dy, dz;
        float  magnitude, radius;
        float  opacity;
        int  num;
        unsigned short  color;
        unsigned char   group;
        unsigned char   type;
} db_star;

typedef  struct  hrec { float  t;  int  num;}  hrec_t;
typedef  struct  mrec { float  mass, x, y, z, vx, vy, vz, rho, temp, sfr, gasmass;
                 int  id, token;}  mrec_t;
  /* end Star Renderer */


#define VDOT( v1, v2 )  ( (v1)->x[0]*(v2)->x[0] + (v1)->x[1]*(v2)->x[1] + (v1)->x[2]*(v2)->x[2] )

#if CAVEMENU
int parti_menuwall, parti_menubox[4];

static MenuHead *pmenu = NULL;
static MenuHead *stubmenu = NULL;

struct ments {
  MenuEnt *alpha, *point, *poly, *label, *psize, *polysize, *fog;
  MenuEnt *step, *fwd, *lumvar, *slum, *colorvar, *where[4];
  MenuEnt *viewall, *setscale, *speed, *fspeed, *every;
  MenuEnt *seedata, *box, *boxlevel[12];
  MenuEnt *annot;
  MenuEnt *tknob[3];
  MenuEnt *menu;
} ment;

struct boxleveler {
  struct stuff *st;
  int level;
};

void set_tknob( float val, MenuEnt *me, void *st );
void set_psize( float psize, MenuEnt *me, void *st );
void set_polysize( float psize, MenuEnt *me, void *st );
void set_alpha( float alpha, MenuEnt *me, void *st );
void set_point( int on, MenuEnt *me, void *st );
void set_poly( int on, MenuEnt *me, void *st );
void set_label( int on, MenuEnt *me, void *st );

void set_step( float time0, MenuEnt *me, void *st );
void set_fwd( int fwd, MenuEnt *me, void *st );
void set_slum( float slum, MenuEnt *me, void *st );
void set_lumvar( int which, MenuEnt *me, void *st );
void set_colorvar( int which, MenuEnt *me, void *st );
void set_seedata( int which, MenuEnt *me, void *st );
void set_speed( float speed, MenuEnt *me, void *st );
void set_fspeed( float fspeed, MenuEnt *me, void *st );
void set_every( float every, MenuEnt *me, void *st );
void set_scale( float logscale, MenuEnt *me, void *st );
void set_viewall( int all, MenuEnt *me, void *st );
void set_allboxes( int on, MenuEnt *me, void *st );
void set_boxlevel( int on, MenuEnt *me, void *bleveler );
void set_where( int unused, MenuEnt *me, void *st );
void set_menu( int on, MenuEnt *me, void *st );
#endif /*CAVEMENU*/

static int defcmap[] = {
  0x11eeee00,
  0x1106ee00, 0x120ea900, 0x1316ce00, 0x1520d500, 0x172ece00, 0x193fcb00,
  0x1c54b400, 0x206da400, 0x24889200, 0x2aa58400, 0x33c07600, 0x3ed96a00,
  0x4eed6100, 0x63fa5a00, 0x7efe5100, 0x99fb4c00, 0xafef4700, 0xc0dc4000,
  0xcbc33900, 0xd4a83000, 0xda8b2900, 0xdf702700, 0xe2572500, 0xe5411700,
  0xe82f1000, 0xea211500, 0xeb171400, 0xed0f1200, 0xee0aee00,
  0xffffff00,
};

int orientboxcolor = PACKRGBA( 0xff, 0xff, 0, 0xff );

void specks_read( struct stuff **stp, char *fname );
void specks_read_boxes( struct stuff *st, char *fname, int timestep );
int  specks_add_box( struct stuff *st, struct AMRbox *box, int timestep );
int  specks_purge( void *vstuff, int nbytes, void *aarena );
int  specks_count( struct specklist *head );
int  specks_gobox( struct stuff *st, int boxno, int argcrest, char *argvrest[] );
int specks_cookcment( struct stuff *st, int cment );
void specks_rgbremap( struct stuff *st );

void specks_reupdate( struct stuff *st, struct specklist *sl );
struct specklist **specks_find_annotation( struct stuff *, struct specklist **);
void specks_set_current_annotation( struct stuff *st, char *annotation );
void specks_add_annotation( struct stuff *st, char *annotation, int timestep );
void specks_timerange( struct stuff *st, double *tminp, double *tmaxp );

extern int getbool( char *str, int defval );
extern double getfloat( char *str, double defval );
extern int getfloats( float *v, int nfloats, int arg0, int argc, char **argv );
extern int specks_set_byvariable( struct stuff *st, char *str, int *val );

static char separ[] = " \t\n";

void strncpyt( char *dst, char *src, int dstsize ) {
  int len = strlen(src);
  if(len >= dstsize) len = dstsize-1;
  memcpy(dst, src, len);
  dst[len] = '\0';
} 

#ifdef sgi
static float defgamma = 1.0;
#else
static float defgamma = 2.5;
#endif

struct stuff *
specks_init( int argc, char *argv[] )
{
  int i;
  struct stuff *st = NewN( struct stuff, 1 );
  char *menuopt;

  i = PACKRGBA(1, 0, 0, 0);
  if(*(char *)&i != 1) {
    msg("specks_init(): trouble: is WORDS_BIGENDIAN mis-set?  Giving up.");
    exit(1);
  }

  memset(st, 0, sizeof(*st));
  st->spacescale = 1.0;
  st->fog = 0;
  st->psize = 1;
  st->alpha = .5;
  st->gamma = defgamma;
  st->rgbgamma[0] = st->rgbgamma[1] = st->rgbgamma[2] = 1.0;
  st->rgbright[0] = st->rgbright[1] = st->rgbright[2] = 1.0;
  specks_rgbremap( st );
  st->alias = NULL;
  st->useme = 1;
  st->usepoly = 0;
  st->usepoint = 1;
  st->usetext = 1;
  st->usetextaxes = 1;
  st->usetextures = 1;
  st->useboxes = 1;
  st->usemeshes = 1;
  st->useellipsoids = 1;
  st->polysizevar = -1;
  st->polyarea = 0;
  st->polyorivar0 = -1;
  st->texturevar = -1;
  st->txscale = .5;
  st->boxlabels = 0;
  st->boxlabelscale = 1.0;
  st->boxlevelmask = ~0;	/* all levels on */
  st->boxaxes = 0;		/* boxes don't show orientation markers */
  st->goboxscale = 1.0;
  st->textsize = .05;
  st->npolygon = 11;
  st->subsample = 1;
  st->everycomp = 1;
  st->maxcomment = sizeof(st->sl->specks->title) - 1;
  st->dyn.enabled = 0;
  st->dyn.data = NULL;
  st->dyn.getspecks = NULL;
  st->dyn.draw = NULL;
  st->dyn.trange = NULL;
  st->dyn.help = NULL;
  st->dyn.ctlcmd = NULL;
  st->dyn.free = NULL;
  st->speckseq = 0;

  st->menudemandfps = 4.0;

  st->pfaint = .05;	/* params for "fast" point-drawing */
  st->plarge = 10;
  st->polymin = .5;	/* don't draw polygons if smaller (pixels) */
  st->polymax = 1e8;	/* don't allow polygons to get bigger than this (pixels) */
  st->polyfademax = 0;
  st->textmin = 2;	/* replace labels with line-segments if smaller (pixels) */
  st->ntextures = 0;
  st->textures = NULL;

  st->fade = F_SPHERICAL;
  st->fadeknee1 = 10.0;
  st->fadeknee2 = 1.0;
  st->knee2steep = 1.0;

  st->gscale = 1.;
  st->gtrans.x[0] = st->gtrans.x[1] = st->gtrans.x[2] = 0;
  st->objTo2w = Tidentity;

  st->ncmap = st->boxncmap = st->textncmap = COUNT(defcmap);
  st->cmap = NewN(struct cment, COUNT(defcmap));
  st->boxcmap = NewN(struct cment, COUNT(defcmap));
  st->textcmap = NewN(struct cment, COUNT(defcmap));
  for(i = 0; i < COUNT(defcmap); i++) {
    st->cmap[i].raw = st->boxcmap[i].raw = st->textcmap[i].raw
	= htonl(defcmap[i]);
    st->cmap[i].cooked = st->boxcmap[i].cooked = st->textcmap[i].cooked
	= specks_cookcment( st, st->cmap[i].raw );
  }
  /* Ensure that textcmap[0] is white by default */
  st->textcmap[0].raw = PACKRGBA( 170, 170, 170, 0 );
  st->textcmap[0].cooked = specks_cookcment( st, st->textcmap[0].raw );

  st->sizedby = 0;
  st->coloredby = 1;
  st->sizeseq = st->colorseq = st->threshseq = 0;
  st->trueradius = 0;
  st->sdbvars = shmstrdup( "mcr" );

  st->useemph = 0;
  selinit( &st->emphsel );
  st->emphfactor = 10;

  selinit( &st->threshsel );
  selinit( &st->seesel );

  st->nghosts = 0;

  st->ntimes = 0;
  st->ndata = 0;
  st->curtime = 0;
  st->curdata = 0;
  st->datatime = 0;
  st->usertrange = 0;
  st->utmin = -HUGE;
  st->utmax = HUGE;
  st->utwrap = 0.1;
  st->sl = NULL;

  st->boxes = NULL;
  st->boxlevels = 0;
  st->boxlinewidth = 0.75;

  st->depthsort = 0;

  st->clk = NewN(SClock, 1);
  clock_init(st->clk);
  clock_set_running(st->clk, 1);

  memset(st->anima, 0, sizeof(st->anima));
  memset(st->annot, 0, sizeof(st->annot));
  memset(st->datafile, 0, sizeof(st->datafile));
  memset(st->fname, 0, sizeof(st->fname));


#if CAVE
  shmrecycler( specks_purge, st );
#endif

  for(i = 1; i < argc; i++)
    specks_read( &st, argv[i] );

#if CAVEMENU
  menu_preinit();

  pmenu = menu_create( .05, .95 );
  menu_setfont( pmenu, font_small );

  ment.annot = menu_addentry( pmenu, " ", NULL, NULL );
  menu_addentry( pmenu, " ", NULL, NULL );

  menuopt = getenv("PARTIMENU");
  if(menuopt == NULL) {
    menuopt = "label";
    st->hidemenu = 1;
  }

  if(strstr(menuopt, "survey")) {	/* SC99DEMO star-formation survey */
    for(i = 0; i < 3; i++) {
	struct boxleveler *tb = NewN( struct boxleveler, 1 );
	tb->st = st;
	tb->level = i;
	ment.tknob[i] = menu_addknob( pmenu, NULL, set_tknob, 0, 0,4, tb );
    }
  }
  menu_addentry( pmenu, " ", NULL, NULL );

  ment.psize = menu_addknob( pmenu, NULL, set_psize, st->psize, 0.1,20, st );
  ment.slum = menu_addknob( pmenu, NULL, set_slum,
		st->vdesc[st->curdata][st->sizedby].lum, .1, 10, st );
  ment.lumvar = menu_addtoggle( pmenu, NULL, set_lumvar, st->sizedby, MAXVAL, st );
  ment.colorvar = menu_addtoggle( pmenu, NULL, set_colorvar, st->coloredby, MAXVAL, st );
  ment.seedata = menu_addtoggle( pmenu, NULL, set_seedata, st->curdata, 2, st );
  if(strstr(menuopt, "label"))
    ment.label = menu_addtoggle( pmenu, NULL, set_label, st->usetext, 2, st );
  if(strstr(menuopt, "poly"))
    ment.poly = menu_addtoggle( pmenu, NULL, set_poly, st->usepoly, 2, st );
  menu_addentry( pmenu, " ", NULL, NULL );	/* spacer */

#ifndef STATIC_SPECKS
  ment.fspeed = menu_addknob( pmenu, NULL, set_fspeed, clock_speed(st->clk), 0, 5, st );
  ment.speed = menu_addknob( pmenu, NULL, set_speed, clock_speed(st->clk), .5, 5, st );
  ment.step = menu_addknob( pmenu, NULL, set_step, clock_time(st->clk), 0, 10, st );
  ment.fwd = menu_addtoggle( pmenu, NULL, set_fwd, clock_fwd(st->clk), 2, st );
#endif

#ifndef NO_WHERE
  ment.where[0] = menu_addentry( pmenu, "", NULL, NULL );
  ment.where[1] = menu_addentry( pmenu, "", NULL, NULL );
  ment.where[2] = menu_addentry( pmenu, "", NULL, NULL );
  ment.where[3] = menu_addentry( pmenu, "", NULL, NULL );
#endif

  ment.every = menu_addknob( pmenu, NULL, set_every, st->subsample, 1, 10, st );

  menu_addentry( pmenu, " ", NULL, NULL );

  ment.setscale = menu_addknob( pmenu, NULL, set_scale, 0., -4., 1., st );
  if(strstr(menuopt, "viewall"))
    ment.viewall = menu_addtoggle( pmenu, NULL, set_viewall, 0, 2, st );

  menu_addentry( pmenu, " ", NULL, NULL );

  ment.box = menu_addtoggle( pmenu, NULL, set_allboxes, st->useboxes, 3, st );
  menu_addentry( pmenu, " ", NULL, NULL );
  {
    struct boxleveler *blev = NewN( struct boxleveler, COUNT(ment.boxlevel) );
    for(i = 0; i < COUNT(ment.boxlevel); i++) {
	blev[i].st = st;
	blev[i].level = i;
	ment.boxlevel[i] = menu_addtoggle( pmenu, NULL, set_boxlevel, 1, 2,
		&blev[i] );
    }
  }

  ment.menu = menu_addtoggle( pmenu, "<Menu>", set_menu, 0, 2, st );
  menu_addentry( pmenu, " ", NULL, NULL );


  stubmenu = menu_create( .05, .15 );
  menu_addtoggle( stubmenu, "<Menu>", set_menu, 0, 2, st );
  menu_sethidden( stubmenu, 1 );
#endif /*NOCAVE*/

  return st;
}

#if CAVEMENU
void specks_refresh_menu(struct stuff *st)
{
  int i;
  int empty = st->ntimes == 0 && st->sl == NULL;

#if MENU_IN_PPR
  menu_check( pmenu, pmenu->cavewall );
  menu_check( stubmenu, stubmenu->cavewall );

#endif
  menu_sethidden( pmenu, st->hidemenu || empty );
  menu_sethidden( stubmenu, st->hidemenu<=0 || empty );

  if(st->vdcmd[0] != '\0') {
    VIDI_queue_commandstr( st->vdcmd );
    st->vdcmd[0] = '\0';
  }
  set_psize( st->psize, ment.psize, st );
  set_lumvar( st->sizedby, ment.lumvar, st );
  set_colorvar( st->coloredby, ment.colorvar, st );
  set_slum( st->vdesc[st->curdata][st->sizedby].lum, ment.slum, st );
  set_every( st->subsample, ment.every, st );
  set_scale( ment.setscale->val, ment.setscale, st );
  if(ment.label) set_label( st->usetext, ment.label, st );
  if(ment.poly) set_poly( st->usepoly, ment.poly, st );
  if(ment.annot)
    menu_settitle( ment.annot, st->annotation ? st->annotation : " " );
  if(st->ntimes > 1) {
    menu_setknobrange( ment.step, 0, st->ntimes-1 );
    set_fwd( st->timefwd, ment.fwd, st );
    set_step( st->time0, ment.step, st );
    set_speed( clock_speed(st->clk), ment.speed, st );
  } else {
    menu_settitle( ment.fwd, "" );
    menu_settitle( ment.step, "" );
    menu_settitle( ment.speed, "" );
    if(ment.viewall)
      menu_settitle( ment.viewall, "" ); /* no need to "view peak", right? */
  }
  if(st->ndata > 1) {
    set_seedata( st->curdata, ment.seedata, st );
  } else {
    menu_settitle( ment.seedata, "" );
  }
  if(st->boxlevels > 0) {
    set_allboxes( st->useboxes, ment.box, st );
  } else {
    menu_settitle( ment.box, "" );
  }

  if(ment.tknob[0])	/* if SC99DEMO, i.e. if PARTIMENU includes "survey" */
    for(i = 0; i < 3; i++)
      set_tknob( ment.tknob[i]->val, ment.tknob[i], ment.tknob[i]->data );


#ifndef NO_WHERE
  for(i = 0; i < COUNT(ment.where) && ment.where[i] != NULL; i++)
    set_where( i, ment.where[i], st );
#endif

  /* Publicize our menu location so other apps -- Matt's vtk AMR code --
   * can avoid interacting if user's wand points there.
   */
  parti_menuwall = pmenu->cavewall;

  parti_menubox[0] = pmenu->x0 * pmenu->wallxpix;
  parti_menubox[1] = ((pmenu->y0 < pmenu->y1) ? pmenu->y0 : pmenu->y1)
			* pmenu->wallypix;

  parti_menubox[2] = pmenu->x1 * pmenu->wallxpix;
  parti_menubox[3] = ((pmenu->y0 > pmenu->y1) ? pmenu->y0 : pmenu->y1)
			* pmenu->wallypix;

}

void specks_evoke_menu( struct stuff *st ) {
  menu_evoke( pmenu );	/* evoke any callbacks */
}

#endif /*CAVEMENU*/

void specks_rethresh( struct stuff *st, struct specklist *sl, int by )
{
  struct valdesc *vd;
  int i;
  int curdata = st->curdata;
  int nel = sl->nspecks;
  float cmin, cmax, normal;
  int ncmap = st->ncmap;
  int index;
  int min, max;
  float threshmin = st->thresh[0], threshmax = st->thresh[1];
  struct speck *p = sl->specks;
  SelMask *sel = sl->sel;

  sl->threshseq = st->threshseq;

  if(sl->text != NULL)	/* specklists with labels shouldn't be thresholded */
    return;

  if(curdata >= st->ndata)
    curdata = 0;

  min = (SMALLSPECKSIZE(by)>=sl->bytesperspeck) ? 0 : st->usethresh&P_THRESHMIN;
  max = (SMALLSPECKSIZE(by)>=sl->bytesperspeck) ? 0 : st->usethresh&P_THRESHMAX;

  for(i = 0, p = sl->specks; i < sl->nspecks; i++, p = NextSpeck( p, sl, 1 )) {
    if((min&&p->val[by]<threshmin) || (max&&p->val[by]>threshmax)) {
	p->rgba |= THRESHBIT;
	SELUNSET( sel[i], &st->threshsel );
    } else {
	p->rgba &= ~THRESHBIT;
	SELSET( sel[i], &st->threshsel );
	st->selcount++;
    }
  }
}

int specks_cookcment( struct stuff *st, int cment )
{
    unsigned char crgba[4];
    memcpy(crgba, &cment, 4);		/* XXX 64-bit bug? */
    return (cment & ~RGBWHITE) |
	   PACKRGBA(
		st->rgbmap[0][crgba[0]],
		st->rgbmap[1][crgba[1]],
		st->rgbmap[2][crgba[2]],
		0 );
}

void specks_rgbremap( struct stuff *st )
{
  int i, k;
  for(i = 0; i < 256; i++) {
    float t = (float)i / 256;
    float v = 255.99f * pow(t, 1/st->rgbgamma[0]);
    k = (int) (st->rgbright[0] * v);
    st->rgbmap[0][i] = (k <= 0) ? 0 : (k > 255) ? 255 : k;
    if(st->rgbgamma[0] != st->rgbgamma[1])
        v = 255.99f * pow(t, 1/st->rgbgamma[1]);
    k = (int) (st->rgbright[1] * v);
    st->rgbmap[1][i] = (k <= 0) ? 0 : (k > 255) ? 255 : k;
    if(st->rgbgamma[1] != st->rgbgamma[2])
        v = 255.99f * pow(t, 1/st->rgbgamma[2]);
    k = (int) (st->rgbright[2] * v);
    st->rgbmap[2][i] = (k <= 0) ? 0 : (k > 255) ? 255 : k;
  }
  /* remake colormap too */
  for(i = 0; i < st->ncmap; i++)
    st->cmap[i].cooked = specks_cookcment( st, st->cmap[i].raw );
  st->colorseq++;
}


void specks_recolor( struct stuff *st, struct specklist *sl, int by )
{
  struct valdesc *vd;
  int i;
  int curdata = st->curdata;
  int nel = sl->nspecks;
  struct speck *sp = sl->specks;
  float cmin, cmax, normal;
  int ncmap = st->ncmap;
  struct cment *cmap = st->cmap;
  int index;
  int rgb565;
  unsigned char (*rgbmap)[256];

  sl->coloredby = by;
  sl->colorseq = st->colorseq;

  if(sl->text != NULL)	/* specklists with labels shouldn't be recolored */
    return;

  if(curdata >= st->ndata)
    curdata = 0;

  if(by == CONSTVAL) {
    /* Hack -- color by given RGB value */
    char crgba[4];
    int r, g, b, rgba;
    vd = &st->vdesc[curdata][CONSTVAL];
    r = (vd->cmin<=0) ? 0 : vd->cmin>=1 ? 255 : (int)(255.99f * vd->cmin);
    g = (vd->cmax<=0) ? 0 : vd->cmax>=1 ? 255 : (int)(255.99f * vd->cmax);
    b = (vd->mean<=0) ? 0 : vd->mean>=1 ? 255 : (int)(255.99f * vd->mean);
    crgba[0] = st->rgbmap[0][r];
    crgba[1] = st->rgbmap[0][g];
    crgba[2] = st->rgbmap[0][b];
    crgba[3] = 0;
    rgba = *(int *)&crgba[0];	/* XXX 64-bit bug? */
    for(i = 0; i < nel; i++, sp = NextSpeck(sp, sl, 1)) {
	sp->rgba = rgba | (sp->rgba & EXTRABITS);
    }
    return; 
  }

  if(by >= MAXVAL || by < 0 || SMALLSPECKSIZE(by) > sl->bytesperspeck)
	by = 0;
	
  vd = &st->vdesc[curdata][by];

  if(vd->vncmap > 0) {
    ncmap = vd->vncmap;
    cmap = vd->vcmap;
  }


  cmin = vd->cmin, cmax = vd->cmax;
  if(cmin == cmax && cmin == 0 || (vd->call && !vd->cexact)) {
    vd->cmin = cmin = vd->min;
    vd->cmax = cmax = vd->max;
  }

  normal = (cmax != cmin && ncmap>2) ? (ncmap-2)/(cmax - cmin) : 0;
  rgb565 = !strcmp(vd->name, "rgb565") || !strcmp(vd->name, "colors565");

  /* cexact field means:
   *   0 (default): scale data range to cmap index 1..ncmap-2;
   *		    use 0 and ncmap-1 for low- and high- out-of-range values.
   *
   *   1 ("exact"): use data value as literal colormap index, 0..ncmap-1.
   */
  rgbmap = &st->rgbmap[0];
  for(i = 0; i < nel; i++, sp = NextSpeck(sp, sl, 1)) {
    if(rgb565) {
	index = sp->val[by];
	sp->rgba = PACKRGBA(
			rgbmap[0][(index&0xF800)>>(11-(8-5))],
			rgbmap[1][(index&0x07E0)>>(5-(8-6))],
			rgbmap[2][(index&0x1F)<<(8-5)],
			0 )
		  | (sp->rgba & EXTRABITS);
    } else {
	index = vd->cexact  ?   sp->val[by] + cmin
			    :  (sp->val[by] - cmin) * normal + 1;

	if(index < 0) index = 0;
	else if(index >= ncmap) index = ncmap-1;

	sp->rgba = cmap[index].cooked | (sp->rgba & EXTRABITS);
    }
  }
}

void specks_resize( struct stuff *st, struct specklist *sl, int by )
{
  int i;
  int nel = sl->nspecks;
  struct speck *sp = sl->specks;
  int curdata = st->curdata;
  struct valdesc *vd;
  float lmin, lmax, normal;

  sl->sizedby = by;
  sl->sizeseq = st->sizeseq;

  if(sl->text != NULL) /* specklists with labels shouldn't be resized */
    return;

  if(curdata < 0 || curdata >= st->ndata)
    curdata = 0;

  if(by == CONSTVAL) {
    vd = &st->vdesc[curdata][CONSTVAL];
    for(i = 0; i < nel; i++, sp = NextSpeck(sp, sl, 1))
	sp->size = vd->lmin;
    return;
  }

  if(by < 0 || by >= MAXVAL)
    by = 0;

  vd = &st->vdesc[curdata][by];

  lmin = vd->lmin, lmax = vd->lmax;
  if(lmin == lmax && lmin == 0 || vd->lall) {
	vd->lmin = lmin = vd->min;
	vd->lmax = lmax = vd->max;
  }

  if(lmax == lmin)
    normal = 1;
  else
    normal = 1 / (lmax - lmin);

  if(st->useemph && st->emphsel.use != SEL_NONE) {
    for(i = 0; i < nel; i++, sp = NextSpeck(sp, sl, 1)) {
	sp->size = (sp->val[by] - lmin) * normal;
	if(SELECTED(sl->sel[i], &st->emphsel)) {
	    sp->size *= st->emphfactor;
	    st->selcount++;
	}
    }
  } else {
    /* normal case */
    for(i = 0; i < nel; i++, sp = NextSpeck(sp, sl, 1))
	sp->size = (sp->val[by] - lmin) * normal;
  }
}

void specks_datawait(struct stuff *st) {
#if USE_IEEEIO
    while(st->fetching && st->fetchpid > 0
	    && st->fetchtime == st->curtime
	    && st->fetchdata == st->curdata)
	usleep(50000);
#endif
}

#ifdef USE_IEEEIO

void specks_ieee_server( void *vst ) {
  struct stuff *st = (struct stuff *)vst;
  struct specklist *sl;

  prctl(PR_TERMCHILD);	/* Die when parent dies */
  /* Await a request */
  for(;;) {
    time_t then;
    while(st->fetching <= 0)
	sginap(5);
    then = time(NULL);
    sl = specks_ieee_read_timestep( st, st->subsample,
				st->fetchdata, st->fetchtime );
    if(then + 5 < time(NULL))
	msg("... got %x (%d particles) from %s", sl, specks_count(sl),
		st->fname[st->fetchdata][st->fetchtime]);
    if(st->curtime == st->fetchtime && st->curdata == st->fetchdata
			&& sl != NULL) {
	st->sl = sl;
	parti_redraw();
    }
    st->fetching = 0;
  }
}

#endif /*USE_IEEEIO*/

void specks_set_speed( struct stuff *st, double newspeed ) {
  clock_set_speed( st->clk, newspeed );
  parti_set_speed( st, newspeed );
}

void specks_set_fspeed( struct stuff *st, double newspeed ) {
  /* ignore it! */
  /* or maybe use fspeed to specify finite clock resolution */
}

void specks_set_timestep( struct stuff *st )
{
  struct specklist *sl = st->sl;
  struct specklist **slp;
  double realtime, tmin, tmax;
  int timestep;

  st->used++;

  tmin = st->clk->tmin;
  tmax = st->clk->tmax;
  specks_timerange( st, &tmin, &tmax );
  if(st->usertrange) {
    if(tmin < st->utmin) tmin = st->utmin;
    if(tmax > st->utmax) tmax = st->utmax;
  }
  clock_set_range( st->clk, tmin, tmax, st->utwrap );
  realtime = clock_time( st->clk );
  timestep = realtime;

  if(st->dyn.enabled && st->dyn.getspecks) {
    struct specklist *sl;

    if(realtime == st->currealtime && st->sl != NULL && st->sl->used >= 0)
	return;
    sl = (*st->dyn.getspecks)(&st->dyn, st, realtime);
    sl->used = 1;
    st->currealtime = realtime;
    st->curtime = timestep;

    slp = specks_timespecksptr( st, 0, 0 );
    st->sl = sl;
    /* Do we need this? */ /* *slp = sl; */ /* Guess not... */
    parti_set_timestep( st, realtime );
    return;
  }

  if(timestep >= st->ntimes) timestep = st->ntimes - 1;
  if(timestep < 0) timestep = 0;

  slp = specks_timespecksptr( st, st->curdata, timestep );
  if(*slp != NULL && (*slp)->subsampled > st->subsample) {
    specks_discard( st, slp );
  }

  sl = *slp;

#ifdef USE_IEEEIO

  if(sl == NULL && (st->fetching == 0 || st->datasync)
		&& st->curdata>=0 && st->curdata<st->ndata
		&& st->datafile[st->curdata] != NULL
		&& timestep>=0 && timestep<st->ntimes
		&& st->datafile[st->curdata][timestep] != NULL) {

    st->fetchdata = st->curdata;
    st->fetchtime = timestep;

    if(getenv("NO_SPROC") || st->datasync) {
	sl = specks_ieee_read_timestep( st, st->subsample,
                                st->fetchdata, st->fetchtime );
        msg("Got %x (%d particles) from %s  (d%d t%d)\n",
		sl, specks_count(sl), st->fname[st->fetchdata][st->fetchtime],
		st->fetchdata, st->fetchtime);

    } else {
	st->fetching = 1;

	if(st->fetchpid <= 0) {
	    st->fetchpid = sproc( specks_ieee_server, PR_SADDR|PR_SFDS, st );
	    if(st->fetchpid < 0)
		perror("sproc");
	}
    }
  }

  if(st->fetchpid > 0) {
    if(kill(st->fetchpid, 0) < 0) {
	perror("specks server vanished: kill -0");
	st->fetching = 0;
	st->fetchpid = 0;
    }
  }
#endif /*USE_IEEEIO*/

  if(sl == NULL && !(st->boxtimes>timestep && st->boxes[timestep]!=NULL))
    return;
  st->sl = sl;   /* st->sl <= anima[][] */
  st->curtime = timestep;
  st->currealtime = timestep;
#if !CAVEMENU
  st->frame_time = st->curtime;	/* if non-CAVE, we have no frame-function */
#endif

  specks_set_current_annotation( st, st->annot[st->curdata][st->curtime]
		? st->annot[st->curdata][st->curtime]->text : NULL );
  parti_set_timestep( st, timestep );

}

void specks_set_current_annotation( struct stuff *st, char *annotation )
{
  st->annotation = annotation;
#ifdef CAVEMENU
  if(ment.annot)
    menu_settitle( ment.annot, st->annotation ? st->annotation : "" );
#endif
}

void specks_add_annotation( struct stuff *st, char *annotation, int timestep )
{
  struct specklist *sl, **slp;
  int curtime = (timestep < 0) ? st->curtime : timestep;

  if(annotation == NULL) annotation = "";

  specks_timespecksptr( st, st->curdata, curtime );
  if(curtime >= st->ntimes) curtime = 0;
  slp = &st->annot[st->curdata][curtime];
  if((sl = *slp) == NULL) {
    sl = *slp = NewN( struct specklist, 1 );
    memset(*slp, 0, sizeof(**slp));
    sl->speckseq = ++st->speckseq;
  } else if(sl->text) {
    Free(sl->text);
  }
  sl->text = NewN( char, strlen(annotation)+1 );
  strcpy( sl->text, annotation );
}


#if !WORDS_BIGENDIAN
void starswap(db_star *st) {
  int i, *wp;
  /* byte-swap x,y,z, dx,dy,dz, magnitude,radius,opacity fields (32-bit float),
   * 		num (32-bit int),
   *		color (16-bit short).
   * group and type fields shouldn't need swapping,
   * assuming the compiler packs bytes into a word in increasing
   * address order.  Seems safe.
   */
  for(i = 0, wp = (int *)st; i < 10; i++)
    wp[i] = htonl(wp[i]);
  st->color = htons(st->color);
}
#endif /*!WORDS_BIGENDIAN*/

void specks_read_sdb( struct stuff *st, char *sdbfname, int timestep )
{
  FILE *inf = fopen(sdbfname, "rb");
  long flen;
  int nspecks, i;
  float min[MAXVAL], max[MAXVAL], sum[MAXVAL];
  struct specklist *sl, **slp;
  register struct speck *sp;
  int dfltvars = (strcmp(st->sdbvars, "mcr") == 0);
  int nvars = strlen(st->sdbvars);
  int needswap = (htonl(1) != 1);

  if(inf == NULL) {
    msg("sdb: %s: cannot open: %s", sdbfname, strerror(errno));
    return;
  }
  /* Just measure file size */
  errno = 0;
  fseek(inf, 0, SEEK_END);
  flen = ftell(inf);
  if(flen == -1 || flen == 0) {
    msg("sdb: %s: can't measure length of file: %s", sdbfname, strerror(errno));

    return;
  }
  nspecks = (flen / sizeof(db_star));

  if(nspecks <= 0) {
    msg("sdb: %s: ignoring empty sdb file", sdbfname);
    return;
  }
  
  sl = NewN(struct specklist, 1);
  memset(sl, 0, sizeof(*sl));
  sl->speckseq = ++st->speckseq;
  
  sl->bytesperspeck = SMALLSPECKSIZE( nvars );
  if(nvars > MAXVAL) nvars = MAXVAL;

  sl->scaledby = st->spacescale;
  sp = NewNSpeck(sl, nspecks);
  sl->specks = sp;
  sl->sel = NewN( SelMask, nspecks );
  sl->nsel = nspecks;
  memset(sl->sel, 0, nspecks*sizeof(SelMask));

  fseek(inf, 0, SEEK_SET);
  for(i = 0; i < nspecks; i++, sp = NextSpeck(sp, sl, 1)) {
    db_star star;
    char *cp;
    int k;
    float *vp;
    if(fread(&star, sizeof(star), 1, inf) <= 0)
	break;
#if !WORDS_BIGENDIAN
    starswap(&star);
#endif
    sp->p.x[0] = star.x * sl->scaledby;
    sp->p.x[1] = star.y * sl->scaledby;
    sp->p.x[2] = star.z * sl->scaledby;
    if(dfltvars) {
	sp->val[0] = exp((-18-star.magnitude)*.921/*log(100)/5*/);
	sp->val[1] = star.color;
	sp->val[2] = star.radius;
    } else {
	for(vp = &sp->val[0], cp = st->sdbvars; *cp; cp++, vp++) {
	    switch(*cp) {
	    case 'm': *vp = exp((-18-star.magnitude)*.921/*log(100)/5*/); break;
	    case 'M': *vp = star.magnitude; break;
	    case 'c': *vp = star.color; break;
	    case 'r': *vp = star.radius; break;
	    case 'o': *vp = star.opacity; break;
	    case 'g': *vp = star.group; break;
	    case 't': *vp = star.type; break;
	    case 'x': *vp = star.dx; break;
	    case 'y': *vp = star.dy; break;
	    case 'z': *vp = star.dz; break;
	    case 'S': *vp = sqrt(star.dx*star.dx + star.dy*star.dy + star.dz*star.dz); break;
	    case 'n': *vp = star.num; break;
	    default: *vp = 1; break;
	    }
	}
    }

    if(i == 0) {
	for(k = 0; k < nvars; k++)
	    sum[k] = min[k] = max[k] = sp->val[k];
    } else {
	for(k = 0; k < nvars; k++) {
	    if(min[k] > sp->val[k]) min[k] = sp->val[k];
	    else if(max[k] < sp->val[k]) max[k] = sp->val[k];
	    sum[k] += sp->val[k];
	}
    }
  }
  sl->nspecks = i;
  sl->sizedby = 0;
  sl->coloredby = 1;

  /* Update statistics */
  if(sl->nspecks > 0) {
    struct valdesc *vdp = &st->vdesc[st->curdata][0];
    for(i = 0; i < nvars; i++, vdp++) {
	if(vdp->min > min[i]) vdp->min = min[i];
	if(vdp->max < max[i]) vdp->max = max[i];
	vdp->nsamples += sl->nspecks;
	vdp->sum += sum[i];
	vdp->mean = vdp->sum / vdp->nsamples;

	if(vdp->name[0] == '\0') {
	    char *name = "unk";
	    switch(st->sdbvars[i]) {
	    case 'm': name = "lumsdb"; break;
	    case 'M': name = "magsdb"; break;
	    case 'c': name = vdp->max > 16384 ? "rgb565" : "colorsdb";
		      vdp->cexact = 1;
		      break;
	    case 'r': name = "radius"; break;
	    case 'o': name = "opacity"; break;
	    case 'g': name = "group"; break;
	    case 't': name = "type"; break;
	    case 'x': name = "dx"; break;
	    case 'y': name = "dy"; break;
	    case 'z': name = "dz"; break;
	    case 'S': name = "speed"; break;
	    case 'n': name = "number"; break;
	    }
	    strcpy(vdp->name, name);
	}
    }
    specks_recolor( st, sl, st->coloredby );
    specks_resize( st, sl, st->sizedby );
  }

  /* Add to running list */
  slp = specks_timespecksptr( st, st->curdata, timestep );
  sl->next = *slp;
  *slp = sl;

  fclose(inf);
    
}

void specks_timerange( struct stuff *st, double *tminp, double *tmaxp )
{
  if(st->dyn.enabled && st->dyn.trange &&
	(*st->dyn.trange)(&st->dyn, st, tminp, tmaxp)) {
    /* OK */
  } else {
    *tminp = 0;
    *tmaxp = (st->ntimes == 0) ? 0 : st->ntimes - 1;
  }
}

int specks_get_datastep( struct stuff *st )
{
  return st->curtime;
}

double specks_get_realtime( struct stuff *st )
{
  return st->currealtime;
}

void set_interest_point( Point *p )
{
#if CAVEMENU
  char cmd[80];
  sprintf(cmd, "interest %g %g %g", p->x[0], p->x[1], p->x[2]);
  VIDI_queue_commandstr( cmd );
#else

  parti_center( p );
#endif
}

/* Only specks_timespecksptr() extends the spans of time or datasets */
struct specklist **
specks_timespecksptr( struct stuff *st, int dataset, int timestep )
{
  int d, needroom;
  struct specklist **na, **nan;
  void **ndf;
  char **nfn;

  if((timestep >= st->ntimes || dataset >= st->ndata)) {
    needroom = st->timeroom;
    if(timestep >= st->timeroom)
	needroom = 2*timestep + 15;

    for(d = 0; d < st->ndata || (dataset < MAXFILES && d <= dataset); d++) {

	if(needroom == st->timeroom && d < st->ndata)
	    continue;

	na = NewN( struct specklist *, needroom );
	nan = NewN( struct specklist *, needroom );
	ndf = NewN( void *, needroom );
	nfn = NewN( char *, needroom );
	memset(na, 0, needroom * sizeof(*na));
	memset(nan, 0, needroom * sizeof(*nan));
	memset(ndf, 0, needroom * sizeof(*ndf));
	memset(nfn, 0, needroom * sizeof(*nfn));
	if(d < st->ndata && st->anima[d])
	    memcpy( na, st->anima[d], st->ntimes * sizeof(*na) );

	if(d < st->ndata && st->annot[d])
	    memcpy( nan, st->annot[d], st->ntimes * sizeof(*nan) );

	if(d < st->ndata && st->datafile[d])
	    memcpy( ndf, st->datafile[d], st->ntimes * sizeof(*ndf) );

	if(d < st->ndata && st->fname[d])
	    memcpy( nfn, st->fname[d], st->ntimes * sizeof(*nfn) );

	/* Don't free old pointers, just in case they're in use. */
	st->anima[d] = na;
	st->annot[d] = nan;
	st->datafile[d] = ndf;
	st->fname[d] = nfn;
    }
    st->timeroom = needroom;

    if(timestep >= st->ntimes)
	st->ntimes = timestep + 1;
    if(dataset >= st->ndata && dataset < MAXFILES)
	st->ndata = dataset + 1;
  }

  return (timestep >= 0 && timestep < st->ntimes &&
			dataset >= 0 && dataset < st->ndata)
	? &st->anima[dataset][timestep] : NULL;
}

struct specklist *
specks_timespecks( struct stuff *st, int dataset, int timestep )
{
  return (timestep >= 0 && timestep < st->ntimes &&
			dataset >= 0 && dataset < st->ndata)
	? st->anima[dataset][timestep] : NULL;
}

void specks_reupdate( struct stuff *st, struct specklist *sl )
{
  struct specklist *tsl;

  if(sl != NULL && sl->threshseq != st->threshseq) {
    for(tsl = sl; tsl != NULL; tsl = tsl->next)
	specks_rethresh( st, tsl, st->threshvar );
  }

  if(sl != NULL && sl->colorseq != st->colorseq) {
    for(tsl = sl; tsl != NULL; tsl = tsl->next)
	specks_recolor( st, tsl, st->coloredby );
  }

  if(sl != NULL && sl->sizeseq != st->sizeseq) {
    for(tsl = sl; tsl != NULL; tsl = tsl->next)
	specks_resize( st, tsl, st->sizedby );
  }
}

void specks_set_time( struct stuff *st, double newtime )
{
  static Point lastinterest;

  clock_set_time( st->clk, newtime );
  specks_set_timestep( st );
  if(st->sl == NULL && st->ntimes > 1 && st->clk->parent == NULL) {
    /* Skip blank time-slots -- keep incrementing until either:
     *  - we find a time-slot that has (or could have) some data, or
     *  - we've run through all time-steps (avoid infinite loops!).
     */
    int nudges, ts;
    for(nudges = 0; nudges < st->ntimes; nudges++) {
	ts = (st->curtime + nudges) % st->ntimes;
	if( specks_timespecks( st, st->curdata, ts ) != NULL ||
	    (st->datafile[st->curdata] != NULL &&
		st->datafile[st->curdata][ts] != NULL) )
	    break;
    }
    clock_set_time( st->clk, ts );
    specks_set_timestep( st );
  }
  /* st->playnext = now + (st->fspeed != 0 ? 1/st->fspeed : 0); */

#if CAVE
  if(ment.tknob[0] == NULL) {	/* if not SC99DEMO */
    struct specklist *sl = st->sl;
    if(sl != NULL && memcmp(&sl->interest, &lastinterest, sizeof(Point))
	    && (sl->interest.x[0]!=0 || sl->interest.x[1]!=0
					    || sl->interest.x[2]!=0)) {
	char str[128];
	float scale = .004;
	lastinterest = sl->interest;
	set_interest_point( &lastinterest );
	sprintf(str, "\002setjump peak %g %g %g 0 0 0 %g",
	    lastinterest.x[0], lastinterest.x[1]-scale*5, lastinterest.x[2] + scale*5,
	    scale);
	VIDI_queue_commandstr( str );
    }
  }
#endif

  specks_reupdate( st, st->sl );

#if CAVE
  specks_refresh_menu( st );
#endif
}

static int specks_freenow( struct specklist **slp, int maxage )
{
  struct specklist *sl, **sprev;
  int any = 0;

  for(sprev = slp; (sl = *sprev) != NULL && sl->used <= maxage; ) {
    if(sl->used <= maxage) {
	*sprev = sl->freelink;
	if(sl->specks != NULL)
	    Free(sl->specks);
	Free(sl);
	any++;
    } else {
	/* too recent - might still be in use */
	sprev = &sl->freelink;
    }
  }
  return any;
}

void specks_discard( struct stuff *st, struct specklist **slp )
{
  struct specklist *sl, *slnext;

  for(sl = *slp; sl != NULL; sl = slnext) {
    slnext = sl->next;
    sl->freelink = st->scrap;
    st->scrap = sl;
  }
  *slp = NULL;
}

int specks_purge( void *vst, int nbytes, void *aarena )
{
  struct stuff *st = (struct stuff *)vst;
  int oldused = st->used;
  int oldtime = -1, oldds = -1;
  int t, ds;
  struct specklist *sl;

#ifdef sgi
  static int first = 1;
  struct mallinfo mi;
  mi = amallinfo( aarena );

  if(first) {
    first = 0;
    msg("Purging %dKbyte shmem arena (currently %dK used in %d blks, %dK free)\n",
	mi.arena>>10, mi.uordblks>>10, mi.ordblks, mi.fordblks>>10);
  }
#endif

  /* Free any known scrap first */

#define OLD_ENOUGH  4

  if(specks_freenow( &st->scrap, st->used - OLD_ENOUGH ) > 0)
    return 1;

  for(t = 0; t < st->ntimes; t++) {
    if(t == st->curtime) continue;
    for(ds = 0; ds < st->ndata; ds++) {
	sl = st->anima[ds][t];
	if(sl != NULL && sl != st->sl && sl->used < oldused) {
	    oldused = st->used;
	    oldtime = t;
	    oldds = ds;
	}
    }
  }
  if(oldtime >= 0) {
    specks_discard( st, &st->anima[oldds][oldtime] );
    specks_freenow( &st->scrap, oldtime );
    return 1;	/* We freed something, so try allocating again */
  } else {
    msg("Ran out of shmem, couldn't find anything more to purge\n");
#ifdef sgi
    msg("%dKbyte shmem arena (currently %dK used in %d blks, %dK free)\n",
	mi.arena>>10, mi.uordblks>>10, mi.ordblks, mi.fordblks>>10);

#endif
    return 0;	/* No progress made -- give up */
  }
}


/*
 * Just stash these values in our frame function so they won't change visibly
 * during a frame.  Each frame function will do this; we'll just hope that
 * they don't change as the various cave-wall processes start.
 */
void specks_current_frame( struct stuff *st, struct specklist *sl )
{
  st->frame_sl = sl;
  st->frame_time = st->curtime;
  st->frame_data = st->curdata;
  st->frame_annotation = st->annotation;
  specks_reupdate( st, sl );
}

#define	MAXXYFAN 16

extern void specks_draw_boxes( struct stuff *st, struct AMRbox *boxes, int levelmask, Matrix Ttext, int oriented );

extern void specks_draw_mesh( struct stuff *st, struct mesh *m, int *texturing );
extern void specks_draw_ellipsoid( struct stuff *st, struct ellipsoid *e );

struct cpoint {
    int rgba;
    Point p;
};

void dumpcpoints( struct cpoint *cp, int n )
{
    while(--n >= 0) {
	glColor4ubv( (GLubyte *)&cp->rgba );
	glVertex3fv( &cp->p.x[0] );
	cp++;
    }
}

static Point depth_fwd;
static float depth_d;

struct order {
  float z;
  struct speck *sp;
  struct specklist *sl;
};

static int depthcmp( const void *a, const void *b )
{
  return ((struct order *)a)->z < ((struct order *)b)->z ?    1
	: ((struct order *)a)->z > ((struct order *)b)->z ? -1 : 0;
}

static int additive_blend;

void sortedpolys( struct stuff *st, struct specklist *slhead, Matrix *Tc2wp, float radperpix, float polysize )
{
  struct speck *sp, *sbase;
  struct order *op, *obase;
  int i, k, total, skip;
  struct specklist *sl;
  int usethresh = st->usethresh&P_USETHRESH ? THRESHBIT : 0;
  int bps = 0;
  int prevrgba = -1;
  int usearea = st->polyarea;
  int sizevar = st->polysizevar;
  int polyorivar = st->polyorivar0;
  int texturevar = st->texturevar;
  int txno;
  int texturing = -1;
  float s;
  float polyminrad = st->polymin * radperpix;
  float polymaxrad = st->polymax * radperpix;
  float mins2d = polyminrad * polyminrad;
  int rgba;
  int alpha = st->alpha * 255;
  int nfan = st->npolygon<MAXXYFAN ? st->npolygon : MAXXYFAN;
  float xyfan[MAXXYFAN][2];
  Point sfan[MAXXYFAN], pfan[MAXXYFAN];
  Matrix Tc2w = *Tc2wp;
  float scl = vlength( (Point *)&Tc2w.m[0*4+0] );
  float fanscale;
  Texture *wanttx;
  int additive = additive_blend;
  int wantblend = additive;

  int useclip = (st->clipbox.level != 0);
  Point clipp0 = st->clipbox.p0;
  Point clipp1 = st->clipbox.p1;

  for(total = 0, sl = slhead; sl != NULL; sl = sl->next) {
    if(sl->text != NULL || sl->nspecks == 0 || sl->special != SPECKS)
	continue;
    skip = st->subsample;
    if(sl->subsampled != 0)	/* if already subsampled */
	skip /= sl->subsampled;
    if(skip <= 0) skip = 1;

    if(bps < sl->bytesperspeck) bps = sl->bytesperspeck;
    if(usethresh) {
	for(i=sl->nspecks, sp=sl->specks; i>0; i-=skip, sp=NextSpeck(sp,sl,skip)) {
	    if((usethresh & sp->rgba) == 0)
		total++;
	}
    } else {
	total += sl->nspecks / skip;
    }
  }

  obase = op = (struct order *)malloc( (total+1) * sizeof(struct order) );
  for(sl = slhead; sl != NULL; sl = sl->next) {
    if(sl->text != NULL || sl->nspecks == 0 || sl->special != SPECKS)
	continue;
    skip = st->subsample;
    if(sl->subsampled != 0)	/* if already subsampled */
	skip /= sl->subsampled;
    if(skip <= 0) skip = 1;
    for(i=sl->nspecks, sp=sl->specks; i > 0; i-=skip, sp=NextSpeck(sp,sl,skip)) {
	float dist;
	if(usethresh & sp->rgba)
	    continue;
	dist = VDOT( &sp->p, &depth_fwd ) + depth_d;
	if(dist < 0)
	    continue;
	if(useclip &&
	  (sp->p.x[0] < clipp0.x[0] ||
	   sp->p.x[0] > clipp1.x[0] ||
	   sp->p.x[1] < clipp0.x[1] ||
	   sp->p.x[1] > clipp1.x[1] ||
	   sp->p.x[2] < clipp0.x[2] ||
	   sp->p.x[2] > clipp1.x[2]))
	    continue;
	op->z = dist;
	op->sp = sp;
	op->sl = sl;
	op++;
    }
  }

  total = op - obase;
  qsort( obase, total, sizeof(*obase), depthcmp );

  prevrgba = 0;

  /* Build prototype fan -- unit disk in screen plane */
  /* Scale the fan big enough that the unit disk is inscribed in our polygon */
  fanscale = 1 / (scl * cos(M_PI/nfan));
  for(i = 0; i < nfan; i++) {
    float theta = 2*M_PI*(i+0.5f)/nfan;
    xyfan[i][0] = cos(theta);
    xyfan[i][1] = sin(theta);
    vcomb( &sfan[i], xyfan[i][0] * fanscale, (Point *)&Tc2w.m[0*4+0],
		     xyfan[i][1] * fanscale, (Point *)&Tc2w.m[1*4+0] );
  }

  if(st->usetextures == 0 || SMALLSPECKSIZE(texturevar) > bps)
    texturevar = -1;
  if(SMALLSPECKSIZE(polyorivar) > bps)
    polyorivar = -1;

  for(i = 0, op = obase; i < total; i++, op++) {
    float dist, size;

    sp = op->sp;
    dist = op->z;
    size = sp->val[sizevar] * polysize;
    if(usearea) {
	if(size < dist * dist * mins2d)
	    continue;
	size = sqrtf(size);
    } else {
	if(size < dist * polyminrad)
	    continue;
    } 
    if(size > dist * polymaxrad)
	size = dist * polymaxrad;

    rgba = sp->rgba & RGBBITS;
    if(rgba != prevrgba) {
	prevrgba = rgba;
	rgba = RGBALPHA( prevrgba, alpha );
	glColor4ubv( (GLubyte *)&rgba );
    }

    if(texturevar >= 0 &&
	    (txno = sp->val[texturevar]) >= 0 &&
	    txno < st->ntextures &&
	    (wanttx = st->textures[txno]) != NULL) {

	txbind( wanttx, &texturing );
	wantblend = (wanttx->flags & TXF_ADD) ? 1 : additive_blend;
    } else if(texturing) {
	glDisable( GL_TEXTURE_2D );
	texturing = 0;
    }

    if(wantblend != additive) {
	additive = wantblend;
	glBlendFunc( GL_SRC_ALPHA, additive ? GL_ONE : GL_ONE_MINUS_SRC_ALPHA );
    }

    if(polyorivar >= 0 && sp->val[polyorivar] < 9) {
	float *xv = &sp->val[polyorivar];
	float *yv = &sp->val[polyorivar+3];

	glBegin( GL_TRIANGLE_FAN );
	if(texturing) {
	    for(k = 0; k < nfan; k++) {
		glTexCoord2fv( &xyfan[k][0] );
		glVertex3f(
		    sp->p.x[0] + size*(xyfan[k][0]*xv[0] + xyfan[k][1]*yv[0]),
		    sp->p.x[1] + size*(xyfan[k][0]*xv[1] + xyfan[k][1]*yv[1]),
		    sp->p.x[2] + size*(xyfan[k][0]*xv[2] + xyfan[k][1]*yv[2]));
	    }
	} else {
	    for(k = 0; k < nfan; k++) {
		glVertex3f(
		    sp->p.x[0] + size*(xyfan[k][0]*xv[0] + xyfan[k][1]*yv[0]),
		    sp->p.x[1] + size*(xyfan[k][0]*xv[1] + xyfan[k][1]*yv[1]),
		    sp->p.x[2] + size*(xyfan[k][0]*xv[2] + xyfan[k][1]*yv[2]));
	    }
	}
	glEnd();

    } else {
	glBegin( GL_TRIANGLE_FAN );
	if(texturing) {
	    for(k = 0; k < nfan; k++) {
		glTexCoord2fv( &xyfan[k][0] );
		glVertex3f(
		    sp->p.x[0] + size*sfan[k].x[0],
		    sp->p.x[1] + size*sfan[k].x[1],
		    sp->p.x[2] + size*sfan[k].x[2] );
	    }
	} else {
	    for(k = 0; k < nfan; k++) {
		glVertex3f(
		    sp->p.x[0] + size*sfan[k].x[0],
		    sp->p.x[1] + size*sfan[k].x[1],
		    sp->p.x[2] + size*sfan[k].x[2] );
	    }
	}
	glEnd();
    }
  }
  free(obase);
  if(texturing > 0)
    txbind( NULL, NULL );
}
  
  
void drawspecks( struct stuff *st )
{
  int i, slno, k;
  int rgba, alpha, prevrgba = 0;
  float prevsize = 0;
  struct specklist *sl, *slhead;
  register struct speck *p;
  Matrix Tw2c, Tc2w, Ttext, Tproj, Ttemp;
  static Point zero = {0,0,0};
  int xywh[4];
  float radperpix;
  Point tp, fan[MAXXYFAN];
  Point eyepoint;
  Point fwd;
  float fwdd;
  float tscale, scl, fanscale;
  int skip;
  static int nxyfan = 0;
  static float xyfan[MAXXYFAN][2];
  static unsigned char randskip[256];
  int randix = 0;
  int fast = st->fast;
  int inpick = st->inpick;
  /*int usethresh = st->usethresh&P_USETHRESH ? THRESHBIT : 0;*/
  SelOp seesel = st->seesel;
  int fixeddist;
  float polyminrad, polymaxrad, polyfaderad;
  float threshmin = st->thresh[0];
  float threshmax = st->thresh[1];
  int useclip = (st->clipbox.level != 0);
  Point clipp0 = st->clipbox.p0;
  Point clipp1 = st->clipbox.p1;

  float plum = st->psize;

  float knee1dist2 = st->fadeknee1 * st->fadeknee1;
  float knee2dist2 = st->fadeknee2 * st->fadeknee2;
  float orthodist2 = st->fadeknee2 * st->fadeknee2;
  float steep2knee2 = st->knee2steep * st->knee2steep / knee2dist2;
  float faderball2 = 1 / (st->fadeknee2 * st->fadeknee2);
  Point fadecen = st->fadecen;
  enum FadeType fademodel = st->fade;

  if(!st->useme)
    return;

  switch(fademodel) {
  case F_CONSTANT:
	if(orthodist2 <= 0)
	    orthodist2 = 1;
	break;
  case F_KNEE12:
	if(st->fadeknee1 >= st->fadeknee2 || st->fadeknee1 <= 0)
	    fademodel = F_KNEE2;	/* and fall into... */
  case F_KNEE2:
	if(st->fadeknee2 <= 0)
	    fademodel = F_SPHERICAL;
	break;
  case F_LREGION:
	if(st->fadeknee2 <= 0) faderball2 = 1;
	break;
  }


  { float r=0,g=0,b=0;	/* Ugh. Allow background to be non-black */
    sscanf(parti_bgcolor(NULL), "%f%f%f", &r,&g,&b);
    additive_blend = (r+g+b == 0);
  }

  if(st->clipbox.level > 0) {
    GLdouble plane[4];

    if(st->clipbox.level == 1) {
	struct AMRbox b[2];
	b[0] = st->clipbox;
	b[0].level = 0;
	b[1].level = -1;
	specks_draw_boxes( st, b, ~0, Tidentity, 1 );
    }
    for(i = 0; i < 3; i++) {
	plane[0] = plane[1] = plane[2] = 0;
	plane[i] = 1;
	plane[3] = -st->clipbox.p0.x[i];
	glClipPlane( GL_CLIP_PLANE0 + i, plane );
	glEnable( GL_CLIP_PLANE0 + i );
	plane[i] = -1;
	plane[3] = st->clipbox.p1.x[i];
	glClipPlane( GL_CLIP_PLANE0 + 3 + i, plane );
	glEnable( GL_CLIP_PLANE0 + 3 + i );
    }
  }


  if(nxyfan != st->npolygon && st->npolygon > 0) {
    if(st->npolygon > MAXXYFAN) st->npolygon = MAXXYFAN;
    nxyfan = st->npolygon;
    fanscale = 1 / cos(M_PI/nxyfan);
    for(i = 0; i < nxyfan; i++) {
	float th = (i+.5f)*2*M_PI / nxyfan;
	xyfan[i][0] = cos(th) * fanscale;
	xyfan[i][1] = sin(th) * fanscale;
    }
    srandom(11);
    for(i = 0; i < 256; i++)
	randskip[i] = random() & 0xFF;
  }


  alpha = st->alpha * 255;
  rgba = RGBALPHA( RGBWHITE, alpha );	/* BIG-ENDIAN (1,1,1,alpha) */

  /* Find displacements which lie in the screen plane, for making
   * billboard-style polygonal patches, and for text.
   */
  glGetFloatv( GL_MODELVIEW_MATRIX, Tw2c.m );
  eucinv( &Tc2w, &Tw2c );
  scl = vlength( (Point *)&Tw2c.m[0] );
  for(i = 0; i < nxyfan; i++) {
    tp.x[0] = scl*st->polysize*xyfan[i][0];
    tp.x[1] = scl*st->polysize*xyfan[i][1];
    tp.x[2] = 0;
    vtfmvector( &fan[i], &tp, &Tc2w );
  }

  /* Find projection matrix and screen (well, viewport) size,
   * so we can convert angular sizes to screen (pixel) sizes,
   * in radians per pixel.
   */
  glGetFloatv( GL_PROJECTION_MATRIX, Tproj.m );
  glGetIntegerv( GL_VIEWPORT, xywh );
  radperpix = 1 / (.5*xywh[2] * Tproj.m[0*4+0]);

  /* Construct a "forward" vector in object coords too, for measuring
   * distance from camera plane.  Note camera looks toward its -Z axis not +Z!
   */
  tp.x[0] = 0, tp.x[1] = 0, tp.x[2] = -1;
  vtfmvector( &fwd, &tp, &Tc2w );
  vunit( &fwd, &fwd );
  /*
   * Actually we want a plane equation, whose value is zero in the
   * eye plane.  Camera-space distance from camera plane = 
   *		vdot( &objectpoint, &fwd ) + fwdd.
   */
  vtfmpoint( &eyepoint, &zero, &Tc2w );
  fwdd = -vdot( &eyepoint, &fwd );

  tscale = scl * st->textsize;
  mcopy( &Ttemp, &Tidentity );
  Ttemp.m[0*4+0] = Ttemp.m[1*4+1] = Ttemp.m[2*4+2] = tscale;
  mmmul( &Ttext, &Ttemp, &Tc2w );
  vsettranslation( &Ttext, &zero );

  glDisable( GL_LIGHTING );

  /* Draw any boxes (even if we have no specks) for this timestep */
  if(st->useboxes && st->boxlevelmask != 0
	&& st->frame_time >= 0 && st->frame_time < st->boxtimes
	&& st->boxes[st->frame_time] != NULL) {
    specks_draw_boxes( st, st->boxes[st->frame_time], st->boxlevelmask, Ttext, 0 );
  }
  if(st->useboxes && st->staticboxes != NULL) {
    specks_draw_boxes( st, st->staticboxes, st->boxlevelmask, Ttext, st->boxaxes );
  }

  slhead = st->frame_sl;	/* st->sl as snapped by specks_ffn() */
  if(slhead == NULL)
    slhead = st->sl;		/* maybe there is no specks_ffn() */

  if(st->dyn.enabled && st->dyn.draw)
    (*st->dyn.draw)( &st->dyn, st, slhead, &Tc2w, radperpix );

  skip = st->subsample;
  if(slhead && slhead->subsampled != 0)	/* if already subsampled */
    skip /= slhead->subsampled;
  if(skip == 0) skip = 1;

  for(sl = slhead; sl != NULL; sl = sl->next)
    sl->used = st->used;

  if((unsigned int)st->sizedby <= MAXVAL
			&& (unsigned int)st->curdata < MAXFILES
			&& st->vdesc[st->curdata][st->sizedby].lum != 0) {
	plum *= st->vdesc[st->curdata][st->sizedby].lum;
  }
  if(st->subsample > 0 && st->everycomp)
      plum *= st->subsample;	/* Compensate for "every" subsampling */


  if(st->alpha >= 1) {
    glDisable(GL_BLEND);
    glEnable(GL_DEPTH_TEST);
  } else {
    glEnable(GL_BLEND);
    glBlendFunc( GL_SRC_ALPHA, additive_blend ? GL_ONE : GL_ONE_MINUS_SRC_ALPHA );
    glEnable(GL_DEPTH_TEST);
    glDepthMask( GL_FALSE );
  }

  if(inpick) glLoadName(0);
  

  if(st->usepoly && st->polysize > 0) {
    int texturing = 0;
    int texturevar = st->usetextures && st->texturevar >= 0
			&& st->texturevar < MAXVAL
		   ? st->texturevar : -1;
    int usearea = st->polyarea;
    int sizevar = st->polysizevar;
    float polysize = st->polysize;
    float mins2d;

    int txno;

    glMatrixMode( GL_TEXTURE );
    glLoadIdentity();
    glTranslatef( .5, .5, 0 );
    glScalef( st->txscale, st->txscale, st->txscale );
    glMatrixMode( GL_MODELVIEW );

    /* Textures seem to be multiplied by some mysterious factor,
     * unless we're in GL_SMOOTH mode.  What gives?
     */
    glShadeModel( texturevar >= 0 ? GL_SMOOTH : GL_FLAT );


    polyminrad = st->polymin * radperpix;
    polymaxrad = st->polymax * radperpix;
    mins2d = polyminrad*polyminrad;

    if(sizevar == -1) {
	/* If polygon size is tied to point size,
	 * then include pointsize scale factors in polygon scaling.
	 */
      if((unsigned int)st->sizedby <= MAXVAL
		&& (unsigned int)st->curdata < MAXFILES
		&& st->vdesc[st->curdata][st->sizedby].lum != 0)
	polysize *= st->vdesc[st->curdata][st->sizedby].lum;
      if(st->subsample > 0 && st->everycomp)
        polysize *= st->subsample; /* Compensate for "every" subsampling */
    }

    if(st->depthsort && !inpick) {
	depth_fwd = fwd;
	depth_d = fwdd;
	sortedpolys( st, slhead, &Tc2w, radperpix, polysize );

    } else {
      int polyalpha = alpha;
      int usepolymax = (st->polymax < 1e8);

#ifdef POLYFADE
      int polyfade = 0;
      if(st->polyfademax > st->polymin) {
	polyfade = 1;
	polyfadesize = 
	glEnable(GL_BLEND);
	glBlendFunc( GL_SRC_ALPHA, additive_blend ? GL_ONE : GL_ONE_MINUS_SRC_ALPHA );
      }
#endif

      for(sl = slhead, slno = 1; sl != NULL; sl = sl->next, slno++) {
	if(sl->text != NULL || sl->special != SPECKS) continue;
	if(inpick) {
	    glLoadName(slno);
	    glPushName(0);
	}
	for(i = 0, p = sl->specks; i < sl->nspecks; i+=skip, p = NextSpeck( p, sl, skip )) {
	    float dist = VDOT( &p->p, &fwd ) + fwdd;
	    float size;

	    if(!SELECTED(sl->sel[i], &seesel))
		continue;

	    size = p->val[sizevar] * polysize;

	    if(dist + size <= 0) continue;
	    if(usearea) {
		if(size < dist * dist * mins2d)
		    continue;
		size = sqrtf(size);
	    } else {
		if(size < dist * polyminrad)
		    continue;
	    } 
	    if(usepolymax && size > dist * polymaxrad && dist > 0)
		size = dist * polymaxrad;

	    rgba = p->rgba & RGBBITS;
	    if(rgba != prevrgba) {
		prevrgba = rgba;
		rgba = RGBALPHA( prevrgba, alpha );
		glColor4ubv( (GLubyte *)&rgba );
	    }
	    if(st->polyorivar0 >= 0 && p->val[st->polyorivar0] < 9) {
		for(k = 0; k < nxyfan; k++) {
		    vcomb( &fan[k],
			size*xyfan[k][0], (Point *)&p->val[st->polyorivar0],
			size*xyfan[k][1], (Point *)&p->val[st->polyorivar0+3] );
		}
	    } else if(p->size != prevsize) {
		float s = scl*size;
		for(k = 0; k < nxyfan; k++) {
		    vcomb( &fan[k], s*xyfan[k][0], (Point *)&Tc2w.m[0*4+0],
				    s*xyfan[k][1], (Point *)&Tc2w.m[1*4+0] );
		}
		prevsize = size;
	    }

#define PFAN(vno, comp)  p->p.x[comp] + fan[vno].x[comp]

	    if(inpick) {
		glLoadName( i );
		glBegin( GL_TRIANGLE_FAN );
		for(k = 0; k < nxyfan; k++) {
		    glVertex3f( PFAN(k,0), PFAN(k,1), PFAN(k,2) );
		}
		glEnd();

	    } else if(texturevar >= 0
		    && (txno = p->val[texturevar]) >= 0
		    && txno < st->ntextures &&
		    st->textures[txno] != NULL) {

		txbind( st->textures[txno], &texturing );

		glBegin( GL_TRIANGLE_FAN );
		for(k = 0; k < nxyfan; k++) {
		    glTexCoord2fv( &xyfan[k][0] );
		    glVertex3f( PFAN(k,0), PFAN(k,1), PFAN(k,2) );
		}
		glEnd();

	    } else {
		if(texturing) {
		    texturing = 0;
		    glDisable( GL_TEXTURE_2D );
		}
		glBegin(GL_TRIANGLE_FAN);
		for(k = 0; k < nxyfan; k++)
		    glVertex3f( PFAN(k,0), PFAN(k,1), PFAN(k,2) );
		glEnd();
	    }
#undef PFAN

	}
	if(inpick) glPopName();
      }
    }
    if(texturing) {
	txbind( NULL, NULL );
	glDisable( GL_TEXTURE_2D );
    }
    glMatrixMode( GL_TEXTURE );
    glLoadIdentity();
    glMatrixMode( GL_MODELVIEW );
  }

  if(st->usepoint && !(st->useboxes == 2)) {

#define MAXPTSIZE 16	/* in half-point units */
#define PERBUCKET 64	/* max points per bucket */

    struct cpoint sized[MAXPTSIZE*2][PERBUCKET];
    int nsized[MAXPTSIZE*2];
    unsigned char invgamma[256];
    float invgam = (st->gamma <= 0) ? 0 : 1/st->gamma;

    for(i = 0; i < 256; i++)
	invgamma[i] = (int) (255.99 * pow( i/255., invgam ));

    if(inpick) {
	for(sl = slhead, slno = 1; sl != NULL; sl = sl->next, slno++) {
	    if(sl->text != NULL || sl->special != SPECKS) continue;
	    glLoadName(slno);
	    glPushName(0);
	    for(i = 0, p = sl->specks; i < sl->nspecks; i+=skip, p = NextSpeck(p, sl, skip)) {
		if(!SELECTED(sl->sel[i], &seesel))
		    continue;

		glLoadName(i);
		glBegin(GL_POINTS);
		glVertex3fv( &p->p.x[0] );
		glEnd();
	    }
	    glPopName();
	}

    } else if(fast) {
	static unsigned char apxsize[MAXPTSIZE*MAXPTSIZE];
	unsigned char faintrand[256];
	int pxsize, oldpxsize;
	int pxmin, pxmax;

	pxmin = 256 * st->pfaint;
	if(st->plarge > MAXPTSIZE) st->plarge = MAXPTSIZE;
	pxmax = 256 * st->plarge * st->plarge;

	if(apxsize[1] == 0) {
	    for(i=0; i<COUNT(apxsize); i++)
		apxsize[i] = (int)ceil(sqrtf(i+1));
	}
	for(i = 0; i < 256; i++)
	    faintrand[i] = randskip[i] * st->pfaint;

	/* Render using fast (non-antialiased) points */
	glDisable( GL_POINT_SMOOTH );
	glEnable( GL_BLEND );
	glBlendFunc( GL_SRC_ALPHA, additive_blend ? GL_ONE : GL_ONE_MINUS_SRC_ALPHA );

	prevrgba = 0;
	prevsize = 0;

	sl = slhead;

	pxsize = oldpxsize = 1;
	glPointSize( pxsize );
	glColor4ubv( (GLubyte *)&rgba );

	for(i = 0; i < MAXPTSIZE*2; i++)
	    nsized[i] = 0;

	glBegin( GL_POINTS );
	for(sl = slhead; sl != NULL; sl = sl->next) {
	    if(sl->text != NULL || sl->special != SPECKS) continue;
	    for(i = 0, p = sl->specks; i < sl->nspecks; i+=skip, p= NextSpeck(p, sl, skip)) {
		int lum, myalpha;
		float dist = VDOT( &p->p, &fwd ) + fwdd;
		if(dist <= 0)	/* Behind eye plane */
		    continue;

		if(!SELECTED(sl->sel[i], &seesel))
		    continue;

		if(useclip &&
		  (p->p.x[0] < clipp0.x[0] ||
		   p->p.x[0] > clipp1.x[0] ||
		   p->p.x[1] < clipp0.x[1] ||
		   p->p.x[1] > clipp1.x[1] ||
		   p->p.x[2] < clipp0.x[2] ||
		   p->p.x[2] > clipp1.x[2]))
		    continue;