library("animation")  #Library to save GIFs
 
#Function to create a blank canvas / scene for drawing objects onto later
createSceneFunc <- function(bottomLeftX, bottomLeftY, width,height,main="",xlab="",ylab="",ann=T,xaxt=NULL,yaxt=NULL){
        plot(c(bottomLeftX, width), c(bottomLeftY,height), type = "n",ann=ann, xaxt=xaxt, yaxt=yaxt,main=main,xlab=xlab,ylab=ylab )
}
 
#Function to draw a box on the scene
createBoxFunc <- function(topLeftX, topLeftY, width, height, fillColour=NA, borderColour="black"){
     polygon(c(topLeftX,topLeftX+width,topLeftX+width,topLeftX),
             c(topLeftY,topLeftY,topLeftY-height,topLeftY-height),
              col = fillColour, border=borderColour)
}
 
#Function to draw a circle on the scene
createCircleFunc <- function(centerX,centerY,radius,fillColour=NA, borderColour="black"){
     symbols(centerX,centerY,circles=radius,inches=F,add=T,fg=borderColour,bg=fillColour)
}
 
#Parameters to control the simulation
simulation.timestep = 0.02
simulation.gravity = 1.8
simulation.numoftimesteps = 2000
 
#Define the size of the scene (used to visualise what is happening in the simulation)
scene.bottomLeftX = 0
scene.bottomLeftY = -1
scene.width = 10
scene.height=10
 
#This is the object the bouncing ball is going to hit
platform.x = scene.bottomLeftX+1
platform.y = 0
platform.width= scene.width - 2
platform.height=0.5
platform.colour = "red"
 
#This is just a box resting on the left end of the platform to practise drawing things
leftbluebox.x = platform.x
leftbluebox.y = platform.y+0.5
leftbluebox.width=0.5
leftbluebox.height=0.5
leftbluebox.colour = "blue"
 
#This is just a box resting on the right end of the platform to practise drawing things
rightbluebox.x = platform.x+platform.width-0.5
rightbluebox.y = platform.y+0.5
rightbluebox.width=0.5
rightbluebox.height=0.5
rightbluebox.colour = "blue"
 
#This is the ball that is going to be bouncing
ball.y=10
ball.x=5
ball.radius = 0.25
ball.yvelocity = 0
ball.yacceleration = 0
ball.coefficientofrestitution = 0.85 #This is a physics term to describe how much velocity as pct is kept after a bounce
ball.colour="purple"
 
#Some variables to store various time series values of the simulation
logger.yposition = rep(NA,simulation.numoftimesteps)
logger.yvelocity = rep(NA,simulation.numoftimesteps)
logger.yacceleration = rep(NA,simulation.numoftimesteps)
 
#Some settings to control the charts used to plot the logged variables
plotcontrol.yposition.ylim = c(0,10)
plotcontrol.yvelocity.ylim = c(-6,6)
plotcontrol.yacceleration.ylim = c(-simulation.gravity,400)
 
runSimulationFunc <- function(){
  #Main simulation loop
  for(i in seq(1,simulation.numoftimesteps)){
 
     #Equations of motion
     ball.yacceleration = -simulation.gravity
     ball.y = ball.y + ball.yvelocity*simulation.timestep
     ball.yvelocity = ball.yvelocity+ball.yacceleration*simulation.timestep
 
     #Logic to check is there has been a collision between the ball and the platform
     if(ball.y-ball.radius <= platform.y){
        #There has been a collision
        newyvelocity = -ball.yvelocity*ball.coefficientofrestitution
        ball.yacceleration = (newyvelocity - ball.yvelocity)/simulation.timestep
        ball.yvelocity = newyvelocity
        ball.y = ball.radius+platform.y
     }
 
     #Log the results of the simulation
     logger.yposition[i] <- ball.y
     logger.yvelocity[i] <- ball.yvelocity
     logger.yacceleration[i] <- ball.yacceleration
 
    #Plot the simulation
    #The layout command arranges the charts
    layout(matrix(c(1,2,1,3,1,4,5,5), 4, 2, byrow = TRUE),heights=c(3,3,3,1))
    par(mar=c(3,4,2,2) + 0.1)
 
    #Create the scene and draw the various objects
    #Create the scene
    createSceneFunc(scene.bottomLeftX,scene.bottomLeftY,scene.width,scene.height,
                 main="Simulation of Bouncing Ball - www.gekkoquant.com",xlab="",ylab="Ball Height",xaxt="n")
 
    #Draw the platform the ball lands on
    createBoxFunc(platform.x,platform.y,platform.width,platform.height,platform.colour)
 
    #Draw a box on the left off the platform
    createBoxFunc(leftbluebox.x,leftbluebox.y,leftbluebox.width,leftbluebox.height,leftbluebox.colour)
 
    #Draw a box on the right of the platform
    createBoxFunc(rightbluebox.x,rightbluebox.y,rightbluebox.width,rightbluebox.height,rightbluebox.colour)
 
    #Draw the ball
    createCircleFunc(ball.x,ball.y,ball.radius,ball.colour,borderColour=NA)
 
    #Plot the logged variables
    plot(logger.yposition, type="l", ylim=plotcontrol.yposition.ylim, ylab="Y POSITION")
    plot(logger.yvelocity, type="l", ylim=plotcontrol.yvelocity.ylim, ylab="Y VELOCITY")
    plot(logger.yacceleration, type="l",ylim=plotcontrol.yacceleration.ylim, ylab="Y ACCELERATION")
 
    #Plot a progress bar
    par(mar=c(2,1,1,1))
    plot(-5, xlim = c(1,simulation.numoftimesteps), ylim = c(0, .3), yaxt = "n", xlab = "", ylab = "", main = "Iteration")
    abline(v=i, lwd=5, col = rgb(0, 0, 255, 255, maxColorValue=255))
 
  }
}
 
oopt = ani.options(ani.width = 1200, ani.height = 800, other.opts = "-pix_fmt yuv420p -b 600k")
saveVideo(runSimulationFunc(),interval = simulation.timestep,ani.options=oopt,video.name="bounce.mp4")
ani.options(oopt)