public Matrix3d(int NumberOfRows, int NumberOfColumns, int NumberOfLayers)
 {
   _data = new DenseMatrix[NumberOfLayers];
   
   for (int i = 0; i < NumberOfLayers; i++)
     _data[i] = new DenseMatrix(NumberOfRows, NumberOfColumns);
 }

        public static DenseMatrix ToMatrix(this object[,] data, int rowStart, int rowEnd, int colStart, int colEnd,
            bool reverseColumns)
        {
            var d = new DenseMatrix(rowEnd - rowStart + 1, colEnd - colStart + 1);

            if (reverseColumns)
            {
                for (int i = rowEnd; i >= rowStart; i--)
                {
                    for (int j = colStart; j <= colEnd; j++)
                    {
                        d[rowEnd - i, j - colStart] = (double) data[i, j];
                    }
                }
            }
            else
            {
                for (int i = rowStart; i <= rowEnd; i++)
                {
                    for (int j = colStart; j <= colEnd; j++)
                    {
                        d[i - rowStart, j - colStart] = (double) data[i, j];
                    }
                }
            }

            return d;
        }

 public SalesmanProblem(DenseMatrix costs, double startRecord, IEnumerable<int> startChain)
     : this(costs)
 {
     Record = startRecord;
     RecordChain = new List<int>(startChain);
     AnswerFound = true;
 }

        /// <summary>
        /// Creates a new Input from double values.
        /// Note that we store each example as a row in the X matrix. While calculating Theta vector, we need to insert the top column of all ones into the X matrix - this will allow us to treat theta0 as just another feature.
        /// </summary>
        internal Input(double[,] x, double[] y, int skip, int take)
        {
            if (take == 0) {
                X = null;
                Y = null;
                return;
            }

            var samples = x.GetLength(0);
            var features = x.GetLength(1);

            //make sure we add first column of ones
            var x1 = new double[take, features + 1];
            var y1 = new double[take];

            for (int sample = 0; sample < samples; sample++) {
                if (sample < skip) {
                    continue;
                }
                for (int feature = 0; feature < features + 1; feature++) {
                    x1[sample - skip, feature] = (feature == 0) ? 1 : x[sample, feature - 1];
                }
                y1[sample - skip] = y[sample];

                take--;
                if (take == 0) {
                    break;
                }
            }

            X = new DenseMatrix(x1);
            Y = new DenseVector(y1).ToColumnMatrix();
        }

        public void Factorize(Matrix meanMatrix)
        {
            int count = 0;
            movieCount = meanMatrix.RowCount;
            userCount = meanMatrix.ColumnCount;
            convergedMovie = DenseMatrix.Build.Dense(K, movieCount, 0);
            convergedUser = DenseMatrix.Build.Dense(userCount, K, 0);

            A = new DenseMatrix(DenseColumnMajorMatrixStorage<double>.OfMatrix(DenseMatrix.Build.Dense(K, movieCount, 0.1).Storage));
            B = new DenseMatrix(DenseColumnMajorMatrixStorage<double>.OfMatrix(DenseMatrix.Build.Dense(userCount, K, 0.1).Storage));

            while (count < 50)
            {
                for (int k = 0; k < K; k++)
                {
                    for (int m = 0; m < movieCount; m++)
                    {
                        for (int u = 0; u < userCount; u++)
                        {
                            train(k, u, m, meanMatrix[m, u]);
                        }
                    }
                }
                Console.WriteLine(count++);
                if (IsConverged())
                    break;
            }
        }

        public void Compute(DenseMatrix x1, int size)
        {
            //Setting Preparation.....................................................................................................
            _estimationLength = x1.Values.Length;

            //Creation of a Z matrix..................................................................................................
            // ReSharper disable once CSharpWarnings::CS0618
            var z = new DenseMatrix(_estimationLength, size, 1.0);
            for (var i = 0; i <= _estimationLength - 1; i++)
            {
                for (var j = 0; j <= size - 1; j++)
                {
                    z[i, j] = Math.Pow(i,j);
                }
            }

            //Computation of Theta matrix.............................................................................................
            var zt = z.Transpose();
            var ztz = zt * z;
            var theta = ztz.Inverse() * zt * x1;

            // ReSharper disable once CSharpWarnings::CS0618
            _output = new DenseMatrix(_estimationLength, 1, 1.0);

            //Output creation.........................................................................................................
            for (var i = 0; i <= _estimationLength - 1; i++)
            {
               for (var j = 0; j <= size - 1; j++)
               {
                   _output[i, 0] += theta[j,0]*Math.Pow(i,j);
               }
            }
            _estimationDone = true;
        }

        /// <summary>
        /// 计算模糊矩阵
        /// </summary>
        private void CalculateFuzzyMatrix()
        {
            MemberShipFun _memeberShipFun = new MemberShipFun();
            for (int i = AHPIndexHierarchyUtil.totalLevelCount - 2; i >= 0; i--)
            {
                List<AHPIndexHierarchy> iLevelAhpIndexs = _ahpIndexUtil.FindbyLevel(i);
                foreach (AHPIndexHierarchy _iLevelAhpIndex in iLevelAhpIndexs)
                {
                    List<string> childrenNames = _iLevelAhpIndex.ChildrenNames;
                    DenseMatrix _iLevelMatrix = new DenseMatrix(childrenNames.Count, MemberShipFun.HealthLevelCount);
                    DenseVector _childrenValue = new DenseVector(childrenNames.Count);
                    for (int j = 0; j < childrenNames.Count; j++)
                    {
                        string name = childrenNames[j];
                        AHPIndexHierarchy _ahpIndex = _ahpIndexUtil.FindbyName(name);
                        if (i == AHPIndexHierarchyUtil.totalLevelCount - 2)//是底层
                        {
                            _ahpIndex.FuzzyValue = _memeberShipFun.TrapezoiMebership(_ahpIndex.IndexValue);
                        }

                        _iLevelMatrix = (DenseMatrix)_iLevelMatrix.InsertRow(j, _ahpIndex.FuzzyValue);
                        _iLevelMatrix = (DenseMatrix)_iLevelMatrix.RemoveRow(j + 1);
                        //_ahpIndex.ChildrenFuzzyMatrix = (DenseMatrix)_iLevelMatrix;
                        _childrenValue[j] = _ahpIndex.IndexValue;
                    }
                    _iLevelAhpIndex.ChildrenFuzzyMatrix = _iLevelMatrix;
                    _iLevelAhpIndex.IndexValue = _iLevelAhpIndex.ChildrenWeightVector * _childrenValue;
                    _iLevelAhpIndex.FuzzyValue = FuzzyOperator.WeightedAverage(_iLevelAhpIndex.ChildrenWeightVector, _iLevelAhpIndex.ChildrenFuzzyMatrix);
                }
            }
        }

        public void Optimize(Dictionary<double, double> values)
        {
            var n = _f.Functions.Count();
              var xs = values.Select(v => v.Key).ToList();
              var ys = values.Select(v => v.Value).ToList();
              var fs = new List<List<double>>(n);

              for (var i = 0; i < n; i++)
              {
            fs[i] = _f.Functions[i].Evaluate(xs);
              }

              var matrix = new DenseMatrix(n, n);
              var vector = new DenseVector(n);
              for (var i = 0; i < n; i++)
              {
            for (var j = 0; j < n; j++)
            {
              matrix[i, j] = fs[i].ScalarProduct(fs[j]);
            }
            vector[i] = ys.ScalarProduct(fs[i]);
              }

              var matrixInverse = matrix.Inverse();

              var result = matrixInverse * vector;

              for (var i = 0; i < n; i++)
              {
            _f.LinearParameters[i].Value = result[i];
              }
        }

    public void Load(string FileName)
    {
      using (StreamReader sr = new StreamReader(FileName))
      {
        string line = sr.ReadLine();
        NumberOfColumns = int.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);
        line = sr.ReadLine();
        NumberOfRows = int.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);
        line = sr.ReadLine();
        XOrigin = double.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);
        line = sr.ReadLine();
        YOrigin = double.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);
        line = sr.ReadLine();
        GridSize = double.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);
        line = sr.ReadLine();
        DeleteValue = double.Parse(line.Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries)[1]);

        Data = new DenseMatrix(NumberOfRows, NumberOfColumns);


        string[] DataRead;
        for (int j = 0; j < NumberOfRows; j++)
        {
          DataRead = sr.ReadLine().Split(new Char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries);
          for (int i = 0; i < NumberOfColumns; i++)
          {
            Data[NumberOfRows - j - 1, i] = double.Parse(DataRead[i]);
          }
        }
      }
    }

        // Returns normalisation matrix : xn = Mx
        // Uses list of points : each point is 3-vector
        // Assusmes that weight of each point is 1
        public static Matrix<double> Normalize2D(List<Vector<double>> points)
        {
            Matrix<double> norm = new DenseMatrix(3, 3);
            int n = points.Count;
            // Compute center of image points
            double xc = 0, yc = 0;
            for(int c = 0; c < n; ++c)
            {
                xc += points[c].At(0);
                yc += points[c].At(1);
            }
            xc /= n;
            yc /= n;
            // Get mean distance of points from center
            double dist = 0;
            for(int c = 0; c < n; ++c)
            {
                dist += Math.Sqrt((points[c].At(0) - xc) * (points[c].At(0) - xc) +
                    (points[c].At(1) - yc) * (points[c].At(1) - yc));
            }
            dist /= n;
            // Normalize in a way that mean dist = sqrt(2)
            double ratio = Math.Sqrt(2) / dist;
            // Noramlize matrix - homonogeus point must be multiplied by it
            norm[0, 0] = ratio;
            norm[1, 1] = ratio;
            norm[0, 2] = -ratio * xc;
            norm[1, 2] = -ratio * yc;
            norm[2, 2] = 1.0;

            return norm;
        }

        // Return convolution A * B
        // Size of returned matrix = size of A, but
        // elements on boundaries ( of length rows(cols) of B / 2 )
        // are not computed and set to 0
        // Size of B must be odd
        public static Matrix Convolve(Matrix A, Matrix B)
        {
            Matrix conv = new DenseMatrix(A.RowCount, A.ColumnCount);

            int row2 = B.RowCount / 2;
            int col2 = B.ColumnCount / 2;
            int xmax = A.ColumnCount - col2;
            int ymax = A.RowCount - row2;

            int y, x, dx, dy;
            double maskSum;

            for ( y = row2; y < ymax; ++y)
                for ( x = col2; x < xmax; ++x)
                {
                    maskSum = 0.0f;
                    for ( dy = -row2; dy <= row2; ++dy)
                        for ( dx = -col2; dx <= col2; ++dx)
                        {
                            maskSum += A[y + dy, x + dx] * B[row2 + dy, col2 + dx];
                        }
                    conv[y, x] = maskSum;
                }

            return conv;
        }

        public static void WriteDataToJson(string[] symbols, DateTime[] dateTimes, DenseMatrix data, string filename)
        {
            StringBuilder jsonStringBuilder = new StringBuilder();

            jsonStringBuilder.Append("[");
            for (int i = 0; i < data.ColumnCount; i++)
            {
                jsonStringBuilder.Append("[");
                jsonStringBuilder.Append(dateTimes[i].Subtract(new DateTime(1970, 1, 1, 0, 0, 0, DateTimeKind.Utc)).TotalMilliseconds.ToString());
                foreach (double d in data.Column(i))
                {
                    jsonStringBuilder.Append(",");
                    jsonStringBuilder.Append((d.Equals(double.NaN) ? "null" : Math.Round(d, 8).ToString()));
                }
                jsonStringBuilder.Append("]");
                if (i != data.ColumnCount - 1) jsonStringBuilder.Append(",\n");
            }
            jsonStringBuilder.Append("]");

            using (var file = new StreamWriter(QSConstants.DEFAULT_DATA_FILEPATH + @filename))
            {
                file.WriteLine(jsonStringBuilder.ToString());
            }

        }

        public static Matrix<double> CreateDCTMatrix(int n)
        {
            if (dctMatrixes.ContainsKey(n))
            {
                return dctMatrixes[n];
            }

            var matrix = new DenseMatrix(n);
            var val = 1d / Math.Sqrt(n);
            for (int i = 0; i < n; i++)
            {
                matrix[0, i] = val;
            }

            var sqrt2 = Math.Sqrt(2d / n);
            for (int x = 0; x < n; x++)
            {
                for (int y = 1; y < n; y++)
                {
                    matrix[x, y] = sqrt2 * Math.Cos((Math.PI / 2 / n) * y * (2 * x + 1));
                }
            }

            dctMatrixes[n] = matrix;
            return matrix;
        }

        DenseMatrix[] sigmas; //covariance of the 3D gaussians.

        #endregion Fields

        #region Constructors

        public CD_HMM(MarkovChain A, DenseVector pi, DenseVector[] mus, DenseMatrix[] sigmas)
        {
            this.A = A;
            this.pi = pi;
            this.mus = mus;
            this.sigmas = sigmas;
        }

 public void FindTransformation(ArrayList kinectCoors, ArrayList projectorCoors)
 {
     PrepareMatrices(kinectCoors, projectorCoors);
     DenseMatrix problem = foundCoordinatesMatrix;
     result = (DenseMatrix) problem.QR().Solve(rightSideMatrix);
     Console.Out.WriteLine(result);
 }

        public MainWindow()
        {
            InitializeComponent();
            Graph = new UserControl1();
            this.Content = Graph;

            var t = Generate();

            //    var s = t.Item1.Svd(true);

            //var nd = t.Item1.Multiply(s.VT().SubMatrix(0, 8, 0, 8));

            var mlp = new MLP(8, 100, 1);

            var r = mlp.Train(t.Item1, t.Item2, null, null, 2500);

            Graph.Set(r.TrainingSquaredError, r.TrainingError);

            DenseMatrix data = new DenseMatrix(new double[,] { { 1, 1 }, { 1, 0 }, { 0, 0 }, { 0, 1 } });
            DenseMatrix labels = new DenseMatrix(new double[,] { { 1, 0 }, { 0, 1 }, { 1, 0 }, { 0, 1 } });

            //var mlp = new MLP(2, 100, 2);

            //var r = mlp.Train(data, labels, data, labels, 5000);

            //Graph.Set(r.TrainingSquaredError, r.TrainingError);
        }

 private static void GetLinRegArgs(IEnumerable<double[]> items, out DenseMatrix xM, out DenseVector yV)
 {
     var original = items;
       var transformed = original.Select(a => new Tuple<double[], double>(Phi(a), a[2])).ToList();
       xM = DenseMatrix.Create(transformed.Count, transformed[0].Item1.Length, (r, c) => transformed[r].Item1[c]);
       yV = DenseVector.OfEnumerable(transformed.Select(t => t.Item2));
 }

        public static Matrix<double> CorrelationMatrix(List<List<double>> x)
        {
            Matrix<double> correlM = new DenseMatrix(x.Count, x.Count);

            //off diagonal elements
            for (int i = 0; i < x.Count; i++)
            {
                var rowx = x[i];
                for (int j = 0; j < x.Count; j++)
                {
                    var rowy = x[j];
                    if (i < j)
                    {
                        correlM[i, j] = Correlation.Pearson(rowx, rowy);
                    }
                    if (i > j)
                    {
                        correlM[i, j] = correlM[j, i];
                    }
                }
            }
            //Diagonal elements
            for (int i = 0; i < x.Count; i++)
            {
                correlM[i, i] = 1;
            }

            return correlM;
        }

 public DenseMatrix ToMatrix()
 {
     var res = new DenseMatrix(3, 1);
     res[0, 0] = this.X;
     res[1, 0] = this.Y;
     res[2, 0] = this.Z;
     return res;
 }

 /// <summary>
 /// ctor with params
 /// </summary>
 /// <param name="N">number of rows of each matrix, point count</param>
 /// <param name="L">number of cols of each matrix, it is height of the tree</param>
 public BasicFuncBoxes(int N,int L)
 {
     n = N;
     l = L;
     x = new DenseMatrix(n, l);
     y = new DenseMatrix(n, l);
     z = new DenseMatrix(n, l);
 }