Prosper Vector Class

' ########################################
' #         Vector Class Module          #
' ########################################
' Auxiliary library to work with Petroleum Experts (TM) OpenServer.
' Part of the SimpleOpenServer package, developed by Patricio Panichelli.
' Some rights reserved (c) www.SimpleOpenServer.com/license
' For documentation, visit www.SimpleOpenServer.com/docs
' SOS Version:              1.0.0
' Developed & tested with:  IPM 9.0
'
Option Explicit
Option Base 0
Private xValues
Enum vDirection
    printDown
    printRight
End Enum
Enum sDirection
    leftToRight
    rightToLeft
End Enum
Const NotFound = -1E+30

' ++++++++++++++++++++++++++++
' Private sub Class_Initialize
' ++++++++++++++++++++++++++++
' Initializes the class with the defaults
Private Sub Class_Initialize()
    ReDim xValues(0 To 1, -1 To -1)
End Sub

' ++++++++++++++++++++
' Public Sub transpose
' ++++++++++++++++++++
' Exchanges the X and Y coordinates of the vector
Public Sub transpose()
    Dim newVec As Variant, i As Integer
    ReDim newVec(0 To 1, 0 To UBound(xValues, 2))
    For i = 0 To UBound(xValues, 2)
        newVec(0, i) = xValues(1, i)
        newVec(1, i) = xValues(0, i)
    Next
    xValues = newVec
End Sub

' +++++++++++++++++
' Public Function Y
' +++++++++++++++++
' Returns the Y value for a given X from a vector.
' Linear interpolation
' Warning: will extrapolate if X is not included in the domain
Public Function Y(ByVal X As Double) As Double
    If Me.points < 2 Then err.Raise 2020, "Y", "Trying to guess Y from X but vector has only " & CStr(Me.points) & " point(s)."
    Dim i As Integer
    For i = 0 To UBound(xValues, 2)
        If xValues(0, i) > X Then Exit For
    Next
    If i = 0 Then i = i + 1
    If i > UBound(xValues, 2) Then i = UBound(xValues, 2)
    If (xValues(0, i) - xValues(0, i - 1)) = 0 Then
        ' Try to avoid division by zero
        Y = xValues(1, i - 1)
    Else
        Y = xValues(1, i - 1) + (X - xValues(0, i - 1)) / (xValues(0, i) - xValues(0, i - 1)) * (xValues(1, i) - xValues(1, i - 1))
    End If
End Function

' ++++++++++++++++++++
' Public Function dYdX
' ++++++++++++++++++++
' Returns the derivative value for a given X from a vector.
' Uses quadratic fit
' CAUTION: will extrapolate if X is not included in the domain
Public Function dYdX(ByVal X As Double) As Double
    Dim numPoints As Integer
    numPoints = Me.points
    If numPoints < 2 Then err.Raise 2020, "dYdX", "Trying to guess dYdX from X but vector has only " & CStr(numPoints) & " point(s)."
    If numPoints = 2 Then
        ' 2 points - Straight line
        dYdX = (xValues(1, 1) - xValues(1, 0)) / (xValues(0, 1) - xValues(0, 0))
    Else
        ' 3 or more points
        Dim i As Long, coeffs
        For i = 0 To numPoints - 1
            If xValues(0, i) > X Then Exit For
        Next
        If i = 0 Then i = i + 1
        If i > numPoints - 1 Then i = i - 1
        If i = numPoints - 1 Then i = i - 1
        coeffs = pFitQuadAtIndex(i)
        dYdX = 2 * coeffs(0) * X + coeffs(1) ' 2 a X + b
    End If
End Function

' ++++++++++++++++++++++
' Public Function XfromY
' ++++++++++++++++++++++
' Returns the X value that corresponds to a given Y from a vector.
' Returns NotFound if the value was not found
Function XfromY(Y As Double, Optional searchDirection As Variant) As Double
    If Me.points < 2 Then err.Raise 2020, "XfromY", "Trying to guess X from Y but vector has only " & CStr(Me.points) & " point(s)."
    If IsMissing(searchDirection) Then searchDirection = sDirection.rightToLeft
    Dim i As Integer
    If searchDirection = sDirection.rightToLeft Then
        For i = UBound(xValues, 2) To 1 Step -1
            If Sgn(xValues(1, i) - Y) <> Sgn(xValues(1, i - 1) - Y) Then Exit For
        Next
        If i = 0 Then
            XfromY = NotFound
            Exit Function
        End If
        XfromY = xValues(0, i - 1) + (Y - xValues(1, i - 1)) / (xValues(1, i) - xValues(1, i - 1)) * (xValues(0, i) - xValues(0, i - 1))
    Else
        For i = 0 To UBound(xValues, 2) - 1 Step 1
            If Sgn(xValues(1, i) - Y) <> Sgn(xValues(1, i + 1) - Y) Then Exit For
        Next
        If i = UBound(xValues, 2) Then
            XfromY = NotFound
            Exit Function
        End If
        XfromY = xValues(0, i) + (xValues(1, i) - Y) / (xValues(1, i) - xValues(1, i + 1)) * (xValues(0, i + 1) - xValues(0, i))
    End If
