00001 /*
00002  * vp_view.c
00003  *
00004  * Routines to compute quantities derived from the view transformation.
00005  *
00006  * Copyright (c) 1994 The Board of Trustees of The Leland Stanford
00007  * Junior University.  All rights reserved.
00008  *
00009  * Permission to use, copy, modify and distribute this software and its
00010  * documentation for any purpose is hereby granted without fee, provided
00011  * that the above copyright notice and this permission notice appear in
00012  * all copies of this software and that you do not sell the software.
00013  * Commercial licensing is available by contacting the author.
00014  * 
00015  * THE SOFTWARE IS PROVIDED "AS IS" AND WITHOUT WARRANTY OF ANY KIND,
00016  * EXPRESS, IMPLIED OR OTHERWISE, INCLUDING WITHOUT LIMITATION, ANY
00017  * WARRANTY OF MERCHANTABILITY OR FITNESS FOR A PARTICULAR PURPOSE.
00018  *
00019  * Author:
00020  *    Phil Lacroute
00021  *    Computer Systems Laboratory
00022  *    Electrical Engineering Dept.
00023  *    Stanford University
00024  */
00025 
00026 /*
00027  * $Date: 2001/12/17 16:16:23 $
00028  * $Revision: 1.1 $
00029  */
00030 
00031 #include "vp_global.h"
00032 
00033 static int FactorAffineView ANSI_ARGS((vpContext *vpc, vpMatrix4 vm));
00034 static int FactorPerspectiveView ANSI_ARGS((vpContext *vpc, vpMatrix4 vm));
00035 static void ComputeAffineOpacityCorrection ANSI_ARGS((vpContext *vpc,
00036     double shear_i, double shear_j, float table[VP_OPACITY_MAX+1]));
00037 static void CheckRenderBuffers ANSI_ARGS((vpContext *vpc));
00038 static void ComputeLightViewTransform ANSI_ARGS((vpContext *vpc,vpMatrix4 vm));
00039 static int FactorLightView ANSI_ARGS((vpContext *vpc, vpMatrix4 vm));
00040 static void CheckShadowBuffer ANSI_ARGS((vpContext *vpc));
00041 
00042 /*
00043  * VPFactorView
00044  *
00045  * Factor the viewing matrix.
00046  */
00047 
00048 vpResult
00049 VPFactorView(vpc)
00050 vpContext *vpc;
00051 {
00052     vpMatrix4 vm;
00053     int retcode;
00054 
00055     if (vpc->factored_view_ready) {
00056         CheckRenderBuffers(vpc);
00057         return(VP_OK);
00058     }
00059 
00060     /* compute the overall view transformation */
00061     VPComputeViewTransform(vpc, vm);
00062 
00063     /* check if transformation is affine and factor it */
00064     if (fabs(vm[3][0]) < VP_EPS && fabs(vm[3][1]) < VP_EPS &&
00065         fabs(vm[3][2]) < VP_EPS && fabs(vm[3][3]-1.) < VP_EPS) {
00066         if ((retcode = FactorAffineView(vpc, vm)) != VP_OK)
00067             return(retcode);
00068         ComputeAffineOpacityCorrection(vpc, vpc->shear_i, vpc->shear_j,
00069                                        vpc->affine_opac_correct);
00070     } else {
00071         FactorPerspectiveView(vpc, vm);
00072     }
00073     CheckRenderBuffers(vpc);
00074 
00075     /* compute viewing transformation from the point of view of the light
00076        source (for calculating shadows) */
00077     if ((retcode = VPCheckShadows(vpc)) != VP_OK)
00078         return(retcode);
00079     if (vpc->enable_shadows) {
00080         ComputeLightViewTransform(vpc, vm);
00081         if ((retcode = FactorLightView(vpc, vm)) != VP_OK)
00082             return(retcode);
00083         ComputeAffineOpacityCorrection(vpc, vpc->shadow_shear_i,
00084                                        vpc->shadow_shear_j,
00085                                        vpc->shadow_opac_correct);
00086         CheckShadowBuffer(vpc);
00087     }
00088     return(VP_OK);
00089 }
00090 
00091 /*
00092  * VPComputeViewTransform
00093  *
00094  * Compute the overall view transformation.
00095  */
00096 
00097 void
00098 VPComputeViewTransform(vpc, vm)
00099 vpContext *vpc;
00100 vpMatrix4 vm;   /* storage for result */
00101 {
00102     vpMatrix4 prematrix;                /* transform volume to unit cube */
00103     vpMatrix4 viewportm;                /* viewport matrix */
00104     vpMatrix4 tmp1m, tmp2m, tmp3m;      /* temporary matrices */
00105     int maxdim;
00106     double scale;
00107 
00108     /* compute prematrix */
00109     vpIdentity4(prematrix);
00110     maxdim = vpc->xlen;
00111     if (vpc->ylen > maxdim)
00112         maxdim = vpc->ylen;
00113     if (vpc->zlen > maxdim)
00114         maxdim = vpc->zlen;
00115     scale = 1. / (double)maxdim;
00116     prematrix[0][0] = scale;
00117     prematrix[1][1] = scale;
00118     prematrix[2][2] = scale;
00119     prematrix[0][3] = -vpc->xlen * scale * 0.5;
00120     prematrix[1][3] = -vpc->ylen * scale * 0.5;
00121     prematrix[2][3] = -vpc->zlen * scale * 0.5;
00122 
00123     /* compute the viewport matrix */
00124     vpIdentity4(viewportm);
00125     viewportm[0][0] = 0.5 * vpc->image_width;
00126     viewportm[0][3] = 0.5 * vpc->image_width;
00127     viewportm[1][1] = 0.5 * vpc->image_height;
00128     viewportm[1][3] = 0.5 * vpc->image_height;
00129     viewportm[2][2] = -0.5;     /* minus sign: switch to left-handed coords. */
00130     viewportm[2][3] = 0.5;
00131 
00132     /* compute the view matrix */
00133     vpMatrixMult4(tmp1m, vpc->transforms[VP_MODEL], prematrix);
00134     vpMatrixMult4(tmp2m, vpc->transforms[VP_VIEW], tmp1m);
00135     vpMatrixMult4(tmp3m, vpc->transforms[VP_PROJECT], tmp2m);
00136     vpMatrixMult4(vm, viewportm, tmp3m);
00137 }
00138 
00139 /*
00140  * FactorAffineView
00141  *
00142  * Factor an affine viewing matrix into two parts:
00143  *  1) A shear and translation which map object coordinates into
00144  *     intermediate image coordinates.
00145  *  2) An affine warp which transforms intermediate image coordinates to
00146  *     image coordinates.
00147  * Return value is a result code.
00148  */
00149 
00150 static int
00151 FactorAffineView(vpc, vm)
00152 vpContext *vpc;
00153 vpMatrix4 vm;
00154 {
00155     vpMatrix4 p;                /* permutation matrix */
00156     vpMatrix4 pvm;              /* permutation of the viewing matrix */
00157     int icount, jcount, kcount; /* dimensions of volume in rotated space */
00158     vpVector4 xobj, yobj, zobj;
00159     vpVector4 xim, yim, zim;
00160     double x, y, z, denom;
00161 
00162     Debug((vpc, VPDEBUG_VIEW, "FactorAffineView\n"));
00163 
00164     vpc->affine_view = 1;
00165 
00166     /*
00167      * Transform the unit x, y and z object-coordinate vectors into image
00168      * space and see which one is most aligned with the view direction
00169      * (which is the z-axis in image coordinates).
00170      */
00171     vpSetVector4(xobj, 1., 0., 0., 0.);
00172     vpSetVector4(yobj, 0., 1., 0., 0.);
00173     vpSetVector4(zobj, 0., 0., 1., 0.);
00174     vpMatrixVectorMult4(xim, vm, xobj);
00175     vpMatrixVectorMult4(yim, vm, yobj);
00176     vpMatrixVectorMult4(zim, vm, zobj);
00177     x = fabs((vpNormalize3(xim) == VPERROR_SINGULAR) ? 0. : xim[2]);
00178     y = fabs((vpNormalize3(yim) == VPERROR_SINGULAR) ? 0. : yim[2]);
00179     z = fabs((vpNormalize3(zim) == VPERROR_SINGULAR) ? 0. : zim[2]);
00180     if (x >= y) {
00181         if (x >= z) {
00182             vpc->best_view_axis = VP_X_AXIS;
00183         } else {
00184             vpc->best_view_axis = VP_Z_AXIS;
00185         }
00186     } else {
00187         if (y >= z) {
00188             vpc->best_view_axis = VP_Y_AXIS;
00189         } else {
00190             vpc->best_view_axis = VP_Z_AXIS;
00191         }
00192     }
00193     switch (vpc->axis_override) {
00194     case VP_X_AXIS:
00195         if (x < VP_EPS)
00196             return(VPSetError(vpc, VPERROR_SINGULAR));
00197         vpc->best_view_axis = VP_X_AXIS;
00198         break;
00199     case VP_Y_AXIS:
00200         if (y < VP_EPS)
00201             return(VPSetError(vpc, VPERROR_SINGULAR));
00202         vpc->best_view_axis = VP_Y_AXIS;
00203         break;
00204     case VP_Z_AXIS:
00205         if (z < VP_EPS)
00206             return(VPSetError(vpc, VPERROR_SINGULAR));
00207         vpc->best_view_axis = VP_Z_AXIS;
00208         break;
00209     default:
00210         break;
00211     }
00212 
00213     /* permute the rows of the viewing matrix so that the third axis is
00214        most parallel to the viewing direction */
00215     bzero(p, sizeof(vpMatrix4));
00216     switch (vpc->best_view_axis) {
00217     case VP_X_AXIS:
00218         p[0][2] = 1.;
00219         p[1][0] = 1.;
00220         p[2][1] = 1.;
00221         p[3][3] = 1.;
00222         icount = vpc->ylen;
00223         jcount = vpc->zlen;
00224         kcount = vpc->xlen;
00225         break;
00226     case VP_Y_AXIS:
00227         p[0][1] = 1.;
00228         p[1][2] = 1.;
00229         p[2][0] = 1.;
00230         p[3][3] = 1.;
00231         icount = vpc->zlen;
00232         jcount = vpc->xlen;
00233         kcount = vpc->ylen;
00234         break;
00235     case VP_Z_AXIS:
00236         p[0][0] = 1.;
00237         p[1][1] = 1.;
00238         p[2][2] = 1.;
00239         p[3][3] = 1.;
00240         icount = vpc->xlen;
00241         jcount = vpc->ylen;
00242         kcount = vpc->zlen;
00243         break;
00244     }
00245     vpMatrixMult4(pvm, vm, p);
00246 
00247     /* compute the shear coefficients */
00248     denom = pvm[0][0]*pvm[1][1] - pvm[0][1]*pvm[1][0];
00249     if (fabs(denom) < VP_EPS)
00250         return(VPSetError(vpc, VPERROR_SINGULAR));
00251     vpc->shear_i = (pvm[0][2]*pvm[1][1] - pvm[0][1]*pvm[1][2]) / denom;
00252     vpc->shear_j = (pvm[0][0]*pvm[1][2] - pvm[1][0]*pvm[0][2]) / denom;
00253     if (pvm[2][0]*vpc->shear_i + pvm[2][1]*vpc->shear_j - pvm[2][2] > 0)
00254         vpc->reverse_slice_order = 0;
00255     else
00256         vpc->reverse_slice_order = 1;
00257 
00258     /* compute the intermediate image size */
00259     vpc->intermediate_width = icount + 1 + (int)ceil((kcount-1)*
00260                               fabs(vpc->shear_i));
00261     vpc->intermediate_height = jcount + 1 + (int)ceil((kcount-1)*
00262                                fabs(vpc->shear_j));
00263 
00264     /* compute the translation coefficients */
00265     if (vpc->shear_i >= 0.)
00266         vpc->trans_i = 1.;
00267     else
00268         vpc->trans_i = 1. - vpc->shear_i * (kcount - 1);
00269     if (vpc->shear_j >= 0.)
00270         vpc->trans_j = 1.;
00271     else
00272         vpc->trans_j = 1. - vpc->shear_j * (kcount - 1);
00273 
00274     /* compute the depth coefficients */
00275     vpc->depth_di = pvm[2][0];
00276     vpc->depth_dj = pvm[2][1];
00277     vpc->depth_dk = pvm[2][2];
00278     vpc->depth_000 = pvm[2][3];
00279 
00280     /* compute the mapping from compositing space to image space */
00281     vpc->warp_2d[0][0] = pvm[0][0];
00282     vpc->warp_2d[0][1] = pvm[0][1];
00283     vpc->warp_2d[0][2] = pvm[0][3] - pvm[0][0]*vpc->trans_i -
00284                          pvm[0][1]*vpc->trans_j;
00285     vpc->warp_2d[1][0] = pvm[1][0];
00286     vpc->warp_2d[1][1] = pvm[1][1];
00287     vpc->warp_2d[1][2] = pvm[1][3] - pvm[1][0]*vpc->trans_i - 
00288                          pvm[1][1]*vpc->trans_j;
00289     vpc->warp_2d[2][0] = 0.;
00290     vpc->warp_2d[2][1] = 0.;
00291     vpc->warp_2d[2][2] = 1.;
00292 
00293     vpc->factored_view_ready = 1;
00294 
00295     Debug((vpc, VPDEBUG_VIEW, "  best_view_axis: %c%c\n", 
00296            vpc->reverse_slice_order ? '-' : '+',
00297            vpc->best_view_axis == VP_X_AXIS ? 'x' :
00298            (vpc->best_view_axis == VP_Y_AXIS ? 'y' : 'z')));
00299     Debug((vpc, VPDEBUG_VIEW, "  shear factors: %g %g\n",
00300            vpc->shear_i, vpc->shear_j));
00301     Debug((vpc, VPDEBUG_VIEW, "  translation: %g %g\n",
00302            vpc->trans_i, vpc->trans_j));
00303     Debug((vpc, VPDEBUG_VIEW, "  depth: d000: %g\n", vpc->depth_000));
00304     Debug((vpc, VPDEBUG_VIEW, "         di:   %g\n", vpc->depth_di));
00305     Debug((vpc, VPDEBUG_VIEW, "         dj:   %g\n", vpc->depth_dj));
00306     Debug((vpc, VPDEBUG_VIEW, "         dk:   %g\n", vpc->depth_dk));
00307     Debug((vpc, VPDEBUG_VIEW, "  intermediate image size: %d %d\n",
00308            vpc->intermediate_width, vpc->intermediate_height));
00309     return(VP_OK);
00310 }
00311 
00312 /*
00313  * FactorPerspectiveView
00314  *
00315  * Factor a perspective view matrix into two parts:
00316  *  1) A shear, translation and scale which map object coordinates into
00317  *     intermediate image coordinates.
00318  *  2) A perspective warp which transforms intermediate image coordinates to
00319  *     image coordinates.
00320  */
00321 
00322 static int
00323 FactorPerspectiveView(vpc, vm)
00324 vpContext *vpc;
00325 vpMatrix4 vm;
00326 {
00327     vpc->affine_view = 0;
00328     return(VP_OK);
00329 
00330 #ifdef notdef
00331     Matrix4 p;          /* permutation matrix */
00332     Matrix4 pvm;        /* permutation of the viewing matrix */
00333     Matrix4 m2d;        /* final warp */
00334     Matrix4 t;
00335     double alpha1, alpha2, alpha3, alpha4;
00336     int icount, jcount, kcount; /* dimensions of volume in rotated space */
00337     Vector4 xobj, yobj, zobj;
00338     double x, y, z, denom;
00339     double i0, j0, i1, j1;
00340     double imin, imax, jmin, jmax;
00341 
00342     Debug((DEBUG_VIEW, "FactorPerspectiveView\n"));
00343 
00344     rbuf->perspective_proj = 1;
00345 
00346     /*
00347      * Transform the unit x, y and z object-coordinate vectors into image
00348      * space and see which one is most aligned with the view direction
00349      * (which is the z-axis in image coordinates).
00350      */
00351     xobj[0] = 1.; xobj[1] = 0.; xobj[2] = 0.; xobj[3] = 0.;
00352     yobj[0] = 0.; yobj[1] = 1.; yobj[2] = 0.; yobj[3] = 0.;
00353     zobj[0] = 0.; zobj[1] = 0.; zobj[2] = 1.; zobj[3] = 0.;
00354     TransformVector4(xobj, view->view_matrix);
00355     TransformVector4(yobj, view->view_matrix);
00356     TransformVector4(zobj, view->view_matrix);
00357     /* normalize each vector to unit length and compare the absolute value
00358        of the z component; note that the w component drops out */
00359     xobj[2] = fabs(xobj[2]);
00360     yobj[2] = fabs(yobj[2]);
00361     zobj[2] = fabs(zobj[2]);
00362     x = (xobj[2] < EPS) ? 0. :
00363         (xobj[2] / sqrt(xobj[0]*xobj[0] + xobj[1]*xobj[1] + xobj[2]*xobj[2]));
00364     y = (yobj[2] < EPS) ? 0. :
00365         (yobj[2] / sqrt(yobj[0]*yobj[0] + yobj[1]*yobj[1] + yobj[2]*yobj[2]));
00366     z = (zobj[2] < EPS) ? 0. :
00367         (zobj[2] / sqrt(zobj[0]*zobj[0] + zobj[1]*zobj[1] + zobj[2]*zobj[2]));
00368     if (x >= y) {
00369         if (x >= z) {
00370             rbuf->best_view_axis = VP_XAXIS;
00371         } else {
00372             rbuf->best_view_axis = VP_ZAXIS;
00373         }
00374     } else {
00375         if (y >= z) {
00376             rbuf->best_view_axis = VP_YAXIS;
00377         } else {
00378             rbuf->best_view_axis = VP_ZAXIS;
00379         }
00380     }
00381 
00382     /* permute the rows of the viewing matrix so that the third axis is
00383        most parallel to the viewing direction */
00384     bzero(p, sizeof(Matrix4));
00385     switch (rbuf->best_view_axis) {
00386     case VP_XAXIS:
00387         p[0][2] = 1.;
00388         p[1][0] = 1.;
00389         p[2][1] = 1.;
00390         p[3][3] = 1.;
00391         icount = ylen;
00392         jcount = zlen;
00393         kcount = xlen;
00394         break;
00395     case VP_YAXIS:
00396         p[0][1] = 1.;
00397         p[1][2] = 1.;
00398         p[2][0] = 1.;
00399         p[3][3] = 1.;
00400         icount = zlen;
00401         jcount = xlen;
00402         kcount = ylen;
00403         break;
00404     case VP_ZAXIS:
00405         p[0][0] = 1.;
00406         p[1][1] = 1.;
00407         p[2][2] = 1.;
00408         p[3][3] = 1.;
00409         icount = xlen;
00410         jcount = ylen;
00411         kcount = zlen;
00412         break;
00413     default:
00414         VPBug("wierd value for best_view_axis in FactorPerspectiveView\n");
00415     }
00416     MatrixMult4(pvm, view->view_matrix, p);
00417 
00418     /* compute the magic alpha coefficients */
00419     alpha1 = pvm[3][1] * (pvm[0][3]*pvm[1][2] - pvm[0][2]*pvm[1][3]) +
00420              pvm[3][2] * (pvm[0][1]*pvm[1][3] - pvm[0][3]*pvm[1][1]) +
00421              pvm[3][3] * (pvm[0][2]*pvm[1][1] - pvm[0][1]*pvm[1][2]);
00422     alpha2 = pvm[3][0] * (pvm[0][2]*pvm[1][3] - pvm[0][3]*pvm[1][2]) +
00423              pvm[3][2] * (pvm[0][3]*pvm[1][0] - pvm[0][0]*pvm[1][3]) +
00424              pvm[3][3] * (pvm[0][0]*pvm[1][2] - pvm[0][2]*pvm[1][0]);
00425     alpha3 = pvm[3][0] * (pvm[0][1]*pvm[1][2] - pvm[0][2]*pvm[1][1]) +
00426              pvm[3][1] * (pvm[0][2]*pvm[1][0] - pvm[0][0]*pvm[1][2]) +
00427              pvm[3][2] * (pvm[0][0]*pvm[1][1] - pvm[0][1]*pvm[1][0]);
00428     alpha4 = pvm[3][0] * (pvm[0][1]*pvm[1][3] - pvm[0][3]*pvm[1][1]) +
00429              pvm[3][1] * (pvm[0][3]*pvm[1][0] - pvm[0][0]*pvm[1][3]) +
00430              pvm[3][3] * (pvm[0][0]*pvm[1][1] - pvm[0][1]*pvm[1][0]);
00431 
00432     /* determine the order of the slices */
00433     if (pvm[2][2] - (pvm[2][0]*alpha1 + pvm[2][1]*alpha2)/(alpha3+alpha4) > 0)
00434         rbuf->reverse_k_order = 1;
00435     else
00436         rbuf->reverse_k_order = 0;
00437 
00438     /* compute the scale coefficients */
00439     rbuf->w_factor = alpha3 / alpha4;
00440     if (rbuf->reverse_k_order)
00441         rbuf->normalize_scale = 1. + rbuf->w_factor*(kcount-1);
00442     else
00443         rbuf->normalize_scale = 1.;
00444 
00445     /* compute the bounding box of the image in compositing space */
00446     denom = 1. / (alpha4 + alpha3*(kcount-1));
00447     i0 = rbuf->normalize_scale*alpha1*(kcount-1) * denom;
00448     j0 = rbuf->normalize_scale*alpha2*(kcount-1) * denom;
00449     i1 = rbuf->normalize_scale*(alpha4*icount + alpha1*(kcount-1)) * denom;
00450     j1 = rbuf->normalize_scale*(alpha4*jcount + alpha2*(kcount-1)) * denom;
00451     imin = MIN(0, i0);
00452     imax = MAX(rbuf->normalize_scale*icount, i1);
00453     jmin = MIN(0, j0);
00454     jmax = MAX(rbuf->normalize_scale*jcount, j1);
00455 
00456     /* compute the size of the intermediate image */
00457     rbuf->intermediate_width = (int)ceil(imax - imin);
00458     rbuf->intermediate_height = (int)ceil(jmax - jmin);
00459 
00460     /* compute the translation and shear coefficients */
00461     rbuf->shear_i = (rbuf->normalize_scale*alpha1 - alpha3*imin) / alpha4;
00462     rbuf->shear_j = (rbuf->normalize_scale*alpha2 - alpha3*jmin) / alpha4;
00463     rbuf->trans_i = -imin;
00464     rbuf->trans_j = -jmin;
00465 
00466     /* compute the depth coefficients */
00467     rbuf->depth_di = pvm[2][0];
00468     rbuf->depth_dj = pvm[2][1];
00469     rbuf->depth_dk = pvm[2][2];
00470     rbuf->depth_000 = pvm[2][3];
00471     rbuf->w_di = pvm[3][0];
00472     rbuf->w_dj = pvm[3][1];
00473     rbuf->w_dk = pvm[3][2];
00474     rbuf->w_000 = pvm[3][3];
00475 
00476     /* compute the mapping from compositing space to image space */
00477     Identity4(t);
00478     t[0][0] = 1. / rbuf->normalize_scale;
00479     t[1][1] = 1. / rbuf->normalize_scale;
00480     t[0][2] = -alpha1 / alpha4;
00481     t[1][2] = -alpha2 / alpha4;
00482     t[3][2] = -alpha3 / alpha4;
00483     t[0][3] = imin / rbuf->normalize_scale;
00484     t[1][3] = jmin / rbuf->normalize_scale;
00485     MatrixMult4(m2d, pvm, t);
00486 
00487     rbuf->warp_2d[0][0] = m2d[0][0];
00488     rbuf->warp_2d[1][0] = m2d[1][0];
00489     rbuf->warp_2d[2][0] = m2d[3][0];
00490     rbuf->warp_2d[0][1] = m2d[0][1];
00491     rbuf->warp_2d[1][1] = m2d[1][1];
00492     rbuf->warp_2d[2][1] = m2d[3][1];
00493     rbuf->warp_2d[0][2] = m2d[0][3];
00494     rbuf->warp_2d[1][2] = m2d[1][3];
00495     rbuf->warp_2d[2][2] = m2d[3][3];
00496 #endif /* notdef */
00497 }
00498 
00499 /*
00500  * ComputeAffineOpacityCorrection
00501  *
00502  * Precompute a lookup table which corrects opacity for an affine viewing
00503  * transformation.  (Opacity correction accounts for variations in the
00504  * apparent thickness of a voxel depending on viewpoint.)
00505  */
00506 
00507 static void
00508 ComputeAffineOpacityCorrection(vpc, shear_i, shear_j, table)
00509 vpContext *vpc;
00510 double shear_i;
00511 double shear_j;
00512 float table[VP_OPACITY_MAX+1];
00513 {
00514     float voxel_size;
00515     int i;
00516 
00517     Debug((vpc, VPDEBUG_OPCCORRECT,
00518            "Computing affine opacity correction table.\n"));
00519     voxel_size = sqrt(1 + shear_i*shear_i + shear_j*shear_j);
00520     for (i = 0; i <= VP_OPACITY_MAX; i++) {
00521 #ifdef NO_OPAC_CORRECT
00522         table[i] = (double)i / (double)VP_OPACITY_MAX;
00523 #else
00524         table[i] = 1.-pow(1.-(double)i/(double)VP_OPACITY_MAX,voxel_size);
00525 #endif
00526     }
00527 }
00528 
00529 /*
00530  * CheckRenderBuffers
00531  *
00532  * Resize the buffers used during rendering, if necessary.
00533  */
00534 
00535 static void
00536 CheckRenderBuffers(vpc)
00537 vpContext *vpc;
00538 {
00539     int new_max_width, new_max_height, new_max_scan;
00540     int resize = 0;
00541 
00542     /* determine if resizing is necessary */
00543     if (vpc->intermediate_width > vpc->max_intermediate_width) {
00544         new_max_width = MAX(vpc->intermediate_width,
00545                             vpc->int_image_width_hint);
00546         resize = 1;
00547     } else {
00548         new_max_width = MAX(vpc->max_intermediate_width,
00549                             vpc->int_image_width_hint);
00550     }
00551     if (vpc->intermediate_height > vpc->max_intermediate_height) {
00552         new_max_height = MAX(vpc->intermediate_height,
00553                              vpc->int_image_height_hint);
00554         resize = 1;
00555     } else {
00556         new_max_height = MAX(vpc->max_intermediate_height,
00557                              vpc->int_image_height_hint);
00558     }
00559     new_max_scan = vpc->xlen;
00560     if (vpc->ylen > new_max_scan)
00561         new_max_scan = vpc->ylen;
00562     if (vpc->zlen > new_max_scan)
00563         new_max_scan = vpc->zlen;
00564     if (new_max_scan > vpc->max_scan_length)
00565         resize = 1;
00566     if (vpc->color_channels != vpc->intermediate_color_channels)
00567         resize = 1;
00568 
00569     /* resize */
00570     if (resize)
00571         VPResizeRenderBuffers(vpc, new_max_width, new_max_height,new_max_scan);
00572 }
00573 
00574 /*
00575  * VPResizeRenderBuffers
00576  *
00577  * Resize the rendering buffers.
00578  */
00579 
00580 void
00581 VPResizeRenderBuffers(vpc, max_width, max_height, max_scan)
00582 vpContext *vpc;
00583 int max_width;  /* new width of the intermediate image */
00584 int max_height; /* new height of the intermediate image */
00585 int max_scan;   /* new max. scanline length */
00586 {
00587     /* free old buffers */
00588     if (vpc->int_image.gray_intim != NULL) {
00589         Debug((vpc, VPDEBUG_RBUF, "Freeing old RenderBuffer(%d,%d,%d,%d)\n",
00590                vpc->max_intermediate_width, vpc->max_intermediate_height,
00591                vpc->max_scan_length, vpc->intermediate_color_channels));
00592         Dealloc(vpc, vpc->int_image.gray_intim);
00593     }
00594 
00595     /* allocate new buffers */
00596     Debug((vpc, VPDEBUG_RBUF, "Allocating RenderBuffer(%d,%d,%d,%d)\n",
00597            max_width, max_height, max_scan, vpc->color_channels));
00598     vpc->max_intermediate_width = max_width;
00599     vpc->max_intermediate_height = max_height;
00600     vpc->max_scan_length = max_scan;
00601     vpc->intermediate_color_channels = vpc->color_channels;
00602     if (max_width > 0) {
00603         if (vpc->color_channels == 1) {
00604             Alloc(vpc, vpc->int_image.gray_intim, GrayIntPixel *,
00605                   max_width * max_height * sizeof(GrayIntPixel), "int_image");
00606         } else {
00607             Alloc(vpc, vpc->int_image.rgb_intim, RGBIntPixel *,
00608                   max_width * max_height * sizeof(RGBIntPixel), "int_image");
00609         }
00610     } else {
00611         vpc->int_image.gray_intim = NULL;
00612     }
00613 }
00614 
00615 /*
00616  * VPResizeDepthCueTable
00617  *
00618  * Resize the depth cueing table.
00619  */
00620 
00621 void
00622 VPResizeDepthCueTable(vpc, entries, copy)
00623 vpContext *vpc;
00624 int entries;    /* new number of table entries */
00625 int copy;       /* if true, copy old entries */
00626 {
00627     float *new_dc_table;
00628 
00629     Debug((vpc, VPDEBUG_DEPTHCUE, "resizing dctable to %d entries (%s)\n",
00630            entries, copy ? "copy" : "nocopy"));
00631     if (entries == 0) {
00632         if (vpc->dc_table != NULL) {
00633             Dealloc(vpc, vpc->dc_table);
00634             vpc->dc_table = NULL;
00635         }
00636         vpc->dc_table_len = 0;
00637     } else {
00638         Alloc(vpc, new_dc_table, float *, entries * sizeof(float), "dc_table");
00639         if (vpc->dc_table != NULL) {
00640             if (copy && vpc->dc_table_len > 0) {
00641                 bcopy(vpc->dc_table, new_dc_table,
00642                       MIN(vpc->dc_table_len, entries) * sizeof(float));
00643             }
00644             Dealloc(vpc, vpc->dc_table);
00645         }
00646         vpc->dc_table = new_dc_table;
00647         vpc->dc_table_len = entries;
00648     }
00649 }
00650 
00651 /*
00652  * VPComputeDepthCueTable
00653  *
00654  * Compute entries in the depth cueing lookup table.  
00655  */
00656 
00657 void
00658 VPComputeDepthCueTable(vpc, first, last)
00659 vpContext *vpc;
00660 int first;      /* first entry to compute */
00661 int last;       /* last entry to compute */
00662 {
00663     int c;
00664     double delta_depth, front_factor, density;
00665 
00666     Debug((vpc, VPDEBUG_DEPTHCUE, "computing dctable entries %d to %d\n",
00667            first, last));
00668     delta_depth = vpc->dc_quantization;
00669     front_factor = vpc->dc_front_factor;
00670     density = vpc->dc_density;
00671     for (c = first; c <= last; c++)
00672         vpc->dc_table[c] = front_factor * exp(-density*(1.0 - c*delta_depth));
00673 }
00674 
00675 /*
00676  * VPSliceDepthCueRatio
00677  *
00678  * Return the ratio of the depth cueing factor for two adjacent slices
00679  * for an affine view.  A constant factor is applied to all voxels in a
00680  * slice, and then a fixup is applied to the pixels of the intermediate
00681  * image.  This produces the correct answer without having to compute
00682  * the depth of each voxel.
00683  */
00684 
00685 float
00686 VPSliceDepthCueRatio(vpc)
00687 vpContext *vpc;
00688 {
00689     float delta_depth;          /* change in depth between adjacent slices */
00690     float slice_dc_ratio;       /* return value */
00691 
00692     if (!vpc->dc_enable)
00693         return(1.);
00694     delta_depth = vpc->depth_dk - vpc->depth_di*vpc->shear_i -
00695                   vpc->depth_dj*vpc->shear_j;
00696     if (vpc->reverse_slice_order)
00697         delta_depth = -delta_depth;
00698     /* slice_dc_ratio = exp(-vpc->dc_density * (-delta_depth)) */
00699     slice_dc_ratio = exp(vpc->dc_density * delta_depth);
00700     Debug((vpc, VPDEBUG_DEPTHCUE, "slice_dc_ratio = %f\n", slice_dc_ratio));
00701     return(slice_dc_ratio);
00702 }
00703 
00704 /*
00705  * VPDepthCueIntImage
00706  *
00707  * Perform depth cueing on the intermediate image.
00708  */
00709 
00710 void
00711 VPDepthCueIntImage(vpc, slicenum)
00712 vpContext *vpc;
00713 int slicenum;   /* slice number corresponding to location of int. image */
00714 {
00715     float pixel_depth_quant;    /* depth of current pixel in image
00716                                    (multiplied by depth_quant) */
00717     int pixel_depth_int;        /* pixel_depth truncated to an integer */
00718     float left_depth;           /* depth of pixel on left edge of current
00719                                    scanline in image */
00720     float left_depth_quant;     /* left_depth * depth_quant */
00721     float *dc_table;            /* depth cueing table */
00722     float depth_di, depth_dj;   /* change in depth for a unit change in each
00723                                    rotated object space coordinate */
00724     float depth_quant;          /* number of quantization levels for depth */
00725     float depth_di_quant;       /* depth_di * depth_quant */
00726     float depth_dj_quant;       /* depth_dj * depth_quant */
00727     float max_depth;            /* maximum (closest) depth in image */
00728     int max_depth_int;          /* maximum quantized depth */
00729     int i, j;                   /* intermediate image coordinates */
00730     float slice_u, slice_v;     /* sheared object space coordinates */
00731     GrayIntPixel *gray_intim;   /* image data (grayscale) */
00732     RGBIntPixel *rgb_intim;     /* image data (RGB) */
00733     int width, height;          /* size of intermediate image */
00734     int c;
00735 #ifdef DEBUG
00736     float pix_depth;
00737 #endif
00738 
00739     Debug((vpc, VPDEBUG_DEPTHCUE, "depth cueing intermediate image\n"));
00740 
00741     /* check the size of the depth cueing table and enlarge if necessary */
00742     width = vpc->intermediate_width;
00743     height = vpc->intermediate_height;
00744     depth_quant = 1.0 / vpc->dc_quantization;
00745     depth_di = vpc->depth_di;
00746     depth_dj = vpc->depth_dj;
00747     slice_u = vpc->shear_i * slicenum + vpc->trans_i;
00748     slice_v = vpc->shear_j * slicenum + vpc->trans_j;
00749     left_depth = vpc->depth_000 + vpc->depth_dk*slicenum -
00750                  slice_u*depth_di - slice_v*depth_dj;
00751     if (depth_di > 0) {
00752         if (depth_dj > 0) {
00753             max_depth = left_depth + depth_di*width + depth_dj*height;
00754         } else {
00755             max_depth = left_depth + depth_di*width;
00756         }
00757     } else {
00758         if (depth_dj > 0) {
00759             max_depth = left_depth + depth_dj*height;
00760         } else {
00761             max_depth = left_depth;
00762         }
00763     }
00764     max_depth_int = max_depth * depth_quant;
00765     if (max_depth_int >= vpc->dc_table_len) {
00766         c = vpc->dc_table_len;
00767         VPResizeDepthCueTable(vpc, max_depth_int+1, 1);
00768         VPComputeDepthCueTable(vpc, c, vpc->dc_table_len-1);
00769     }
00770     dc_table = vpc->dc_table;
00771     depth_di_quant = depth_di * depth_quant;
00772     depth_dj_quant = depth_dj * depth_quant;
00773     left_depth_quant = left_depth * depth_quant;
00774 
00775 #ifdef DEBUG
00776     Debug((vpc, VPDEBUG_DEPTHCUE, "depth cueing at image corners:\n"));
00777     pix_depth = left_depth + 0*depth_di + 0*depth_dj;
00778     pixel_depth_int = (int)(pix_depth * depth_quant);
00779     if (pixel_depth_int < 0) pixel_depth_int = 0;
00780     Debug((vpc, VPDEBUG_DEPTHCUE,
00781            "  %3d %3d: depth = %10.6f, factor = %10.6f, table[%d] = %10.6f\n",
00782            0, 0, pix_depth,
00783            vpc->dc_front_factor * exp(-vpc->dc_density * (1.0 - pix_depth)),
00784            pixel_depth_int, dc_table[pixel_depth_int]));
00785     pix_depth = left_depth + width*depth_di + 0*depth_dj;
00786     pixel_depth_int = (int)(pix_depth * depth_quant);
00787     if (pixel_depth_int < 0) pixel_depth_int = 0;
00788     Debug((vpc, VPDEBUG_DEPTHCUE,
00789            "  %3d %3d: depth = %10.6f, factor = %10.6f, table[%d] = %10.6f\n",
00790            width, 0, pix_depth,
00791            vpc->dc_front_factor * exp(-vpc->dc_density * (1.0 - pix_depth)),
00792            pixel_depth_int, dc_table[pixel_depth_int]));
00793     pix_depth = left_depth + width*depth_di + height*depth_dj;
00794     pixel_depth_int = (int)(pix_depth * depth_quant);
00795     if (pixel_depth_int < 0) pixel_depth_int = 0;
00796     Debug((vpc, VPDEBUG_DEPTHCUE,
00797            "  %3d %3d: depth = %10.6f, factor = %10.6f, table[%d] = %10.6f\n",
00798            width, height, pix_depth,
00799            vpc->dc_front_factor * exp(-vpc->dc_density * (1.0 - pix_depth)),
00800            pixel_depth_int, dc_table[pixel_depth_int]));
00801     pix_depth = left_depth + 0*depth_di + height*depth_dj;
00802     pixel_depth_int = (int)(pix_depth * depth_quant);
00803     if (pixel_depth_int < 0) pixel_depth_int = 0;
00804     Debug((vpc, VPDEBUG_DEPTHCUE,
00805            "  %3d %3d: depth = %10.6f, factor = %10.6f, table[%d] = %10.6f\n",
00806            0, height, pix_depth,
00807            vpc->dc_front_factor * exp(-vpc->dc_density * (1.0 - pix_depth)),
00808            pixel_depth_int, dc_table[pixel_depth_int]));
00809 #endif /* DEBUG */
00810 
00811     /* foreach pixel, compute depth and scale color by dc factor */
00812     if (vpc->color_channels == 1) {
00813         gray_intim = vpc->int_image.gray_intim;
00814         for (j = height; j > 0; j--) {
00815             pixel_depth_quant = left_depth_quant;
00816             left_depth_quant += depth_dj_quant;
00817             for (i = width; i > 0; i--) {
00818                 pixel_depth_int = pixel_depth_quant;
00819                 pixel_depth_quant += depth_di_quant;
00820                 if (pixel_depth_int < 0)
00821                     pixel_depth_int = 0;
00822                 if (pixel_depth_int >= vpc->dc_table_len) {
00823                     VPBug("VPDepthCueIntImage: depth too large (%d >= %d)",
00824                           pixel_depth_int, vpc->dc_table_len);
00825                 }
00826                 gray_intim->clrflt *= dc_table[pixel_depth_int];
00827                 gray_intim++;
00828             } /* for i */
00829         } /* for j */
00830     } else {
00831         rgb_intim = vpc->int_image.rgb_intim;
00832         for (j = height; j > 0; j--) {
00833             pixel_depth_quant = left_depth_quant;
00834             left_depth_quant += depth_dj_quant;
00835             for (i = width; i > 0; i--) {
00836                 pixel_depth_int = pixel_depth_quant;
00837                 pixel_depth_quant += depth_di_quant;
00838                 if (pixel_depth_int < 0)
00839                     pixel_depth_int = 0;
00840                 if (pixel_depth_int >= vpc->dc_table_len) {
00841                     VPBug("VPDepthCueIntImage: depth too large (%d >= %d)",
00842                           pixel_depth_int, vpc->dc_table_len);
00843                 }
00844                 rgb_intim->rclrflt *= dc_table[pixel_depth_int];
00845                 rgb_intim->gclrflt *= dc_table[pixel_depth_int];
00846                 rgb_intim->bclrflt *= dc_table[pixel_depth_int];
00847                 rgb_intim++;
00848             } /* for i */
00849         } /* for j */
00850     }
00851 }
00852 
00853 /*
00854  * ComputeLightViewTransform
00855  *
00856  * Compute the view transformation from the point of view of the one
00857  * light source that produces shadows.
00858  */
00859 
00860 static void
00861 ComputeLightViewTransform(vpc, vm)
00862 vpContext *vpc;
00863 vpMatrix4 vm;   /* storage for result */
00864 {
00865     vpMatrix4 prematrix;        /* transform volume to unit cube */
00866     vpMatrix4 viewportm;                /* viewport matrix */
00867     vpVector3 vpn;              /* view plane normal */
00868     vpVector3 vup;              /* view up vector */
00869     vpVector3 tmp1v, tmp2v;     /* temporary vectors */
00870     vpMatrix4 view;             /* transform world coordinates to
00871                                    eye coordinates, with view
00872                                    direction equal to light vector */
00873     vpMatrix4 tmp1m, tmp2m;     /* temporary matrices */
00874     double lx, ly, lz;          /* components of light vector */
00875     int light_index;
00876     int maxdim;
00877     double scale;
00878 
00879     /* check for errors */
00880     ASSERT(vpc->shadow_light_num >= VP_LIGHT0 &&
00881            vpc->shadow_light_num <= VP_LIGHT5);
00882     ASSERT(vpc->light_enable[vpc->shadow_light_num - VP_LIGHT0]);
00883 
00884     /* compute prematrix */
00885     vpIdentity4(prematrix);
00886     maxdim = vpc->xlen;
00887     if (vpc->ylen > maxdim)
00888         maxdim = vpc->ylen;
00889     if (vpc->zlen > maxdim)
00890         maxdim = vpc->zlen;
00891     scale = 1. / (double)maxdim;
00892     prematrix[0][0] = scale;
00893     prematrix[1][1] = scale;
00894     prematrix[2][2] = scale;
00895     prematrix[0][3] = -vpc->xlen * scale * 0.5;
00896     prematrix[1][3] = -vpc->ylen * scale * 0.5;
00897     prematrix[2][3] = -vpc->zlen * scale * 0.5;
00898 
00899     /* compute the world-to-eye coordinate transformation */
00900     light_index = vpc->shadow_light_num - VP_LIGHT0;
00901     lx = vpc->light_vector[light_index][0];
00902     ly = vpc->light_vector[light_index][1];
00903     lz = vpc->light_vector[light_index][2];
00904     vpSetVector3(vpn, lx, ly, lz);
00905     if (fabs(lx) < fabs(ly)) {
00906         if (fabs(lx) < fabs(lz)) {
00907             vpSetVector3(vup, 1.0, 0.0, 0.0);
00908         } else {
00909             vpSetVector3(vup, 0.0, 0.0, 1.0);
00910         }
00911     } else {
00912         if (fabs(ly) < fabs(lz)) {
00913             vpSetVector3(vup, 0.0, 1.0, 0.0);
00914         } else {
00915             vpSetVector3(vup, 0.0, 0.0, 1.0);
00916         }
00917     }
00918     vpCrossProduct(tmp1v, vup, vpn);
00919     vpCrossProduct(tmp2v, vpn, tmp1v);
00920 
00921     vpIdentity4(view);
00922     view[0][0] = tmp1v[0];
00923     view[0][1] = tmp1v[1];
00924     view[0][2] = tmp1v[2];
00925     view[1][0] = tmp2v[0];
00926     view[1][1] = tmp2v[1];
00927     view[1][2] = tmp2v[2];
00928     view[2][0] = vpn[0];
00929     view[2][1] = vpn[1];
00930     view[2][2] = vpn[2];
00931 
00932     /* initialize matrix to switch to left-handed coords. */
00933     vpIdentity4(viewportm);
00934     viewportm[2][2] = -1.0;
00935 
00936     /* compute the view matrix */
00937     vpMatrixMult4(tmp1m, vpc->transforms[VP_MODEL], prematrix);
00938     vpMatrixMult4(tmp2m, view, tmp1m);
00939     vpMatrixMult4(vm, viewportm, tmp2m);
00940 }
00941 
00942 /*
00943  * FactorLightView
00944  *
00945  * Factor an affine viewing matrix that specifies the view seen by
00946  * a light source.  Most of the parameters of the factorization are
00947  * taken from the factorization of the normal viewing matrix.
00948  */
00949 
00950 static int
00951 FactorLightView(vpc, vm)
00952 vpContext *vpc;
00953 vpMatrix4 vm;
00954 {
00955     vpMatrix4 p;                /* permutation matrix */
00956     vpMatrix4 pvm;              /* permutation of the viewing matrix */
00957     int icount, jcount, kcount; /* dimensions of volume in rotated space */
00958     double denom;
00959 
00960     Debug((vpc, VPDEBUG_SHADOW, "FactorLightView\n"));
00961 
00962     /* permute the rows of the viewing matrix according to the best viewing
00963        axis for the viewing direction */
00964     bzero(p, sizeof(vpMatrix4));
00965     switch (vpc->best_view_axis) {
00966     case VP_X_AXIS:
00967         p[0][2] = 1.;
00968         p[1][0] = 1.;
00969         p[2][1] = 1.;
00970         p[3][3] = 1.;
00971         icount = vpc->ylen;
00972         jcount = vpc->zlen;
00973         kcount = vpc->xlen;
00974         break;
00975     case VP_Y_AXIS:
00976         p[0][1] = 1.;
00977         p[1][2] = 1.;
00978         p[2][0] = 1.;
00979         p[3][3] = 1.;
00980         icount = vpc->zlen;
00981         jcount = vpc->xlen;
00982         kcount = vpc->ylen;
00983         break;
00984     case VP_Z_AXIS:
00985         p[0][0] = 1.;
00986         p[1][1] = 1.;
00987         p[2][2] = 1.;
00988         p[3][3] = 1.;
00989         icount = vpc->xlen;
00990         jcount = vpc->ylen;
00991         kcount = vpc->zlen;
00992         break;
00993     }
00994     vpMatrixMult4(pvm, vm, p);
00995 
00996     /* compute the shear coefficients */
00997     denom = pvm[0][0]*pvm[1][1] - pvm[0][1]*pvm[1][0];
00998     if (fabs(denom) < VP_EPS)
00999         return(VPSetError(vpc, VPERROR_SINGULAR));
01000     vpc->shadow_shear_i = (pvm[0][2]*pvm[1][1] - pvm[0][1]*pvm[1][2]) / denom;
01001     vpc->shadow_shear_j = (pvm[0][0]*pvm[1][2] - pvm[1][0]*pvm[0][2]) / denom;
01002 
01003     /* check that light direction is compatible with compositing direction */
01004     if (pvm[2][0]*vpc->shadow_shear_i + pvm[2][1]*vpc->shadow_shear_j -
01005         pvm[2][2] > 0) {
01006         if (vpc->reverse_slice_order != 0)
01007             return(VPERROR_BAD_SHADOW);
01008     } else {
01009         if (vpc->reverse_slice_order != 1)
01010             return(VPERROR_BAD_SHADOW);
01011     }
01012 
01013     /* compute the shadow buffer image size */
01014     vpc->shadow_width = (int)ceil((kcount-1)*fabs(vpc->shadow_shear_i)) +
01015         icount + 1;
01016     vpc->shadow_height = (int)ceil((kcount-1)*fabs(vpc->shadow_shear_j)) +
01017         jcount + 1;
01018 
01019     /* compute the translation coefficients */
01020     if (vpc->shadow_shear_i >= 0.)
01021         vpc->shadow_trans_i = 1.;
01022     else
01023         vpc->shadow_trans_i = 1. - vpc->shadow_shear_i * (kcount - 1);
01024     if (vpc->shadow_shear_j >= 0.)
01025         vpc->shadow_trans_j = 1.;
01026     else
01027         vpc->shadow_trans_j = 1. - vpc->shadow_shear_j * (kcount - 1);
01028 
01029     Debug((vpc, VPDEBUG_SHADOW, "  shadow shear factors: %g %g\n",
01030            vpc->shadow_shear_i, vpc->shadow_shear_j));
01031     Debug((vpc, VPDEBUG_SHADOW, "  shadow translation: %g %g\n",
01032            vpc->shadow_trans_i, vpc->shadow_trans_j));
01033 
01034     return(VP_OK);
01035 }
01036 
01037 /*
01038  * CheckShadowBuffer
01039  *
01040  * Resize the shadow buffer, if necessary.
01041  */
01042 
01043 static void
01044 CheckShadowBuffer(vpc)
01045 vpContext *vpc;
01046 {
01047     int new_max_width, new_max_height;
01048     int resize = 0;
01049 
01050     /* determine if resizing is necessary */
01051     if (vpc->shadow_width > vpc->max_shadow_width) {
01052         new_max_width = MAX(vpc->shadow_width, vpc->shadow_width_hint);
01053         resize = 1;
01054     } else {
01055         new_max_width = MAX(vpc->max_shadow_width, vpc->shadow_width_hint);
01056     }
01057     if (vpc->shadow_height > vpc->max_shadow_height) {
01058         new_max_height = MAX(vpc->shadow_height, vpc->shadow_height_hint);
01059         resize = 1;
01060     } else {
01061         new_max_height = MAX(vpc->max_shadow_height, vpc->shadow_height_hint);
01062     }
01063 
01064     /* resize */
01065     if (resize)
01066         VPResizeShadowBuffer(vpc, new_max_width, new_max_height);
01067 }
01068 
01069 /*
01070  * VPResizeShadowBuffer
01071  *
01072  * Resize the shadow buffer.
01073  */
01074 
01075 void
01076 VPResizeShadowBuffer(vpc, max_width, max_height)
01077 vpContext *vpc;
01078 int max_width;  /* new width of the intermediate image */
01079 int max_height; /* new height of the intermediate image */
01080 {
01081     /* free old buffers */
01082     if (vpc->shadow_buffer != NULL) {
01083         Dealloc(vpc, vpc->shadow_buffer);
01084     }
01085 
01086     /* allocate new buffers */
01087     vpc->max_shadow_width = max_width;
01088     vpc->max_shadow_height = max_height;
01089     if (max_width > 0) {
01090         Alloc(vpc, vpc->shadow_buffer, GrayIntPixel *,
01091               max_width * max_height * sizeof(GrayIntPixel), "shadow_buffer");
01092     } else {
01093         vpc->shadow_buffer = NULL;
01094     }
01095 }