/* 
Title:     Rolling Racers, a classic physics demonstration for moment of inertia
Version:   1.0
Author:    Lucas V. Barbosa
Date:      December 23, 2012
Coded in:  POV-Ray 3.6
License:   Public Domain


Updated 2021-06-16


Info
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
Author user page: http://en.wikipedia.org/wiki/User:Kieff

The final animation at Wikimedia Commons:
     OGG Theora: http://commons.wikimedia.org/wiki/File:Rolling_Racers_-_Moment_of_inertia (HD).ogv
                 http://commons.wikimedia.org/wiki/File:Rolling_Racers_-_Moment_of_inertia.ogv
   Animated GIF: http://commons.wikimedia.org/wiki/File:Rolling_Racers_-_Moment_of_inertia.gif
                              
POV-Ray command line options: +A0.2 +AM2 +R4 -j +w1920 +h1080 +kff252 File_Gamma=sRGB
                              +A0.2           = anti-aliasing
                              +AM2            = anti-aliasing mode
                              +R4             = anti-aliasing quality
                              -j              = no jittering. Jittering sucks.
                              +w1920          = 1920 px of width
                              +h1080          = 1080 px of height (so, widescreen, Full HD)
                              +kff252         = render 252 frames (4.2 seconds * 60 frames per second)
                              File_Gamma=sRGB = use sRGB colour space
                              
Aspect Ratio: 16:9 (widescreen)
 
 
Notes
¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
I have not written a full-blown physics simulator here!

All objects are assumed to have the same mass and radius (which, in reality, is irrelevant), and to roll without slipping.
The idea is to illustrate how important the geometry of the object is in defining its rotational properties.

It shouldn't be too hard to add objects with different properties.
*/

/* ----------------------- BASIC DIRECTIVES ----------------------- */

#version 3.7;
#include "rad_def.inc"

// In the beginning, there were global settings
#local p_start		=	64/image_width;
#local p_end_final	=	1/image_width;

global_settings {
  assumed_gamma 1.0
  radiosity {
    pretrace_start p_start
    pretrace_end   p_end_final
    count 60
    nearest_count 20
    error_bound 0.4
    recursion_limit 3
    minimum_reuse 0.005
  }
}

// ... Lights! ...
// An area light (creates soft shadows)
light_source {
  0*x                 // light's position (translated below)
  color rgb 1.0       // light's color
  area_light
  <5, 0, 0> <0, 0, 5> // lights spread out across this distance (x * z)
  10, 10                // total number of lights in grid (4x*4z = 16 lights)
  adaptive 1          // 0,1,2,3...
  circular            // make the shape of the light circular
  orient              // orient light
  translate <10, 50, -40>   // <x y z> position of light
}



// ... Camera! ...
camera {
  location  <-0.1, 2.5, -4.0>*1
  right     x*image_width/image_height
  look_at   <0.0, 0.0,  0.0>
  direction 2.2*z
  rotate -30*y
  translate -1.25*x
}



/* ----------------------- MACROS ----------------------- */

// Generates a ramp using a prism and the given parameters
// Downwards is towards a positive x
// Angle = incline angle, in degrees
// Length, Width = self-explanatory, in metres
#macro Ramp(Angle, Length, Width)
  #local Height = RampLength * tan(radians(Angle));
    prism {
      linear_sweep
      linear_spline
      0,
      Width,
      4,
      <0,0>, <-Length,0>, <-Length,Height>, <0,0>
      rotate -90*x
      translate Width/2*z
    }
#end



// Generates a red or green flag
// Color of flag is based on position of the associated object
#macro Flag(ObjectPosition)
    union {
        union {
            cylinder {
                0, y*0.3, 0.0075
            }
            sphere {
                0.3*y, 0.0075*1.5
            }
            pigment { color rgb 0.9 }
            finish {
                specular 0.5
                roughness 0.2
                reflection 0.2
            }
        }
        prism {
          linear_sweep
          linear_spline
          0,
          1/100,
          4,
          <0,0>, <0,1>, <1,0.75>, <0,0>
          rotate 90*x
          scale 0.08*y
          scale 0.15*x
          translate y*0.29
            pigment {
                color rgb
                // If the center of mass of the object has positive x, then it's past the finish line
                #if (ObjectPosition.x > 0)
                    <0,1,0>*1.22
                #else
                    <1,0,0>*1.2
                #end
            }
        }
    }
#end



/* ----------------------- SCENE SETUP ----------------------- */

// Ramp settings
#declare RampAngle = 5; // degrees
#declare RampLength = 3; // metres
#declare RampWidth = 1.25;  // metres

// Physics settings
#declare Gravity = 9.8;  // m/s²

// User settings
#declare PhysicsAwesome = true; // VERY important!
#declare Contrast = 1.6; // multiplier for pattern color for objects. 1 = solid color = lame