End Function

' +++++++++++++++++++
' Public Property min
' +++++++++++++++++++
' Returns the Minimum value of a vector
' By default is the Y value, but pass an "X" as
' parameter and you will get the X corresponding
' to the maximum.
Public Property Get min(Optional coordinate)
    Dim i As Integer, mV As Double, returnX As Boolean
    If UCase(CStr(coordinate)) = "X" Then returnX = True
    mV = 1E+30
    For i = 0 To UBound(xValues, 2)
        If CDbl(xValues(1, i)) < mV Then
            mV = xValues(1, i)
            If returnX Then
                min = xValues(0, i)
            Else
                min = xValues(1, i)
            End If
        End If
    Next
End Property

' +++++++++++++++++++
' Public Property max
' +++++++++++++++++++
' Returns the Maximum value of a vector
' By default is the Y value, but pass an "X" as
' parameter and you will get the X corresponding
' to the maximum.
Public Property Get max(Optional coordinate)
    Dim i As Integer, mV As Double, returnX As Boolean
    If UCase(CStr(coordinate)) = "X" Then returnX = True
    mV = -1E+30
    For i = 0 To UBound(xValues, 2)
        If CDbl(xValues(1, i)) > mV Then
            mV = xValues(1, i)
            If returnX Then
                max = xValues(0, i)
            Else
                max = xValues(1, i)
            End If
        End If
    Next
End Property

' ++++++++++++++++++++
' Public Property last
' ++++++++++++++++++++
' Returns the last value of a vector
' By default is the Y value, but pass an "X" as
' parameter and you will get the X.
Public Property Get last(Optional coordinate)
    If UCase(CStr(coordinate)) = "X" Then
        last = xValues(0, UBound(xValues, 2))
    Else
        last = xValues(1, UBound(xValues, 2))
    End If
End Property

' +++++++++++++++++++++
' Public Property first
' +++++++++++++++++++++
' Returns the first value of a vector.
' By default is the Y value, but pass an "X" as
' parameter and you will get the X.
Public Property Get first(Optional coordinate)
    If UCase(CStr(coordinate)) = "X" Then
        first = xValues(0, 0)
    Else
        first = xValues(1, 0)
    End If
End Property

' ++++++++++++++++++++++
' Public Property values
' ++++++++++++++++++++++
' Returns the values as a variant variable
Public Property Get values() As Variant
    values = xValues
End Property

' ++++++++++++++++++++++++
' Public Sub newFromValues
' ++++++++++++++++++++++++
' Sets the vector equal to a variant object.
' Useful to initialize with values
Public Sub newFromValues(ByVal values As Variant)
    xValues = values
    pReverseIfNecessary
End Sub

' +++++++++++++++++++++++
' Public Sub newFromRange
' +++++++++++++++++++++++
' Sets the vector equal to a range.
' Must be a 2-column vertical range, with X values in the left and Y values in the right
Public Sub newFromRange(ByVal rng)
    If TypeName(rng) = "String" Then rng = range(rng)
    Dim newVec, i As Integer, numRows As Integer
    If TypeName(rng) = "Range" Then
        numRows = rng.rows.Count - 1
    Else
        numRows = UBound(rng, 1) - 1
    End If
    ReDim newVec(0 To 1, 0 To numRows)
    For i = 0 To UBound(newVec, 2)
        If rng(i + 1, 1) = "" Or rng(i + 1, 2) = "" Then Exit For
        newVec(0, i) = CDbl(rng(i + 1, 1))
        newVec(1, i) = CDbl(rng(i + 1, 2))
    Next
    ReDim Preserve newVec(0 To 1, 0 To i - 1)
    xValues = newVec
    pReverseIfNecessary
End Sub

' +++++++++++++++++++
' Public Sub addPoint
' +++++++++++++++++++
' Adds the given point (x, y) to the vector
' Useful to create point-based vectors
Public Sub addPoint(X As Double, Y As Double)
    Dim numRows As Integer
    numRows = UBound(xValues, 2) + 1
    ReDim Preserve xValues(0 To 1, 0 To numRows)
    xValues(0, numRows) = X
    xValues(1, numRows) = Y
