/brlcad/tags/rel-7-16-10/src/rt/do.c

/*                            D O . C
 * BRL-CAD
 *
 * Copyright (c) 1987-2010 United States Government as represented by
 * the U.S. Army Research Laboratory.
 *
 * This program is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Lesser General Public License
 * version 2.1 as published by the Free Software Foundation.
 *
 * This program is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public
 * License along with this file; see the file named COPYING for more
 * information.
 */
/** @file do.c
 *
 * Routines that process the various commands, and manage the overall
 * process of running the raytracing process.
 *
 */

#include "common.h"

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

#ifdef HAVE_SYS_TYPES_H
#  include <sys/types.h>
#endif
#ifdef HAVE_SYS_STAT_H
#  include <sys/stat.h>
#endif

#include "bio.h"
#include "vmath.h"
#include "raytrace.h"
#include "fb.h"
#include "rtprivate.h"
#include "bu.h"

/***** Variables shared with viewing model *** */
extern FILE *outfp;			/* optional pixel output file */
extern double azimuth, elevation;
extern mat_t view2model;
extern mat_t model2view;
/***** end of sharing with viewing model *****/

extern void grid_setup(void);
extern void worker(int cpu, genptr_t arg);

/***** variables shared with worker() ******/
extern struct application ap;
extern int hypersample;		/* number of extra rays to fire */
extern fastf_t aspect;		/* view aspect ratio X/Y */
extern fastf_t cell_width;	/* model space grid cell width */
extern fastf_t cell_height;	/* model space grid cell height */
extern point_t eye_model;	/* model-space location of eye */
extern fastf_t eye_backoff;	/* dist of eye from center */
extern fastf_t rt_perspective;	/* persp (degrees X) 0 => ortho */
extern size_t width;		/* # of pixels in X */
extern size_t height;		/* # of lines in Y */
extern mat_t Viewrotscale;	/* view orientation quaternion */
extern fastf_t viewsize;
extern int incr_mode;		/* !0 for incremental resolution */
extern int incr_level;		/* current incremental level */
extern int incr_nlevel;		/* number of levels */
extern int npsw;
extern struct resource resource[];
/***** end variables shared with worker() */

/***** variables shared with rt.c *****/
extern char *string_pix_start;	/* string spec of starting pixel */
extern char *string_pix_end;	/* string spec of ending pixel */
extern int pix_start;		/* pixel to start at */
extern int pix_end;		/* pixel to end at */
extern int nobjs;		/* Number of cmd-line treetops */
extern char **objtab;		/* array of treetop strings */
extern int matflag;		/* read matrix from stdin */
extern int desiredframe;	/* frame to start at */
extern int finalframe;		/* frame to halt at */
extern int curframe;		/* current frame number */
extern char *outputfile;	/* name of base of output file */
extern int interactive;		/* human is watching results */
extern int save_overlaps;	/* flag for setting rti_save_overlaps */
/***** end variables shared with rt.c *****/

/***** variables shared with viewg3.c *****/
struct bu_vls ray_data_file;  /* file name for ray data output */
/***** end variables shared with viewg3.c *****/

/***** variables for frame buffer black pixel rendering *****/
unsigned char *pixmap = NULL; /* Pixel Map for rerendering of black pixels */


void def_tree(register struct rt_i *rtip);
void do_ae(double azim, double elev);
void res_pr(void);
void memory_summary(void);

extern struct bu_image_file *bif;

/**
 * O L D _ F R A M E
 *
 * Acquire particulars about a frame, in the old format.  Returns -1
 * if unable to acquire info, 0 if successful.
 */
int
old_frame(FILE *fp)
{
    register int i;
    char number[128];

    /* Visible part is from -1 to +1 in view space */
    if (fscanf(fp, "%128s", number) != 1) return -1;
    viewsize = atof(number);
    if (fscanf(fp, "%128s", number) != 1) return -1;
    eye_model[X] = atof(number);
    if (fscanf(fp, "%128s", number) != 1) return -1;
    eye_model[Y] = atof(number);
    if (fscanf(fp, "%128s", number) != 1) return -1;
    eye_model[Z] = atof(number);
    for (i=0; i < 16; i++) {
	if (fscanf(fp, "%128s", number) != 1)
	    return -1;
	Viewrotscale[i] = atof(number);
    }
    return 0;		/* OK */
}


/**
 * O L D _ W A Y
 *
 * Determine if input file is old or new format, and if old format,
 * handle process.  Returns 0 if new way, 1 if old way (and all done).
 * Note that the rewind() will fail on ttys, pipes, and sockets
 * (sigh).
 */
