6.16. Криві та інтегрування

Математичні сторінки охопили операції, що приймають масив і повертають масив (або скаляр): арифметику, редукції, трансляцію. Ця сторінка охоплює інший клас операцій: ті, що розглядають масив як дискретизовану функцію і ставлять запитання про саму функцію. Інтерполяція між відліками, апроксимація кривої по них, інтегрування під ними, згортка з іншим буфером.

Усі вони приймають вхідні дані ndarray і повертають або скаляр з плаваючою крапкою, або ndarray з плаваючою крапкою.

6.16.1. Інтерполяція

interp() виконує одновимірну лінійну інтерполяцію. xp — монотонно зростаючий одновимірний масив незалежних значень; fp — відповідні залежні значення; x — точки, в яких потрібно обчислити інтерполянт:

xp = np.array([0.0, 1.0, 2.0, 3.0])
fp = np.array([10.0, 20.0, 25.0, 30.0])
x  = np.array([0.5, 1.5, 2.5])

np.interp(x, xp, fp)
# array([15.0, 22.5, 27.5])

Поза діапазоном [xp[0], xp[-1]] результат обрізається до fp[0] і fp[-1] відповідно; ключові слова left= і right= замінюють ці крайні значення.

Використовується на камері для перемаштабування таблиці калібрування до довільних позицій вибірки – таблиця температура-напруга термістора, нелінійна крива відгуку пікселя, гамма-таблиця пошуку для кожного каналу. Один виклик бібліотеки, один прохід по вхідних даних, без циклу Python.

6.16.2. Поліноміальна апроксимація та обчислення

polyfit() апроксимує поліном степеня deg до точок даних (x, y) методом найменших квадратів і повертає коефіцієнти (починаючи з найвищого степеня):

x = np.array([0.0, 1.0, 2.0, 3.0, 4.0])
y = np.array([0.0, 1.1, 3.9, 9.1, 15.8])

coeffs = np.polyfit(x, y, 2)
# array([1.0, 0.0, 0.0])  (approximately)

Якщо x не вказано, використовується range(len(y)) – зручно для швидкої поліноміальної апроксимації рівномірно дискретизованого буфера без осі x:

np.polyfit(y, 2)

polyval() обчислює значення полінома з коефіцієнтами p у точці x. Вхідний x може бути скаляром (повертає float) або ndarray (повертає ndarray):

p = np.polyfit(x, y, 2)
fitted = np.polyval(p, x)

Природне поєднання – викликати polyfit один раз під час калібрування, зберегти коефіцієнти і викликати polyval для обчислення отриманої кривої на кожному кадрі. Крок обчислення полінома — це декілька операцій з плаваючою крапкою на відлік, що є недорогим навіть на найменших камерах.

6.16.3. Згортка

convolve() повертає повнодовжню дискретну лінійну згортку двох одновимірних масивів. Реалізовано лише режим full; довжина виводу – len(a) + len(v) - 1. Зрізайте результат для отримання ефекту, аналогічного режимам same і valid у настільному numpy

a = np.array([1.0, 2.0, 3.0])
v = np.array([0.5, 0.5])

np.convolve(a, v)
# array([0.5, 1.5, 2.5, 1.5])

Корисно для коротких FIR-фільтрів та ядер згладжування (прямокутне, трикутне, гаусівське), де налаштування ланцюжка SOS надлишкове. Час виконання пропорційний добутку довжин двох масивів – прийнятно для коротких ядер, але швидко стає дорожчим за FFT-згортку для довгих.

6.16.4. Трапецієподібне інтегрування

trapz() інтегрує дискретизовану функцію за складеним трапецієподібним правилом:

x = np.linspace(0, np.pi, num=128)
y = np.sin(x)
np.trapz(y, x)             # ~ 2.0

Передавайте dx=, коли крок між відліками рівномірний і важливий лише розмір кроку; передавайте x=, коли відліки нерівномірно розподілені. Правильний виклик для інтегрування вже захоплених даних датчика, де аналітична форма недоступна.

Для вибіркових даних, що були смугово обмежені (аудіобуфер після антиаліасингового фільтра, наприклад), трапецієподібне правило збігається як квадрат кількості відліків, тобто подвоєння довжини буфера зменшує похибку вчетверо.

Якщо підінтегральна функція є не буфером відліків, а функцією Python, яку застосунок може обчислювати у довільних точках, то правильним вибором є інший клас розв’язувачів.