// Number of objects to use, so they can be equally spaced automatically
#declare NumObjects = 4; 




/* ----------------------- OBJECTS ----------------------- */


// Radius and mass of objects, assumed to be equal for all of them
#declare Object_Radius = 0.1; // metres
#declare Object_Mass = 0.5; // kg

// Base moment of inertia
// This bit is constant for all objects with circular symmetry along the rotation axis, so it can be reused
#declare MR2 = Object_Mass * pow(Object_Radius,2);



// Solid Sphere (spoiler: our winner!)
#declare I_SolidSphere = 2/5 * MR2;
#declare Object_SolidSphere = sphere {
    0, Object_Radius
    pigment {
        checker
        color rgb <1, 0.7, 0.3>,
        color rgb <1, 0.7, 0.3>*Contrast
        translate <1,0,1>*0.5
        scale Object_Radius/3.3*y
    }
    finish {
        specular 0.3
        roughness 0.01
    }
}



// Hollow Sphere
#declare I_HollowSphere = 2/3 * MR2;
#declare Object_HollowSphere = difference {
    sphere {
        0, Object_Radius
    }
    sphere {
        0, Object_Radius * 0.95
    }
    pigment {
        checker
        // Added some transparency to the hollow sphere so it looks thin-walled. Seems to be a bit more convincing.
        color rgb <1, 0.2, 0.1> filter 0.3 transmit 0.5,
        color rgb <1, 0.2, 0.1>*Contrast filter 0.3 transmit 0
        translate <1,0,1>*0.5
        scale Object_Radius/3.3*y
    }
    finish {
        specular 0.3
        roughness 0.01
    }
}



// Solid Cylinder
#declare I_SolidCylinder = 1/2 * MR2;
#declare Object_SolidCylinder = cylinder {
    -Object_Radius*0.75*z, Object_Radius*0.75*z, Object_Radius
    pigment {
        checker
        color rgb <0.3, 0.7, 1>,
        color rgb <0.3, 0.7, 1>*Contrast
    }
    finish {
        specular 0.3
        roughness 0.01
    }
}



// Hollow Cylinder (never bet on this one!)
#declare I_HollowCylinder = MR2;
#declare Object_HollowCylinder = difference {
    cylinder {
        -Object_Radius*0.75*z, Object_Radius*0.75*z, Object_Radius
    }
    cylinder {
        -Object_Radius*2*z, Object_Radius*2*z, Object_Radius*0.95
    }
    pigment {
        checker
        color rgb <0.7, 1, 0.3>,
        color rgb <0.7, 1, 0.3>*Contrast
    }
    finish {
        specular 0.3
        roughness 0.01
    }
}




/* ----------------------- INTERNAL VARIABLES ----------------------- */
// No need to change these!

#declare InitialPosition = RampLength * 0.9; // how far to the left objects start

#declare InitialHeight = InitialPosition * tan(radians(RampAngle));

#declare VerticalOffset = Object_Radius / cos(radians(RampAngle)); // the vertical offset places the objects tangential to the plane
                                                                   // Therefore, I can treat positions of the center of mass without worrying about contact angles

#declare ObjectSpacing = RampWidth / NumObjects; // distributes available ramp width across objects

#declare BasePosition = < -InitialPosition , InitialHeight + VerticalOffset , -RampWidth/2 + ObjectSpacing/2 >; // base position


/* ----------------------- PHYSICS ----------------------- */

// Duration of simulation
#declare Duration = 4.2;          // seconds - 4.2 seconds is just a good time for when everything has left the frame
#declare time = clock * Duration; // time, in seconds, at the current frame


// ... Action!


// This macro returns the angular acceleration alpha (in radians/second²) for a given body of radius R, mass M and moment of inertia I.
// As you can see, and it is somewhat surprising, the mass and radius are irrelevant for objects of the same shape.
// In other words, two solid cylinders with different masses and radii will roll at the same rate!
// I could have simplified the code, but it would become a bit cryptic to those studying its inner workings and the physics behind it.
// See the file's description page (check header at the top of this source code) for the derivation of this formula
#macro Alpha(Radius, Mass, I)
  (Gravity * sin(radians(RampAngle))) / ((1 + I/(Mass*pow(Radius,2))) * Radius)
#end


// half of t^2. Integration term for t
#declare t2 = pow(time,2) / 2;

// Angular accelerations (alphas), in radians/s²)
#declare Alpha_SolidCylinder  = Alpha(Object_Radius, Object_Mass, I_SolidCylinder);
#declare Alpha_HollowCylinder = Alpha(Object_Radius, Object_Mass, I_HollowCylinder);
#declare Alpha_SolidSphere    = Alpha(Object_Radius, Object_Mass, I_SolidSphere);
#declare Alpha_HollowSphere   = Alpha(Object_Radius, Object_Mass, I_HollowSphere);

