// Persistence of Vision Ray Tracer Scene Description File
// File: Hypercube.pov
// Vers: 1.0
// Desc: Render 4d objects
// Date: mm/dd/yy
// Auth: J. Harshbarger (Hypercubed, http://www.hypercubed.com/)          
// Thanks to Steven Richard Hollasch and his Thesis on Four-Space Visualization of 4D Objects (August 1991)
//              http://research.microsoft.com/~hollasch/thesis/default.htm

#version 3.1;

#include "colors.inc"


global_settings
{
  assumed_gamma 1.0
}

// Scene Settins
camera {
  perspective
  //orthographic
  //fisheye
  location <-1.2, -.6, .6>
  up <0, 1, 0>
  sky <0, 0, 1>
  look_at <0, 0, 0>
}
 
 
background { color White }


light_source
{
    0*x
    color red 1.0  green 1.0  blue 1.0

    shadowless
    translate <-50, -150, 100>
}

// 4-d conversion functions. 

 #macro v4to3 (A)
        #local V = v4dotM(A-From4, TM);  //Translate to eye coordinates
        #if (Per4 = 1)
                #local R = (V.t*Tp);  //4d Perspective
        #else
                #local R = v4len(From4);
        #end        
        <V.x/R,V.y/R,V.z/R>  //3d vector
 #end
     
 #macro v4len(V)  //Length of 4 vector
     sqrt(v4dot(V, V))
 #end
 
 #macro Norm4(V)  //Normalize 4 vector
     V/v4len(V)
 #end
 
#macro Cross4 (U, V, W)  // 4d cross product

    // Calculate intermediate values.

    #local A = (V.x * W.y) - (V.y * W.x);
    #local B = (V.x * W.z) - (V.z * W.x);
    #local C = (V.x * W.t) - (V.t * W.x);
    #local D = (V.y * W.z) - (V.z * W.y);
    #local E = (V.y * W.t) - (V.t * W.y);
    #local F = (V.z * W.t) - (V.t * W.z);

    // Calculate the result-vector components.

    #local X = ((U.y * F) - (U.z * E) + (U.t * D));
    #local Y = ((U.x * F) + (U.z * C) - (U.t * B));
    #local Z = ((U.x * E) - (U.y * C) + (U.t * A));
    #local T = (-(U.x * D) + (U.y * B) - (U.z * A));
    
    <X, Y, Z, T>
#end
    
     
 #macro v4dot (A, B)  //4d dot product
       A.x*B.x+A.y*B.y+A.z*B.z+A.t*B.t
 #end
 
 #macro v4dotM (A, M)  //4d vector dot 4x4 matrix
       <v4dot(A, M[0]), v4dot(A, M[1]), v4dot(A, M[2]), v4dot(A, M[3])>
 #end

// 4-d Objects
 
 #macro point3d (P, pradius, pcolor)
         sphere
        {
                P
                pradius
                pigment { pcolor }
        }
 #end
  
 #macro point4d (P, pradius, pcolor)
        point3d (v4to3(P), pradius, pcolor)
 #end 
 
 #macro vector3d (S, E, vradius, vcolor)
        #if (!((S.x = E.x) & (S.y = E.y) & (S.z = E.z)))
        cone
        {
                S, vradius,
                E, vradius
                pigment { vcolor }
        }
        #end 
 #end
 
 #macro vector4d (S, E, vradius, vcolor)
        vector3d (v4to3(S), v4to3(E), segradius, vcolor)
 #end

 #macro Seg(S, E, segcolor)  
        union
        {
                point3d (S, segradius, segcolor)
                #if (!((S.x = E.x) & (S.y = E.y) & (S.z = E.z)))
                        point3d (E, segradius, segcolor)
                        vector3d (S, E, segradius, segcolor)
                #end                     
        }
#end

#macro Seg4(S, E, segcolor)   //Simply a cone with end points
        Seg(v4to3(S), v4to3(E), segcolor)
#end




// Environment Varaibles

#declare T = pi/2*clock;
#declare Tp = tan((2*pi/360)*45);

#declare From4 = 4*<cos(T), 0, -sin(T), 0>;
#declare To4 = <0, 0, 0, 0>;
#declare Up4    = <0, 1, 0, 0>; 
#declare Over4  = <0, 0, 1, 0>; 
#declare Per4 = 1;


#declare D = Norm4(To4 - From4);
#declare A = Norm4(Cross4(Up4, Over4, D));
#declare B = Norm4(Cross4(Over4, D, A));
#declare C = Cross4(D, A, B);

#declare TM =   //Translation Matrix
array[4]
{ A, B, C, D }


#declare segcolor = Red;
#declare vcolor = Yellow;
#declare pcolor = Red;
#declare segradius = .015;

// Hypercube

#macro corner(C)
        Seg4(C, <-C.x, C.y, C.z, C.t>, segcolor)
        Seg4(C, <C.x, -C.y, C.z, C.t>, segcolor)
        Seg4(C, <C.x, C.y, -C.z, C.t>, segcolor)
        Seg4(C, <C.x, C.y, C.z, -C.t>, segcolor)
#end

corner(<1, 1, 1, 1>)
corner(<-1, 1, 1, 1>)
corner(<1, -1, 1, 1>)
corner(<-1, -1, 1, 1>)
corner(<1, 1, -1, 1>)
corner(<-1, 1, -1, 1>)
corner(<1, -1, -1, 1>) 
corner(<-1, -1, -1, 1>)

corner(<1, 1, 1, -1>)
corner(<-1, 1, 1, -1>)
corner(<1, -1, 1, -1>)
corner(<-1, -1, 1, -1>)
corner(<1, 1, -1, -1>)
corner(<-1, 1, -1, -1>)
corner(<1, -1, -1, -1>) 
corner(<-1, -1, -1, -1>)