最小二乘法的曲线拟合

最小二乘法的曲线拟合

BOOL CDataDistillView::LeastDoubleMultiplication(LONG *pX, LONG *pY, long M, long N, double *result, double &warp){//data--  X,Y两轴的坐标//M--  结果变量组数//N--  采样数目//result--  结果参数//warp--偏差平方和register long i, j, k;double Z, D1, D2, C, P, G, Q;double *B, *T, *S;B = new double[N];T = new double[N];S = new double[N];if(M > N){M = N;}for(i = 0; i < M; i++){result[i] = 0;}Z = 0;B[0] = 1;D1 = N;P = 0;C = 0;for(i = 0; i < N; i++){P = P + pX[i] - Z;C = C + pY[i];}C = C / D1;P = P / D1;result[0] = C * B[0];if(M > 1){T[1] = 1;T[0] = -P;D2 = 0;C = 0;G = 0;for(i = 0; i < N; i++){Q = pX[i] - Z - P;D2 = D2 + Q * Q;C = pY[i] * Q + C;G = (pX[i] - Z) * Q * Q + G;}C = C / D2;P = G / D2;Q = D2 / D1;D1 = D2;result[1] = C * T[1];result[0] = C * T[0] + result[0];}for(j = 2; j < M; j++){S[j] = T[j - 1];S[j - 1] = -P * T[j - 1] + T[j - 2];if(j >= 3){for(k=j-2;k>=1;k--){S[k] = -P * T[k] + T[k - 1] - Q * B[k];}}S[0] = -P * T[0] - Q * B[0];D2 = 0;C = 0;G = 0;for(i = 0; i < N; i++){Q=S[j];for(k = j - 1; k >= 0; k--){Q = Q * (pX[i] - Z) + S[k];}D2 = D2 + Q * Q;C = pY[i] * Q + C;G = (pX[i] - Z) * Q * Q + G;}C = C / D2;P = G / D2;Q = D2 / D1;D1 = D2;result[j] = C * S[j];T[j] = S[j];for(k = j - 1;k >= 0;k--){result[k] = C * S[k] + result[k];B[k] = T[k];T[k] = S[k];}}warp = 0;double dt; for(i = 0; i < N; i++)   {Q = result[M - 1];   for(k = M-1; k>=1; k--)   {Q = Q * (pX[i] - Z) + result[k - 1];   }dt = Q - pY[i];   warp = warp + dt * dt;   }delete [] B;delete [] T;delete [] S;return TRUE;}