'---------------------------------------------------------------------------------------
'
' P U R P O S E
'
' The function lin_xyzt uses linear interpolation (or extrapolation) to find the value
' z(x=x_value,y=y_value,t=t_value) from multiple tables (x, y, z) of different t values.
' Each table is arranged so that the 1st column and the 1st row contain x and y values,
' respectively.  The z values are at the intersection of its x and y values.
' The program assumes x and y in each table are in ascending order.
'
' I N P U T
'
' Table     = a range of cells in the worksheet that defines the multiple tables.
' row_col   = a range of cells in the following format (one line for each table in Table):
'             row_1, row_2, col_1, col_2, t
'             where
'             row_1  = starting row number of table, relative to the first cell in Table.
'             row_2  = ending row number of table.
'             col_1  = starting column number of table, relatve to the first cell in Table.
'             col_2  = ending column number of table
'             t      = table values, arranged in ascending order.
' x_value   = the x value where the z value is to be determined.
' y_value   = the y value where the z value is to be determined.
' t_value   = the t value where the z value is to be determined.
' extrapolate (optional)
'           = if omitted, will return #N/A if x_value, y_value or t_value is out-of-bound.
'           = "bd", will return table bounded value if x_value, y_value or t_value is out-of-bound.
'           = "et", will linearly extrapolate if x_value, y_value or t_value is out-of-bound.
'
' O U T P U T
'
' lin_xyzt  = the z value at (x=x_value, y=y_value, t=t_value).
'             For multiple x_values or y_values, the first pair will be returned.
'             If x_value, y_value or t_value is out-of-bounds, the argument extrapolate
'             will determine the returned value.
'
'---------------------------------------------------------------------------------------
Function lin_xyzt(Table As Object, row_col As Object, x_value, y_value, t_value, Optional extrapolate)

Dim i, j As Integer              ' counter.
Dim i1, i2 As Integer
Dim z1, z2 As Single
Dim active_all(1) As Boolean
Dim extr

Dim row_1, row_2, col_1, col_2 As Integer

Dim T As Object
Dim Table_1 As Range
Dim Table_2 As Range

If IsMissing(extrapolate) Then
   extr = ""
Else
   extr = LCase(extrapolate)
End If

Select Case row_col.Rows.count
Case 1                           ' single table.  t_value is ignored.
   row_1 = row_col(1, 1)
   row_2 = row_col(1, 2)
   col_1 = row_col(1, 3)
   col_2 = row_col(1, 4)
      
   Set Table_1 = Range(Table.Cells(row_1, col_1), Table.Cells(row_2, col_2))
   
   lin_xyzt = lin_xyz(Table_1, x_value, y_value, extrapolate)
