Biegral Blog

using System;
using System.Collections;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.IO;
namespace CsTest
{   //可以确定的是n个插值点的拉格朗日多项式的最大次数为n
    public class Ploy   //多项式
    {
        public double[] a;//系数,指数即索引值
        public Ploy(int n = 1)
        {
            a = new double[n]; //当未给数组中的每一项赋值时，其值为0
        }
        public static Ploy operator +(Ploy a1, Ploy b1)//多项式a1和b1的项数相同
        {
            Ploy p = new Ploy(a1.a.Length);
            for (int i = 0; i < a1.a.Length; i++)
            {
                p.a[i] = a1.a[i] + b1.a[i];//系数相加
            }
            return p;
        }
        public static Ploy operator *(Ploy a1, Ploy b1)
        {
            Ploy p = new Ploy(a1.a.Length);
            for (int i = 0; i < a1.a.Length; i++)
            {
                if (a1.a[i] == 0)
                {
                    continue;  //此项系数为0
                }
                for (int j = 0; j < a1.a.Length; j++)
                {
                    if (b1.a[j] == 0)//此项系数为0
                    {
                        continue;
                    }
                    if (i + j > a1.a.Length)
                    {
                        Console.WriteLine("运算出现错误");
                        break;
                    }
                    p.a[i + j] += a1.a[i] * b1.a[j];
                }
            }
            return p;
        }
    }
    class Program
    {
        static double NewtonDiff(int[] X, int[] Y, int m)//当m=1时差商为对应Y值
        {
            //double result=0;
            if (m == 2)
            {
                return (Y[m - 1] - Y[m - 2]) / (X[m - 1] - X[m - 2]);
            }
            else
            {
                int[] xfore = new int[m - 1];
                int[] yfore = new int[m - 1];
                int[] xlast = new int[m - 1];
                int[] ylast = new int[m - 1];
                for (int i = 0; i < m - 1; i++)
                {
                    xfore[i] = X[i];
                    yfore[i] = Y[i];
                    xlast[i] = X[i + 1];
                    ylast[i] = Y[i + 1];
                }
                return (NewtonDiff(xlast, ylast, m - 1) - NewtonDiff(xfore, yfore, m - 1)) / (X[m - 1] - X[0]);
            }
        }

        static void Newton(int[] X, int[] Y, out Ploy Newt)
        {

            Ploy tm = new Ploy(X.Length);
            tm.a[0] = Y[0];  //第一个差商即函数值本身，是一个常数项，即对应x零次方的系数
            Newt = tm;
            for (int i = X.Length; i != 1; i--)
            {
                Ploy multy = new Ploy(X.Length);
                multy.a[0] = NewtonDiff(X, Y, i);//对应第i个差商
                for (int j = 1; j < i; j++)
                {
                    Ploy tm1 = new Ploy(X.Length);
                    tm1.a[1] = 1;
                    tm1.a[0] = -X[j - 1];
                    multy = multy * tm1;
                }
                Newt = Newt + multy;
            }
        }
        static void Main(string[] args)
        {
            int n = 4;//没选取n个节点进行一次牛顿运算
            string Inpath = @"E://360壁纸//color.txt";//原数据文件路径
            string Outpath = @"F://Rout.txt";//输出数据文件
            StreamReader sr = new StreamReader(Inpath, Encoding.Default);
            StreamWriter sw = new StreamWriter(Outpath, true);
            string temps;
            // ArrayList aL = new ArrayList();//将节点存入数组
            while ((temps = sr.ReadLine()) != null)
            {
                int Length = 0;
                string[] s = temps.Split(new char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);//根据空格区分数据
                Length = s.Length; //获取一行有多少节点

                int[] Z = new int[Length * 2 - 1];//存放当前行所做牛顿插值分段运算的值
                for (int i = 0; i < Length - n + 1; i++)//对当前行进行分段牛顿插值运算
                {
                    int[] X = new int[n];
                    for (int j = 0; j < n; j++)
                    {
                        X[j] = j;
                    }
                    int[] Y = new int[n];
                    for (int k = 0; k < n; k++)
                    {
                        Y[k] = Convert.ToInt32(s[i + k]);
                    }
                    Ploy Newt;
                    Program.Newton(X, Y, out Newt);

                    int[] r = new int[n * 2 - 1];
                    for (int m = 0; m < n; m++)//带入横坐标求值
                    {
                        double xhalf = Convert.ToDouble(m) + 0.5;
                        double result = 0;
                        for (int t = n - 1; t >= 0; t--)//求对应于result的牛顿插值
                        {
                            result += Newt.a[t] * Math.Pow(xhalf, t);
                        }
                        r[m * 2] = Y[m];
                        if (result > 255)
                            result = 255;
                        if (m != n - 1)
                        {
                            r[m * 2 + 1] = Convert.ToInt32(result);
                        }

                    }
                    if (i == 0) //说明是第一个指定n个节点所做的牛顿运算
                    {
                        for (int m = 0; m < r.Length; m++)
                        {
                            Z[i + m] = r[m];
                        }
                    }
                    else    //倒数第2个之前的值与Z中对应的值求平均
                    {
                        for (int m = 0; m < r.Length - 2; m++)
                        {
                            Z[i * 2 + m] = (Z[i * 2 + m] + r[m]) / 2;
                        }
                        Z[i * 2 + r.Length - 2] = r[r.Length - 2];
                        Z[i * 2 + r.Length - 1] = r[r.Length - 1];
                    }

                    //将Z中的值输出到文件中
                    for (int m = 0; m < Z.Length; m++)
                    {
                        sw.Write(Z[m].ToString() + " ");
                    }
                    sw.WriteLine();
                }
            }
            sw.Flush();
            sw.Close();
            sr.Close();
        }
    }
}