配置Eigen矩阵运算库
后续计算需要利用矩阵运算来求解拟合系数,用到Eigen矩阵运算库,配置方法可看VS参考安装连接。
记得前两个画圈的地方改为所有配置和所有平台
拟合原理
以二次曲线拟合为例,拟合曲线应满足以下方程:
a ∗ X 2 + b ∗ X + c = y
如要要拟合的4个点为:
(3, 6)
(5, 8)
(6, 12)
(7, 4)
准备数据
由于Python的Matplotlib库可以很好地可视化数据,所以选用Python来生成数据文件——points.txt。
代码如下,程序中设置的系数为a = 2.2, b = 1.4, c = -1.3:
import numpy as np
import matplotlib.pyplot as plt
# 设置点数
num = 100
# 设置噪声系数
rate = 0
# 生成数据
x = np.linspace(-10, 10, num)
y = []
a = 2.2
b = 1.4
c = -1.3
for i in x:
noise = np.random.random()
temp = i + noise * rate
y.append(a * pow(temp, 2) + b * temp + c)
# 绘制曲线
plt.plot(x, y)
plt.show()
# 写入文件
with open("./points.txt", 'w') as f:
for i in range(num):
f.write(str(x[i]))
f.write("\t")
f.write(str(y[i]))
f.write("\n")
曲线为:
程序会在当前目录下生成points.txt文件。
C++最小二乘拟合代码
以下代码理论上可以完成N阶的曲线拟合,需要对参数N进行设置,并提供正确的数据点。
#include <Eigen/Core>
#include <Eigen/Dense>
#include <Eigen/Geometry>
#include <Eigen/Eigenvalues>
#include <iostream>
#include <fstream>
#include <string>
#include <vector>
#include <cmath>
using namespace std;
int main() {
// txt点数据文件路径
const string fileName = "./points.txt";
// 设置是几次拟合
const int N = 2;
// 创建两个vector
vector<float> x, y;
// 读取文件
ifstream f(fileName);
if (!f) {
cout << "数据文件打开失败" << endl;
exit(EXIT_FAILURE);
}
float tempx, tempy;
while (f >> tempx >> tempy) {
x.push_back(tempx);
y.push_back(tempy);
}
if (x.size() != y.size()) {
cout << "数据文件内容有误" << endl;
exit(EXIT_FAILURE);
}
// 创建A矩阵
Eigen::MatrixXd A(x.size(), N + 1);
for (unsigned int i = 0; i < x.size(); ++i) { // 遍历所有点
for (int n = N, dex = 0; n >= 1; --n, ++dex) { // 遍历N到1阶
A(i, dex) = pow(x[i], n);
}
A(i, N) = 1; //
}
// 创建B矩阵
Eigen::MatrixXd B(y.size(), 1);
for (unsigned int i = 0; i < y.size(); ++i) {
B(i, 0) = y[i];
}
// 创建矩阵W
Eigen::MatrixXd W;
W = (A.transpose() * A).inverse() * A.transpose() * B;
// 打印W结果
cout << W << endl;
}