myfunction.vb

vb.net编辑插值与拟合软件

Imports System.Math
Public Class MyFunction
    Function ZhiXian_ChaZhi(ByVal x() As Single, ByVal y() As Single, ByVal myx As Single) As Single '直线差值。
        Dim ChaZhi As Single = 0

        For i As Integer = 0 To x.GetUpperBound(0) - 1

            If myx > x(i) And myx < x(1 + i) Then
                ChaZhi = y(i) + (y(i + 1) - y(i)) * (myx - x(i)) / (x(i + 1) - x(i))
                Return ChaZhi
            ElseIf myx = x(i) Then
                ChaZhi = y(i)
                Return ChaZhi
            ElseIf myx = x(i + 1) Then

                ChaZhi = y(i + 1)
                Return ChaZhi
                Exit For

            End If
        Next
    End Function


    Public Function spline(ByVal x() As Single, ByVal y() As Single, ByVal myx As Single, ByVal JingDu As Integer) As Double

        Dim ChaZhi_spline As Double '存储三次样条差值的结果
        Dim u1(100000) As Single, v1(100000) As Single
        Dim num As Single

        Dim n As Integer = x.GetUpperBound(0)
        Dim ch As Integer = 2
        Dim tx1 As Single = 0
        Dim tx2 As Single = 0
        Dim ty1 As Single = 0
        Dim ty2 As Single = 0

        Dim a(10000) As Single, b(10000) As Single, c(10000) As Single, dx(10000) As Single, dy(10000) As Single
        Dim qx(10000) As Single, qy(10000) As Single
        Dim tt As Single, bx3 As Single, bx4 As Single, by3 As Single, by4 As Single
        Dim cx As Single, cy As Single, t(10000) As Single, px(10000) As Single, py(10000) As Single
        Dim u(3) As Single, v(3) As Single, i As Single
        num = 0
        For i = 1 To n
            t(i) = (x(i) - x(i - 1)) ^ 2 + (y(i) - y(i - 1)) ^ 2
            t(i) = Sqrt(t(i))
        Next i
        Select Case ch
            Case 0     '抛物条件 有错误啊
                u(0) = (x(1) - x(0)) / t(1) : u(1) = (x(2) - x(1)) / t(2)
                u(2) = (u(1) - u(0)) / (t(2) + t(1))
                tx1 = u(0) - u(2) * t(1)
                u(0) = (y(1) - y(0)) / t(1) : u(1) = (y(2) - y(1)) / t(2)
                u(2) = (u(1) - u(0)) / (t(2) + t(1))
                ty1 = u(0) - u(2) * t(1)
                u(0) = (x(n) - x(n - 1)) / t(n) : u(1) = (x(n - 1) - x(n - 2)) / t(n - 1)
                u(2) = (u(0) - u(1)) / (t(n) + t(n - 1))
                tx2 = u(0) + u(2) * t(n)
                u(0) = (y(n) - y(n - 1)) / t(n) : u(1) = (y(n - 1) - y(n - 2)) / t(n - 1)
                u(2) = (u(0) - u(1)) / (t(n) + t(n - 1))
                ty2 = u(0) + u(2) * t(n)

                a(0) = 1 : c(0) = 0
                a(n) = 1 : b(n) = 0
            Case 1     '夹持条件
                a(0) = 1 : c(0) = 0 : dx(0) = tx1 : dy(0) = ty1
                a(n) = 1 : b(n) = 0 : dx(n) = tx2 : dy(n) = ty2


            Case 2     '自由条件
                a(0) = 2 : c(0) = 1
                dx(0) = 3 * (x(1) - x(0)) / t(1) : dy(0) = 3 * (y(1) - y(0)) / t(1)
                a(n) = 2 : b(n) = 1
                dx(n) = 3 * (x(n) - x(n - 1)) / t(n) : dy(n) = 3 * (y(n) - y(n - 1)) / t(n)


            Case 3     '循环条件
                a(0) = 2 : c(0) = 1
                dx(0) = 3 * (x(1) - x(0)) / t(1) - (t(1) * (x(2) - x(1)) / t(2) - x(1) + x(0)) / (t(1) + t(2))
                dy(0) = 3 * (y(1) - y(0)) / t(1) - (t(1) * (y(2) - y(1)) / t(2) - y(1) + y(0)) / (t(1) + t(2))
                a(n) = 2 : b(n) = 1
                dx(n) = 3 * (x(n) - x(n - 1)) / t(n)
                dx(n) = dx(n) + (x(n) - x(n - 1) - t(n) * (x(n - 1) - x(n - 2)) / t(n - 1)) / (t(n) + t(n - 1))
                dy(n) = 3 * (y(n) - y(n - 1)) / t(n)
                dy(n) = dy(n) + (y(n) - y(n - 1) - t(n) * (y(n - 1) - y(n - 2)) / t(n - 1)) / (t(n) + t(n - 1))


        End Select

        '计算方程组系数阵和常数阵
        For i = 1 To n - 1
            a(i) = 2 * (t(i) + t(i + 1)) : b(i) = t(i + 1) : c(i) = t(i)
            dx(i) = 3 * (t(i) * (x(i + 1) - x(i)) / t(i + 1) + t(i + 1) * (x(i) - x(i - 1)) / t(i))
            dy(i) = 3 * (t(i) * (y(i + 1) - y(i)) / t(i + 1) + t(i + 1) * (y(i) - y(i - 1)) / t(i))
        Next i

        '采用追赶法解方程组
        If ch = 0 Then
            c(0) = 0
        Else
            c(0) = c(0) / a(0)
        End If

        For i = 1 To n - 1
            a(i) = a(i) - b(i) * c(i - 1) : c(i) = c(i) / a(i)
        Next i
        a(n) = a(n) - b(n) * c(i - 1)

        If ch = 0 Then
            qx(0) = 0 : qy(0) = 0
        Else
            qx(0) = dx(0) / a(0) : qy(0) = dy(0) / a(0)
        End If


        For i = 1 To n
            qx(i) = (dx(i) - b(i) * qx(i - 1)) / a(i)
            qy(i) = (dy(i) - b(i) * qy(i - 1)) / a(i)
        Next i
        px(n) = qx(n) : py(n) = qy(n)
        For i = n - 1 To 0 Step -1
            px(i) = qx(i) - c(i) * px(i + 1)
            py(i) = qy(i) - c(i) * py(i + 1)
        Next i
        '计算曲线上点的坐标
        For i = 0 To n - 1
            bx3 = (3 * (x(i + 1) - x(i)) / t(i + 1) - 2 * px(i) - px(i + 1)) / t(i + 1)
            bx4 = ((2 * (x(i) - x(i + 1)) / t(i + 1) + px(i) + px(i + 1)) / t(i + 1)) / t(i + 1)
            by3 = (3 * (y(i + 1) - y(i)) / t(i + 1) - 2 * py(i) - py(i + 1)) / t(i + 1)
            by4 = ((2 * (y(i) - y(i + 1)) / t(i + 1) + py(i) + py(i + 1)) / t(i + 1)) / t(i + 1)
            tt = 0
            While (tt <= t(i + 1))
                cx = x(i) + (px(i) + (bx3 + bx4 * tt) * tt) * tt
                cy = y(i) + (py(i) + (by3 + by4 * tt) * tt) * tt
                u1(num) = cx : v1(num) = cy : num = num + 1 : tt = tt + (t(n) - t(0)) / (JingDu * n) '插入些节点用来控制曲线。
            End While
            u1(num) = x(i + 1) : v1(num) = y(i + 1) : num = num + 1
        Next i



        u1(0) = x(0) : v1(0) = y(0)

        '利用直线差值计算出我们要求的Y值。
        For i = 0 To num - 1
            If myx = u1(i) And i = num - 1 Then

                ChaZhi_spline = v1(i)
                Exit For

            ElseIf myx >= u1(i) And myx < u1(1 + i) Then

                ChaZhi_spline = v1(i) + (v1(i + 1) - v1(i)) * (myx - u1(i)) / (u1(i + 1) - u1(i))
                Exit For

            End If
        Next
        Return ChaZhi_spline

    End Function
End Class