Excel. Долгая дорога оцифровки. Часть 6. Кусочная интерполяция⁠⁠

Ну теперь пора перейти именно к интерполяции исходных данных. Итак, я напомню - у нас был лист с распечатанным графиком, мы его отсканировали, получили набор точек ХY и... имеем вот такую (в лучшем случае) картину (см.первый скрин). Т.е. по данному набору точек невозможно сделать корректную апроксимацию полиномом.

Выходом из данной ситуации является разбиение данных на несколько частей, в данном случае 2 с общей точкой Do=428, создание кусучнозаданной функции (ЕСЛИ меньше 428 одна функция, если больше - вторая функция). Но так в данном случае, а если надо сделать два, три...и больше разбиений? Кропотливая работа. Но зачем, если можно заставить Excel в автоматическом режиме выбирать малое количество точек и проводить через них интерполяционную функцию.

Отчасти кусочную интерполяцию показывал в прошлом посте серии ( Excel. Долгая дорога оцифровки. Часть 5. Создание пользовательской функции для двух аргументов. Ручной вариант ) при поиске решения между заданных критериев.

Как простые варианты рассмотрим кусочную интерполяцию по двум и по четырём точкам для заданных ниже данных.

Допустим нужно определить значение Y при X = 2.5.

При кусочной интерполяции по двум точкам используются две ближайшие заданные точки к X = 2.5. , т.е. 2 и 3, через данные точки провидится линия,и по ней находится Y при X = 2.5..

При кусочной интерполяции по четырём точкам используются две ближайшие заданные точки к X = 2.5. справа и две две ближайшие заданные точки к X = 2.5. слева, т.е. 1 и 2 и 3 и 4, через данные точки провидится кривая (полином 3-й степени),и по ней находится Y при X = 2.5.

Это справедливо для данных за 2-й и перед предпоследней известной точкой. Для данных отрезков и для экстраполяции использую линейную интерполяцию (по 2-м точкам).

Как видно использование кусочной интерполяции по 4-м точкам (голубая линия)немного сглаживает итоговую функцию, что позволяет снимать при оцифровке чуть меньше точек :) .

Вообще правило такое, в зависимости от вида графика:

Ну и собственно с помощью чего сие выполняется:

Макрос кусочной интерполяции при использовании данных с листа

======

' Интерполяция по 2-м, 3-м или 4-м ближайшим (до и после) к заданной (Xisk) точке

' В подпрограмму передаются все заданные точки

' Интерполяция происходит косочно, по количеству заданных точек с учётом расположения заданного Х

' Данные из листа Excel.

Public Function kus_interp_Ex(Xt As Range, Yt As Range, Xisk As Single, Optional ByVal toch As Integer = 2) As Variant

Dim i As Long

Dim xd() As Double

Dim yd() As Double

Dim cd() As Double

' toch - указание поиска решения с использованием количества точек (2, 3, 4).

Select Case toch

Case 2 ' Уравнение а·х+b

kus_interp_Ex = linterp(Xt.Rows(Xt.Count - 1), Xt.Rows(Xt.Count), Yt.Rows(Xt.Count - 1), Yt.Rows(Xt.Count), Xisk)

For i = 1 To Xt.Count - 1

If Xisk < Xt.Rows(i + 1) Then

kus_interp_Ex = linterp(Xt.Rows(i), Xt.Rows(i + 1), Yt.Rows(i), Yt.Rows(i + 1), Xisk)

Exit For

End If

Next i

Case 3 ' Уравнение а·х^2+b·x+c Интерполяция по принципу х1 Х х2 х3

kus_interp_Ex = kubterp(Xt.Rows(Xt.Count - 2), Xt.Rows(Xt.Count - 1), Xt.Rows(Xt.Count), _

Yt.Rows(Xt.Count - 2), Yt.Rows(Xt.Count - 1), Yt.Rows(Xt.Count), Xisk)

For i = 1 To Xt.Count - 2

If Xisk < Xt.Rows(i + 1) Then

kus_interp_Ex = kubterp(Xt.Rows(i), Xt.Rows(i + 1), Xt.Rows(i + 2), _

Yt.Rows(i), Yt.Rows(i + 1), Yt.Rows(i + 2), Xisk)

Exit For

End If

Next i

Case 4 ' Уравнение а·х^3+b·x^2+c·x+d Интерполяция по принципу х1 х2 X х3 x4

ReDim xd(1 To 4) As Double

ReDim yd(1 To 4) As Double

