using MathNet.Numerics;
using MathNet.Numerics.LinearAlgebra.Double;
using NPOI.SS.Formula.Functions;
using System;
using System.Collections.Generic;
using System.Drawing;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
using System.Windows.Forms.DataVisualization.Charting;

namespace Hydro.CommonBase
{
    public class CurveFitHelper
    {


        public DenseVector coefficients;
        int degree;
        double xMax;
        double xMin;
        List<PointF> data;
        bool isFitted = false;

        public CurveFitHelper(List<PointF> data, int degree = 2)
        {
            this.degree = degree;
            this.data = data;

            double[] xlist = data.Select(p => (double)p.X).ToArray();
            if (data.Count>2)
            {
                coefficients = Fit.Polynomial(xlist, data.Select(p => (double)p.Y).ToArray(), 2);
                isFitted = true;
            }
            else if(data.Count==2)
            {
                //使用两个点的线性拟合
                var tuple = Fit.Line(xlist, data.Select(p => (double)p.Y).ToArray());
                coefficients= new DenseVector(new double[] { tuple.Item1, tuple.Item2 });
                isFitted = true;
            }
            else if (data.Count==1)
            {
                coefficients = new DenseVector(new double[] { data[0].Y });
                isFitted = true;
            }
            else
            {
                coefficients = new DenseVector(new double[] { 0});
                isFitted = true;
            }
            
            var list = xlist.ToList();
            xMax = list.Max();
            xMin = list.Min();
        }

        public double Evaluate(double x)
        {
            if (!isFitted) return double.NaN;
            double y = 0;
            for (int i = 0; i < coefficients.Count; i++)
            {
                y += (float)(coefficients[i] * Math.Pow(x, i));
            }
            return y;
        }


        public List<PointF> GetFitCurve(int count)
        {
            List<PointF> curve = new List<PointF>();
            
            for (double x = xMin; x <= xMax && curve.Count<count; x += (xMax - xMin) / count)
            {
                double y = Evaluate(x);
                curve.Add(new PointF((float)x, (float)y));
            }
            return curve;

        }