int
old_way(FILE *fp)
{
    int c;

    viewsize = -42.0;

    /* Sneek a peek at the first character, and then put it back */
    if ((c = fgetc(fp)) == EOF) {
	/* Claim old way, all (ie, nothing) done */
	return 1;
    }
    if (ungetc(c, fp) != c)
	bu_exit(EXIT_FAILURE, "do.c:old_way() ungetc failure\n");

    /*
     * Old format files start immediately with a %.9e format, so the
     * very first character should be a digit or '-'.
     */
    if ((c < '0' || c > '9') && c != '-') {
	return 0;		/* Not old way */
    }

    if (old_frame(fp) < 0 || viewsize <= 0.0) {
	rewind(fp);
	return 0;		/* Not old way */
    }
    bu_log("Interpreting command stream in old format\n");

    def_tree(ap.a_rt_i);	/* Load the default trees */

    curframe = 0;
    do {
	if (finalframe >= 0 && curframe > finalframe)
	    return 1;
	if (curframe >= desiredframe)
	    do_frame(curframe);
	curframe++;
    }  while (old_frame(fp) >= 0 && viewsize > 0.0);
    return 1;			/* old way, all done */
}


/**
 * C M _ S T A R T
 *
 * Process "start" command in new format input stream
 */
int cm_start(int argc, char **argv)
{
    char *buf = (char *)NULL;
    int frame;

    frame = atoi(argv[1]);
    if (finalframe >= 0 && frame > finalframe)
	return -1;	/* Indicate EOF -- user declared a halt */
    if (frame >= desiredframe) {
	curframe = frame;
	return 0;
    }

    /* Skip over unwanted frames -- find next frame start */
    while ((buf = rt_read_cmd(stdin)) != (char *)0) {
	register char *cp;

	cp = buf;
	while (*cp && isspace(*cp)) cp++;	/* skip spaces */
	if (strncmp(cp, "start", 5) != 0) continue;
	while (*cp && !isspace(*cp)) cp++;	/* skip keyword */
	while (*cp && isspace(*cp)) cp++;	/* skip spaces */
	frame = atoi(cp);
	bu_free(buf, "rt_read_cmd command buffer (skipping frames)");
	buf = (char *)0;
	if (finalframe >= 0 && frame > finalframe)
	    return -1;			/* "EOF" */
	if (frame >= desiredframe) {
	    curframe = frame;
	    return 0;
	}
    }
    return -1;		/* EOF */
}


int cm_vsize(int argc, char **argv)
{
    viewsize = atof(argv[1]);
    return 0;
}


int cm_eyept(int argc, char **argv)
{
    register int i;

    for (i=0; i<3; i++)
	eye_model[i] = atof(argv[i+1]);
    return 0;
}


int cm_lookat_pt(int argc, char **argv)
{
    point_t pt;
    vect_t dir;
    int yflip = 0;
    quat_t quat;

    if (argc < 4)
	return -1;
    pt[X] = atof(argv[1]);
    pt[Y] = atof(argv[2]);
    pt[Z] = atof(argv[3]);
    if (argc > 4)
	yflip = atoi(argv[4]);

    /*
     * eye_pt should have been specified before here and different
     * from the lookat point or the lookat point will be from the
     * "front"
     */
    if (VAPPROXEQUAL(pt, eye_model, VDIVIDE_TOL)) {
	VSETALLN(quat, 0.5, 4);
	quat_quat2mat(Viewrotscale, quat); /* front */
    } else {
	VSUB2(dir, pt, eye_model);
	VUNITIZE(dir);
	bn_mat_lookat(Viewrotscale, dir, yflip);
    }

    return 0;
}


int cm_vrot(int argc, char **argv)
{
    register int i;

    for (i=0; i<16; i++)
	Viewrotscale[i] = atof(argv[i+1]);
    return 0;
}


int cm_orientation(int argc, char **argv)
{
    register int i;
    quat_t quat;

    for (i=0; i<4; i++)
	quat[i] = atof(argv[i+1]);
    quat_quat2mat(Viewrotscale, quat);
    return 0;
}


int cm_end(int argc, char **argv)
{
    struct rt_i *rtip = ap.a_rt_i;

    if (rtip && BU_LIST_IS_EMPTY(&rtip->HeadRegion)) {
	def_tree(rtip);		/* Load the default trees */
    }

    /* If no matrix or az/el specified yet, use params from cmd line */
    if (Viewrotscale[15] <= 0.0)
	do_ae(azimuth, elevation);

    if (do_frame(curframe) < 0) return -1;
    return 0;
}