Case Else
   Set T = Range(row_col.Cells(1, 5), row_col.Cells(row_col.Rows.count, 5))
      
   active_all(1) = True
   Call bracket(T, t_value, i1, i2, active_all())
   
   row_1 = row_col(i1, 1)
   row_2 = row_col(i1, 2)
   col_1 = row_col(i1, 3)
   col_2 = row_col(i1, 4)
   Set Table_1 = Range(Table.Cells(row_1, col_1), Table.Cells(row_2, col_2))
         
   row_1 = row_col(i2, 1)
   row_2 = row_col(i2, 2)
   col_1 = row_col(i2, 3)
   col_2 = row_col(i2, 4)
   Set Table_2 = Range(Table.Cells(row_1, col_1), Table.Cells(row_2, col_2))
      
   z1 = lin_xyz(Table_1, x_value, y_value, extrapolate)
   z2 = lin_xyz(Table_2, x_value, y_value, extrapolate)
   
   If t_value >= T(i1) And t_value <= T(i2) Then
      lin_xyzt = lin_2xy(T(i1), z1, T(i2), z2, t_value)
   Else
      Select Case extr
      Case "bd"
         If t_value < T(i1) Then
            lin_xyzt = z1
         ElseIf t_value > T(i2) Then
            lin_xyzt = z2
         Else
            lin_xyzt = [#N/A]
         End If
      Case "et"
         lin_xyzt = lin_2xy(T(i1), z1, T(i2), z2, t_value)
      Case Else
         lin_xyzt = [#N/A]
      End Select
   End If
   
End Select

End Function

'---------------------------------------------------------------------------------------
'
' P U R P O S E
'
' The function lin_2xy uses linear interpolation to find the value y(x=x_value) from
' a pair (x1,y1), (x2,y2).
'
' I N P U T
'
' x1        = 1st x value.
' y1        = 1st y value.
' x2        = 2nd x value.
' y2        = 2nd y value.
' x_value   = the x value where the y value is to be determined.
'
' O U T P U T
'
' lin_2xy   = the y value at x_value.
'
'---------------------------------------------------------------------------------------
Function lin_2xy(x1, y1, x2, y2, x_value)

If x1 = x2 Then
   lin_2xy = [#N/A]
Else
   lin_2xy = y1 + (y2 - y1) / (x2 - x1) * (x_value - x1)
End If

End Function

'---------------------------------------------------------------------------------------
'
' P U R P O S E
'
' The function bracket find x1=x(i1) and x2=x(i2) from a table x
' so that x1 <= x_value <= x2.
'
' I N P U T
'
' x         = a range of cells representing x.
' x_value   = the x value to be bracketed.
' active    = an array of x.count elements.
'           = 0 if the corresponding x is to be ignored.
'
' O U T P U T
'
' i1        = the x index which will enclose x_value.
' i2        = the x index which will enclose x_value.
'           x(i1) <= x_value <= x(i2)
'           If x_value is out-of-bound, the two closest x values to x_value
'           will be returned.
'           If there are identical values of x(i1) or x(i2) in x, the first one
'           encountered will be used.
'
'---------------------------------------------------------------------------------------
Sub bracket(x As Object, x_value, i1, i2, active() As Boolean)

Dim i As Integer           ' counter.
Dim x1, xi, x2 As Single   ' = x(i1), x(i), x(i2), respectively.
Dim nactive As Integer
Dim usable As Boolean

nactive = UBound(active)

' If dx of first pair of x times number of points matches the last x,
' the program assumes x is equally spaced.

If x.count >= 2 And nactive <= 1 Then
   If Abs((x(x.count) - x(1)) / (x.count - 1) - (x(2) - x(1))) <= 0.01 * Abs(x(2) - x(1)) Then
      If x(x.count) >= x(1) Then
         i1 = Int((x_value - x(1)) / (x(2) - x(1))) + 1
         If x_value < x(i1) Then          ' Minor adjustment in case of round-off error.
            i1 = i1 - 1
         End If
         If i1 < 1 Then
            i1 = 1
         ElseIf i1 > x.count - 1 Then
            i1 = x.count - 1
         End If
         i2 = i1 + 1
      Else
         i1 = x.count - Int((x_value - x(x.count)) / (x(1) - x(2)))
         If x_value < x(i1) Then          ' Minor adjustment in case of round-off error.
            i1 = i1 + 1
         End If
         If i1 < 2 Then
            i1 = 2
         ElseIf i1 > x.count Then
            i1 = x.count
         End If
         i2 = i1 - 1
      End If
      If x_value >= x(i1) And x_value <= x(i2) Then
         Exit Sub
      ElseIf x_value < x(i1) And (i1 = 1 Or i1 = x.count) Then
         Exit Sub
      ElseIf x_value > x(i2) And (i2 = 1 Or i2 = x.count) Then
         Exit Sub
      End If
   End If
End If

' If the above fails, a more thorough search is performed.

i1 = 0
i2 = 0
For i = 1 To x.count
   xi = x(i)
   If nactive <= 1 Then
      usable = True
   Else
      usable = active(i)
   End If
   
   If Not usable Or Not IsNumeric(xi) Or VarType(xi) = vbEmpty Then
      
      ' Do nothing.
      
   ElseIf i1 = i2 Then
      If i1 = 0 Then
         i1 = i
         i2 = i
         x1 = xi
         x2 = xi
      ElseIf xi > x2 Then
         i2 = i
         x2 = xi
      ElseIf xi < x1 Then
         i1 = i
         x1 = xi
      End If
   ElseIf x1 <= x_value And x2 >= x_value Then
      If xi <= x1 Or xi >= x2 Then
         ' Do nothing.
      ElseIf xi < x_value Then
         i1 = i
         x1 = xi
      Else
         i2 = i
         x2 = xi
      End If
   ElseIf x1 >= x_value Then
      If xi < x1 Then
         i2 = i1
         x2 = x1
         i1 = i
         x1 = xi
      ElseIf xi <> x1 And xi < x2 Then
         i2 = i
         x2 = xi
      End If
   ElseIf x2 <= x_value Then
      If xi > x2 Then
         i1 = i2
         x1 = x2
         i2 = i
         x2 = xi
      ElseIf xi <> x2 And xi > x1 Then
         i1 = i
         x1 = xi
      End If
   End If
Next i

End Sub