rusefi-1/firmware/util/math/interpolation.cpp

159 lines
4.0 KiB
C++
Raw Normal View History

2015-07-10 06:01:56 -07:00
/**
* @file interpolation.cpp
* @brief Linear interpolation algorithms
*
2018-08-19 07:53:25 -07:00
* See test_interpolation_3d.cpp
*
*
2015-07-10 06:01:56 -07:00
* @date Oct 17, 2013
2020-01-13 18:57:43 -08:00
* @author Andrey Belomutskiy, (c) 2012-2020
2015-07-10 06:01:56 -07:00
* @author Dmitry Sidin, (c) 2015
*/
#include "pch.h"
2015-07-10 06:01:56 -07:00
#include "interpolation.h"
#define BINARY_PERF true
2016-06-29 17:02:00 -07:00
#if BINARY_PERF && ! EFI_UNIT_TEST
2015-07-10 06:01:56 -07:00
#define COUNT 10000
float array16[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16 };
static void testBinary() {
2015-07-10 06:01:56 -07:00
const int size16 = 16;
uint32_t totalOld = 0;
uint32_t totalNew = 0;
for (int v = 0; v <= 16; v++) {
uint32_t timeOld;
{
2019-05-07 16:32:08 -07:00
uint32_t start = getTimeNowLowerNt();
2015-07-10 06:01:56 -07:00
int temp = 0;
for (int i = 0; i < COUNT; i++) {
temp += findIndex(array16, size16, v);
}
2019-05-07 16:32:08 -07:00
timeOld = getTimeNowLowerNt() - start;
2015-07-10 06:01:56 -07:00
}
uint32_t timeNew;
{
2019-05-07 16:32:08 -07:00
uint32_t start = getTimeNowLowerNt();
2015-07-10 06:01:56 -07:00
int temp = 0;
for (int i = 0; i < COUNT; i++) {
temp += findIndex2(array16, size16, v);
}
2019-05-07 16:32:08 -07:00
timeNew = getTimeNowLowerNt() - start;
2015-07-10 06:01:56 -07:00
}
efiPrintf("for v=%d old=%d ticks", v, timeOld);
efiPrintf("for v=%d new=%d ticks", v, timeNew);
2015-07-10 06:01:56 -07:00
totalOld += timeOld;
totalNew += timeNew;
}
efiPrintf("totalOld=%d ticks", totalOld);
efiPrintf("totalNew=%d ticks", totalNew);
2015-07-10 06:01:56 -07:00
}
#endif
/** @brief Linear interpolation by two points
*
* @param x1 key of the first point
* @param y1 value of the first point
* @param x2 key of the second point
* @param y2 value of the second point
* @param X key to be interpolated
*
2018-06-12 02:45:11 -07:00
* @note For example, "interpolateMsg("", engineConfiguration.tpsMin, 0, engineConfiguration.tpsMax, 100, adc);"
2021-11-02 20:34:44 -07:00
* @see interpolateClamped
2015-07-10 06:01:56 -07:00
*/
2016-01-09 12:01:41 -08:00
float interpolateMsg(const char *msg, float x1, float y1, float x2, float y2, float x) {
2019-06-27 19:23:18 -07:00
if (cisnan(x1) || cisnan(x2) || cisnan(y1) || cisnan(y2)) {
warning(CUSTOM_INTEPOLATE_ERROR, "interpolate%s: why param", msg);
return NAN;
}
if (cisnan(x)) {
warning(CUSTOM_INTEPOLATE_ERROR, "interpolate%s: why X", msg);
return NAN;
}
2015-07-10 06:01:56 -07:00
// todo: double comparison using EPS
if (x1 == x2) {
/**
* we could end up here for example while resetting bins while changing engine type
*/
2018-01-23 09:05:14 -08:00
warning(CUSTOM_INTEPOLATE_ERROR, "interpolate%s: Same x1 and x2 in interpolate: %.2f/%.2f", msg, x1, x2);
2015-07-10 06:01:56 -07:00
return NAN;
}
// a*x1 + b = y1
// a*x2 + b = y2
2018-07-25 20:03:04 -07:00
// efiAssertVoid(CUSTOM_ERR_ASSERT_VOID, x1 != x2, "no way we can interpolate");
2015-07-10 06:01:56 -07:00
float a = INTERPOLATION_A(x1, y1, x2, y2);
2019-06-27 19:23:18 -07:00
if (cisnan(a)) {
warning(CUSTOM_INTEPOLATE_ERROR, "interpolate%s: why a", msg);
return NAN;
}
2015-07-10 06:01:56 -07:00
float b = y1 - a * x1;
float result = a * x + b;
#if DEBUG_FUEL
2018-01-23 09:05:14 -08:00
printf("x1=%.2f y1=%.2f x2=%.2f y2=%.2f\r\n", x1, y1, x2, y2);
printf("a=%.2f b=%.2f result=%.2f\r\n", a, b, result);
2015-07-10 06:01:56 -07:00
#endif
return result;
}
2021-11-02 20:34:44 -07:00
/**
* @see interpolateMsg
*/
float interpolateClamped(float x1, float y1, float x2, float y2, float x) {
if (x <= x1)
return y1;
if (x >= x2)
return y2;
float a = INTERPOLATION_A(x1, y1, x2, y2);
float b = y1 - a * x1;
return a * x + b;
}
2017-05-28 21:02:22 -07:00
/**
* Another implementation, which one is faster?
*/
2015-07-10 06:01:56 -07:00
int findIndex2(const float array[], unsigned size, float value) {
2018-07-25 20:30:00 -07:00
efiAssert(CUSTOM_ERR_ASSERT, !cisnan(value), "NaN in findIndex2", 0);
efiAssert(CUSTOM_ERR_ASSERT, size > 1, "size in findIndex", 0);
2015-07-10 06:01:56 -07:00
// if (size <= 1)
// return size && *array <= value ? 0 : -1;
signed i = 0;
//unsigned b = 1 << int(log(float(size) - 1) / 0.69314718055994530942);
unsigned b = size >> 1; // in our case size is always a power of 2
2018-07-25 20:30:00 -07:00
efiAssert(CUSTOM_ERR_ASSERT, b + b == size, "Size not power of 2", -1);
2015-07-10 06:01:56 -07:00
for (; b; b >>= 1) {
unsigned j = i | b;
/**
* it should be
* "if (j < size && array[j] <= value)"
* but in our case size is always power of 2 thus size is always more then j
*/
2018-07-25 20:30:00 -07:00
// efiAssert(CUSTOM_ERR_ASSERT, j < size, "size", 0);
2015-07-10 06:01:56 -07:00
if (array[j] <= value)
i = j;
}
return i || *array <= value ? i : -1;
}
2017-05-28 21:02:22 -07:00
int findIndex(const float array[], int size, float value) {
return findIndexMsg("", array, size, value);
}
void initInterpolation() {
2015-07-10 06:01:56 -07:00
#if BINARY_PERF && ! EFI_UNIT_TEST
addConsoleAction("binarytest", testBinary);
#endif
}