Skip to content

Commit c375a5e

Browse files
committed
fix(quadratic_fit): 拟合过程使用双精度浮点,避免精度丢失
1 parent f10f786 commit c375a5e

File tree

1 file changed

+10
-10
lines changed

1 file changed

+10
-10
lines changed

src/algorithms/polyfit.c

Lines changed: 10 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -130,7 +130,7 @@ void find_max_y_x(double coef[4], double x_min, double x_max, double *x_max_y, d
130130
}
131131
void linear_curve_fit(float *x, float *y, int size, float *slope, float *intercept)
132132
{
133-
float sum_x = 0, sum_y = 0, sum_xy = 0, sum_x2 = 0;
133+
double sum_x = 0, sum_y = 0, sum_xy = 0, sum_x2 = 0;
134134

135135
for (int i = 0; i < size; i++)
136136
{
@@ -146,13 +146,13 @@ void linear_curve_fit(float *x, float *y, int size, float *slope, float *interce
146146
// 二次拟合函数
147147
void quadratic_fit(float *x, float *y, int n, float *coeff)
148148
{
149-
float sum_x = 0, sum_x2 = 0, sum_x3 = 0, sum_x4 = 0;
150-
float sum_y = 0, sum_xy = 0, sum_x2y = 0;
149+
double sum_x = 0, sum_x2 = 0, sum_x3 = 0, sum_x4 = 0;
150+
double sum_y = 0, sum_xy = 0, sum_x2y = 0;
151151

152152
for (int i = 0; i < n; i++)
153153
{
154-
float xi = x[i];
155-
float yi = y[i];
154+
double xi = x[i];
155+
double yi = y[i];
156156
sum_x += xi;
157157
sum_x2 += xi * xi;
158158
sum_x3 += xi * xi * xi;
@@ -162,20 +162,20 @@ void quadratic_fit(float *x, float *y, int n, float *coeff)
162162
sum_x2y += xi * xi * yi;
163163
}
164164

165-
float A[3][3] = {
165+
double A[3][3] = {
166166
{ n, sum_x, sum_x2},
167167
{ sum_x, sum_x2, sum_x3},
168168
{sum_x2, sum_x3, sum_x4}
169169
};
170170

171-
float B[3] = {sum_y, sum_xy, sum_x2y};
171+
double B[3] = {sum_y, sum_xy, sum_x2y};
172172

173173
// 解线性方程组
174174
for (int i = 0; i < 3; i++)
175175
{
176176
for (int j = i + 1; j < 3; j++)
177177
{
178-
float ratio = A[j][i] / A[i][i];
178+
double ratio = A[j][i] / A[i][i];
179179
for (int k = 0; k < 3; k++) { A[j][k] -= ratio * A[i][k]; }
180180
B[j] -= ratio * B[i];
181181
}
@@ -230,8 +230,8 @@ int calculateBoundedExtremum(const float coeff[3], float x_min, float x_max, flo
230230
else
231231
{
232232
// 极值点在范围外,计算两个端点的值
233-
float y_min = coeff[0] + coeff[1] * x_min + coeff[2] * x_min * x_min;
234-
float y_max = coeff[0] + coeff[1] * x_max + coeff[2] * x_max * x_max;
233+
double y_min = coeff[0] + coeff[1] * x_min + coeff[2] * x_min * x_min;
234+
double y_max = coeff[0] + coeff[1] * x_max + coeff[2] * x_max * x_max;
235235

236236
if (coeff[2] > 0)
237237
{

0 commit comments

Comments
 (0)