UNCHAR3   **   CContourJudge::ImgInterp(int   k,   float   imgfactor,   float   mode,UNCHAR3   **resimg)  
  {  
  int   i,j;  
  int   pos;  
  //UNCHAR3   **   reimg;  
  long   x,y;  
  long   newimgw,newimgh;  
  float   xScale,   yScale,   fX,   fY;  
  MYVOXEL   voxnode;  
  CMyPoint   tmpt,rgb1,rgb2,rgb3,rgb4,test;  
  newimgw   =   (long)(pDoc->xdim   *   imgfactor);  
  newimgh   =   (long)(pDoc->ydim   *   imgfactor);  
   
  xScale   =   (float)(pDoc->xdim)   /   (float)newimgw;  
  yScale   =   (float)(pDoc->ydim)   /   (float)newimgh;  
   
  for(i=0;   i<newimgh;   i++)  
  {  
  for(j=0;   j<newimgw;   j++)  
  {  
  resimg[i][j][0]   =   0;  
  resimg[i][j][1]   =   0;  
  resimg[i][j][2]   =   0;  
  }  
  }  
   
  //CContourListProc   m_cclistobj;  
   
  switch((long)mode)  
  {  
  case   1:     //nearest   neighbour  
  for(y=0;   y<newimgh;   y++)  
  {  
  fY   =   y   *   yScale;  
  for(x=0;   x<newimgw;   x++)  
  {  
  fX   =   x   *   xScale;  
   
  test.x   =   (double)fX;  
  test.y   =   (double)fY;  
  test.z   =   (double)k;  
  float   LayHeight=   (float)pDoc->layerthickness*k+pDoc->minz;  
   
   
  pos   =   pDoc->CompPos((long)fX,(long)fY,k);  
  voxnode   =   pDoc->myvoxarray[pos];  
  if(voxnode.flag   !=   255)  
  {  
  resimg[y][x][0]   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  resimg[y][x][1]   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  resimg[y][x][2]   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);    
  }  
  else  
  {  
  resimg[y][x][0]   =   0;  
  resimg[y][x][1]   =   0;  
  resimg[y][x][2]   =   0;    
  }  

   
  }  
  }  
  break;  
  case   2:   //   bicubic   interpolation   by   Blake   L.   Carlson   <blake-carlson(at)uiowa(dot)edu  
  float   f_x,   f_y,   a,   b,   rr,   gg,   bb,   r1,   r2;  
  int       i_x,   i_y,   xx,   yy;  
  for(y=0;   y<newimgh;   y++)  
  {  
  f_y   =   (float)   y   *   yScale;  
  i_y   =   (int)   floor(f_y);  
  a       =   f_y   -   (float)floor(f_y);  
  for(x=0;   x<newimgw;   x++)  
  {  
  f_x   =   (float)   x   *   xScale;  
  i_x   =   (int)   floor(f_x);  
  b       =   f_x   -   (float)floor(f_x);  
   
  rr   =   gg   =   bb   =   0.0F;  
  for(int   m=-1;   m<3;   m++)    
  {  
  r1   =   b3spline((float)   m   -   a);  
  for(int   n=-1;   n<3;   n++)    
  {  
  r2   =   b3spline(-1.0F*((float)n   -   b));    
  xx   =   i_x+n+2;  
  yy   =   i_y+m+2;  
  if   (xx<0)   xx=0;  
  if   (yy<0)   yy=0;  
  if   (xx>=pDoc->xdim)   xx=pDoc->xdim   -   1;  
  if   (yy>=pDoc->ydim)   yy=pDoc->ydim   -   1;  
   
  pos   =   pDoc->CompPos(xx,yy,k);  
   
  voxnode   =   pDoc->myvoxarray[pos];  
  if(voxnode.flag   !=   255)  
  {  
  tmpt.x   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  tmpt.y   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  tmpt.z   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);  
   
  rr   +=   (float)tmpt.x   *   r1   *   r2;  
  gg   +=   (float)tmpt.y   *   r1   *   r2;  
  bb   +=   (float)tmpt.z   *   r1   *   r2;  
  }  
  }  
  }  
  resimg[y][x][0]   =   (unsigned   char)   rr;  
  resimg[y][x][1]   =   (unsigned   char)   gg;  
  resimg[y][x][2]   =   (unsigned   char)   bb;  
  }  
  }  
  break;  
  default:   //bilinear   interpolation  
  {  
  long   ifX,   ifY,   ifX1,   ifY1,   xmax,   ymax;  
  float   ir1,   ir2,   ig1,   ig2,   ib1,   ib2,   dx,   dy;  
  BYTE   r,g,b;  
   
  xmax   =   pDoc->xdim-1;  
  ymax   =   pDoc->ydim-1;  
  for(y=0;   y<newimgh-1;   y++){  
  fY   =   y   *   yScale;  
  ifY   =   (int)fY;  
  ifY1   =   min(ymax,   ifY+1);  
  dy   =   fY   -   ifY;  
  for(x=0;   x<newimgw-1;   x++){  
  fX   =   x   *   xScale;  
  ifX   =   (int)fX;  
  ifX1   =   min(xmax,   ifX+1);  
  dx   =   fX   -   ifX;  
   
  //   Interpolate   using   the   four   nearest   pixels   in   the   source  
  pos   =   pDoc->CompPos(ifX1,ifY1,k);  
  voxnode   =   pDoc->myvoxarray[pos];  
  rgb1.x   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  rgb1.y   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  rgb1.z   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);  
   
  pos   =   pDoc->CompPos(ifX1+1,ifY1,k);  
  voxnode   =   pDoc->myvoxarray[pos];  
  rgb2.x   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  rgb2.y   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  rgb2.z   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);  
   
  pos   =   pDoc->CompPos(ifX1,ifY1+1,k);  
  voxnode   =   pDoc->myvoxarray[pos];  
  rgb3.x   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  rgb3.y   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  rgb3.z   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);  
   
  pos   =   pDoc->CompPos(ifX1+1,ifY1+1,k);  
  voxnode   =   pDoc->myvoxarray[pos];  
  rgb4.x   =   (unsigned   char)   (255   *   voxnode.matvalue[0]);  
  rgb4.y   =   (unsigned   char)   (255   *   voxnode.matvalue[1]);  
  rgb4.z   =   (unsigned   char)   (255   *   voxnode.matvalue[2]);  
   
  //   Interplate   in   x   direction:  
  ir1   =   (float)rgb1.x     *   (1   -   dy)   +   (float)rgb3.x   *   dy;  
  ig1   =   (float)rgb1.y     *   (1   -   dy)   +   (float)rgb3.y   *   dy;  
  ib1   =   (float)rgb1.z     *   (1   -   dy)   +   (float)rgb3.z   *   dy;  
  ir2   =   (float)rgb2.x     *   (1   -   dy)   +   (float)rgb4.x   *   dy;  
  ig2   =   (float)rgb2.y     *   (1   -   dy)   +   (float)rgb4.y   *   dy;  
  ib2   =   (float)rgb2.z     *   (1   -   dy)   +   (float)rgb4.z   *   dy;  
  //   Interpolate   in   y:  
  r   =   (BYTE)(ir1   *   (1   -   dx)   +   ir2   *   dx);  
  g   =   (BYTE)(ig1   *   (1   -   dx)   +   ig2   *   dx);  
  b   =   (BYTE)(ib1   *   (1   -   dx)   +   ib2   *   dx);  
  //   Set   output  
  resimg[y][x][0]   =   r;  
  resimg[y][x][1]   =   g;  
  resimg[y][x][2]   =   b;  
   
  }  
  }    
  }    
  break;  
  }  
   
  return   resimg;  
  }