====== Nötige Methode für die Bezier-Klasse ======
* x(t): Liefert die Position zum Parameter $t$.
* v(t): Liefert die Geschwindigkeit zum Parameter $t$ ("exakt").
* len(tmin, tmax): Liefert die Länge der Kurve zwischen den Parametern tmin, tmax (angenähert).
* vnorm(t): Geschwindkeitsvektor mit Länge 1 zum Parameter $t$ ("exakt").
* adl(t): Die Ableitung von vnorm nach der **Länge** der Kurve
*
===== Todo =====
* forward: Interpolation der verbleibenden Länge, wenn über die Kurve hinaus gegangen wird
* adl(t):
* a(t,v,g): Effektiver Beschleunigungsvektor beim Punkt $p(t)$, wenn der Geschwindigkeitsbetrag v und Inverse Erdbeschleunigung gegeben ist.
* koordsyst(t): Orthonormales Koordinatensystem (Liste von 3 Vektoren) mit $x$ parallel zu $v$, $y$ parallel zu $a$ und $z = x \times y$.
# -*- coding: utf-8 -*-
from vector import Vector
class Bezier:
def __init__(self,pts):
if len(pts)==4:
self.pts = [p.copy() for p in pts]
if len(pts)==3:
self.pts = [pts[0].copy(), 1/3*pts[0]+2/3*pts[1], 1/3*pts[2]*2/3*pts[1], pts[2].copy()]
if len(pts)==2:
self.pts = [pts[0].copy(), 2/3*pts[0]+1/3*pts[1], 1/3*pts[0]+2/3*pts[1], pts[1].copy()]
# Punkt für Parameter t (in [0,1])
def x(self,t):
return (1.0-t)**3*self.pts[0]+\
3*(1.0-t)**2*t*self.pts[1]+\
3*(1.0-t)*t**2*self.pts[2]+\
t**3*self.pts[3]
# Geschwindigkeit (dx/dt) in t-Parametrierung
def v(self,t):
return -3*(1.0-t)**2*self.pts[0]+\
3*(1.0-4*t+3*t**2)*self.pts[1]+\
3*(2*t-3*t**2)*self.pts[2]+\
3*t**2*self.pts[3]
# Normalisierte Geschwindigkeit
def vnorm(self,t):
return self.v(t).normalize()
def adl(self,t, dt=0.000001): # Ableitung der normalisierten Geschwindigkeit nach der Länge
dl = self.v(t).length*dt # (self.x(t+dt)-self.x(t-dt)).length
return (self.vnorm(t+dt/2)-self.vnorm(t-dt/2))/dl
# Effektive Beschleunigung im Punkt zum Parameter t, mit
# gegebenem Geschwindigkeitsbetrag und Gravitationsvektor
def aeff(self, t, v, g):
return self.adl(t)*v*v-g
# Orthonnormales Bahn-Koordinatensystem im Punkt zum Parameter t, mit
# gegebenem Geschwindigkeitsbetrag und Gravitationsvektor.
# Resultat (v,a, v cross a)
def koordsyst(self, t, v, g):
x = self.vnorm(t)
# Effektive Beschleunigung
aeff = self.aeff(t,v,g)
# Projektion von aeff auf x (d.h. normalisierte Geschwindigkeit)
av = x.dot(aeff)/(x.dot(x))*x
# Effektive Beschleunigung in Normalebene zur Bahn.
# (ist mathematisch rechtwinklig!)
aa = (aeff-av).normalize()
z = x.cross(aa)
return (x,aa,z)
# Laenge der Bezier-Kurve
def length(self, tmin=0.0, tmax=1.0, steps=1000):
l = 0.0
# Zeitschritt
dt = (tmax-tmin)/steps
for i in range(steps):
# Den i-ten Zeitpunkt von steps Zeitpunkten
ti = (i+0.5)*dt+tmin
l += self.v(ti).length*dt
return l
# Input: Startzeitpunkt t, Streckenlaenge l
# Output: Endzeitpunkt t oder -Restlaenge
def forward(self, l, startt, dt=0.0001):
lastl = l
lastt = startt
while (startt<1.0 and l>0):
dl = self.v(startt).length*dt
lastl = l
lastt = startt
l -= dl
startt += dt
if startt<1.0: # Noch auf der Bezierkurve
return (startt-lastt)*lastl/(lastl+abs(l))+lastt
if l>0: # Es ist noch Laenge uebrig, also auf naechster Kurve
return -(lastl - (self.x(1)-self.x(lastt)).length)
return 1.0 # Sonst genau auf dem Ende der Kurve
# Zeichnet den Spline
# Achtung: Funktioniert nur, wenn diese Datei ausgeführt wird...
def draw(self, steps=1000):
move(self.pts[0][0], self.pts[0][1]);
for i in range(1,steps+1):
p = self.x(1.0*i/steps)
lineTo((p[0],p[1]))
if __name__== "__main__":
from gpanel import *
makeGPanel(Size(600,600))
window(0, 599, 0, 599)
b = Bezier([Vector([100,500,0]), Vector([200,0,0]), Vector([300,300,0]), Vector([500,300,0])])
b.draw()
print(b.length(0,1,10000))
t = 0.0
while (t>=0.0):
p = b.x(t)
a = b.adl(t)
pa = p+a*10000
# move(p.x, p.y)
# circle(5)
line(p.x,p.y, pa.x, pa.y)
t = b.forward(10,t)
print(t)
====== Bahn-Koordinatensystem ======
===== Bewegung mit konstantem Geschwindigkeitsbetrag 1 =====
Kreisbewegung mit Geschwindigkeitsbetrag 1:
\begin{align*}
x(t) & = r \cdot \textrm{e}^{\textrm{i}t \cdot \frac{1}{r}}\\
v(t) & = x'(t) = \textrm{i}\textrm{e}^{\textrm{i}t \cdot \frac{1}{r}}\\
a(t) & = v'(t) = -\frac{1}{r} \cdot \textrm{e}^{\textrm{i}t \cdot \frac{1}{r}}\\
\end{align*}
Beträgt jetzt der Betrag der Geschwindigkeit $\lambda$ anstatt 1, ändern sich die Gleichungen wie folgt:
\begin{align*}
x(t) & = r \cdot \textrm{e}^{\textrm{i} \lambda t \cdot \frac{1}{r}}\\
v(t) & = x'(t) = \textrm{i}\lambda \cdot \textrm{e}^{\textrm{i} \lambda t \cdot \frac{1}{r}}\\
a(t) & = v'(t) = -\lambda^2 \cdot \frac{1}{r} \cdot \textrm{e}^{\textrm{i} \lambda t \cdot \frac{1}{r}}\\
\end{align*}
D.h., wenn der Betrag der Geschwindigkeit von $1$ auf $\lambda$ erhöht wird, wird der Betrag der Beschleunigung mit $\lambda^2$ multipliziert.
===== Beschleunigung bei konstantem Geschwindigkeisbetrag 1 =====
Sei $\ell(t)$ die Bogenlänge der Bezierkurve, also $\ell(t) = \int_0^t |p'(s)| \textrm{d}s$. Insbesondere gilt dann
$\ell'(t) = |p'(t)|$. (Die Länge ändert sich so stark, wie der Punkt in $t$ schnell ist).
Sei $t(\ell)$ die Umkehrfunktion, d.h. der $t$-Parameter für eine bestimmte Länge. Die Ableitung ist dann
(nach der Formel $(f^{-1})' = \frac{1}{f'(f^{-1})}$ mit $f^{-1}=t$ und $f=\ell$):
$$
\frac{\textrm{d}t(\ell)}{\textrm{d}\ell} = \frac{1}{\ell'(t(\ell))} = \frac{1}{|p'(t(\ell))|}.
$$
(D.h. der Parameter $t$ ändert sich mit der Länge indirekt proportional zur $t$-Geschwindigkeit).
Anstatt der "normalen" Parametrierung in $t \in [0,1]$, stellen wir uns eine Parametrierung $p_{\ell}$ vor, so dass die Bogenlänge der Bezierkurve von $p_{\ell}(0)$ bis $p_{\ell}(\ell)$ genau $\ell$ beträgt.
Es gilt $p_{\ell}(\ell) = p(t(\ell))$ und
\[
v_n(t(\ell)) := \frac{\textrm{d}p_{\ell}(\ell)}{\textrm{d} \ell} = p'(t(\ell)) \cdot \frac{\textrm{d}t(\ell)}{\textrm{d}\ell} =
p'(t(\ell)) \cdot \frac{1}{|p'(t(\ell))|}
\]
d.h. man erhält den auf den Betrag 1 normalisierten Geschwindgkeitsvektor (was ja genau der Sinn der $p_{\ell}$ Parametrierung ist).
Die zweite Ableitung von $p_{\ell}$ nach $\ell$ ergibt die Beschleunigung $a_n(\ell)$, die rechtwinklig auf $v_n$ steht (sonst würde sich der Betrag von $v_n$ ändern).
Diese zweite Ableitung berechnen wir nummerisch durch Ableiten der ersten:
\[
a_n(t(\ell)) := \frac{\mathrm{d}v_n(\ell)}{\mathrm{d}\ell} \approx \frac{v_n(t(\ell)+\Delta t)-v_n(t(\ell)-\Delta t)}{\ell(t+\Delta t)-\ell(t-\Delta t)} \approx \frac{v_n(t(\ell)+\Delta t)-v_n(t(\ell)-\Delta t)}{|p(t+\Delta t)-p(t-\Delta t)|}
\]
Algebraisch erhält man folgendes:
\[
\frac{\mathrm{d}}{\mathrm{d}\ell} \left(\frac{p'}{|p'|}\right) =
\frac{p'' - e \cdot \left( e \cdot p''\right)}{|p'|^2}
\]
mit $e=\frac{p'}{|p'|}$. Das ist bis auf den Faktor $|p'|^2$ das Gram-Schmidt Verfahren. Sachen gibts...
===== Effektive Beschleunigung und Komponente in Bahnnormalebene =====
Sei $v_{\text{eff}}(t) \in \mathbb{R}$ der Betrag der effektiven Bahngeschwindigkeit im Punkt zum entsprechenden $t$-Parameter.
Die effektive Beschleunigung ist also
\[
a_{\text{eff}}(t) = (v_{\text{eff}}(t))^2 \cdot a_n(t) -g
\]
wobei $g$ die Gravitationsbeschleunigung ist.
Diese Beschleunigung zerlegen wir in zwei Komponenten, eine tangential zur Bahn $a_t$, eine rechtwinklig dazu $a_r$. Die Komponente $a_t$ sorgt für die Beschleunigung des Zugs, die Komponente $a_r$ ist für die Bahnneigung relevant. Diese sollte nämlich rechtwinklig zu $a_r$ sein.
$a_t$ ist die Projektion von $a_{\text{eff}}$ auf $v_n$. Es gilt (mit $|v_n|=1$):
\[
a_t = (\vec a_{\text{eff}} \cdot \vec v_n) \cdot \vec v_n \qquad \text{Skalarprodukt in der Klammer!}
\]
Und damit
\[
a_r = a_{\text{eff}}-a_t
\]
Damit bilden wir ein Koordinatensystem $K(t)$ mit Ursprung $p(t)$ und Einheitsvektoren
\[
e_1 = v_n(t), \qquad e_2 = \frac{a_r(t)}{|a_r(t)|}, \qquad e_3 = v_n(t) \times a_r(t) \cdot \frac{1}{|a_r(t)|}
\]
===== Bau der Bahn =====
In geeigneter Schrittgrösse der Kurve entlang gehen, an jedem Punkt mit $K(t)$ Schienen-Punkte definieren.
Eventuell Stützen definieren.
Blender Code analog zur Generierung der glatten Kurve.
===== Simulation der Bewegung =====
Zustand: $t$ (Ort auf der Bahn), $v_{\text{eff}}$ (aktueller Betrag der Geschwindigkeit)
Schritt: Zeit um 1/framerate vorrücken, der Bahn folgen (z.B. um die Strecke, die mit $v_{\text{eff}}$ in dieser Zeit zurückgelegt würde, oder genauere schrittweise Simulation). Aus der Höhendifferenz und eventuell Reibung die neue Geschwindigkeit berechnen.
Kamera entsprechend positionieren und Keyframe setzen:
cam = bpy.data.objects['Camera']
frame = 0
cam.animation_data_clear()
cam.matrix_world = ( (y.x,y.y,y.z,1), (-an.x, -an.y, -an.z, 1), (-vv.x,-vv.y,-vv.z,1), (pp.x, pp.y, pp.z, 0))
cam.keyframe_insert(data_path="rotation_euler", frame=frame)
cam.keyframe_insert(data_path="location", frame=frame)
frame+=1
Siehe auch https://blender.stackexchange.com/questions/108938/how-to-interpret-the-camera-world-matrix
D.h. die erste Koordinatenrichtung ist rechts, die zweite oben und die dritte ist entgegen der Blickrichtung.
===== Blender =====
Bezier Klasse laden in Blender:
# Nimmt die Bezierkurven aus myspline und erzeugt
# die Bahn und die Kamera-Animation
# Import in Blender 2.8 (see https://devtalk.blender.org/t/2-80-using-multiple-internal-scripts-breaking-change/6980 )
Bezier = bpy.data.texts["bezier.py"].as_module().Bezier
obj = bpy.data.objects['mySpline']
# Kurvenpunkte auslesen
mypoints=[]
if obj.type == 'CURVE':
for subcurve in obj.data.splines:
curvetype = subcurve.type
if curvetype == 'BEZIER':
for bezpoint in subcurve.bezier_points:
mypoints.append(bezpoint.handle_left)
mypoints.append(bezpoint.co)
mypoints.append(bezpoint.handle_right)
# Sammlung von Bezierkurven erzeugen
mySplines = []
numpoints = len(mypoints)
totalLength = 0
for i in range(numpoints//3):
mySplines.append(Bezier((mypoints[i*3+1],
mypoints[i*3+2],
mypoints[(i*3+3)%numpoints],
mypoints[(i*3+4)%numpoints])))
totalLength+=mySplines[-1].length()
# Bahn erzeugen
try:
bpy.ops.collection.objects_remove(bpy.data.collections['Rails'])
except:
pass
railsCol = bpy.data.collections.new('Rails')
linksCol = bpy.data.collections.new('RailLinks')
railsCol.children.link(linksCol)
bpy.context.scene.collection.children.link(railsCol)
abstand = 0.2 # Bahnpunkte
ldone = 0 # Erledigte Bahnstrecke
i=0 # Aktuelle Bezierkurve
t = 0 # Aktuelle t-Parameter
g = Vector(0,0,-9.81) # Gravitationbeschleunigung
hmax = 40 # Hoehe fuer v=0
# Bahnpunkte: Ctrl-Links, Knoten, Ctrl-Rechts
railspts=[[],[],[]] # Bahnpunkte, Schiene L, Schiene R, Träger
while(ldone