This week only. Pastebin PRO Accounts Christmas Special! Don't miss out!Want more features on Pastebin? Sign Up, it's FREE!
Guest

Interception of objects with linear motion

By: a guest on Mar 31st, 2013  |  syntax: Python  |  size: 3.17 KB  |  views: 29  |  expires: Never
download  |  raw  |  embed  |  report abuse  |  print
Text below is selected. Please press Ctrl+C to copy to your clipboard. (⌘+C on Mac)
  1. # PYTHON 2.7
  2.  
  3. # examples:
  4. def test():
  5.         a = {
  6.                 'p': (-1.0,4.0),
  7.                 'v': (0.0,-2.0),
  8.         }
  9.        
  10.         b = {
  11.                 'p': (10.0, 8.0),
  12.                 'v': (-5.0, 5.0),
  13.         }
  14.        
  15.         print InterceptRelative(a,b,1)
  16.         print InterceptRelative(a,b,10)
  17.         print InterceptRelative(a,b,100)
  18.        
  19.         print InterceptAbsolute(a,b,1)
  20.         print InterceptAbsolute(a,b,10)
  21.         print InterceptAbsolute(a,b,100)
  22.  
  23. # suppose source will shoot at target
  24. # the projectile's speed relative to source will equal speed
  25. def InterceptRelative(source,target,speed):
  26.         # we will work in 'source-relative-space'
  27.         # it is as if source is stationary at the origin
  28.         # and we just observe target
  29.         delta_position = Vsub(target['p'], source['p']) # t-s gives the vector FROM s TO t
  30.         delta_velocity = Vsub(target['v'], source['v'])
  31.        
  32.         return Intercept(delta_position, delta_velocity, speed)
  33.  
  34. # suppose source shoots at target
  35. # but the projectile's speed relative to the universe equal's speed
  36. def InterceptAbsolute(source,target,speed):
  37.         # we will work in 'source-relative-space'
  38.         # BUT we pretend source's velocity equals zero
  39.         delta_position = Vsub(target['p'], source['p'])
  40.         delta_velocity = target['v']
  41.        
  42.         return Intercept(delta_position, delta_velocity, speed)
  43.  
  44. # the actual interception function
  45. # there is a sphere (circle in 2D) at the origin expanding at rate speed
  46. # when does it touch the target?
  47. def Intercept(target_pos, target_vel, speed):
  48.         # there's some math going on here
  49.         # roughly speaking, we imagine a sphere expanding from the origin at rate speed
  50.         # and compute when that sphere will intersect the point
  51.         a = Vsqmag(target_vel) - (speed**2.0)
  52.         b = 2.0*Vdot(target_pos, target_vel)
  53.         c = Vsqmag(target_pos)
  54.        
  55.         # use our solver
  56.         times = Quadratic(a,b,c)
  57.        
  58.         if times is None:
  59.                 # there are no real solutions - therefore it's not possible to hit the target
  60.                 return None
  61.        
  62.         # choose the non-negative times
  63.         times = [t for t in times if t >= 0]
  64.        
  65.         if len(times) == 0:
  66.                 # we can only hit it in the past! - so not really at all :(
  67.                 return None
  68.        
  69.         # choose the soonest time when we can hit it
  70.         t = min(times)
  71.         collisionpt = Vadd(target_pos, Vscale(target_vel, t))
  72.         print t, collisionpt
  73.         return Vscale(Vnormalize(collisionpt), speed)
  74.        
  75.  
  76. # some vector operations
  77. # note that they are dimensionless!
  78. from itertools import izip
  79. from math import sqrt, fsum
  80.  
  81. def Vadd(v1, v2):
  82.         return tuple(x+y for (x,y) in izip(v1,v2))
  83.  
  84. def Vsub(v1, v2):
  85.         return tuple(x-y for (x,y) in izip(v1,v2))
  86.  
  87. def Vscale(v, s):
  88.         return tuple(x*s for x in v)
  89.  
  90. def Vdot(v1,v2):
  91.         return fsum(x*y for (x,y) in izip(v1,v2))
  92.  
  93. def Vsqmag(v):
  94.         return Vdot(v,v)
  95.  
  96. def Vmag(v):
  97.         return sqrt(Vsqmag(v))
  98.  
  99. def Vnormalize(v):
  100.         return Vscale(v, 1.0/Vmag(v))
  101.  
  102. # a simple quadratic solver for the form
  103. # Ax^2 + Bx + C = 0
  104. # returns results as a tuple (low x, high x)
  105. # returns None if the result is complex
  106. def Quadratic(A,B,C):
  107.         if A == 0.0:
  108.                 # this is not quadratic
  109.                 if B == 0.0:
  110.                         # this is not even linear!
  111.                         return None
  112.                 else:
  113.                         return (-C/B,)
  114.                
  115.         radical = B*B - 4.0*A*C
  116.         if(radical < 0):
  117.                 return None
  118.        
  119.         part1 = -B / (2.0*A)
  120.         part2 = sqrt(radical)/(2.0*A)
  121.        
  122.         return (part1 - part2, part1 + part2)
  123.  
  124. if(__name__ == "__main__"):
  125.         test()
clone this paste RAW Paste Data