int cm_tree(int argc, const char **argv)
{
    register struct rt_i *rtip = ap.a_rt_i;
    struct bu_vls times;

    if (argc <= 1) {
	def_tree(rtip);		/* Load the default trees */
	return 0;
    }
    bu_vls_init(&times);

    rt_prep_timer();
    if (rt_gettrees(rtip, argc-1, &argv[1], npsw) < 0)
	bu_log("rt_gettrees(%s) FAILED\n", argv[0]);
    (void)rt_get_timer(&times, NULL);

    if (rt_verbosity & VERBOSE_STATS)
	bu_log("GETTREE: %s\n", bu_vls_addr(&times));
    bu_vls_free(&times);
    return 0;
}


int cm_multiview(int argc, char **argv)
{
    register struct rt_i *rtip = ap.a_rt_i;
    size_t i;
    static int a[] = {
	35,   0,
	0,  90, 135, 180, 225, 270, 315,
	0,  90, 135, 180, 225, 270, 315
    };
    static int e[] = {
	25, 90,
	30, 30, 30, 30, 30, 30, 30,
	60, 60, 60, 60, 60, 60, 60
    };

    if (rtip && BU_LIST_IS_EMPTY(&rtip->HeadRegion)) {
	def_tree(rtip);		/* Load the default trees */
    }
    for (i=0; i<(sizeof(a)/sizeof(a[0])); i++) {
	do_ae((double)a[i], (double)e[i]);
	(void)do_frame(curframe++);
    }
    return -1;	/* end RT by returning an error */
}


/**
 * C M _ A N I M
 *
 * Experimental animation code
 *
 * Usage:  anim <path> <type> args
 */
int cm_anim(int argc, const char **argv)
{

    if (db_parse_anim(ap.a_rt_i->rti_dbip, argc, argv) < 0) {
	bu_log("cm_anim:  %s %s failed\n", argv[1], argv[2]);
	return -1;		/* BAD */
    }
    return 0;
}


/**
 * C M _ C L E A N
 *
 * Clean out results of last rt_prep(), and start anew.
 */
int cm_clean(int argc, char **argv)
{
    /* Allow lighting model clean up (e.g. lights, materials, etc) */
    view_cleanup(ap.a_rt_i);

    rt_clean(ap.a_rt_i);

    if (R_DEBUG&RDEBUG_RTMEM_END)
	bu_prmem("After cm_clean");
    return 0;
}


/**
 * C M _ C L O S E D B
 *
 * To be invoked after a "clean" command, to close out the ".g"
 * database.  Intended for memory debugging, to help chase down memory
 * "leaks".  This terminates the program, as there is no longer a
 * database.
 */
int cm_closedb(int argc, char **argv)
{
    db_close(ap.a_rt_i->rti_dbip);
    ap.a_rt_i->rti_dbip = DBI_NULL;

    bu_free((genptr_t)ap.a_rt_i, "struct rt_i");
    ap.a_rt_i = RTI_NULL;

    bu_prmem("After _closedb");
    bu_exit(0, NULL);

    return 1;	/* for compiler */
}


/* viewing module specific variables */
extern struct bu_structparse view_parse[];

struct bu_structparse set_parse[] = {
    {"%d",	1, "width",	bu_byteoffset(width),		BU_STRUCTPARSE_FUNC_NULL },
    {"%d",	1, "height",	bu_byteoffset(height),		BU_STRUCTPARSE_FUNC_NULL },
    {"%d",	1, "save_overlaps", bu_byteoffset(save_overlaps),	BU_STRUCTPARSE_FUNC_NULL },
    {"%f",	1, "perspective", bu_byteoffset(rt_perspective),	BU_STRUCTPARSE_FUNC_NULL },
    {"%f",	1, "angle",	bu_byteoffset(rt_perspective),	BU_STRUCTPARSE_FUNC_NULL },
#if !defined(_WIN32) || defined(__CYGWIN__)
    {"%d",  1, "rt_bot_minpieces", bu_byteoffset(rt_bot_minpieces), BU_STRUCTPARSE_FUNC_NULL },
    {"%d",  1, "rt_bot_tri_per_piece", bu_byteoffset(rt_bot_tri_per_piece), BU_STRUCTPARSE_FUNC_NULL },
    {"%f",  1, "rt_cline_radius", bu_byteoffset(rt_cline_radius), BU_STRUCTPARSE_FUNC_NULL },
#endif
    {"%V",  1, "ray_data_file", bu_byteoffset(ray_data_file), BU_STRUCTPARSE_FUNC_NULL },
    {"i", bu_byteoffset(view_parse[0]), "View_Module-Specific Parameters", 0, BU_STRUCTPARSE_FUNC_NULL },
    {"",	0, (char *)0,	0,				BU_STRUCTPARSE_FUNC_NULL }
};


