3、最小二乘法拟合：使用最小二乘法计算校准参数、矩阵运算与求解、误差分析与精度评估

好，咱们接着聊。上一节我们把触摸数据采回来了，一堆 (LCD_x, LCD_y, ADC_x, ADC_y) 的配对点。现在问题来了：怎么从这些离散的点里，算出最靠谱的校准参数？

你可能会想，直接拿两组点解方程不就行了？嗯，理论上可以。但实际采样总有噪声，手指按下去也有偏差。用两组点算出来的参数，换几个点测可能就偏了。我早期做一款 7 寸屏时就吃过这个亏——用 3 个点算参数，结果屏幕左上角准了，右下角差了十几个像素。

所以，我们需要一种方法：让所有采样点都“参与投票”，找出一个整体误差最小的解。这就是最小二乘法干的事。

3.1 最小二乘法的核心思想

说白了，最小二乘法就是找一组参数，让所有采样点的“预测值”和“实际值”的差的平方和最小。为什么用平方？因为平方可以放大偏差，让大的误差更显眼，优化时会优先把它们压下去。

对于触摸校准，我们通常用线性模型：

LCD_x = A * ADC_x + B * ADC_y + C
LCD_y = D * ADC_x + E * ADC_y + F

这里有 6 个未知数（A, B, C, D, E, F）。理论上，3 组不共线的点就能解。但用 N 组点（N > 3），就变成了一个超定方程组。最小二乘法就是解这个超定方程组的最佳工具。

3.2 矩阵运算与求解步骤

别被矩阵吓到。实际做的时候，我们一步步来。我个人习惯用 C 语言手写矩阵运算，因为嵌入式环境不一定有现成的库。

假设我们采集了 N 个点，每个点有 (adc_x, adc_y, lcd_x, lcd_y)。

第一步：构造矩阵 X 和向量 Y

// 对于 X 矩阵，每行是 [adc_x, adc_y, 1]
// 对于 Y 向量，先算 X 方向
// 矩阵大小：X 是 N×3，Y 是 N×1

float X[N][3];
float Yx[N];  // 用于计算 X 方向的参数
float Yy[N];  // 用于计算 Y 方向的参数

for (int i = 0; i < N; i++) {
    X[i][0] = adc_x[i];
    X[i][1] = adc_y[i];
    X[i][2] = 1.0f;
    Yx[i] = lcd_x[i];
    Yy[i] = lcd_y[i];
}

第二步：计算 X^T * X 和 X^T * Y

这一步是矩阵乘法。X^T * X 是一个 3×3 的对称矩阵，X^T * Y 是一个 3×1 的向量。

float XTX[3][3] = {0};
float XTYx[3] = {0};
float XTYy[3] = {0};

for (int i = 0; i < N; i++) {
    // XT * X
    XTX[0][0] += X[i][0] * X[i][0];
    XTX[0][1] += X[i][0] * X[i][1];
    XTX[0][2] += X[i][0] * X[i][2];
    XTX[1][0] += X[i][1] * X[i][0];
    XTX[1][1] += X[i][1] * X[i][1];
    XTX[1][2] += X[i][1] * X[i][2];
    XTX[2][0] += X[i][2] * X[i][0];
    XTX[2][1] += X[i][2] * X[i][1];
    XTX[2][2] += X[i][2] * X[i][2];

    // XT * Yx
    XTYx[0] += X[i][0] * Yx[i];
    XTYx[1] += X[i][1] * Yx[i];
    XTYx[2] += X[i][2] * Yx[i];

    // XT * Yy
    XTYy[0] += X[i][0] * Yy[i];
    XTYy[1] += X[i][1] * Yy[i];
    XTYy[2] += X[i][2] * Yy[i];
}

第三步：求解线性方程组

现在我们有 X^TX * K = X^TY。解这个 3 元一次方程组。我常用高斯消元法，稳定且容易实现。

// 高斯消元求解 3x3 方程组
// 输入：A[3][3] 系数矩阵，b[3] 右侧向量
// 输出：x[3] 解向量

int gauss_3x3(float A[3][3], float b[3], float x[3]) {
    float M[3][4];
    for (int i = 0; i < 3; i++) {
        for (int j = 0; j < 3; j++) M[i][j] = A[i][j];
        M[i][3] = b[i];
    }

    // 前向消元
    for (int col = 0; col < 3; col++) {
        // 选主元（这里简化，实际建议加选主元逻辑）
        float pivot = M[col][col];
        if (fabs(pivot) < 1e-10) return -1;  // 奇异矩阵

        for (int row = col + 1; row < 3; row++) {
            float factor = M[row][col] / pivot;
            for (int k = col; k <= 3; k++) {
                M[row][k] -= factor * M[col][k];
            }
        }
    }

    // 回代
    for (int i = 2; i >= 0; i--) {
        x[i] = M[i][3];
        for (int j = i + 1; j < 3; j++) {
            x[i] -= M[i][j] * x[j];
        }
        x[i] /= M[i][i];
    }
    return 0;
}

// 调用两次，分别求 X 方向和 Y 方向的参数
float Kx[3], Ky[3];
gauss_3x3(XTX, XTYx, Kx);  // Kx = [A, B, C]
gauss_3x3(XTX, XTYy, Ky);  // Ky = [D, E, F]

3.3 误差分析与精度评估

参数算出来了，但效果怎么样？不能光靠感觉。我们需要量化评估。

残差计算

对每个采样点，用算出的参数反算预测的 LCD 坐标，然后和实际值比较：

float total_error_x = 0;
float total_error_y = 0;
float max_error_x = 0;
float max_error_y = 0;

for (int i = 0; i < N; i++) {
    float pred_x = Kx[0] * adc_x[i] + Kx[1] * adc_y[i] + Kx[2];
    float pred_y = Ky[0] * adc_x[i] + Ky[1] * adc_y[i] + Ky[2];

    float err_x = fabs(pred_x - lcd_x[i]);
    float err_y = fabs(pred_y - lcd_y[i]);

    total_error_x += err_x;
    total_error_y += err_y;

    if (err_x > max_error_x) max_error_x = err_x;
    if (err_y > max_error_y) max_error_y = err_y;
}

float rmse_x = sqrt(total_error_x / N);
float rmse_y = sqrt(total_error_y / N);

评估指标

指标	含义	可接受范围
RMSE（均方根误差）	整体偏差的“平均大小”	< 3 像素
最大误差	最差的那个点偏差多少	< 5 像素
平均误差	简单平均偏差	< 2 像素

3.4 知识体系总览

下面这张图把最小二乘法校准的完整流程串起来了，你可以对照着看：

3.5 实战中的几个坑

讲几个我踩过的坑，你遇到了可以少走弯路：

• 采样点数量：我曾经只用 4 个点，结果参数对采样点拟合很好，但采样点之间的区域偏差很大。后来我改成 9 个点（3×3 网格），效果稳定多了。
• 浮点精度：在 Cortex-M3 上做浮点运算，如果 ADC 值很大（比如 12 位 ADC，范围 0~4095），X^TX 矩阵里的数值会很大。我遇到过因为 float 精度不够导致矩阵求逆失败的情况。解决方案：先把 ADC 值归一化到 0~1 范围再计算。
• 验证集：不要用参与拟合的点来评估精度。我习惯留 2~3 个点不参与计算，专门用来验证。这样测出来的误差才是真实的。