www.pudn.com > Polyfit.rar > Polyfit.cs, change:2012-12-26,size:5716b

using System; 
using System.Collections.Generic; 
using System.Text; 
 
namespace PolyFit 
{ 
    class Polyfit 
    { 
        private static double[,] MatrixA; 
        private static double[] MatrixB; 
        private static double[,] MatrixMix; 
 
        /// <summary> 
        /// 一元N次多项式系数 
        /// </summary> 
        public static double[] MatrixAk; 
        //private static int m; // 数据点 个数 
        private static int N; 
 
        private static void getReady(int n, double[] dArrXSrc, double[] dArrYSrc) 
        { 
            N = n; 
            MatrixA = new double[N + 1, N + 1]; 
            MatrixMix = new double[N + 1, N + 1]; 
            MatrixB = new double[N + 1]; 
            MatrixAk = new double[N + 1]; 
            indx = new int[N + 1]; 
 
            for (int iRow = 0; iRow <= N; iRow++) 
            { 
                for (int iCol = 0; iCol <= N; iCol++) 
                { 
                    double SumOfNPowOfXi = 0; 
                    for (int i = 0; i < dArrXSrc.Length; i++) 
                    { 
                        SumOfNPowOfXi += Math.Pow(dArrXSrc[i], iCol + iRow); 
                    } 
                    MatrixA[iRow, iCol] = SumOfNPowOfXi; 
                } 
                double SumOfNPowOfXiYi = 0; 
                for (int i = 0; i < dArrYSrc.Length; i++) 
                { 
                    SumOfNPowOfXiYi += (Math.Pow(dArrXSrc[i], iRow) * dArrYSrc[i]); 
                } 
                MatrixB[iRow] = SumOfNPowOfXiYi; 
            } 
 
            //if (dArrYSrc.Length > 0) 
            //    m = dArrXSrc.Length; 
        } 
 
        private static int[] indx; 
        private static int parity = 1; 
 
        private static void LUCmp() 
        { 
            double tiny = 1.0e-20; 
            int iMax = 0; 
            double big, dum, temp; 
            double[] vv = new double[N + 1]; 
 
            for (int iRow = 0; iRow <= N; iRow++) 
                for (int iCol = 0; iCol <= N; iCol++) 
                    MatrixMix[iRow, iCol] = MatrixA[iRow, iCol]; 
 
            for (int i = 0; i <= N; i++) 
            { 
                big = 0.0; 
                for (int j = 0; j <= N; j++) 
                { 
                    if ((temp = Math.Abs(MatrixMix[i, j])) > big) 
                        big = temp; 
                } 
                if (big == 0.0); 
                vv[i] = 1.0 / big; 
            } 
 
            for (int j = 0; j <= N; j++) 
            { 
                for (int i = 0; i < j; i++) 
                { 
                    double dSumOfAikBkj = MatrixMix[i, j]; 
                    for (int k = 0; k < i; k++) 
                    { 
                        dSumOfAikBkj -= MatrixMix[i, k] * MatrixMix[k, j]; 
                    } 
                    MatrixMix[i, j] = dSumOfAikBkj; 
                } 
                big = 0.0; 
                for (int i = j; i <= N; i++) 
                { 
                    double dSumOfAikBkj = MatrixMix[i, j]; 
                    for (int k = 0; k < j; k++) 
                    { 
                        dSumOfAikBkj -= MatrixMix[i, k] * MatrixMix[k, j]; 
                    } 
                    MatrixMix[i, j] = dSumOfAikBkj; 
                    if ((dum = vv[i] * Math.Abs(dSumOfAikBkj)) >= big) 
                    { 
                        big = dum; 
                        iMax = i; 
                    } 
                } 
                if (j != iMax) 
                { 
                    for (int k = 0; k <= N; k++) 
                    { 
                        dum = MatrixMix[iMax, k]; 
                        MatrixMix[iMax, k] = MatrixMix[j, k]; 
                        MatrixMix[j, k] = dum; 
                    } 
                    parity = -parity; 
                    dum = vv[iMax]; 
                    vv[iMax] = vv[j]; 
                    vv[j] = dum; 
                } 
                indx[j] = iMax; 
                if (MatrixMix[j, j] == 0.0) 
                    MatrixMix[j, j] = tiny; 
                if (j != N) 
                { 
                    dum = 1.0 / MatrixMix[j, j]; 
                    for (int i = j + 1; i <= N; i++) 
                        MatrixMix[i, j] *= dum; 
                } 
            } 
        } 
 
        private static void LUbksb() 
        { 
            for (int i = 0; i < MatrixB.Length; i++) 
                MatrixAk[i] = MatrixB[i]; 
 
            double sum; 
            int ii = 0; 
            for (int i = 0; i <= N; i++) 
            { 
                int ip = indx[i]; 
                sum = MatrixAk[ip]; 
                MatrixAk[ip] = MatrixAk[i]; 
                if (ii != 0) 
                    for (int j = ii - 1; j < i; j++) 
                        sum -= MatrixMix[i, j] * MatrixAk[j]; 
                else if (sum != 0.0) 
                    ii = i + 1; 
                MatrixAk[i] = sum; 
            } 
            for (int i = N; i >= 0; i--) 
            { 
                sum = MatrixAk[i]; 
                for (int j = i + 1; j <= N; j++) 
                    sum -= MatrixMix[i, j] * MatrixAk[j]; 
                MatrixAk[i] = sum / MatrixMix[i, i]; 
            } 
        } 
 
        /// <summary> 
        /// 基于最小二乘法的多项式拟合 
        /// </summary> 
        /// <param name="n">一元n次多项式</param> 
        /// <param name="dArrX">测试点x坐标值</param> 
        /// <param name="dArrY">测试点y坐标值</param> 
        public static void fittingOfAPolynomial(int n, double[] dArrX, double[]dArrY) 
        { 
            getReady(n, dArrX, dArrY); 
            LUCmp(); 
            LUbksb(); 
        } 
    } 
}