/**
 * C M _ S E T
 *
 * Allow variable values to be set or examined.
 */
int cm_set(int argc, char **argv)
{
    struct bu_vls str;

    if (argc <= 1) {
	bu_struct_print("Generic and Application-Specific Parameter Values",
			set_parse, (char *)0);
	return 0;
    }
    bu_vls_init(&str);
    bu_vls_from_argv(&str, argc-1, (const char **)argv+1);
    if (bu_struct_parse(&str, set_parse, (char *)0) < 0) {
	bu_vls_free(&str);
	return -1;
    }
    bu_vls_free(&str);
    return 0;
}


/**
 * C M _ A E
 */
int cm_ae(int argc, char **argv)
{
    azimuth = atof(argv[1]);	/* set elevation and azimuth */
    elevation = atof(argv[2]);
    do_ae(azimuth, elevation);

    return 0;
}


/**
 * C M _ O P T
 */
int cm_opt(int argc, char **argv)
{
    int old_bu_optind=bu_optind;	/* need to restore this value after calling get_args() */

    if (get_args(argc, argv) <= 0) {
	bu_optind = old_bu_optind;
	return -1;
    }
    bu_optind = old_bu_optind;
    return 0;
}


/**
 * D E F _ T R E E
 *
 * Load default tree list, from command line.
 */
void
def_tree(register struct rt_i *rtip)
{
    struct bu_vls times;

    RT_CK_RTI(rtip);

    bu_vls_init(&times);
    rt_prep_timer();
    if (rt_gettrees(rtip, nobjs, (const char **)objtab, npsw) < 0)
	bu_log("rt_gettrees(%s) FAILED\n", objtab[0]);
    (void)rt_get_timer(&times, NULL);

    if (rt_verbosity & VERBOSE_STATS)
	bu_log("GETTREE: %s\n", bu_vls_addr(&times));
    bu_vls_free(&times);
    memory_summary();
}


/**
 * D O _ P R E P
 *
 * This is a separate function primarily as a service to REMRT.
 */
void
do_prep(struct rt_i *rtip)
{
    struct bu_vls times;

    RT_CHECK_RTI(rtip);
    if (rtip->needprep) {
	/* Allow lighting model to set up (e.g. lights, materials, etc) */
	view_setup(rtip);

	/* Allow RT library to prepare itself */
	bu_vls_init(&times);
	rt_prep_timer();
	rt_prep_parallel(rtip, npsw);

	(void)rt_get_timer(&times, NULL);
	if (rt_verbosity & VERBOSE_STATS)
	    bu_log("PREP: %s\n", bu_vls_addr(&times));
	bu_vls_free(&times);
    }
    memory_summary();
    if (rt_verbosity & VERBOSE_STATS) {
	bu_log("%s: %d nu, %d cut, %d box (%d empty)\n",
	       rtip->rti_space_partition == RT_PART_NUGRID ?
	       "NUGrid" : "NUBSP",
	       rtip->rti_ncut_by_type[CUT_NUGRIDNODE],
	       rtip->rti_ncut_by_type[CUT_CUTNODE],
	       rtip->rti_ncut_by_type[CUT_BOXNODE],
	       rtip->nempty_cells);
    }
}


/**
 * D O _ F R A M E
 *
 * Do all the actual work to run a frame.
 *
 * Returns -1 on error, 0 if OK.
 */
