/****************************************************************************/
/*                                                                          */
/*  VolVis is a volume visualization system for investigating, manipulating */
/*  and rendering geometric and volumetric data.                            */
/*                                                                          */
/*  Copyright (C) 1993 by the Research Foundation of the State University   */
/*                            of New York                                   */
/*                                                                          */
/*  This program is free software; you can redistribute it and/or modify    */
/*  it under the terms of the GNU General Public License as published by    */
/*  the Free Software Foundation; either version 1, or (at your option)     */
/*  any later version.                                                      */
/*                                                                          */
/*  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 General Public License for more details.                            */
/*                                                                          */
/*  You should have received a copy of the GNU General Public License       */
/*  along with this program; if not, write to the Free Software             */
/*  Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.               */
/*                                                                          */
/*  For information on VolVis, contact us at:                               */
/*                                                                          */
/*                volvis@cs.sunysb.edu                         (email)      */
/*                                                                          */
/*                Lisa Sobierajski & Ricardo Avila             (US Mail)    */
/*                Department of Computer Science                            */
/*                State University of New York at Stony Brook               */
/*                Stony Brook, New York  11794-4400                         */
/*                                                                          */
/****************************************************************************/



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

# include "C_volvis.h"
# include "C_raytrace.h"

# define debug(x)


extern C_IntersectionList *C_bounding_intersection_list();


	
int C_trace_ray( rt_info, world, view, intersect_info )
C_RayTraceInfo			*rt_info;
C_World				*world;
C_View				*view;
C_IntersectInfo			*intersect_info;
{
	extern void C_rt_volume_simple_range();
	extern C_transform_fposition(); 
	extern C_transform_fvector(); 
	extern void C_shade_point();

	int			exiting_volume;
	int			exited_volume;
	int			vloop;
	int			loop;
	C_Volume		*volume;
	C_IntersectionList	*intersections,
				*iptr;
	C_Voxel_8bit		*dptr,
				low,
				high;
	C_StepType		prev_step_type;
	float			dot_product;	
	float			t, distance;
	C_FPosition		new_ray_origin;
	C_FPosition		point_one, point_two;
	C_FPosition		midpoint;
	C_FPosition		intersection_point;
	C_FPosition		junk_point;
	C_FVector		new_ray_direction;
	C_FVector		normal;
	C_FVector		world_normal;
	float			closest_distance;
	float			current_intersection;
	float			num_steps;
	int			xstep, ystep, zstep,
				found_intersection,
				voxel_x, voxel_y, voxel_z;
	float			tmax_x, tmax_y, tmax_z,
				tdelta_x, tdelta_y, tdelta_z;
	int			tstep_x, tstep_y, tstep_z;
	int			sign_diffx, sign_diffy, sign_diffz;
	float			diffx, diffy, diffz;

	register float		A, B, C, D, E, F, G, H;
	float
				x_weight1, y_weight1, z_weight1,
				x_weight2, y_weight2, z_weight2,
				mid_x_weight, mid_y_weight, mid_z_weight,
				midpointvalue,
				value1, value2;
	C_ScalarData_8bit	*scan_data;


	float			quad_eq_a,
				quad_eq_b,
				quad_eq_c,
				t1, t2,
				radical_term;
	C_GeoData		*table_entry;

	C_Polygon		*polygon;

	C_GeoData		*geo_data;	

	C_FPosition		*vertex_ptr1, *vertex_ptr2;
	int			cross_count;

	C_MultiPoly		*multi_poly;
	int			multi_loop;

	C_ShadeInfo		shade_info;

	int			start_in_range;

	float			delta_t;
	float			scale_distance;
	float			accum_opacity;
	float			accum_color_red;
	float			accum_color_green;
	float			accum_color_blue;
	float			curr_opacity;

	C_FPosition		old_origin;


/***************************************************************************/
/****                                                                   ****/
/****  Create the intersections list...                                 ****/
/****                                                                   ****/
/****  C_bounding_intersection_list (located in C_bounding_intersect.c) ****/
/****  is called to compute the intersection of this ray with each      ****/
/****  volume in the world.  The returned list contains the following   ****/
/****  information for each intersection:				****/
/****                                                                   ****/
/****                                                                   ****/
/****           position - an (x,y,z) location (in the local coordinates****/
/****                      of the intersected volume, in voxels)        ****/
/****                      indicating the position of the intersection. ****/
/****                                                                   ****/
/****           distance - the distance at which the intersection       ****/
/****                      occurred (in units, relative to the world    ****/
/****                      coordinate system).                          ****/
/****                                                                   ****/
/****           volume   - the index of the intersected volume.         ****/
/****                                                                   ****/
/****           start_step - indicates which face of the volume was     ****/
/**** 			   intersected by this ray                      ****/
/****                                                                   ****/
/****           max_distance - the maximum distance (within this 	****/
/****			   volume) that the true intersection can have. ****/
/****			   Used for limiting intersections due to       ****/
/****			   light sources (so a shadow ray does not      ****/
/****			   intersect a shadowing object after passing   ****/
/****			   the light source) and for cut geometries.    ****/
/****                                                                   ****/
/****           next_intersection - link to the next structure.         ****/
/****                                                                   ****/
/***************************************************************************/

	if ( rt_info->ray_type != C_PARC_RAY )
	{
	    if ( rt_info->ray_type != C_PRIMARY_RAY )
	    {
		old_origin.x = rt_info->ray_origin.x;
		old_origin.y = rt_info->ray_origin.y;
		old_origin.z = rt_info->ray_origin.z;

		rt_info->ray_origin.x += 0.1*rt_info->ray_direction.x;
		rt_info->ray_origin.y += 0.1*rt_info->ray_direction.y;
		rt_info->ray_origin.z += 0.1*rt_info->ray_direction.z;
	    }
	    intersections = C_bounding_intersection_list(
				&(rt_info->ray_origin), 
				&(rt_info->ray_direction), 
			        world, rt_info->max_distance );
	}
	else
	{
	    intersections = C_New( C_IntersectionList );

	    intersections->start_cut = NULL;
	    intersections->end_cut = NULL;
	    
	    intersections->volume = intersect_info->intersect_volume;

	    point_one.x = rt_info->ray_origin.x +  
		rt_info->ray_direction.x * 
		( intersect_info->intersect_distance - 0.0 );
	    point_one.y = rt_info->ray_origin.y +  
		rt_info->ray_direction.y * 
		( intersect_info->intersect_distance - 0.0 );
	    point_one.z = rt_info->ray_origin.z +  
		rt_info->ray_direction.z * 
		( intersect_info->intersect_distance - 0.0 );

	    C_transform_fposition( &(point_one), 
				   &(intersections->position),
		   &(world->volume[intersections->volume]->wtol_voxels ) );
				    
				   
	    intersections->distance = intersect_info->intersect_distance;
	    intersections->max_distance = rt_info->max_distance;
	    intersections->start_step = C_UNDEFINED_STEP;
	    intersections->next_intersection = NULL;
	}

/*-------------------------------------------------------------------------*/


/***************************************************************************/
/****									****/
/****  Find the closest true intersection...				****/
/****									****/
/****  In this section of the code, the closest true intersection is	****/
/****  determined.  This is done using the following algorithm.		****/
/****									****/
/****		1.  Let iptr point to the first intersection in		****/
/****		    intersections.					****/
/****									****/
/****		2.  Let closest_distance = C_BIG_FLOAT 			****/
/****               (the largest number that can be put into a float).	****/
/****									****/
/****		3.  While (iptr != NULL && closest_distance > 		****/
/****			   iptr->distance) do				****/
/****									****/
/****			3a.  Traverse the ray through the volume	****/
/****			     iptr->volume to determine if the ray	****/
/**** 			     intersects an object in this volume.	****/
/****									****/
/****			3b.  If an intersection does occur, find the	****/
/****			     distance from the ray_origin to this	****/
/****			     intersection point and store this		****/
/****			     in current_intersection.			****/
/****									****/
/****			3c.  if (current_intersection < 		****/
/****				 closest_distance) then			****/
/****			     closest_distance =		 		****/
/****			         current_intersection and		****/
/****			     store all other intersection information	****/
/****			     including the volume, intersection point   ****/
/****			     in local and world coordinates.		****/
/****									****/
/****			3d.  iptr = iptr->next_intersection 		****/
/****									****/
/****									****/
/****  At the end of this, closest_distance contains the closest 	****/
/****  intersection along the ray.					****/
/****									****/
/***************************************************************************/

	iptr = intersections;

	closest_distance = C_BIG_FLOAT;

	intersect_info->intersect_color.red = 
	intersect_info->intersect_color.green = 
	intersect_info->intersect_color.blue = 255;

	while ( iptr != NULL  && 
		iptr->distance < closest_distance )
	{
	  switch (world->volume[iptr->volume]->data_type)
	  {
	    case C_SCALAR_DATA_8BIT:

	      switch (world->volume[iptr->volume]->
			data.scalar_data_8bit->seg_type)
	      {
		case C_ISO_SURFACE:

		  C_rt_volume_iso_surface( world, view, rt_info, iptr, 
					    &closest_distance, &shade_info );

		  break;

		case C_SCALAR_OPACITY:

		  C_rt_volume_scalar_opacity( world, view, rt_info, iptr, 
					      &closest_distance, &shade_info );

		  break;

		default:
			
		      C_error_message("Unsupported segmentation type\n");
		      return C_ERROR;

		      break;
	      }

	      break;

	    case C_GEOMETRIC_DATA:

		switch ( world->volume[iptr->volume]->data.geo_data->geo_type )
		{
		    case C_SPHERE:
		    case C_HOLLOW_SPHERE:

			C_rt_geo_sphere( world, rt_info, iptr, 
					 &closest_distance, &shade_info );

			break;

		    case C_CYLINDER:
		    case C_HOLLOW_CYLINDER:

			C_rt_geo_cylinder( world, rt_info, iptr, 
					   &closest_distance, &shade_info );

		     	break;

		    case C_TORUS:
		    case C_HOLLOW_TORUS:

			C_error_message("How did a torus get in here?\n");

			break;

		    case C_POLYGON:

			C_rt_geo_polygon( world, rt_info, iptr, 
					   &closest_distance, &shade_info );

		     	break;

		    case C_MULTIPOLY:
			
			C_rt_geo_multi_polygon( world, rt_info, iptr, 
					   &closest_distance, &shade_info );

		     	break;

		    default:

			C_error_message("Unknown Geo Object Type!!!\n");
			return C_ERROR;

			break;

		  }

		  break;
		  	

	    default:
		
		C_error_message("Error: Unknow data type.\n");
		return C_ERROR;

		break;
	  }

	  volume = world->volume[iptr->volume];

	   if ( closest_distance < C_BIG_FLOAT && 
		rt_info->ray_type == C_SHADOW_RAY )
	     if ( volume->shade_type == C_SIMPLE_SHADE )
		if ( volume->local_shade.simple_shade->transmit_coef > 0.0 )
	        {
		  if (volume->color_type == C_SINGLE_COLOR)
		  {
			closest_distance = C_BIG_FLOAT - 1;
			found_intersection = FALSE;
			intersect_info->intersect_color.red *= 
			    volume->local_shade.simple_shade->transmit_coef *
			    (float)(
			    volume->color_data.single_color->red)/(255.0);
			intersect_info->intersect_color.green *= 
			    volume->local_shade.simple_shade->transmit_coef *
			    (float)(
			    volume->color_data.single_color->green)/(255.0);
			intersect_info->intersect_color.blue *= 
			    volume->local_shade.simple_shade->transmit_coef *
			    (float)(
			    volume->color_data.single_color->blue)/(255.0);
		  }
	        }
	        else
		{
			intersect_info->intersect_color.red = 
			intersect_info->intersect_color.green = 
			intersect_info->intersect_color.blue = 0;
		}
	  
	   iptr = iptr->next_intersection;

	} 


	/* Now that we have closest_intersection ..... */

	intersect_info->intersect_distance = closest_distance;


	if ( closest_distance < C_BIG_FLOAT )
	{
	  intersect_info->intersect_volume = shade_info.volume_num;
	  intersect_info->intersect_local_position.x = 
					shade_info.local_point.x;
	  intersect_info->intersect_local_position.y = 
					shade_info.local_point.y;
	  intersect_info->intersect_local_position.z = 
					shade_info.local_point.z;
	  intersect_info->intersect_world_position.x = 
					shade_info.world_point.x;
	  intersect_info->intersect_world_position.y = 
					shade_info.world_point.y;
	  intersect_info->intersect_world_position.z = 
					shade_info.world_point.z;

	  switch ( rt_info->ray_type )
	  {
	    case C_PRIMARY_RAY:
	    case C_REFLECT_RAY:
	    case C_TRANSMIT_RAY:
	  	shade_info.ray_type = rt_info->ray_type;
	  	shade_info.ray_level = rt_info->ray_level;
	
	  	shade_info.ray_direction.x = rt_info->ray_direction.x;
	  	shade_info.ray_direction.y = rt_info->ray_direction.y;
	  	shade_info.ray_direction.z = rt_info->ray_direction.z;
	  	C_shade_point( world, view, &shade_info, 
				&(intersect_info->intersect_color) );
		break;
	
	    case C_SHADOW_RAY:

		C_Set_Color( (intersect_info->intersect_color), 0, 0, 0 );
		break;

	    case C_PROBE_RAY:

		C_rt_normal( world->volume[shade_info.volume_num],
			     &shade_info );

		break;

	    case C_PARC_RAY:
	  	shade_info.ray_type = rt_info->ray_type;
	  	shade_info.ray_level = rt_info->ray_level;
	
	  	shade_info.ray_direction.x = rt_info->ray_direction.x;
	  	shade_info.ray_direction.y = rt_info->ray_direction.y;
	  	shade_info.ray_direction.z = rt_info->ray_direction.z;
	  	C_shade_point( world, view, &shade_info, 
				&(intersect_info->intersect_color) );
		break;
	  }

	  intersect_info->intersect_world_normal.x =
					shade_info.world_normal.x;
	  intersect_info->intersect_world_normal.y =
					shade_info.world_normal.y;
	  intersect_info->intersect_world_normal.z =
					shade_info.world_normal.z;

	}
	else
	{
	  intersect_info->intersect_volume = -1;

	  switch ( rt_info->ray_type )
	  {
	    case C_PRIMARY_RAY:
	    case C_TRANSMIT_RAY:
		

		C_Set_Color( (intersect_info->intersect_color),
			        world->world_shade.bkgrnd_color.red,
			        world->world_shade.bkgrnd_color.green,
			        world->world_shade.bkgrnd_color.blue );

		break;

	    case C_REFLECT_RAY:

		C_Set_Color( (intersect_info->intersect_color), 0, 0, 0 );

		break;

	    case C_SHADOW_RAY:

		break;

	    case C_PROBE_RAY:

		break;

	    case C_PARC_RAY:

		C_Set_Color( (intersect_info->intersect_color),
			        world->world_shade.bkgrnd_color.red,
			        world->world_shade.bkgrnd_color.green,
			        world->world_shade.bkgrnd_color.blue );
		break;
	  }
	}
		
	iptr = intersections;
	loop = 0;

	while (iptr)
	{
		loop++;
		iptr = iptr->next_intersection;
		free(intersections);
		intersections = iptr;
	}
 
	if ( rt_info->ray_type != C_PARC_RAY )
	{
	    if ( rt_info->ray_type != C_PRIMARY_RAY )
	    {
		rt_info->ray_origin.x = old_origin.x;
		rt_info->ray_origin.y = old_origin.y;
		rt_info->ray_origin.z = old_origin.z;
	    }
	}
/****
printf("%d --- ray ending:  (%f,%f,%f), level %d\n",
	C_Parallel_Get_Processor_Id(), rt_info->ray_origin.x,
	rt_info->ray_origin.y, rt_info->ray_origin.z, rt_info->ray_level);
****/
}