End Sub

' +++++++++++++++++++
' Public Sub multiply
' +++++++++++++++++++
' Multiplies two vectors' Y parameters (v1 + v2)
' Or will multiply a vector by a number
Public Sub multiply(ByVal vectorOrValue As Variant)
    Dim typN As String
    typN = TypeName(vectorOrValue)
    If typN = "Integer" Or typN = "Double" Or typN = "Long" Then
        pMultByNumber CDbl(vectorOrValue)
    ElseIf typN = "Vector" Then
        pMultByVector vectorOrValue.values
    Else
        err.Raise 2020, "add", "Trying to add vector or value, but no value or vector was passed"
    End If
End Sub

' +++++++++++++++++++++++++
' Private Sub pMultByNumber
' +++++++++++++++++++++++++
' Multiplies the Y-coordinate by  a double
Private Sub pMultByNumber(value As Double)
    Dim i As Integer
    For i = 0 To UBound(xValues, 2)
        xValues(1, i) = xValues(1, i) * value
    Next
End Sub

' +++++++++++++++++++++++++
' Private Sub pMultByVector
' +++++++++++++++++++++++++
' Multiplies the Y-coordinates of two vectors
' Both vectors need to be of the same length and the same x-coords
Private Sub pMultByVector(v2 As Variant)
    If UBound(xValues, 2) = UBound(v2, 2) Then
        Dim i As Integer
        For i = 0 To UBound(xValues, 2)
            If xValues(0, i) <> v2(0, i) Then err.Raise 2020, "multiply", "When multiplying two vectors, the X-Coordinates were not compatible."
            xValues(1, i) = xValues(1, i) * v2(1, i)
        Next
    Else
        err.Raise 2020, "add", "Tried to multiply two vectors, but sizes are not compatible."
    End If
End Sub

' ++++++++++++++
' Public Sub add
' ++++++++++++++
' Adds two vectors' Y parameters (v1 + v2)
' Or will add a vector to a number
Public Sub add(ByVal vectorOrValue As Variant)
    Dim typN As String
    typN = TypeName(vectorOrValue)
    If typN = "Integer" Or typN = "Double" Or typN = "Long" Then
        pAddNumber CDbl(vectorOrValue)
    ElseIf typN = "Vector" Then
        pAddVector vectorOrValue.values
    Else
        err.Raise 2020, "add", "Trying to add vector or value, but no value or vector was passed"
    End If
End Sub

' ++++++++++++++++++++
' Public Sub subtract
' ++++++++++++++++++++
' Subtracts two vectors' Y parameters (v1 - v2)
' Or will subtract a number from a vector
Public Sub subtract(ByVal vectorOrValue As Variant)
    Dim typN As String
    typN = TypeName(vectorOrValue)
    If typN = "Integer" Or typN = "Double" Or typN = "Long" Then
        pAddNumber -CDbl(vectorOrValue)
    ElseIf typN = "Vector" Then
        vectorOrValue.multiply -1
        pAddVector vectorOrValue.values
        vectorOrValue.multiply -1
    Else
        err.Raise 2020, "subtract", "Trying to subtract vector or value, but no value or vector was passed"
    End If
End Sub

' ++++++++++++++++++++++
' Private Sub pAddNumber
' ++++++++++++++++++++++
' Adds a double to the Y-coordinate
Private Sub pAddNumber(value As Double)
    Dim i As Integer
    For i = 0 To UBound(xValues, 2)
        xValues(1, i) = xValues(1, i) + value
    Next
End Sub

' ++++++++++++++++++++++
' Private Sub pAddVector
' ++++++++++++++++++++++
' Adds a Y-coordinates of two vectors
' Both vectors need to be of the same length and the same x-coords
Private Sub pAddVector(v2 As Variant)
    If UBound(xValues, 2) = UBound(v2, 2) Then
        Dim i As Integer
        For i = 0 To UBound(xValues, 2)
            If xValues(0, i) <> v2(0, i) Then err.Raise 2020, "add", "When adding two vectors, the X-Coordinates were not compatible."
            xValues(1, i) = xValues(1, i) + v2(1, i)
        Next
    Else
        err.Raise 2020, "add", "Tried to add two vectors, but sizes are not compatible."
    End If
End Sub