int
do_frame(int framenumber)
{
    struct bu_vls times;
    char framename[128] = {0};		/* File name to hold current frame */
    struct rt_i *rtip = ap.a_rt_i;
    double utime = 0.0;			/* CPU time used */
    double nutime = 0.0;		/* CPU time used, normalized by ncpu */
    double wallclock = 0.0;		/* # seconds of wall clock time */
    int npix = 0;			/* # of pixel values to be done */
    int lim = 0;
    vect_t work, temp;
    quat_t quat;

    if (rt_verbosity & VERBOSE_FRAMENUMBER)
	bu_log("\n...................Frame %5d...................\n",
	       framenumber);

    /* Compute model RPP, etc */
    do_prep(rtip);

    if (rt_verbosity & VERBOSE_VIEWDETAIL)
	bu_log("Tree: %d solids in %d regions\n", rtip->nsolids, rtip->nregions);

    if (Query_one_pixel) {
	query_rdebug = R_DEBUG;
	query_debug = RT_G_DEBUG;
	rt_g.debug = rdebug = 0;
    }

    if (rtip->nsolids <= 0)
	bu_exit(3, "ERROR: No primitives remaining.\n");

    if (rtip->nregions <= 0)
	bu_exit(3, "ERROR: No regions remaining.\n");

    if (rt_verbosity & VERBOSE_VIEWDETAIL)
	bu_log("Model: X(%g, %g), Y(%g, %g), Z(%g, %g)\n",
	       rtip->mdl_min[X], rtip->mdl_max[X],
	       rtip->mdl_min[Y], rtip->mdl_max[Y],
	       rtip->mdl_min[Z], rtip->mdl_max[Z]);

    /*
     * Perform Grid setup.
     * This may alter cell size or width/height.
     */
    grid_setup();
    /* az/el 0, 0 is when screen +Z is model +X */
    VSET(work, 0, 0, 1);
    MAT3X3VEC(temp, view2model, work);
    bn_ae_vec(&azimuth, &elevation, temp);

    if (rt_verbosity & VERBOSE_VIEWDETAIL)
	bu_log(
	    "View: %g azimuth, %g elevation off of front view\n",
	    azimuth, elevation);
    quat_mat2quat(quat, model2view);
    if (rt_verbosity & VERBOSE_VIEWDETAIL) {
	bu_log("Orientation: %g, %g, %g, %g\n", V4ARGS(quat));
	bu_log("Eye_pos: %g, %g, %g\n", V3ARGS(eye_model));
	bu_log("Size: %gmm\n", viewsize);
#if 0
	/*
	 * This code shows how the model2view matrix can be
	 * reconstructed using the information from the Orientation,
	 * Eye_pos, and Size messages.
	 */
	{
	    mat_t rotscale, xlate;
	    mat_t new;
	    quat_t newquat;

	    bn_mat_print("model2view", model2view);
	    quat_quat2mat(rotscale, quat);
	    rotscale[15] = 0.5 * viewsize;
	    MAT_IDN(xlate);
	    MAT_DELTAS_VEC_NEG(xlate, eye_model);
	    bn_mat_mul(new, rotscale, xlate);
	    bn_mat_print("reconstructed m2v", new);
	    quat_mat2quat(newquat, new);
	    HPRINT("reconstructed orientation:", newquat);
	}
#endif
	bu_log("Grid: (%g, %g) mm, (%d, %d) pixels\n",
	       cell_width, cell_height,
	       width, height);
	bu_log("Beam: radius=%g mm, divergence=%g mm/1mm\n",
	       ap.a_rbeam, ap.a_diverge);
    }

    /* Process -b and ??? options now, for this frame */
    if (pix_start == -1) {
	pix_start = 0;
	pix_end = height * width - 1;
    }
    if (string_pix_start) {
	int xx, yy;
	register char *cp = string_pix_start;

	xx = atoi(cp);
	while (*cp >= '0' && *cp <= '9') cp++;
	while (*cp && (*cp < '0' || *cp > '9')) cp++;
	yy = atoi(cp);
	bu_log("only pixel %d %d\n", xx, yy);
	if (xx * yy >= 0) {
	    pix_start = yy * width + xx;
	    pix_end = pix_start;
	}
    }
    if (string_pix_end) {
	int xx, yy;
	register char *cp = string_pix_end;

	xx = atoi(cp);
	while (*cp >= '0' && *cp <= '9') cp++;
	while (*cp && (*cp < '0' || *cp > '9')) cp++;
	yy = atoi(cp);
	bu_log("ending pixel %d %d\n", xx, yy);
	if (xx * yy >= 0) {
	    pix_end = yy * width + xx;
	}
    }

    /*
     * After the parameters for this calculation have been
     * established, deal with CPU limits and priorities, where
     * appropriate.  Because limits exist, they better be adequate.
     * We assume that the Cray can produce MINRATE pixels/sec on
     * images with extreme amounts of glass & mirrors.
     */
#define MINRATE 65
    npix = width*height*(hypersample+1);
    if ((lim = bu_cpulimit_get()) > 0) {
	bu_cpulimit_set(lim + npix / MINRATE + 100);
    }

    /* Allocate data for pixel map for rerendering of black pixels */
    if (pixmap == NULL) {
	pixmap = (unsigned char*)bu_calloc(sizeof(RGBpixel), width*height, "pixmap allocate");
    }

    /*
     * If this image is unlikely to be for debugging,
     * be gentle to the machine.
     */
    if (!interactive) {
	if (npix > 256*256)
	    bu_nice_set(10);
	else if (npix > 512*512)
	    bu_nice_set(14);
    }

    /*
     * Determine output file name
     * On UNIX only, check to see if this is a "restart".
     */
    if (outputfile != (char *)0) {
	if (framenumber <= 0) {
	    snprintf(framename, 128, "%s", outputfile);
	} else {
	    snprintf(framename, 128, "%s.%d", outputfile, framenumber);
	}

#ifdef HAVE_SYS_STAT_H
	/*
	 * This code allows the computation of a particular frame to a
	 * disk file to be resumed automatically.  This is worthwhile
	 * crash protection.  This use of stat() and fseek() is
	 * UNIX-specific.
	 *
	 * It is not appropriate for the RT "top part" to assume
	 * anything about the data that the view module will be
	 * storing.  Therefore, it is the responsibility of
	 * view_2init() to also detect that some existing data is in
	 * the file, and take appropriate measures (like reading it
	 * in).
	 *
	 * view_2init() can depend on the file being open for both
	 * reading and writing, but must do it's own positioning.
	 */
	{
	    struct stat sb;
	    if (stat(framename, &sb) >= 0 &&
		sb.st_size > 0 &&
		sb.st_size < width*height*sizeof(RGBpixel)) {
		/* File exists, with partial results */
		register int fd;
		if ((fd = open(framename, 2)) < 0 ||
		    (outfp = fdopen(fd, "r+")) == NULL) {
		    perror(framename);
		    if (matflag) return 0;	/* OK */
		    return -1;			/* Bad */
		}
		/* Read existing pix data into the frame buffer */
		if (sb.st_size > 0) {
		    (void)fread(pixmap, 1, (size_t)sb.st_size, outfp);
		}
	    }
	}
#endif

	/* Ordinary case for creating output file */
	if (outfp == NULL) {
	    bif = bu_image_save_open(framename, BU_IMAGE_AUTO_NO_PIX, width, height, 3);
	    if (bif == NULL && (outfp = fopen(framename, "w+b")) == NULL) {
		perror(framename);
		if (matflag) return 0;	/* OK */
		return -1;			/* Bad */
	    }
	}

	if (rt_verbosity & VERBOSE_OUTPUTFILE)
	    bu_log("Output file is '%s' %dx%d pixels\n",
		   framename, width, height);
    }

    /* initialize lighting, may update pix_start */
    view_2init(&ap, framename);

    /* Just while doing the ray-tracing */
    if (R_DEBUG&RDEBUG_RTMEM)
	bu_debug |= (BU_DEBUG_MEM_CHECK|BU_DEBUG_MEM_LOG);

    rtip->nshots = 0;
    rtip->nmiss_model = 0;
    rtip->nmiss_tree = 0;
    rtip->nmiss_solid = 0;
    rtip->nmiss = 0;
    rtip->nhits = 0;
    rtip->rti_nrays = 0;

    if (rt_verbosity & (VERBOSE_LIGHTINFO|VERBOSE_STATS))
	bu_log("\n");
    fflush(stdout);
    fflush(stderr);

    /*
     * Compute the image
     * It may prove desirable to do this in chunks
     */
    rt_prep_timer();

    if (incr_mode) {
	for (incr_level = 1; incr_level <= incr_nlevel; incr_level++) {
	    if (incr_level > 1)
		view_2init(&ap, framename);

	    do_run(0, (1<<incr_level)*(1<<incr_level)-1);
	}
    } else {
	do_run(pix_start, pix_end);

	/* Reset values to full size, for next frame (if any) */
	pix_start = 0;
	pix_end = height*width - 1;
    }
    bu_vls_init(&times);
    utime = rt_get_timer(&times, &wallclock);

    /*
     * End of application.  Done outside of timing section.
     * Typically, writes any remaining results out.
     */
    view_end(&ap);

    /* Stop memory debug printing until next frame, leave full checking on */
    if (R_DEBUG&RDEBUG_RTMEM)
	bu_debug &= ~BU_DEBUG_MEM_LOG;

    /*
     * Certain parallel systems (eg, Alliant) count the entire
     * multi-processor complex as one computer, and charge only once.
     * This matches the desired behavior here.  Other vendors (eg,
     * SGI) count each processor separately, and charge for all of
     * them.  These results need to be normalized.  Otherwise, all we
     * would know is that a given workload takes about the same amount
     * of CPU time, regardless of the number of CPUs.
     */
#if !defined(alliant)
    if (npsw > 1) {
	int avail_cpus;
	int ncpus;

	avail_cpus = bu_avail_cpus();
	if (npsw > avail_cpus) {
	    ncpus = avail_cpus;
	} else {
	    ncpus = npsw;
	}
	nutime = utime / ncpus;			/* compensate */
    } else
#endif
	nutime = utime;

    /* prevent a bogus near-zero time to prevent infinate and
     * near-infinate results without relying on IEEE floating point
     * zero comparison.
     */
    if (NEAR_ZERO(nutime, VDIVIDE_TOL)) {
	bu_log("WARNING:  Raytrace timings are likely to be meaningless\n");
	nutime = VDIVIDE_TOL;
    }

    /*
     * All done.  Display run statistics.
     */
    if (rt_verbosity & VERBOSE_STATS)
	bu_log("SHOT: %s\n", bu_vls_addr(&times));
    bu_vls_free(&times);
    memory_summary();
    if (rt_verbosity & VERBOSE_STATS) {
	bu_log("%ld solid/ray intersections: %ld hits + %ld miss\n",
	       rtip->nshots, rtip->nhits, rtip->nmiss);
	bu_log("pruned %.1f%%:  %ld model RPP, %ld dups skipped, %ld solid RPP\n",
	       rtip->nshots>0?((double)rtip->nhits*100.0)/rtip->nshots:100.0,
	       rtip->nmiss_model, rtip->ndup, rtip->nmiss_solid);
	bu_log("Frame %2d: %10d pixels in %9.2f sec = %12.2f pixels/sec\n",
	       framenumber,
	       width*height, nutime, ((double)(width*height))/nutime);
	bu_log("Frame %2d: %10d rays   in %9.2f sec = %12.2f rays/sec (RTFM)\n",
	       framenumber,
	       rtip->rti_nrays, nutime, ((double)(rtip->rti_nrays))/nutime);
	bu_log("Frame %2d: %10d rays   in %9.2f sec = %12.2f rays/CPU_sec\n",
	       framenumber,
	       rtip->rti_nrays, utime, ((double)(rtip->rti_nrays))/utime);
	bu_log("Frame %2d: %10d rays   in %9.2f sec = %12.2f rays/sec (wallclock)\n",
	       framenumber,
	       rtip->rti_nrays,
	       wallclock, ((double)(rtip->rti_nrays))/wallclock);
    }
    if (bif != NULL)
	bu_image_save_close(bif);
    bif = NULL;
    if (outfp != NULL) {
	/* Protect finished product */
	if (outputfile != (char *)0)
	    (void)bu_fchmod(outfp, 0444);

	(void)fclose(outfp);
	outfp = NULL;
    }

    if (R_DEBUG&RDEBUG_STATS) {
	/* Print additional statistics */
	res_pr();
    }

    bu_log("\n");
    bu_free(pixmap, "pixmap allocate");
    pixmap = (unsigned char *)NULL;
    return 0;		/* OK */
}


