Interpolation cubique (courbe bézier...)

[comtois]

C'est ça que tu cherches à faire ?

http://dumas.ccsd.cnrs.fr/docs/00/63/68/24/PDF/Yol.pdf

[blendman]

C'est très intéressant, mais ce n'est pas tout à fait ça. Là, c'est assez complexe à mettre en place.
Je souhaite quelque chose de plus simple.
La technique de G-rom me semble très intéressante pour ça .

En fait, ça ressemble à l'explication qu'a donné G-rom.
L'idée est la suivante :
- pour le moment, lorsque je dessine avec mon outil, si je dessine trop vite, cela créé des lignes droites, même si mon dessin était avec lignes plus arrondis
- je voudrais donc arrondir un peu ces lignes droites, un peu avec une courbe de bézier, ou comme avec la technique de G-rom.

[G-Rom]

un peu avec une courbe de bézier, ou comme avec la technique de G-rom.

C'est la même chose ^^




avec des points de contrôle , c'est plus parlant ( en bleu )

[blendman]

C'est la même chose ^^

oui, j'avais bien compris, mais ton explication est bien classe car vraiment simple à mettre en oeuvre

[G-Rom]

Voila une petite fonction intéressante qui calcule le point du milieu entre deux vecteurs :

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 

[blendman]

ah oui, c'est intéressant ça, merci !

[falsam]

Tout simplement géniale cette fonction. Merci G-Rom. Je pense que je vais l'adapter dans mon projet MindMap.

[Mesa]

J'ai essayé un truc mais ça ne marche pas comme je l'espérais mais ça peut peut-être te donner une idée.

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

Mesa.

[blendman]

salut

une petit boulette ici :

Code : Tout sélectionner

Global pourentage.f

pourcentage=0.2

If Abs(pente)>Abs(pente0)+pourcentage*pente0 ; si la pente est sup de 20%

Mais ça ne change pas grand chose en corrigeant.
C'est sensé faire quoi normalement ?

[comtois]

les Splines c'est pas mal aussi pour faire des tracés , ou pour adoucir un chemin calculé dans un pathfinding, dans l'exemple qui suit j'ai fait un mixte (tracé d'une courbe et sur le tracé je déplace un cercle)

C'est une adaptation rapide du code SimpleSpline d'ogre; donc ça demande encore pas mal de boulot pour le mettre en forme, le simplifier, etc. Mais comme je n'aurai pas forcément le temps de le faire voici une démo de ce que ça donne :

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)

[comtois]

Une petite modif pour définir les points avec la souris.

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)

[Mesa]

J'ai une erreur, index de tableau hors limite ligne 159.

Mesa.

[comtois]

C'est exact, le problème apparait avec le débogueur activé. J'avais fait la correction chez moi et pas ici ! C'est corrigé, tu peux tester à nouveau