        public List<List<PointF>> GetConfidenceCurve(int count)
        {
            List<PointF> up_list = new List<PointF>();
            List<PointF> down_list = new List<PointF>();
            List<PointF> avg_list = new List<PointF>();

            double X =0;
            double Y=0;
            for (double x = xMin; x < xMax; x += (xMax - xMin) / count)
            {
                var list = data.FindAll(p => x - (xMax - xMin) / count / 2 <= p.X && p.X < x + (xMax - xMin) / count/2).Select(p => (double)p.Y).ToList();
                if (list.Count >= 3)
                {
                    var Bounds = GetConfidenceRange(list, 95);
                    var x_avg = x;
                    if (Bounds[0] == double.NaN || Bounds[1] == double.NaN) continue;
                    down_list.Add(new PointF((float)x_avg, (float)Bounds[0]));
                    up_list.Add(new PointF((float)x_avg, (float)Bounds[1]));
                    X = x;
                    //平均值
                    Y = (Bounds[0] + Bounds[1] )/ 2;
                    avg_list.Add(new PointF((float)x_avg, (float)Y));
                }
                else if (list.Count == 2)
                {
                    var x_avg = x;
                    down_list.Add(new PointF((float)x_avg, (float)list.Min()));
                    up_list.Add(new PointF((float)x_avg, (float)list.Max()));
                    X = x;
                    Y = (list.Min() + list.Max()) / 2;
                    avg_list.Add(new PointF((float)x_avg, (float)Y));
                }
                else if (list.Count == 1)
                                       {
                    var x_avg = x;
                    down_list.Add(new PointF((float)x_avg, (float)list[0]));
                    up_list.Add(new PointF((float)x_avg, (float)list[0]));
                    X = x; 
                    Y = list[0];
                    avg_list.Add(new PointF((float)x_avg, (float)Y));
                }
                else if(list.Count == 0)
                {
                    //使用前一个点的值,X和Y都是前一个点的值
                    down_list.Add(new PointF((float)X, (float)Y));
                    up_list.Add(new PointF((float)X, (float)Y));
                    avg_list.Add(new PointF((float)X, (float)Y));
                }

            }

            //CurveFitHelper h0 = new CurveFitHelper(down_list);
            //var curve0 = h0.GetFitCurve(count);
            //CurveFitHelper h1 = new CurveFitHelper(up_list);
            //var curve1 = h1.GetFitCurve(count);
            return new List<List<PointF>> { down_list, up_list, avg_list };
        }
        public static double[] GetConfidenceRange(List<double> data, double confidenceLevel = 90, bool mode = true)//true是中位数，置信区间是false
        {
            //如果model==true,这计算中位数,比例为confidenceLevel
            if (mode)
            {

                data.Sort();

                double lowerBound = 0;
                double upperBound = 0;
                //取百分位confidenceLevel/100和1-confidenceLevel/100.0的数
                int upperIndex = (int)Math.Round(data.Count * confidenceLevel / 100.0);
                int lowerIndex = (int)Math.Round(data.Count * (1 - confidenceLevel / 100.0));
                if (upperIndex == data.Count)
                {
                    upperIndex = data.Count - 1;
                }
                if (lowerIndex < 0)
                {
                    lowerIndex = 0;
                }
                return new double[] { data[lowerIndex], data[upperIndex] };


            }

            else
            {
                // 使用循环计算数据的最小值和最大值
                double minValue = double.MaxValue;
                double maxValue = double.MinValue;
                foreach (double value in data)
                {
                    if (value < minValue)
                    {
                        minValue = value;
                    }
                    if (value > maxValue)
                    {
                        maxValue = value;
                    }
                }

                // 使用循环计算样本均值
                double sum = 0;
                foreach (double value in data)
                {
                    sum += value;
                }
                double mean = sum / data.Count;

                // 使用循环计算标准差
                double sumOfSquares = 0;
                foreach (double value in data)
                {
                    sumOfSquares += Math.Pow(value - mean, 2);
                }
                double stdDev = Math.Sqrt(sumOfSquares / (data.Count - 1));

                double zValue = 0;
                double alpha = 1 - confidenceLevel / 100.0;

                double zValueOneTail = GetZValueOneTail(alpha / 2);
                double marginOfError = zValueOneTail * stdDev;
                double lowerBound = mean - marginOfError;
                double upperBound = mean + marginOfError;

                return new double[] { lowerBound, upperBound };
            }



        }
        static double[,] zTable = {
      {0.0000, 0.0039, 0.0078, 0.0117, 0.0156, 0.0195, 0.0234, 0.0274, 0.0314, 0.0353},
      {0.0392, 0.0432, 0.0471, 0.0510, 0.0550, 0.0589, 0.0628, 0.0668, 0.0708, 0.0749},
      {0.0789, 0.0830, 0.0871, 0.0912, 0.0953, 0.0994, 0.1035, 0.1076, 0.1117, 0.1158},
      {0.1190, 0.1231, 0.1272, 0.1314, 0.1357, 0.1398, 0.1441, 0.1483, 0.1525, 0.1568},
      {0.1611, 0.1654, 0.1697, 0.1740, 0.1783, 0.1826, 0.1869, 0.1912, 0.1955, 0.1999},
      {0.2042, 0.2085, 0.2129, 0.2172, 0.2216, 0.2260, 0.2304, 0.2348, 0.2392, 0.2436},
      {0.2480, 0.2524, 0.2568, 0.2612, 0.2658, 0.2703, 0.2747, 0.2793, 0.2838, 0.2883},
      {0.2929, 0.2974, 0.3020, 0.3066, 0.3112, 0.3159, 0.3205, 0.3251, 0.3299, 0.3346},
      {0.3393, 0.3441, 0.3487, 0.3535, 0.3583, 0.3632, 0.3679, 0.3727, 0.3776, 0.3823},
      {0.3872, 0.3920, 0.3968, 0.4017, 0.4066, 0.4115, 0.4164, 0.4214, 0.4264, 0.4313}
    };
        static double GetZValueOneTail(double alpha)
        {
            if (alpha == 0.5)
            {
                return 0;
            }
            else if (alpha > 0.5)
            {
                // Use symmetry property of the normal distribution
                return -GetZValueOneTail(1 - alpha);
            }
            double[] alphaValues = { 0.1, 0.05, 0.025, 0.01, 0.005, 0.001 };
            double[] zValues = { 1.28, 1.645, 1.96, 2.33, 2.575, 3.09 };
            int n = alphaValues.Length;
            int i = 0;
            while (i < n && alpha < alphaValues[i])
            {
                i++;
            }

            if (i == 0)
            {
                return 0;
            }
            else if (i == n)
            {
                return zValues[n - 1];
            }

            double x0 = alphaValues[i - 1];
            double x1 = alphaValues[i];
            double y0 = zValues[i - 1];
            double y1 = zValues[i];

            double zValue = ((y1 - y0) / (x1 - x0)) * (alpha - x0) + y0;

            return zValue;
        }