/**
 * D O _ A E
 *
 * Compute the rotation specified by the azimuth and elevation
 * parameters.  First, note that these are specified relative to the
 * GIFT "front view", ie, model (X, Y, Z) is view (Z, X, Y): looking
 * down X axis, Y right, Z up.
 *
 * A positive azimuth represents rotating the *eye* around the
 * Y axis, or, rotating the *model* in -Y.
 * A positive elevation represents rotating the *eye* around the
 * X axis, or, rotating the *model* in -X.
 */
void
do_ae(double azim, double elev)
{
    vect_t temp;
    vect_t diag;
    mat_t toEye;
    struct rt_i *rtip = ap.a_rt_i;

    if (rtip->nsolids <= 0)
	bu_exit(EXIT_FAILURE, "ERROR: no primitives active\n");

    if (rtip->nregions <= 0)
	bu_exit(EXIT_FAILURE, "ERROR: no regions active\n");

    if (rtip->mdl_max[X] >= INFINITY) {
	bu_log("do_ae: infinite model bounds? setting a unit minimum\n");
	VSETALL(rtip->mdl_min, -1);
    }
    if (rtip->mdl_max[X] <= -INFINITY) {
	bu_log("do_ae: infinite model bounds? setting a unit maximum\n");
	VSETALL(rtip->mdl_max, 1);
    }

    /*
     * Enlarge the model RPP just slightly, to avoid nasty effects
     * with a solid's face being exactly on the edge NOTE: This code
     * is duplicated out of librt/tree.c/rt_prep(), and has to appear
     * here to enable the viewsize calculation to match the final RPP.
     */
    rtip->mdl_min[X] = floor(rtip->mdl_min[X]);
    rtip->mdl_min[Y] = floor(rtip->mdl_min[Y]);
    rtip->mdl_min[Z] = floor(rtip->mdl_min[Z]);
    rtip->mdl_max[X] = ceil(rtip->mdl_max[X]);
    rtip->mdl_max[Y] = ceil(rtip->mdl_max[Y]);
    rtip->mdl_max[Z] = ceil(rtip->mdl_max[Z]);

    MAT_IDN(Viewrotscale);
    bn_mat_angles(Viewrotscale, 270.0+elev, 0.0, 270.0-azim);

    /* Look at the center of the model */
    MAT_IDN(toEye);
    toEye[MDX] = -((rtip->mdl_max[X]+rtip->mdl_min[X])/2.0);
    toEye[MDY] = -((rtip->mdl_max[Y]+rtip->mdl_min[Y])/2.0);
    toEye[MDZ] = -((rtip->mdl_max[Z]+rtip->mdl_min[Z])/2.0);

    /* Fit a sphere to the model RPP, diameter is viewsize, unless
     * viewsize command used to override.
     */
    if (viewsize <= 0) {
	VSUB2(diag, rtip->mdl_max, rtip->mdl_min);
	viewsize = MAGNITUDE(diag);
	if (aspect > 1) {
	    /* don't clip any of the image when autoscaling */
	    viewsize *= aspect;
	}
    }
    Viewrotscale[15] = 0.5*viewsize;	/* Viewscale */
    bn_mat_mul(model2view, Viewrotscale, toEye);
    bn_mat_inv(view2model, model2view);
    VSET(temp, 0, 0, eye_backoff);
    MAT4X3PNT(eye_model, view2model, temp);
}


