Re: Interpolation cubique (courbe bézier...)
Publié : mar. 17/juil./2012 21:19
Forums PureBasic - Français
https://www.purebasic.fr/french/
C'est la même chose ^^un peu avec une courbe de bézier, ou comme avec la technique de G-rom.
oui, j'avais bien compris, mais ton explication est bien classe car vraiment simple à mettre en oeuvreG-Rom a écrit :C'est la même chose ^^
Code : Tout sélectionner
Structure vector2f
x.f
y.f
EndStructure
Procedure.i newVector2f(x.f = 0 ,y.f = 0)
*v.vector2f = AllocateMemory(SizeOf(vector2f))
If *v
*v\x = x
*v\y = y
ProcedureReturn *v
Else
ProcedureReturn #Null
EndIf
EndProcedure
Procedure freeVector2f(*v.vector2f)
If *v <> #Null
FreeMemory(*v)
EndIf
EndProcedure
Procedure.i midPoint(*a.vector2f, *b.vector2f)
*v.vector2f = AllocateMemory(SizeOf(vector2f))
If *v
*v\x = ( *a\x + *b\x ) * 0.5
*v\y = ( *a\y + *b\y ) * 0.5
ProcedureReturn *v
Else
ProcedureReturn #Null
EndIf
EndProcedure
; vive les macros...
Macro position(v)
v\x,v\y
EndMacro
*A.vector2f = newVector2f(100 , 100)
*B.vector2f = newVector2f(500 , 450)
*MidPoint.vector2f = midPoint(*A, *B)
If OpenWindow(0,0,0,640,640,"")
Repeat
event = WindowEvent()
; mouvement
*A\x = 100 - 80 * Cos( (ElapsedMilliseconds() / 10) * #PI /180)
*A\y = 100 - 50 * Sin( (ElapsedMilliseconds() / 10) * #PI /180)
*B\x = 500 + 100 * Cos( (ElapsedMilliseconds() / 10) * #PI /180)
*B\y = 450 + 150 * Sin( (ElapsedMilliseconds() / 10) * #PI /180)
freeVector2f(*MidPoint)
*MidPoint = midPoint(*A, *B)
If timer < ElapsedMilliseconds()
timer = ElapsedMilliseconds() + 10
StartDrawing(WindowOutput(0))
Box(0,0,640,640,$FFFFFF)
LineXY(position(*A),position(*B), $BBBBBB)
Circle(position(*A),2,0)
Circle(position(*B),2,0)
Circle(position(*MidPoint), 4, $FF)
StopDrawing()
EndIf
Delay(1)
Until event = #PB_Event_CloseWindow
freeVector2f(*A)
freeVector2f(*B)
freeVector2f(*MidPoint)
End
EndIf
Code : Tout sélectionner
;http://www.purebasic.fr/english/viewtopic.php?f=12&t=47084
;IncludeFile "Maths Matrix Calculation.pb" ; Guimauve
#IMAGE_Content =0
Enumeration
#GADGET_Canvas
#GADGET_Brush
#GADGET_Clear
EndEnumeration
Global Pas = 2,MouseX_Old,MouseY_Old
; Structure Point
; x.i
; y.i
; EndStructure
;Global pointxy.point
Global NewList listexy.POINT()
Global NewList listexy2.POINT()
Global poucentage.f
Global matrice_de_points.s
pourcentage=0.2
; Procedure lagrange(matrice_de_points.s) ;pemet de trouver l'equation d'une courbe à partir de ses points
; ;String_To_Matrix(PointXY.Matrix, "[1,2; 4,5; 3,3; 5,1; 7,2]")
; String_To_Matrix(PointXY.Matrix, matrice_de_points)
; LagrangePolynomialMatrix(GetMatrixLine(PointXY), Lagrange2D.Matrix, LagrangeInverse2D.Matrix)
; ProductMatrix(Equation2D.Matrix, LagrangeInverse2D, PointXY)
; DebugMatrix(Equation2D, 3)
; Debug ""
; DebugPolynomialMatrix(Equation2D, 3)
; EndProcedure
; Draw the mouse action result depending on the currently selected mode and event type
Macro point_distance(x1,y1,x2,y2)
Int(Sqr((x2-x1)*(x2-x1) + (y2-y1)*(y2-y1)) )
EndMacro
Macro point_direction(x1,y1,x2,y2)
ATan2((y2- y1),(x2- x1))
EndMacro
Procedure DrawBrush(x1,y1)
x2 = MouseX_Old
y2 = MouseY_Old
flech = point_distance(x1,y1,x2,y2)
d.d = point_direction(x1,y1,x2,y2)
deg.d=Degree(d)
;Debug deg
AddElement(listexy()) ; on stock les points
listexy()\x = x1
listexy()\y = y1
For i = 1 To flech
i+pas
x_result = x1 + Sin(d) *i
y_result = y1 + Cos(d) *i
Circle(x_result,y_result,3,0)
;DrawImage(ImageID(#IMAGE_Brush),x_result-brush\size/2,y_result-brush\size/2)
Next i
EndProcedure
Procedure DrawAction(x, y, EventType)
If StartDrawing(CanvasOutput(#GADGET_Canvas))
If EventType = #PB_EventType_LeftButtonDown Or EventType = #PB_EventType_MouseMove
DrawBrush(x,y)
EndIf
StopDrawing()
EndIf
EndProcedure
CreateImage(#IMAGE_Content, 380, 380, 24)
If OpenWindow(0, 0, 0, 800, 600, "Interpolation", #PB_Window_SystemMenu | #PB_Window_ScreenCentered)
CanvasGadget(#GADGET_Canvas, 10, 10, WindowWidth(0)-80, WindowHeight(0)-20, #PB_Canvas_ClipMouse)
ButtonGadget(#GADGET_Clear, WindowWidth(0)-60, 100, 50, 25, "Clear")
SetGadgetAttribute(#GADGET_Canvas, #PB_Canvas_Cursor, #PB_Cursor_Cross)
Repeat
Event = WaitWindowEvent()
If Event = #PB_Event_Gadget
Select EventGadget()
Case #GADGET_Canvas
X = GetGadgetAttribute(#GADGET_Canvas, #PB_Canvas_MouseX)
Y = GetGadgetAttribute(#GADGET_Canvas, #PB_Canvas_MouseY)
Type = EventType()
Select EventType()
Case #PB_EventType_LeftButtonDown
If StartDrawing(ImageOutput(#IMAGE_Content))
DrawImage(GetGadgetAttribute(#GADGET_Canvas, #PB_Canvas_Image), 0, 0)
StopDrawing()
EndIf
MouseY_Old = y
MouseX_Old = x
DrawAction(X, Y, EventType())
Case #PB_EventType_LeftButtonUp
DrawAction(X, Y, EventType())
;Pour chaque point de la courbe
FirstElement(Listexy())
;point0
x0=listexy()\x
y0=listexy()\y
AddElement(listexy2())
;matrice_de_points="["+Str(x0)+","+Str(y0)+";" ;"[1,2; 4,5; 3,3; 5,1; 7,2]")
;point1
NextElement(Listexy())
pente0.f=(listexy()\y - y0)/(listexy()\x - x0) ; calcul de la pente
x0=listexy()\x
y0=listexy()\y
;matrice_de_points+Str(x0)+","+Str(y0)+";"
;point2 et suivant
For n=2 To ListSize(Listexy())
NextElement(Listexy())
pente.f=(listexy()\y - y0)/(listexy()\x - x0)
If Abs(pente)>Abs(pente0)+poucentage*pente0 ; si la pente est sup de 20%
;Il faut changer cet élément
listexy()\x=x0
listexy()\y=Int(pente*(listexy()\x-x0)+y0) ; bug ici
EndIf
x0=listexy()\x
y0=listexy()\y
;matrice_de_points+Str(x0)+","+Str(y0)+";"
pente0=pente
Next n
;matrice_de_points+"]"
;Debug matrice_de_points
FirstElement(Listexy())
StartDrawing(CanvasOutput(#GADGET_Canvas))
For i=1 To ListSize(Listexy())
;Debug listexy()\x
;Debug listexy()\y
Circle(listexy()\x,listexy()\y,5,$ff0000)
NextElement(Listexy())
Next i
StopDrawing()
Case #PB_EventType_MouseMove
If GetGadgetAttribute(#GADGET_Canvas, #PB_Canvas_Buttons) & #PB_Canvas_LeftButton
DrawAction(X, Y, EventType())
MouseY_Old = y
MouseX_Old = x
EndIf
EndSelect
Case #GADGET_Clear
If StartDrawing(CanvasOutput(#GADGET_Canvas))
Box(0, 0, GadgetWidth(#GADGET_Canvas), GadgetHeight(#GADGET_Canvas), $FFFFFF)
StopDrawing()
EndIf
EndSelect
EndIf
Until Event = #PB_Event_CloseWindow
EndIf
Code : Tout sélectionner
Global pourentage.f
pourcentage=0.2
If Abs(pente)>Abs(pente0)+pourcentage*pente0 ; si la pente est sup de 20%
Code : Tout sélectionner
; This source file is part of OGRE
; (Object-oriented Graphics Rendering Engine)
; For the latest info, see http://www.ogre3d.org/
;
; Copyright (c) 2000-2009 Torus Knot Software Ltd
InitSprite()
InitKeyboard()
OpenScreen(1024, 768, 32, "spline")
Structure Vector3
x.f
y.f
z.f
EndStructure
Macro CopyVector3(a, b)
a\x = b\x
a\y = b\y
a\z = b\z
EndMacro
Global Dim mCoeffs.f(3,3)
Global mAutoCalc = #True
Global NewList lPoints.Vector3()
Global Dim mPoints.Vector3(0)
Global Dim mTangents.Vector3(0)
Declare recalcTangents()
Declare interpolate2(*pos.Vector3, fromIndex, t.f)
; Set up matrix
; Hermite polynomial
mCoeffs(0,0) = 2
mCoeffs(0,1) = -2
mCoeffs(0,2) = 1
mCoeffs(0,3) = 1
mCoeffs(1,0) = -3
mCoeffs(1,1) = 3
mCoeffs(1,2) = -2
mCoeffs(1,3) = -1
mCoeffs(2,0) = 0
mCoeffs(2,1) = 0
mCoeffs(2,2) = 1
mCoeffs(2,3) = 0
mCoeffs(3,0) = 1
mCoeffs(3,1) = 0
mCoeffs(3,2) = 0
mCoeffs(3,3) = 0
Procedure concatenate(Array r.f(2), Array m.f(2), Array m2.f(2))
r(0,0) = m(0,0) * m2(0,0) + m(0,1) * m2(1,0) + m(0,2) * m2(2,0) + m(0,3) * m2(3,0);
r(0,1) = m(0,0) * m2(0,1) + m(0,1) * m2(1,1) + m(0,2) * m2(2,1) + m(0,3) * m2(3,1);
r(0,2) = m(0,0) * m2(0,2) + m(0,1) * m2(1,2) + m(0,2) * m2(2,2) + m(0,3) * m2(3,2);
r(0,3) = m(0,0) * m2(0,3) + m(0,1) * m2(1,3) + m(0,2) * m2(2,3) + m(0,3) * m2(3,3);
r(1,0) = m(1,0) * m2(0,0) + m(1,1) * m2(1,0) + m(1,2) * m2(2,0) + m(1,3) * m2(3,0);
r(1,1) = m(1,0) * m2(0,1) + m(1,1) * m2(1,1) + m(1,2) * m2(2,1) + m(1,3) * m2(3,1);
r(1,2) = m(1,0) * m2(0,2) + m(1,1) * m2(1,2) + m(1,2) * m2(2,2) + m(1,3) * m2(3,2);
r(1,3) = m(1,0) * m2(0,3) + m(1,1) * m2(1,3) + m(1,2) * m2(2,3) + m(1,3) * m2(3,3);
r(2,0) = m(2,0) * m2(0,0) + m(2,1) * m2(1,0) + m(2,2) * m2(2,0) + m(2,3) * m2(3,0);
r(2,1) = m(2,0) * m2(0,1) + m(2,1) * m2(1,1) + m(2,2) * m2(2,1) + m(2,3) * m2(3,1);
r(2,2) = m(2,0) * m2(0,2) + m(2,1) * m2(1,2) + m(2,2) * m2(2,2) + m(2,3) * m2(3,2);
r(2,3) = m(2,0) * m2(0,3) + m(2,1) * m2(1,3) + m(2,2) * m2(2,3) + m(2,3) * m2(3,3);
r(3,0) = m(3,0) * m2(0,0) + m(3,1) * m2(1,0) + m(3,2) * m2(2,0) + m(3,3) * m2(3,0);
r(3,1) = m(3,0) * m2(0,1) + m(3,1) * m2(1,1) + m(3,2) * m2(2,1) + m(3,3) * m2(3,1);
r(3,2) = m(3,0) * m2(0,2) + m(3,1) * m2(1,2) + m(3,2) * m2(2,2) + m(3,3) * m2(3,2);
r(3,3) = m(3,0) * m2(0,3) + m(3,1) * m2(1,3) + m(3,2) * m2(2,3) + m(3,3) * m2(3,3);
EndProcedure
Procedure Multiplication(Array r.f(1), Array v.f(1), Array m.f(2))
r(0) = m(0,0) * v(0) + m(1,0) * v(1) + m(2,0) * v(2) + m(3,0) * v(3)
r(1) = m(0,1) * v(0) + m(1,1) * v(1) + m(2,1) * v(2) + m(3,1) * v(3)
r(2) = m(0,2) * v(0) + m(1,2) * v(1) + m(2,2) * v(2) + m(3,2) * v(3)
r(3) = m(0,3) * v(0) + m(1,3) * v(1) + m(2,3) * v(2) + m(3,3) * v(3)
EndProcedure
Procedure addPoint(*p.Vector3)
AddElement(lPoints())
lPoints()\x = *p\x
lPoints()\y = *p\y
lPoints()\z = *p\z
Dim mPoints(ListSize(lPoints()))
i = 0
ForEach lPoints()
mPoints(i) = lPoints()
i + 1
Next
If (mAutoCalc)
recalcTangents()
EndIf
EndProcedure
Procedure interpolate(*pos.Vector3, t.f)
;Currently assumes points are evenly spaced, will cause velocity
;change where this is Not the Case
;TODO: base on arclength?
; Work out which segment this is in
fSeg.f = t * (ListSize(lPoints()) - 1)
segIdx = Int(fSeg)
; Apportion t
t = fSeg - segIdx
interpolate2(*pos, segIdx, t)
EndProcedure
Procedure interpolate2(*pos.Vector3, fromIndex, t.f)
;// Bounds check
;assert (fromIndex < mPoints.size() && "fromIndex out of bounds")
If (fromIndex + 1) = ListSize(lPoints())
;// Duff request, cannot blend To nothing
;// Just Return source
CopyVector3(*pos, mPoints(fromIndex))
ProcedureReturn
EndIf
;// Fast special cases
If t = 0.0
CopyVector3(*pos, mPoints(fromIndex))
ProcedureReturn
ElseIf t = 1.0
CopyVector3(*pos, mPoints(fromIndex + 1))
ProcedureReturn
EndIf
;// Real interpolation
;// Form a vector of powers of t
Define.f t2, t3
t2 = t * t
t3 = t2 * t
Dim powers.f(3)
powers(0) = t3
powers(1) = t2
powers(2) = t
powers(3) = 1
;// Algorithm is ret = powers * mCoeffs * Matrix4(point1, point2, tangent1, tangent2)
Define.Vector3 point1, point2, tan1, tan2
CopyVector3(point1, mPoints(fromIndex))
CopyVector3(point2, mPoints(fromIndex+1))
CopyVector3(tan1 , mTangents(fromIndex))
CopyVector3(tan2 , mTangents(fromIndex+1))
Dim pt.f(3,3)
pt(0,0) = point1\x
pt(0,1) = point1\y
pt(0,2) = point1\z
pt(0,3) = 1.0
pt(1,0) = point2\x
pt(1,1) = point2\y
pt(1,2) = point2\z
pt(1,3) = 1.0
pt(2,0) = tan1\x
pt(2,1) = tan1\y
pt(2,2) = tan1\z
pt(2,3) = 1.0
pt(3,0) = tan2\x
pt(3,1) = tan2\y
pt(3,2) = tan2\z
pt(3,3) = 1.0
Dim ret.f(3)
Dim r.f(3,3)
;ret = powers * mCoeffs * pt;
concatenate(r(), mCoeffs(), pt())
Multiplication(ret(), powers(), r())
*pos\x = ret(0)
*pos\y = ret(1)
*pos\z = ret(2)
;ProcedureReturn Vector3(ret.x, ret.y, ret.z);
EndProcedure
Procedure recalcTangents()
; // Catmull-Rom approach
; //
; // tangent[i] = 0.5 * (point[i+1] - point[i-1])
; //
; // Assume endpoint tangents are parallel With line With neighbour
Define.i i, numPoints, isClosed
numPoints = ListSize(lPoints())
If numPoints < 2
; Can't do anything yet
ProcedureReturn
EndIf
; Closed Or open?
If mPoints(0) = mPoints(numPoints-1)
isClosed = #True
Else
isClosed = #False
EndIf
Dim mTangents.Vector3(numPoints)
For i = 0 To numPoints-1
If i =0
; Special Case start
If isClosed
; Use numPoints-2 since numPoints-1 is the last point And == [0]
mTangents(i)\x = 0.5 * (mPoints(1)\x - mPoints(numPoints-2)\x)
mTangents(i)\y = 0.5 * (mPoints(1)\y - mPoints(numPoints-2)\y)
mTangents(i)\z = 0.5 * (mPoints(1)\z - mPoints(numPoints-2)\z)
Else
mTangents(i)\x = 0.5 * (mPoints(1)\x - mPoints(0)\x)
mTangents(i)\y = 0.5 * (mPoints(1)\y - mPoints(0)\y)
mTangents(i)\z = 0.5 * (mPoints(1)\z - mPoints(0)\z)
EndIf
ElseIf i = numPoints-1
; Special Case End
If isClosed
; Use same tangent As already calculated For [0]
mTangents(i)\x = mTangents(0)\x
mTangents(i)\y = mTangents(0)\y
mTangents(i)\z = mTangents(0)\z
Else
mTangents(i)\x = 0.5 * (mPoints(i)\x - mPoints(i-1)\x)
mTangents(i)\y = 0.5 * (mPoints(i)\y - mPoints(i-1)\y)
mTangents(i)\z = 0.5 * (mPoints(i)\z - mPoints(i-1)\z)
EndIf
Else
mTangents(i)\x = 0.5 * (mPoints(i+1)\x - mPoints(i-1)\x)
mTangents(i)\y = 0.5 * (mPoints(i+1)\y - mPoints(i-1)\y)
mTangents(i)\z = 0.5 * (mPoints(i+1)\z - mPoints(i-1)\z)
EndIf
Next
EndProcedure
;- Ajoute quelques points et ensuite on trace la spline
p.Vector3
pos.Vector3
p\x = 100 : p\y = 0 : p\z = 100
addPoint(@p)
p\x = 300 : p\y = 0 : p\z = 450
addPoint(@p)
p\x = 570 : p\y = 0 : p\z = 360
addPoint(@p)
p\x = 730 : p\y = 0 : p\z = 570
addPoint(@p)
p\x = 330 : p\y = 0 : p\z = 590
addPoint(@p)
p\x = 30 : p\y = 0 : p\z = 90
addPoint(@p)
xd.f = 100
yd.f = 100
pas.f = 1
time.f = 0
Repeat
ClearScreen(0)
StartDrawing(ScreenOutput())
For t=0 To 100
interpolate(@pos, t/100.0)
If t > 0
LineXY(xd, yd, pos\x, pos\z, $AABBCC)
EndIf
xd = pos\x
yd = pos\z
Next t
interpolate(@pos, time)
Circle(pos\x, pos\z, 10, $FFAAAA)
time + pas * 0.005
If time > 1
Time = 1
pas = - pas
ElseIf time < 0
Time = 0
pas = - pas
EndIf
For i=0 To ListSize(lPoints())-1
Circle(mPoints(i)\x, mPoints(i)\z, 2, $FF)
Next
StopDrawing()
ExamineKeyboard()
FlipBuffers()
Until KeyboardPushed(#PB_Key_Escape)
Code : Tout sélectionner
; This source file is part of OGRE
; (Object-oriented Graphics Rendering Engine)
; For the latest info, see http://www.ogre3d.org/
;
; Copyright (c) 2000-2009 Torus Knot Software Ltd
InitSprite()
InitKeyboard()
InitMouse()
ExamineDesktops()
Scx = DesktopWidth(0)
Scy = DesktopHeight(0)
OpenScreen(Scx, Scy, DesktopDepth(0), "spline")
Structure Vector3
x.f
y.f
z.f
EndStructure
Macro CopyVector3(a, b)
a\x = b\x
a\y = b\y
a\z = b\z
EndMacro
Global Dim mCoeffs.f(3,3)
Global mAutoCalc = #True
Global NewList lPoints.Vector3()
Global Dim mPoints.Vector3(0)
Global Dim mTangents.Vector3(0)
Declare recalcTangents()
Declare interpolate2(*pos.Vector3, fromIndex, t.f)
; Set up matrix
; Hermite polynomial
mCoeffs(0,0) = 2
mCoeffs(0,1) = -2
mCoeffs(0,2) = 1
mCoeffs(0,3) = 1
mCoeffs(1,0) = -3
mCoeffs(1,1) = 3
mCoeffs(1,2) = -2
mCoeffs(1,3) = -1
mCoeffs(2,0) = 0
mCoeffs(2,1) = 0
mCoeffs(2,2) = 1
mCoeffs(2,3) = 0
mCoeffs(3,0) = 1
mCoeffs(3,1) = 0
mCoeffs(3,2) = 0
mCoeffs(3,3) = 0
Procedure concatenate(Array r.f(2), Array m.f(2), Array m2.f(2))
r(0,0) = m(0,0) * m2(0,0) + m(0,1) * m2(1,0) + m(0,2) * m2(2,0) + m(0,3) * m2(3,0);
r(0,1) = m(0,0) * m2(0,1) + m(0,1) * m2(1,1) + m(0,2) * m2(2,1) + m(0,3) * m2(3,1);
r(0,2) = m(0,0) * m2(0,2) + m(0,1) * m2(1,2) + m(0,2) * m2(2,2) + m(0,3) * m2(3,2);
r(0,3) = m(0,0) * m2(0,3) + m(0,1) * m2(1,3) + m(0,2) * m2(2,3) + m(0,3) * m2(3,3);
r(1,0) = m(1,0) * m2(0,0) + m(1,1) * m2(1,0) + m(1,2) * m2(2,0) + m(1,3) * m2(3,0);
r(1,1) = m(1,0) * m2(0,1) + m(1,1) * m2(1,1) + m(1,2) * m2(2,1) + m(1,3) * m2(3,1);
r(1,2) = m(1,0) * m2(0,2) + m(1,1) * m2(1,2) + m(1,2) * m2(2,2) + m(1,3) * m2(3,2);
r(1,3) = m(1,0) * m2(0,3) + m(1,1) * m2(1,3) + m(1,2) * m2(2,3) + m(1,3) * m2(3,3);
r(2,0) = m(2,0) * m2(0,0) + m(2,1) * m2(1,0) + m(2,2) * m2(2,0) + m(2,3) * m2(3,0);
r(2,1) = m(2,0) * m2(0,1) + m(2,1) * m2(1,1) + m(2,2) * m2(2,1) + m(2,3) * m2(3,1);
r(2,2) = m(2,0) * m2(0,2) + m(2,1) * m2(1,2) + m(2,2) * m2(2,2) + m(2,3) * m2(3,2);
r(2,3) = m(2,0) * m2(0,3) + m(2,1) * m2(1,3) + m(2,2) * m2(2,3) + m(2,3) * m2(3,3);
r(3,0) = m(3,0) * m2(0,0) + m(3,1) * m2(1,0) + m(3,2) * m2(2,0) + m(3,3) * m2(3,0);
r(3,1) = m(3,0) * m2(0,1) + m(3,1) * m2(1,1) + m(3,2) * m2(2,1) + m(3,3) * m2(3,1);
r(3,2) = m(3,0) * m2(0,2) + m(3,1) * m2(1,2) + m(3,2) * m2(2,2) + m(3,3) * m2(3,2);
r(3,3) = m(3,0) * m2(0,3) + m(3,1) * m2(1,3) + m(3,2) * m2(2,3) + m(3,3) * m2(3,3);
EndProcedure
Procedure Multiplication(Array r.f(1), Array v.f(1), Array m.f(2))
r(0) = m(0,0) * v(0) + m(1,0) * v(1) + m(2,0) * v(2) + m(3,0) * v(3)
r(1) = m(0,1) * v(0) + m(1,1) * v(1) + m(2,1) * v(2) + m(3,1) * v(3)
r(2) = m(0,2) * v(0) + m(1,2) * v(1) + m(2,2) * v(2) + m(3,2) * v(3)
r(3) = m(0,3) * v(0) + m(1,3) * v(1) + m(2,3) * v(2) + m(3,3) * v(3)
EndProcedure
Procedure addPoint(*p.Vector3)
AddElement(lPoints())
lPoints()\x = *p\x
lPoints()\y = *p\y
lPoints()\z = *p\z
Dim mPoints(ListSize(lPoints()))
i = 0
ForEach lPoints()
mPoints(i) = lPoints()
i + 1
Next
If (mAutoCalc)
recalcTangents()
EndIf
EndProcedure
Procedure interpolate(*pos.Vector3, t.f)
;Currently assumes points are evenly spaced, will cause velocity
;change where this is Not the Case
;TODO: base on arclength?
; Work out which segment this is in
fSeg.f = t * (ListSize(lPoints()) - 1)
segIdx = Int(fSeg)
; Apportion t
t = fSeg - segIdx
interpolate2(*pos, segIdx, t)
EndProcedure
Procedure interpolate2(*pos.Vector3, fromIndex, t.f)
;// Bounds check
;assert (fromIndex < mPoints.size() && "fromIndex out of bounds")
If ListSize(lPoints()) = 0 Or fromIndex >= ListSize(lPoints())
ProcedureReturn
EndIf
If (fromIndex + 1) = ListSize(lPoints())
;// Duff request, cannot blend To nothing
;// Just Return source
CopyVector3(*pos, mPoints(fromIndex))
ProcedureReturn
EndIf
;// Fast special cases
If t = 0.0
CopyVector3(*pos, mPoints(fromIndex))
ProcedureReturn
ElseIf t = 1.0
CopyVector3(*pos, mPoints(fromIndex + 1))
ProcedureReturn
EndIf
;// Real interpolation
;// Form a vector of powers of t
Define.f t2, t3
t2 = t * t
t3 = t2 * t
Dim powers.f(3)
powers(0) = t3
powers(1) = t2
powers(2) = t
powers(3) = 1
;// Algorithm is ret = powers * mCoeffs * Matrix4(point1, point2, tangent1, tangent2)
Define.Vector3 point1, point2, tan1, tan2
CopyVector3(point1, mPoints(fromIndex))
CopyVector3(point2, mPoints(fromIndex+1))
CopyVector3(tan1 , mTangents(fromIndex))
CopyVector3(tan2 , mTangents(fromIndex+1))
Dim pt.f(3,3)
pt(0,0) = point1\x
pt(0,1) = point1\y
pt(0,2) = point1\z
pt(0,3) = 1.0
pt(1,0) = point2\x
pt(1,1) = point2\y
pt(1,2) = point2\z
pt(1,3) = 1.0
pt(2,0) = tan1\x
pt(2,1) = tan1\y
pt(2,2) = tan1\z
pt(2,3) = 1.0
pt(3,0) = tan2\x
pt(3,1) = tan2\y
pt(3,2) = tan2\z
pt(3,3) = 1.0
Dim ret.f(3)
Dim r.f(3,3)
;ret = powers * mCoeffs * pt;
concatenate(r(), mCoeffs(), pt())
Multiplication(ret(), powers(), r())
*pos\x = ret(0)
*pos\y = ret(1)
*pos\z = ret(2)
;ProcedureReturn Vector3(ret.x, ret.y, ret.z);
EndProcedure
Procedure recalcTangents()
; // Catmull-Rom approach
; //
; // tangent[i] = 0.5 * (point[i+1] - point[i-1])
; //
; // Assume endpoint tangents are parallel With line With neighbour
Define.i i, numPoints, isClosed
numPoints = ListSize(lPoints())
If numPoints < 2
; Can't do anything yet
ProcedureReturn
EndIf
; Closed Or open?
If mPoints(0) = mPoints(numPoints-1)
isClosed = #True
Else
isClosed = #False
EndIf
Dim mTangents.Vector3(numPoints)
For i = 0 To numPoints-1
If i =0
; Special Case start
If isClosed
; Use numPoints-2 since numPoints-1 is the last point And == [0]
mTangents(i)\x = 0.5 * (mPoints(1)\x - mPoints(numPoints-2)\x)
mTangents(i)\y = 0.5 * (mPoints(1)\y - mPoints(numPoints-2)\y)
mTangents(i)\z = 0.5 * (mPoints(1)\z - mPoints(numPoints-2)\z)
Else
mTangents(i)\x = 0.5 * (mPoints(1)\x - mPoints(0)\x)
mTangents(i)\y = 0.5 * (mPoints(1)\y - mPoints(0)\y)
mTangents(i)\z = 0.5 * (mPoints(1)\z - mPoints(0)\z)
EndIf
ElseIf i = numPoints-1
; Special Case End
If isClosed
; Use same tangent As already calculated For [0]
mTangents(i)\x = mTangents(0)\x
mTangents(i)\y = mTangents(0)\y
mTangents(i)\z = mTangents(0)\z
Else
mTangents(i)\x = 0.5 * (mPoints(i)\x - mPoints(i-1)\x)
mTangents(i)\y = 0.5 * (mPoints(i)\y - mPoints(i-1)\y)
mTangents(i)\z = 0.5 * (mPoints(i)\z - mPoints(i-1)\z)
EndIf
Else
mTangents(i)\x = 0.5 * (mPoints(i+1)\x - mPoints(i-1)\x)
mTangents(i)\y = 0.5 * (mPoints(i+1)\y - mPoints(i-1)\y)
mTangents(i)\z = 0.5 * (mPoints(i+1)\z - mPoints(i-1)\z)
EndIf
Next
EndProcedure
;- Ajoute quelques points et ensuite on trace la spline
p.Vector3
pos.Vector3
xd.f = 0
yd.f = 0
pas.f = 1
time.f = 0
Repeat
ClearScreen(0)
ExamineMouse()
If MouseButton(1)
If flag = 0
Flag = 1
p\x = MouseX()
p\z = MouseY()
addPoint(@p)
If xd = 0 And yd = 0
xd = p\x
yd = p\z
EndIf
EndIf
Else
flag = 0
EndIf
StartDrawing(ScreenOutput())
For t=0 To 300
interpolate(@pos, t/300.0)
If t > 0
LineXY(xd, yd, pos\x, pos\z, $AABBCC)
EndIf
xd = pos\x
yd = pos\z
Next t
interpolate(@pos, time)
Circle(pos\x, pos\z, 10, $FFAAAA)
time + pas * 0.02 / ListSize(lPoints())
If time > 1
Time = 1
pas = - pas
ElseIf time < 0
Time = 0
pas = - pas
EndIf
For i=0 To ListSize(lPoints())-1
Circle(mPoints(i)\x, mPoints(i)\z, 2, $FF)
Next
LineXY(0, MouseY(), Scx, MouseY(), $AAAAAA)
LineXY(MouseX(), 0, MouseX(), Scy, $AAAAAA)
StopDrawing()
ExamineKeyboard()
FlipBuffers()
Until KeyboardPushed(#PB_Key_Escape)