If Xisk < Xt.Rows(2) Then ' Экстраполяция ДО и интерполяция ДО второй известной точки - линейна

kus_interp_Ex = linterp(Xt.Rows(1), Xt.Rows(2), Yt.Rows(1), Yt.Rows(2), Xisk)

Else

If Xisk >= Xt.Rows(Xt.Count - 1) Then ' Экстраполяция ЗА и интерполяция ПОСЛЕ второй известной точки - линейна

kus_interp_Ex = linterp(Xt.Rows(Xt.Count - 1), Xt.Rows(Xt.Count), Yt.Rows(Xt.Count - 1), Yt.Rows(Xt.Count), Xisk)

Else ' Между ними считаю по интерполяции полиномом с расположением заданного икса между двух пар точек

For i = 3 To Xt.Count - 1

If Xisk < Xt.Rows(i) Then

xd(1) = Xt.Rows(i - 2): xd(2) = Xt.Rows(i - 1): xd(3) = Xt.Rows(i): xd(4) = Xt.Rows(i + 1)

yd(1) = Yt.Rows(i - 2): yd(2) = Yt.Rows(i - 1): yd(3) = Yt.Rows(i): yd(4) = Yt.Rows(i + 1)

Linia_trenda yd, xd, 3, cd

kus_interp_Ex = cd(1) * Xisk ^ 3 + cd(2) * Xisk ^ 2 + cd(3) * Xisk ^ 1 + cd(4)

Exit For

End If

Next i

End If

End If

Case Else

End Select

End Function

Входные данные:

Xt - Столбец исходных Х

Yt  - Столбец исходных Y

Xisk  - X при котором требуется определить Y

toch - количество используемых точек при интерполяции.

Наблюдательный заметит, что в макросе присутствует и интерполяция с использованием 3-х точек. (0_о)

Дополнительная функция, требуемая для макроса кусочной интерполяции - определение коэффициентов полинома линии тренда:

======

' Проведение интерполяции с использованием функционала Excel

' На выходе - коэффициенты полинома. Число точек должно быть минимум на одну больше, чем степень полинома.

' Данные берутся из программы

Public Sub Linia_trenda(ByRef Y() As Double, ByRef x() As Double, ByVal PolyStep As Integer, ByRef c() As Double, Optional ByRef r2 As Double)

Dim stepen As Long

' Ввожу проверку не превышения степени массива

If (UBound(Y) - LBound(Y) - 1) < PolyStep Then

stepen = UBound(Y) - LBound(Y)

Else

stepen = PolyStep

End If

' Объявляю переменные, создаю матрицы под размер данных и степень полинома.

Dim X1() As Double, Y1() As Double

ReDim X1(LBound(Y) To UBound(Y), 1 To stepen) As Double

ReDim Y1(LBound(Y) To UBound(Y), 1 To 1) As Double

ReDim c(1 To stepen + 1) As Double

' Заполню массив Х в соответствии со степенью уравнения.

For i = LBound(x) To UBound(x)

Y1(i, 1) = Y(i)

X1(i, 1) = x(i)

For N = 2 To stepen

X1(i, N) = X1(i, 1) ^ N

Next N

Next i

' Нахожу уравнение.

Dim Coefs As Variant

Coefs = WorksheetFunction.LinEst(Y1, X1, True, True)

' Вытаскиваю коэффициенты полинома.

For i = 1 To stepen + 1

c(i) = Coefs(1, i)

Next i

' Вытаскиваю величину достоверности апроксимации.

r2 = Coefs(3, 1)

End Sub

Макрос linterp был представлен в прошлый раз.

Одним из замечательных применений кусочной интерполяции является возможность автоматического создания макросов функций без проблем с невозможностью достоверной апроксимации исходных данных вот в таком виде:

' Поправки Сербия Панчево Страница 39 из 77

Public Function ТЭХ_ПТ80_Рис3(ByRef Go As Single) As Variant

Dim Xt As Variant

Dim Yt As Variant

Xt = Array(13.0042194092827, 13.4767932489451, 14.0675105485232)

Yt = Array(-6.38888888888889E-02, -6.38888888888889E-02, -0.06875)

ТЭХ_ПТ80_Рис3 = kus_interp(Xt, Yt, Go, 4, 2)

End Function

Учтите что kus_interp при использовании данных с листа и из макроса отличаются...

Но об этом в следующий раз.

===========================

Планы на будущее

1. Часть 7. Автоматическое создание макроса функции с использованием кусочной интерполяции.

2. Построим поиск решения.

3. Строим график функции.