Index: src/rt/worker.c
===================================================================
--- src/rt/worker.c	(revision 63934)
+++ src/rt/worker.c	(working copy)
@@ -153,12 +153,10 @@
 
 
 void
-do_pixel(int cpu, int pat_num, int pixelnum)
+do_pixel(int cpu, int pat_num, int pixelnum, struct xray *rays)
 {
     struct application a;
     struct pixel_ext pe;
-    vect_t stereo_point;		/* Ref point on eye or view plane */
-    vect_t point;		/* Ref point on eye or view plane */
     vect_t colorsum = {(fastf_t)0.0, (fastf_t)0.0, (fastf_t)0.0};
     int samplenum = 0;
     static const double one_over_255 = 1.0 / 255.0;
@@ -240,8 +238,7 @@
 	}
     }
 
-    /* our starting point, used for non-jitter */
-    VJOIN2 (point, viewbase_model, a.a_x, dx_model, a.a_y, dy_model);
+    a.a_ray = rays[a.a_y * width + a.a_x];
 
     /* not tracing the corners of a prism by default */
     a.a_pixelext=(struct pixel_ext *)NULL;
@@ -261,7 +258,7 @@
 	/****************/
 
 	if (jitter & JITTER_CELL) {
-	    jitter_start_pt(point, &a, samplenum, pat_num);
+	    jitter_start_pt(rt_perspective > 0.0 ? a.a_ray.r_dir : a.a_ray.r_pt, &a, samplenum, pat_num);
 	}
 
 	if (a.a_rt_i->rti_prismtrace) {
@@ -275,9 +272,6 @@
 	}
 
 	if (rt_perspective > 0.0) {
-	    VSUB2(a.a_ray.r_dir, point, eye_model);
-	    VUNITIZE(a.a_ray.r_dir);
-	    VMOVE(a.a_ray.r_pt, eye_model);
 	    if (a.a_rt_i->rti_prismtrace) {
 		VSUB2(pe.corner[0].r_dir, pe.corner[0].r_pt, eye_model);
 		VSUB2(pe.corner[1].r_dir, pe.corner[1].r_pt, eye_model);
@@ -285,9 +279,6 @@
 		VSUB2(pe.corner[3].r_dir, pe.corner[3].r_pt, eye_model);
 	    }
 	} else {
-	    VMOVE(a.a_ray.r_pt, point);
-	    VMOVE(a.a_ray.r_dir, APP.a_ray.r_dir);
-
 	    if (a.a_rt_i->rti_prismtrace) {
 		VMOVE(pe.corner[0].r_dir, a.a_ray.r_dir);
 		VMOVE(pe.corner[1].r_dir, a.a_ray.r_dir);
@@ -309,13 +300,12 @@
 
 	    right = CRT_BLEND(a.a_color);
 
-	    VSUB2(stereo_point, point, left_eye_delta);
 	    if (rt_perspective > 0.0) {
-		VSUB2(a.a_ray.r_dir, stereo_point, eye_model);
+		VSUB2(a.a_ray.r_dir, a.a_ray.r_dir, left_eye_delta);
 		VUNITIZE(a.a_ray.r_dir);
-		VADD2(a.a_ray.r_pt, eye_model, left_eye_delta);
+		VADD2(a.a_ray.r_pt, a.a_ray.r_pt, left_eye_delta);
 	    } else {
-		VMOVE(a.a_ray.r_pt, stereo_point);
+		VSUB2(a.a_ray.r_pt, a.a_ray.r_pt, left_eye_delta);
 	    }
 	    a.a_level = 0;		/* recursion level */
 	    a.a_purpose = "left eye ray";
@@ -341,7 +331,7 @@
 	    /**********************/
 
 	    if (jitter & JITTER_CELL) {
-		jitter_start_pt(point, &a, samplenum, pat_num);
+		jitter_start_pt(rt_perspective > 0.0 ? a.a_ray.r_dir : a.a_ray.r_pt, &a, samplenum, pat_num);
 	    }
 
 	    if (a.a_rt_i->rti_prismtrace) {
@@ -355,9 +345,6 @@
 	    }
 
 	    if (rt_perspective > 0.0) {
-		VSUB2(a.a_ray.r_dir, point, eye_model);
-		VUNITIZE(a.a_ray.r_dir);
-		VMOVE(a.a_ray.r_pt, eye_model);
 		if (a.a_rt_i->rti_prismtrace) {
 		    VSUB2(pe.corner[0].r_dir, pe.corner[0].r_pt, eye_model);
 		    VSUB2(pe.corner[1].r_dir, pe.corner[1].r_pt, eye_model);
@@ -365,9 +352,6 @@
 		    VSUB2(pe.corner[3].r_dir, pe.corner[3].r_pt, eye_model);
 		}
 	    } else {
-		VMOVE(a.a_ray.r_pt, point);
-		VMOVE(a.a_ray.r_dir, APP.a_ray.r_dir);
-
 		if (a.a_rt_i->rti_prismtrace) {
 		    VMOVE(pe.corner[0].r_dir, a.a_ray.r_dir);
 		    VMOVE(pe.corner[1].r_dir, a.a_ray.r_dir);
@@ -389,13 +373,12 @@
 
 		right = CRT_BLEND(a.a_color);
 
-		VSUB2(stereo_point, point, left_eye_delta);
 		if (rt_perspective > 0.0) {
-		    VSUB2(a.a_ray.r_dir, stereo_point, eye_model);
+		    VSUB2(a.a_ray.r_dir, a.a_ray.r_dir, left_eye_delta);
 		    VUNITIZE(a.a_ray.r_dir);
-		    VADD2(a.a_ray.r_pt, eye_model, left_eye_delta);
+		    VADD2(a.a_ray.r_pt, a.a_ray.r_pt, left_eye_delta);
 		} else {
-		    VMOVE(a.a_ray.r_pt, stereo_point);
+		    VSUB2(a.a_ray.r_pt, a.a_ray.r_pt, left_eye_delta);
 		}
 		a.a_level = 0;		/* recursion level */
 		a.a_purpose = "left eye ray";
@@ -460,11 +443,12 @@
  * For a general-purpose version, see LIBRT rt_shoot_many_rays()
  */
 void
-worker(int cpu, void *UNUSED(arg))
+worker(int cpu, void *arg)
 {
     int pixel_start;
     int pixelnum;
     int pat_num = -1;
+    struct xray *rays = (struct xray *) arg;
 
     /* The more CPUs at work, the bigger the bites we take */
     if (per_processor_chunk <= 0) per_processor_chunk = npsw;
@@ -510,7 +494,7 @@
 		if (pixelnum < last_pixel)
 		    return;
 
-		do_pixel(cpu, pat_num, pixelnum);
+		do_pixel(cpu, pat_num, pixelnum, rays);
 	    }
 	}
     } else if (random_mode) {
@@ -521,7 +505,7 @@
 	       with multi-threaded RNG */
 	    pixelnum = rand()*1.0/RAND_MAX*(last_pixel + 1);
 	    if (pixelnum >= last_pixel) pixelnum = last_pixel;
-	    do_pixel(cpu, pat_num, pixelnum);
+	    do_pixel(cpu, pat_num, pixelnum, rays);
 	}
 
     } else {
@@ -539,7 +523,7 @@
 		if (pixelnum > last_pixel)
 		    return;
 
-		do_pixel(cpu, pat_num, pixelnum);
+		do_pixel(cpu, pat_num, pixelnum, rays);
 	    }
 	}
     }
