rendered paste body/**************************************************************************** * spheres.cpp * * This module implements the sphere primitive. * * from Persistence of Vision(tm) Ray Tracer version 3.6. * Copyright 1991-2003 Persistence of Vision Team * Copyright 2003-2004 Persistence of Vision Raytracer Pty. Ltd. *--------------------------------------------------------------------------- * NOTICE: This source code file is provided so that users may experiment * with enhancements to POV-Ray and to port the software to platforms other * than those supported by the POV-Ray developers. There are strict rules * regarding how you are permitted to use this file. These rules are contained * in the distribution and derivative versions licenses which should have been * provided with this file. * * These licences may be found online, linked from the end-user license * agreement that is located at http://www.povray.org/povlegal.html *--------------------------------------------------------------------------- * This program is based on the popular DKB raytracer version 2.12. * DKBTrace was originally written by David K. Buck. * DKBTrace Ver 2.0-2.12 were written by David K. Buck & Aaron A. Collins. *--------------------------------------------------------------------------- * $File: //depot/povray/3.6-release/source/spheres.cpp $ * $Revision: #2 $ * $Change: 2939 $ * $DateTime: 2004/07/05 03:43:26 $ * $Author: root $ * $Log$ *****************************************************************************/#include "frame.h"#include "povray.h"#include "vector.h"#include "povproto.h"#include "bbox.h"#include "matrices.h"#include "objects.h"#include "spheres.h"BEGIN_POV_NAMESPACE/****************************************************************************** Local preprocessor defines******************************************************************************/const DBL DEPTH_TOLERANCE = 1.0e-6;/****************************************************************************** Static functions******************************************************************************/static int All_Sphere_Intersections (OBJECT *Object, RAY *Ray, ISTACK *Depth_Stack);static int All_Ellipsoid_Intersections (OBJECT *Object, RAY *Ray, ISTACK *Depth_Stack);static int Inside_Sphere (VECTOR IPoint, OBJECT *Object);static int Inside_Ellipsoid (VECTOR IPoint, OBJECT *Object);static void Sphere_Normal (VECTOR Result, OBJECT *Object, INTERSECTION *Inter);static void Sphere_UVCoord (UV_VECT Result, OBJECT *Object, INTERSECTION *Inter);static void Ellipsoid_Normal (VECTOR Result, OBJECT *Object, INTERSECTION *Inter);static void Translate_Sphere (OBJECT *Object, VECTOR Vector, TRANSFORM *Trans);static void Rotate_Sphere (OBJECT *Object, VECTOR Vector, TRANSFORM *Trans);static void Scale_Sphere (OBJECT *Object, VECTOR Vector, TRANSFORM *Trans);static void Invert_Sphere (OBJECT *Object);/****************************************************************************** Local variables******************************************************************************/METHODS Sphere_Methods ={ All_Sphere_Intersections, Inside_Sphere, Sphere_Normal, Sphere_UVCoord, (COPY_METHOD)Copy_Sphere, Translate_Sphere, Rotate_Sphere, Scale_Sphere, Transform_Sphere, Invert_Sphere, Destroy_Sphere};METHODS Ellipsoid_Methods ={ All_Ellipsoid_Intersections, Inside_Ellipsoid, Ellipsoid_Normal, Sphere_UVCoord, (COPY_METHOD)Copy_Sphere, Translate_Sphere, Rotate_Sphere, Scale_Sphere, Transform_Sphere, Invert_Sphere, Destroy_Sphere};/******************************************************************************* FUNCTION** All_Sphere_Intersection** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static int All_Sphere_Intersections(OBJECT *Object, RAY *Ray, ISTACK *Depth_Stack){ register int Intersection_Found; DBL Depth1, Depth2; VECTOR IPoint; SPHERE *Sphere = (SPHERE *)Object; Intersection_Found = false; if (Intersect_Sphere(Ray, Sphere->Center, Sqr(Sphere->Radius), &Depth1, &Depth2)) { if ((Depth1 > DEPTH_TOLERANCE) && (Depth1 < Max_Distance)) { VEvaluateRay(IPoint, Ray->Initial, Depth1, Ray->Direction); if (Point_In_Clip(IPoint, Object->Clip)) { push_entry(Depth1, IPoint, Object, Depth_Stack); Intersection_Found = true; } } if ((Depth2 > DEPTH_TOLERANCE) && (Depth2 < Max_Distance)) { VEvaluateRay(IPoint, Ray->Initial, Depth2, Ray->Direction); if (Point_In_Clip(IPoint, Object->Clip)) { push_entry(Depth2, IPoint, Object, Depth_Stack); Intersection_Found = true; } } } return(Intersection_Found);}/******************************************************************************* FUNCTION** All_Ellipsoid_Intersection** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static int All_Ellipsoid_Intersections(OBJECT *Object, RAY *Ray, ISTACK *Depth_Stack){ register int Intersection_Found; DBL Depth1, Depth2, len; VECTOR IPoint; RAY New_Ray; SPHERE *Sphere = (SPHERE *)Object; /* Transform the ray into the ellipsoid's space */ MInvTransPoint(New_Ray.Initial, Ray->Initial, ((SPHERE *)Object)->Trans); MInvTransDirection(New_Ray.Direction, Ray->Direction, ((SPHERE *)Object)->Trans); VLength(len, New_Ray.Direction); VInverseScaleEq(New_Ray.Direction, len); Intersection_Found = false; if (Intersect_Sphere(&New_Ray, Sphere->Center, Sqr(Sphere->Radius), &Depth1, &Depth2)) { if ((Depth1 > DEPTH_TOLERANCE) && (Depth1 < Max_Distance)) { VEvaluateRay(IPoint, New_Ray.Initial, Depth1, New_Ray.Direction); MTransPoint(IPoint, IPoint, ((SPHERE *)Object)->Trans); if (Point_In_Clip(IPoint, Object->Clip)) { push_entry(Depth1 / len, IPoint, Object, Depth_Stack); Intersection_Found = true; } } if ((Depth2 > DEPTH_TOLERANCE) && (Depth2 < Max_Distance)) { VEvaluateRay(IPoint, New_Ray.Initial, Depth2, New_Ray.Direction); MTransPoint(IPoint, IPoint, ((SPHERE *)Object)->Trans); if (Point_In_Clip(IPoint, Object->Clip)) { push_entry(Depth2 / len, IPoint, Object, Depth_Stack); Intersection_Found = true; } } } return(Intersection_Found);}/******************************************************************************* FUNCTION** Intersect_Sphere** INPUT** Ray - Ray to test intersection with* Center - Center of the sphere* Radius2 - Squared radius of the sphere* Depth1 - Lower intersection distance* Depth2 - Upper intersection distance* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/int Intersect_Sphere(RAY *Ray, VECTOR Center, DBL Radius2, DBL *Depth1, DBL *Depth2){ DBL OCSquared, t_Closest_Approach, Half_Chord, t_Half_Chord_Squared; VECTOR Origin_To_Center; Increase_Counter(stats[Ray_Sphere_Tests]); VSub(Origin_To_Center, Center, Ray->Initial); VDot(OCSquared, Origin_To_Center, Origin_To_Center); VDot(t_Closest_Approach, Origin_To_Center, Ray->Direction); if ((OCSquared >= Radius2) && (t_Closest_Approach < EPSILON)) { return(false); } t_Half_Chord_Squared = Radius2 - OCSquared + Sqr(t_Closest_Approach); if (t_Half_Chord_Squared > EPSILON) { Half_Chord = sqrt(t_Half_Chord_Squared); *Depth1 = t_Closest_Approach - Half_Chord; *Depth2 = t_Closest_Approach + Half_Chord; Increase_Counter(stats[Ray_Sphere_Tests_Succeeded]); return(true); } return(false);}/******************************************************************************* FUNCTION** Inside_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static int Inside_Sphere(VECTOR IPoint, OBJECT *Object){ DBL OCSquared; VECTOR Origin_To_Center; VSub(Origin_To_Center, ((SPHERE *)Object)->Center, IPoint); VDot(OCSquared, Origin_To_Center, Origin_To_Center); if (Test_Flag(Object, INVERTED_FLAG)) { return(OCSquared > Sqr(((SPHERE *)Object)->Radius)); } else { return(OCSquared < Sqr(((SPHERE *)Object)->Radius)); }}/******************************************************************************* FUNCTION** Inside_Ellipsoid** INPUT** OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static int Inside_Ellipsoid(VECTOR IPoint, OBJECT *Object){ DBL OCSquared; VECTOR Origin_To_Center, New_Point; /* Transform the point into the sphere's space */ MInvTransPoint(New_Point, IPoint, ((SPHERE *)Object)->Trans); VSub(Origin_To_Center, ((SPHERE *)Object)->Center, New_Point); VDot(OCSquared, Origin_To_Center, Origin_To_Center); if (Test_Flag(Object, INVERTED_FLAG)) { return(OCSquared > Sqr(((SPHERE *)Object)->Radius)); } else { return(OCSquared < Sqr(((SPHERE *)Object)->Radius)); }}/******************************************************************************* FUNCTION** Sphere_Normal** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Sphere_Normal(VECTOR Result, OBJECT *Object, INTERSECTION *Inter){ VSub(Result, Inter->IPoint, ((SPHERE *)Object)->Center); VInverseScaleEq(Result, ((SPHERE *)Object)->Radius);}/******************************************************************************* FUNCTION** Ellipsoid_Normal** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Ellipsoid_Normal(VECTOR Result, OBJECT *Object, INTERSECTION *Inter){ VECTOR New_Point; /* Transform the point into the sphere's space */ MInvTransPoint(New_Point, Inter->IPoint, ((SPHERE *)Object)->Trans); VSub(Result, New_Point, ((SPHERE *)Object)->Center); MTransNormal(Result, Result, ((SPHERE *)Object)->Trans); VNormalize(Result, Result);}/******************************************************************************* FUNCTION** Copy_Shere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/SPHERE *Copy_Sphere(OBJECT *Object){ SPHERE *New; New = Create_Sphere(); *New = *((SPHERE *)Object); New->Trans = Copy_Transform(((SPHERE *)Object)->Trans); return(New);}/******************************************************************************* FUNCTION** Translate_Sphere** INPUT** OUTPUT** RETURNS** AUTHOR** ?** DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Translate_Sphere(OBJECT *Object, VECTOR Vector, TRANSFORM *Trans){ SPHERE *Sphere = (SPHERE *) Object; if (Sphere->Trans == NULL) { VAddEq(Sphere->Center, Vector); Compute_Sphere_BBox(Sphere); } else { Transform_Sphere(Object, Trans); }}/******************************************************************************* FUNCTION** Rotate_Sphere** INPUT** OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Rotate_Sphere(OBJECT *Object, VECTOR, TRANSFORM *Trans){ SPHERE *Sphere = (SPHERE *) Object; if (Sphere->Trans == NULL) { MTransPoint(Sphere->Center, Sphere->Center, Trans); Compute_Sphere_BBox(Sphere); } else { Transform_Sphere(Object, Trans); }}/******************************************************************************* FUNCTION** Scale_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Scale_Sphere(OBJECT *Object, VECTOR Vector, TRANSFORM *Trans){ SPHERE *Sphere = (SPHERE *) Object; if ((Vector[X] != Vector[Y]) || (Vector[X] != Vector[Z])) { if (Sphere->Trans == NULL) { Sphere->Methods = &Ellipsoid_Methods; Sphere->Trans = Create_Transform(); } } if (Sphere->Trans == NULL) { VScaleEq(Sphere->Center, Vector[X]); Sphere->Radius *= fabs(Vector[X]); Compute_Sphere_BBox(Sphere); } else { Transform_Sphere(Object, Trans); }}/******************************************************************************* FUNCTION** Invert_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/static void Invert_Sphere(OBJECT *Object){ Invert_Flag(Object, INVERTED_FLAG);}/******************************************************************************* FUNCTION** Create_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/SPHERE *Create_Sphere(){ SPHERE *New; New = (SPHERE *)POV_MALLOC(sizeof(SPHERE), "sphere"); INIT_OBJECT_FIELDS(New, SPHERE_OBJECT, &Sphere_Methods) Make_Vector(New->Center, 0.0, 0.0, 0.0); New->Radius = 1.0; New->Trans = NULL; return(New);}/******************************************************************************* FUNCTION** Transform_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/void Transform_Sphere(OBJECT *Object, TRANSFORM *Trans){ SPHERE *Sphere = (SPHERE *)Object; if (Sphere->Trans == NULL) { Sphere->Methods = &Ellipsoid_Methods; Sphere->Trans = Create_Transform(); } Compose_Transforms(Sphere->Trans, Trans); Compute_Sphere_BBox(Sphere);}/******************************************************************************* FUNCTION** Destroy_Sphere** INPUT* * OUTPUT* * RETURNS* * AUTHOR** ?* * DESCRIPTION** -** CHANGES** -*******************************************************************************/void Destroy_Sphere(OBJECT *Object){#if(DUMP_OBJECT_DATA == 1) Debug_Info("{ // SPHERE \n"); DUMP_OBJECT_FIELDS(Object); Debug_Info("\t{ %f, %f, %f }, // Center\n", \ (DBL)((SPHERE *)Object)->Center[X], \ (DBL)((SPHERE *)Object)->Center[Y], \ (DBL)((SPHERE *)Object)->Center[Z]); \ Debug_Info("\t%f // Radius\n", (DBL)((SPHERE *)Object)->Radius); Debug_Info("}\n");#endif Destroy_Transform(((SPHERE *)Object)->Trans); POV_FREE (Object);}/******************************************************************************* FUNCTION** Compute_Sphere_BBox** INPUT** Sphere - Sphere* * OUTPUT** Sphere* * RETURNS* * AUTHOR** Dieter Bayer* * DESCRIPTION** Calculate the bounding box of a sphere.** CHANGES** Aug 1994 : Creation.*******************************************************************************/void Compute_Sphere_BBox(SPHERE *Sphere){ Make_BBox(Sphere->BBox, Sphere->Center[X] - Sphere->Radius, Sphere->Center[Y] - Sphere->Radius, Sphere->Center[Z] - Sphere->Radius, 2.0 * Sphere->Radius, 2.0 * Sphere->Radius, 2.0 * Sphere->Radius); if (Sphere->Trans != NULL) { Recompute_BBox(&Sphere->BBox, Sphere->Trans); }}/******************************************************************************* FUNCTION** Sphere_UVCoord** INPUT** OUTPUT** RETURNS** AUTHOR** Nathan Kopp (adapted from spherical_image_map)** DESCRIPTION** Find the UV coordinates of a sphere.* Map a point (x, y, z) on a sphere of radius 1 to a 2-d image. (Or is it the* other way around?)** CHANGES** -*******************************************************************************/static void Sphere_UVCoord(UV_VECT Result, OBJECT *Object, INTERSECTION *Inter){ DBL len, phi, theta; DBL x,y,z; VECTOR New_Point, New_Center; SPHERE *Sphere = (SPHERE *)Object; /* Transform the point into the sphere's space */ if (Object->UV_Trans != NULL) { MInvTransPoint(New_Point, Inter->IPoint, Object->UV_Trans); if (Sphere->Trans != NULL) MTransPoint(New_Center, Sphere->Center, Sphere->Trans); else Assign_Vector(New_Center, Sphere->Center); MInvTransPoint(New_Center, New_Center, Object->UV_Trans); } else { Assign_Vector(New_Point, Inter->IPoint); Assign_Vector(New_Center, Sphere->Center); } x = New_Point[X]-New_Center[X]; y = New_Point[Y]-New_Center[Y]; z = New_Point[Z]-New_Center[Z]; len = sqrt(x * x + y * y + z * z); if (len == 0.0) return; else { x /= len; y /= len; z /= len; } /* Determine its angle from the x-z plane. */ phi = 0.5 + asin(y) / M_PI; /* This will be from 0 to 1 */ /* Determine its angle from the point (1, 0, 0) in the x-z plane. */ len = x * x + z * z; if (len > EPSILON) { len = sqrt(len); if (z == 0.0) { if (x > 0) theta = 0.0; else theta = M_PI; } else { theta = acos(x / len); if (z < 0.0) theta = TWO_M_PI - theta; } theta /= TWO_M_PI; /* This will be from 0 to 1 */ } else /* This point is at one of the poles. Any value of xcoord will be ok... */ theta = 0; Result[U] = theta; Result[V] = phi;}END_POV_NAMESPACE