/**
 * R E S _ P R
 */
void
res_pr(void)
{
    register struct resource *res;
    register int i;

    fprintf(stderr, "\nResource use summary, by processor:\n");
    res = &resource[0];
    for (i=0; i<npsw; i++, res++) {
	fprintf(stderr, "---CPU %d:\n", i);
	if (res->re_magic != RESOURCE_MAGIC) {
	    fprintf(stderr, "Bad magic number!!\n");
	    continue;
	}
	fprintf(stderr, "seg       len=%10ld get=%10ld free=%10ld\n",
		res->re_seglen, res->re_segget, res->re_segfree);
	fprintf(stderr, "partition len=%10ld get=%10ld free=%10ld\n",
		res->re_partlen, res->re_partget, res->re_partfree);
#if 0
	fprintf(stderr, "bitv_elem len=%10ld get=%10ld free=%10ld\n",
		res->re_bitvlen, res->re_bitvget, res->re_bitvfree);
#endif
	fprintf(stderr, "boolstack len=%10ld\n",
		res->re_boolslen);
    }
}


/**
 * Command table for RT control script language
 */
struct command_tab rt_cmdtab[] = {
    {"start", "frame number", "start a new frame",
     cm_start,	2, 2},
    {"viewsize", "size in mm", "set view size",
     cm_vsize,	2, 2},
    {"eye_pt", "xyz of eye", "set eye point",
     cm_eyept,	4, 4},
    {"lookat_pt", "x y z [yflip]", "set eye look direction, in X-Y plane",
     cm_lookat_pt,	4, 5},
    {"viewrot", "4x4 matrix", "set view direction from matrix",
     cm_vrot,	17, 17},
    {"orientation", "quaturnion", "set view direction from quaturnion",
     cm_orientation,	5, 5},
    {"end", 	"", "end of frame setup, begin raytrace",
     cm_end,		1, 1},
    {"multiview", "", "produce stock set of views",
     cm_multiview,	1, 1},
    {"anim", 	"path type args", "specify articulation animation",
     cm_anim,	4, 999},
    {"tree", 	"treetop(s)", "specify alternate list of tree tops",
     cm_tree,	1, 999},
    {"clean", "", "clean articulation from previous frame",
     cm_clean,	1, 1},
    {"_closedb", "", "Close .g database, (for memory debugging)",
     cm_closedb,	1, 1},
    {"set", 	"", "show or set parameters",
     cm_set,		1, 999},
    {"ae", "azim elev", "specify view as azim and elev, in degrees",
     cm_ae,		3, 3},
    {"opt", "-flags", "set flags, like on command line",
     cm_opt,		2, 999},
    {(char *)0, (char *)0, (char *)0,
     0,		0, 0	/* END */}
};


/*
 * Local Variables:
 * mode: C
 * tab-width: 8
 * indent-tabs-mode: t
 * c-file-style: "stroustrup"
 * End:
 * ex: shiftwidth=4 tabstop=8
 */