@@ -700,6 +684,55 @@
 {
     int cpu;
 
+    struct xray center_ray;
+    struct xrays xrays;
+    struct xrays *ray;
+    int i = 0;
+
+    struct xray *rays = (struct xray *)bu_calloc(1, sizeof(struct xray), "ray bundle");
+
+    BU_LIST_INIT(&(xrays.l));
+
+    if (rt_perspective > 0.0) {
+	fastf_t pw, ph;
+	VMOVE(center_ray.r_pt, eye_model);
+	/* The center ray is slightly off center in order to
+	 * accommodate rt_gen_frustum()'s fully-closed set behavior
+	 * (see above). */
+	VJOIN2(center_ray.r_dir, viewbase_model, width / 2 - 0.5, dx_model, height / 2 - 0.5, dy_model);
+	VSUB2(center_ray.r_dir, center_ray.r_dir, eye_model);
+	/* Calculating the perspective angle to be used is somewhat
+	 * complicated due to rts half-open set behavior.  This
+	 * re-derives the perspective from the width and the magnitude
+	 * of center_ray, making the width one cell narrower, so that
+	 * it encloses (with the appropriate center point) the
+	 * half-open set of rt as a fully closed set.
+	 */
+	pw = atan((width - 1) / 2.0 * cell_width / MAGNITUDE(center_ray.r_dir));
+	ph = aspect * pw;
+	rt_gen_frustum(&xrays, &center_ray, dx_model, dy_model, pw, ph, width, height);
+    } else {
+	vect_t w, h;
+	VSCALE(w, dx_model, width / 2 - 0.5);
+	VSCALE(h, dy_model, height / 2 - 0.5);
+	VJOIN2(center_ray.r_pt, viewbase_model, width / 2 - 0.5, dx_model, height / 2 - 0.5, dy_model);
+	VMOVE(center_ray.r_dir, APP.a_ray.r_dir);
+	rt_gen_rect(&xrays, &center_ray, w, h, MAGNITUDE(dx_model), MAGNITUDE(dy_model));
+    }
+
+    i = bu_list_len(&(xrays.l));
+    if (i == 0) {
+	/* No rays to shoot, so just exit */
+	return;
+    }
+    rays = (struct xray *)bu_calloc(i, sizeof(struct xray), "ray bundle");
+    i--;
+    while (BU_LIST_WHILE(ray, xrays, &(xrays.l))) {
+	BU_LIST_DEQUEUE(&(ray->l));
+	rays[i--] = ray->ray;
+	bu_free(ray, "bundled ray");
+    }
+
 #ifdef USE_FORKED_THREADS
     int pid, wpid;
     int waitret;
@@ -723,7 +756,7 @@
 	 * SERIAL case -- one CPU does all the work.
 	 */
 	npsw = 1;
-	worker(0, NULL);
+	worker(0, (void *)rays);
     } else {
 	/*
 	 * Parallel case.
@@ -747,7 +780,7 @@
 	} else if (pid == 0) {
 #endif
 
-	    bu_parallel(worker, npsw, NULL);
+	    bu_parallel(worker, npsw, (void *)rays);
 
 #ifdef USE_FORKED_THREADS
 	    /* send raytrace instance data back to the parent */
Index: src/librt/mkbundle.c
===================================================================
--- src/librt/mkbundle.c	(revision 63934)
+++ src/librt/mkbundle.c	(working copy)
@@ -226,9 +231,10 @@
 
     fastf_t x, y;
 
-    fastf_t a_length = tan(a_theta);
-    fastf_t b_length = tan(b_theta);
+    fastf_t a_length = tan(a_theta) * MAGNITUDE(center_ray->r_dir);
+    fastf_t b_length = tan(b_theta) * MAGNITUDE(center_ray->r_dir);
     fastf_t a_inc = 2 * a_length / (a_num - 1);
+    fastf_t b_inc = 2 * b_length / (b_num - 1);
 
     vect_t a_dir, b_dir;
 
@@ -245,7 +251,7 @@
     /* This adds BN_TOL_DIST to the *_length variables in the
      * condition because in some cases, floating-point problems can
      * make extremely close numbers compare incorrectly. */
-    for (y = -b_length; y <= b_length + BN_TOL_DIST;) {
+    for (y = -b_length; y <= b_length + BN_TOL_DIST; y += b_inc) {
 	for (x = -a_length; x <= a_length + BN_TOL_DIST; x += a_inc) {
 	    BU_ALLOC(xrayp, struct xrays);
 	    VMOVE(xrayp->ray.r_pt, start);
@@ -288,8 +294,8 @@
     VMOVE(b_dir, b_vec);
     VUNITIZE(b_dir);
 
-    for (y = -b_length; y <= b_length; y += db) {
-	for (x = -a_length; x <= a_length; x += da) {
+    for (y = -b_length; y <= b_length + BN_TOL_DIST; y += db) {
+	for (x = -a_length; x <= a_length +  BN_TOL_DIST; x += da) {
 	    BU_ALLOC(xrayp, struct xrays);
 	    VJOIN2(xrayp->ray.r_pt, orig_start, x, a_dir, y, b_dir);
 	    VMOVE(xrayp->ray.r_dir, dir);