' ++++++++++++++++++++++
' Public Sub printInCell
' ++++++++++++++++++++++
' prints the values of a vector in the worksheet, starting in the range specified
' Optional values > direction: [printDown] | printRight
'                 > printX:    [True]      | False
Public Sub printInCell(inRange As range, Optional direction As Variant, Optional printX As Variant)
    If IsMissing(printX) Then printX = True
    printX = CBool(printX)
    If IsMissing(direction) Then direction = vDirection.printDown

    If printX Then
        If direction = vDirection.printRight Then
            inRange.Resize(2, UBound(xValues, 2) + 1) = xValues
        Else
            inRange.Resize(UBound(xValues, 2) + 1, 2) = Application.transpose(xValues)
        End If
    Else
        If direction = vDirection.printRight Then
            inRange.Resize(1, UBound(xValues, 2) + 1) = pGetYValues()
        Else
            inRange.Resize(UBound(xValues, 2) + 1, 1) = Application.transpose(pGetYValues())
        End If
    End If
End Sub

' +++++++++++++++++++++++++++++++
' Private Sub pReverseIfNecessary
' +++++++++++++++++++++++++++++++
' If x Values are in descending order, it reverses the vector.
' Will only check for first two values to see if it is reversed
Private Sub pReverseIfNecessary()
    Dim nPoints As Long
    nPoints = Me.points
    If nPoints <= 1 Then Exit Sub
    If xValues(0, 0) > xValues(0, 1) Then
        Dim newValues, i As Long
        ReDim newValues(0 To 1, 0 To nPoints - 1)
        For i = 0 To nPoints - 1
            newValues(0, i) = xValues(0, nPoints - 1 - i)
            newValues(1, i) = xValues(1, nPoints - 1 - i)
        Next
        xValues = newValues
    End If
End Sub

' ++++++++++++++++++++++
' Public Property points
' ++++++++++++++++++++++
' Returns the amount of points in a vector.
Public Property Get points() As Integer
    points = UBound(xValues, 2) + 1
End Property

' +++++++++++++++++++
' Public Function MSE
' +++++++++++++++++++
' Returns the mean square error between two vectors
' calculated as SUM( dif^2 ) / nPoints
' Used in the match function
Public Function MSE(ByVal points As Vector) As Double
    Dim i As Integer, sqD As Double, PointValues, numPoints As Integer
    numPoints = points.points
    PointValues = points.values
    For i = 0 To numPoints - 1
        sqD = sqD + (Me.Y(CDbl(PointValues(0, i))) - PointValues(1, i)) ^ 2
    Next
    MSE = sqD / numPoints
End Function

' ++++++++++++++++++++++++++++++++
' Private function pFitQuadAtIndex
' ++++++++++++++++++++++++++++++++
' fits a quadratic around point at index i
' using a quadratic function ax2 + bx + c
' returns [ a, b, c ]
Private Function pFitQuadAtIndex(i As Long)
    pFitQuadAtIndex = pFitQuad(CDbl(xValues(0, i - 1)), CDbl(xValues(1, i - 1)), _
                               CDbl(xValues(0, i)), CDbl(xValues(1, i)), _
                               CDbl(xValues(0, i + 1)), CDbl(xValues(1, i + 1)))
End Function

' +++++++++++++++++++++++++
' Private function pFitQuad
' +++++++++++++++++++++++++
' given three points 1,2,3: fits them using a quadratic function ax2 + bx + c
' returns [ a, b, c ]
Private Function pFitQuad(x1 As Double, y1 As Double, x2 As Double, y2 As Double, x3 As Double, y3 As Double) As Variant
    Dim result(2), X(1 To 3, 1 To 3) As Double, Y(1 To 3, 1 To 1) As Double, c
    X(1, 1) = 1
    X(1, 2) = x1
    X(1, 3) = x1 * x1
    X(2, 1) = 1
    X(2, 2) = x2
    X(2, 3) = x2 * x2
    X(3, 1) = 1
    X(3, 2) = x3
    X(3, 3) = x3 * x3
    Y(1, 1) = y1
    Y(2, 1) = y2
    Y(3, 1) = y3
    c = Application.MMult(Application.MInverse(X()), Y())
    result(0) = c(3, 1)
    result(1) = c(2, 1)
    result(2) = c(1, 1)
    pFitQuad = result
End Function


' ++++++++++++++++++++++++++++
' Private Function pGetYValues
' ++++++++++++++++++++++++++++
' Returns a 1-dimension (real vector) of the Y values of our vector
Private Function pGetYValues() As Variant
    Dim returnedArr As Variant, i As Integer
    ReDim returnedArr(0 To UBound(xValues, 2))
    For i = 0 To UBound(xValues, 2)
        returnedArr(i) = xValues(1, i)
    Next
    pGetYValues = returnedArr
End Function