        //数据频率分布图的英文怎么说,
        /// <summary>
        /// 获取数据频率的分布图
        /// </summary>
        /// <param name="data"></param>
        /// <param name="numIntervals"></param>
        /// <param name="filterNum">以最小值作为过滤的系数</param>
        /// <returns></returns>
        public static List<DataPoint> GetFrequencyDistribution(List<Result> data, out List<DataPoint>[] datapoints, int numIntervals = 10,double filterNum=5 )
        {

            List<DataPoint> resultPoints = new List<DataPoint>();
            datapoints= new List<DataPoint>[numIntervals];
            for(int i=0;i<numIntervals;i++)
            {
                datapoints[i]= new List<DataPoint>();
            }
            // 设置置信水平


            // 设置区间数量


            ////在data中过滤和最大的5%的数据
            //data = data.Where(d => d <= data.Min() * filterNum).ToList();



            // 计算数据的最小值和最大值
            double minValue = data.Min(x=>x.ObjFunctionValue);
            double maxValue = data.Max(x => x.ObjFunctionValue);

            // 计算区间宽度
            double intervalWidth = (maxValue - minValue) / numIntervals;

            // 初始化频率分布数组
            int[] freqDist = new int[numIntervals];

            // 计算各个区间的频率
            foreach (var xd in data)
            {
                var x = xd.ObjFunctionValue;
                int index = (int)((x - minValue) / intervalWidth);
                if (index < 0) // 数据小于最小值
                {
                    freqDist[0]++;
                    datapoints[0].Add(new DataPoint() {XValue= x,Tag=xd });
                }
                else if (index >= numIntervals) // 数据大于最大值
                {
                    freqDist[numIntervals - 1]++;
                    datapoints[numIntervals - 1].Add(new DataPoint() { XValue = x, Tag = xd });
                }
                else // 数据在区间内
                {
                    freqDist[index]++;
                    datapoints[index].Add(new DataPoint() { XValue = x, Tag = xd });
                }
            }
            double[] freqDistPersent =new double[numIntervals];
            double sum = freqDist.Sum();
            //计算freqDist的频率百分比
            for (int i = 0; i < numIntervals; i++)
            {
                freqDistPersent[i] = (freqDist[i] / (double)sum );
            }

            //// 计算样本均值和标准差
            //double mean = data.Average();
            //double stdDev = Math.Sqrt(data.Select(x => (x - mean) * (x - mean)).Sum() / (data.Length - 1));

            //double zValue = 0;
            //double alpha = 1 - confidenceLevel / 100.0;

            //double zValueOneTail = GetZValueOneTail(alpha / 2);

            ////zValue = mean + zValueOneTail * stdDev / Math.Sqrt(data.Length);
            //// 计算置信区间的范围
            ////double z = GetZValue(confidenceLevel); //chart1.DataManipulator.Statistics.InverseNormalDistribution(1 - confidenceLevel / 200.0); // 根据置信水平查找标准正态分布的临界值
            //double marginOfError = zValueOneTail * stdDev;/// Math.Sqrt(data.Length);
            //double lowerBound = mean - marginOfError;
            //double upperBound = mean + marginOfError;
            for (int i = 0; i < numIntervals; i++)
            {
                double xValue = minValue + i * intervalWidth + intervalWidth / 2;
                double yValue = freqDistPersent[i];
                if (isValid(xValue) && isValid(yValue))
                    resultPoints.Add(new DataPoint(xValue, yValue));
            }


            return resultPoints;



        }
        
        static bool isValid(double value)
        {
            if (value == (int)value && value < 0) return false;
            return !Double.IsNaN(value) && !Double.IsInfinity(value);
        }
       

    }

    //public static class ConfidenceRangeHelper
    //{


    //}
}