// Rotation amount (thetas), in radians, from integrating for constant acceleration
#declare Theta_SolidCylinder  = -(Alpha_SolidCylinder  * t2);
#declare Theta_HollowCylinder = -(Alpha_HollowCylinder * t2);
#declare Theta_SolidSphere    = -(Alpha_SolidSphere    * t2);
#declare Theta_HollowSphere   = -(Alpha_HollowSphere   * t2);


// Gets position for an object given its theta and radius
#macro GetPosition(Theta, Radius)

  // Distance travelled considering rotation without slipping
  // This is neither vertical or horizontal, but total length
  #local DistanceTraveled = -Theta * Radius;
  
  #local dist_x = + DistanceTraveled * cos(radians(RampAngle)); // total motion in the x direction
  #local dist_y = - DistanceTraveled * sin(radians(RampAngle)); // total motion in the y direction
  
  
  // Position so far
  #declare pos = BasePosition
               + x * dist_x
               + y * dist_y
  ;
  
  // If the object has sunk into the ground so far, then it's actually off the ramp
  #if (pos.y <= Radius)
      // Therefore, we disregard the movemevent in the y direction after this point
      // and put it all in the x direction. The math is simplified because the ramp ends at the origin
      #local r = pos.x / cos(radians(RampAngle));
      #declare pos =
               + x * r                    // move it the correct amount horizontally
               + y * Radius                // stick it to the ground
               + z * BasePosition.z       // keep z position
    ;
  #end

  pos
#end


// Defines the new position vectors
#declare Pos_SolidSphere    = GetPosition(Theta_SolidSphere, Object_Radius)    + z * ( ObjectSpacing*2 );
#declare Pos_SolidCylinder  = GetPosition(Theta_SolidCylinder, Object_Radius)  + z * ( ObjectSpacing*0 );
#declare Pos_HollowCylinder = GetPosition(Theta_HollowCylinder, Object_Radius) + z * ( ObjectSpacing*1 );
#declare Pos_HollowSphere   = GetPosition(Theta_HollowSphere, Object_Radius)   + z * ( ObjectSpacing*3 );




/* ----------------------- SCENE DECLARATIONS ----------------------- */


// Ground! Ha! I wonder if it will be friends with me?
plane {
  y, 0
  pigment { color rgb <1,1,1>*1.5 }
     finish {
        specular 0.5
        roughness 0.1
        //reflection 0.3     // reflection seemed superfluous and added nothing to it
    }
}


// Place ramp in the scene
union {
    object {
        Ramp(RampAngle, RampLength, RampWidth)
        // Add some color. The pigment creates the stripes for each "track"
        pigment {
          checker
          color rgb <0.8,1,1>
          color rgb <0.8,1,1>*0.8
          scale (RampWidth/4)*z
          scale RampLength*x
        }
        // Ooh, shiny!
        finish {
            specular 0.5
            roughness 0.1
            reflection 0.2
        }
    }
    // Side of the ramp, just so it stands out better
    object {
        Ramp(RampAngle, RampLength, 0.01)
        pigment {
            color rgb <0.8,0.8,1>*0.9
        }
        translate -(RampWidth/2 + 0.01/2)*z
    }
}


// Finish line
box {
    <0,0,-RampWidth/2>, <0.15,1/10000,RampWidth/2>
    pigment {
        checker
        color rgb <1,1,1>*1.7
        color rgb <1,1,0.9>*0.2
        scale RampWidth/(NumObjects*4)
    }
}


// The four little flags
object {
    Flag(Pos_HollowSphere)
    translate z*(RampWidth/2*z - ObjectSpacing*0)
}
object {
    Flag(Pos_SolidSphere)
    translate z*(RampWidth/2*z - ObjectSpacing*1)
}
object {
    Flag(Pos_HollowCylinder)
    translate z*(RampWidth/2*z - ObjectSpacing*2)
}
object {
    Flag(Pos_SolidCylinder)
    translate z*(RampWidth/2*z - ObjectSpacing*3)
}



// Finally, our objects
object {
    Object_SolidCylinder
    rotate z*degrees(Theta_SolidCylinder)
    translate Pos_SolidCylinder
}

object {
    Object_HollowCylinder
    rotate z*degrees(Theta_HollowCylinder)
    translate Pos_HollowCylinder
}

object {
    Object_SolidSphere
    rotate z*degrees(Theta_SolidSphere)
    translate Pos_SolidSphere
}

object {
    Object_HollowSphere
    rotate z*degrees(Theta_HollowSphere)
    translate Pos_HollowSphere
}


// ... and cut!

/* END OF FILE */