C# Matrix.At方法代码示例

 // Converts 3x3 rotation matrix to XYZ euler angles (assumes vector is correctly allocated)
 public static void MatrixToEuler(Vector<double> euler, Matrix<double> matrix)
 {
     if(matrix.At(0, 2) < 1.0 - 1e-9)
     {
         if(matrix.At(0, 2) > -1.0 + 1e-9)
         {
             //      t he t aY = a s i n(r 0 2);
             //     t he t aX = a t a n 2(−r12, r 2 2);
             //    t h e t aZ = a t a n 2(−r01, r 0 0);
             euler.At(1, Math.Asin(matrix.At(0, 2)));
             euler.At(0, Math.Atan2(-matrix.At(1, 2), matrix.At(2, 2)));
             euler.At(2, Math.Atan2(-matrix.At(0, 1), matrix.At(0, 0)));
         }
         else // r 0 2 = −1
         {
             // Not a u n i q u e s o l u t i o n : t h e t aZ − t he t aX = a t a n 2 ( r10 , r 1 1 )
             //  t he t aY = −PI / 2;
             // t he t aX = −a t a n 2(r10, r 1 1);
             // t h e t aZ = 0;
             euler.At(1, -Math.PI * 0.5);
             euler.At(0, -Math.Atan2(-matrix.At(1, 0), matrix.At(1, 1)));
             euler.At(2, 0.0);
         }
     }
     else // r 0 2 = +1
     {
         // Not a u n i q u e s o l u t i o n : t h e t aZ + t he t aX = a t a n 2 ( r10 , r 1 1 )
         // t he t aY = +PI / 2;
         // t he t aX = a t a n 2(r10, r 1 1);
         // t h e t aZ = 0;
         euler.At(1, Math.PI * 0.5);
         euler.At(0, -Math.Atan2(-matrix.At(1, 0), matrix.At(1, 1)));
         euler.At(2, 0.0);
     }
 }

 public override void DisturbMatrix(Matrix<double> matToBeDisturbed)
 {
     double[] samples = new double[matToBeDisturbed.RowCount * matToBeDisturbed.ColumnCount];
     _gauss.Samples(samples);
     int i = 0;
     for(int c = 0; c < matToBeDisturbed.ColumnCount; ++c)
     {
         for(int r = 0; r < matToBeDisturbed.RowCount; ++r)
         {
             matToBeDisturbed.At(r, c,
                matToBeDisturbed.At(r, c) + samples[i]);
             ++i;
         }
     }
 }

        // Computes coefficients of real polynomial (or estimates if values are noised) using Svd
        // Each row of matrix should contain [x, P(x)], at least rank+1 rows
        // Supplied x-es best have magintude in range [1-2], so resulting coefficient matrix is well conditioned
        public static Polynomial EstimatePolynomial(Matrix<float> values, int rank)
        {
            // 1) Create equation Xa = b
            // | x1^n x1^n-1 ... x1 1 | | a0 |   | P(x1) |
            // |                      | | ...| = |       |
            // | xk^n xk^n-1 ... xk 1 | | an |   | P(xk) |
            Matrix<float> X = new DenseMatrix(values.RowCount, rank + 1);
            Vector<float> P = new DenseVector(values.RowCount);

            for(int i = 0; i < values.RowCount; ++i)
            {
                X[i, rank] = 1.0f;
                for(int c = rank - 1; c >= 0; --c)
                {
                    X[i, c] = X[i, c + 1] * values.At(i, 0);
                }
                P[i] = values[i, 1];
            }

            return new Polynomial()
            {
                Coefficents = SvdSolver.Solve(X, P),
                Rank = rank
            };
        }

        // Returns normalisation matrix : xn = Mx
        // Uses list of points in matrix : each point is column 3-vector (so xi is [0,i])
        // Assusmes that weight of each point is 1
        public static Matrix<double> Normalize2D(Matrix<double> points)
        {
            Matrix<double> norm = new DenseMatrix(3, 3);
            int n = points.ColumnCount;
            // Compute center of image points
            double xc = 0, yc = 0;
            for(int c = 0; c < n; ++c)
            {
                xc += points.At(0, c);
                yc += points.At(1, c);
            }
            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.At(0, c) - xc) * (points.At(0, c) - xc) +
                    (points.At(1, c) - yc) * (points.At(1, c) - 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;
        }

        public static float predictionQuality(Matrix<float> matrix)
        {
            int batchSize = matrix.RowCount;
            int numOfPositive = batchSize / 2;

            int corretClassified = 0;
            int cc = matrix.ColumnCount - 1;

            int i = 0;

            for (; i < numOfPositive; ++i)
            {
                if (matrix.At(i, cc) > 0.5f) ++corretClassified;
                Console.WriteLine("label " + i + ": " + matrix.At(i, cc));
            }

            for (; i < batchSize; ++i)
            {
                if (matrix.At(i, cc) < 0.5f) ++corretClassified;
                Console.WriteLine("label " + i + ": " + matrix.At(i, cc));
            }

            return (float)(corretClassified) / (float)(batchSize);
        }

        public static float reconstructionError(Matrix<float> original, Matrix<float> reconstruction)
        {
            float error = 0.0f;

            for (int row = 0; row < original.RowCount; ++row)
            {
                for (int column = 0; column < original.ColumnCount; ++column)
                {
                    error += Math.Abs(original.At(row, column) - reconstruction.At(row, column));
                }
            }

            error /= original.RowCount * original.ColumnCount;

            return error;
        }

 // Converts XYZ euler angles to 3x3 rotation matrix (assumes matrix is correctly allocated)
 public static void EulerToMatrix(double[] euler, Matrix<double> matrix)
 {
     // | cycz           −cysz           sy    |
     // | czsxsy + cxsz  cxcz − sxsysz   −cysx |
     // | −cxczsy + sxsz czsx + cxsysz   cxcy  |
     //
     double sx = Math.Sin(euler[0]);
     double cx = Math.Cos(euler[0]);
     double sy = Math.Sin(euler[1]);
     double cy = Math.Cos(euler[1]);
     double sz = Math.Sin(euler[2]);
     double cz = Math.Cos(euler[2]);
     matrix.At(0, 0, cy * cz);
     matrix.At(0, 1, -cy * sz);
     matrix.At(0, 2, sy);
     matrix.At(1, 0, cz * sx * sy + cx * sz);
     matrix.At(1, 1, cx * cz - sx * sy * sz);
     matrix.At(1, 2, -cy * sx);
     matrix.At(2, 0, -cx * cz * sy + sz * sz);
     matrix.At(2, 1, cz * sx + cx * sy * sz);
     matrix.At(2, 2, cx * cy);
 }

        /// <summary>
        /// Multiplies this matrix with another matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        protected override void DoMultiply(Matrix<double> other, Matrix<double> result)
        {
            var denseOther = other as DenseMatrix;
            var denseResult = result as DenseMatrix;
            if (denseOther != null && denseResult != null)
            {
                Control.LinearAlgebraProvider.MatrixMultiply(
                    _values,
                    _rowCount,
                    _columnCount,
                    denseOther._values,
                    denseOther._rowCount,
                    denseOther._columnCount,
                    denseResult._values);
                return;
            }

            var diagonalOther = other.Storage as DiagonalMatrixStorage<double>;
            if (diagonalOther != null)
            {
                var diagonal = diagonalOther.Data;
                var d = Math.Min(ColumnCount, other.ColumnCount);
                if (d < other.ColumnCount)
                {
                    result.ClearSubMatrix(0, RowCount, ColumnCount, other.ColumnCount - ColumnCount);
                }
                int index = 0;
                for (int j = 0; j < d; j++)
                {
                    for (int i = 0; i < RowCount; i++)
                    {
                        result.At(i, j, _values[index]*diagonal[j]);
                        index++;
                    }
                }
                return;
            }

            base.DoMultiply(other, result);
        }

        /// <summary>
        /// Subtracts another matrix from this matrix.
        /// </summary>
        /// <param name="other">The matrix to subtract.</param>
        /// <param name="result">The matrix to store the result of the subtraction.</param>
        protected override void DoSubtract(Matrix<double> other, Matrix<double> result)
        {
            // dense + dense = dense
            var denseOther = other.Storage as DenseColumnMajorMatrixStorage<double>;
            var denseResult = result.Storage as DenseColumnMajorMatrixStorage<double>;
            if (denseOther != null && denseResult != null)
            {
                Control.LinearAlgebraProvider.SubtractArrays(_values, denseOther.Data, denseResult.Data);
                return;
            }

            // dense + diagonal = matrix
            var diagonalOther = other.Storage as DiagonalMatrixStorage<double>;
            if (diagonalOther != null)
            {
                CopyTo(result);
                var diagonal = diagonalOther.Data;
                for (int i = 0; i < diagonal.Length; i++)
                {
                    result.At(i, i, result.At(i, i) - diagonal[i]);
                }
                return;
            }

            base.DoSubtract(other, result);
        }

        /// <summary>
        /// Adds another matrix to this matrix.
        /// </summary>
        /// <param name="other">The matrix to add to this matrix.</param>
        /// <param name="result">The matrix to store the result of the addition.</param>
        /// <exception cref="ArgumentNullException">If the other matrix is <see langword="null"/>.</exception>
        /// <exception cref="ArgumentOutOfRangeException">If the two matrices don't have the same dimensions.</exception>
        protected override void DoAdd(Matrix<double> other, Matrix<double> result)
        {
            var sparseOther = other as SparseMatrix;
            var sparseResult = result as SparseMatrix;
            if (sparseOther == null || sparseResult == null)
            {
                base.DoAdd(other, result);
                return;
            }

            if (ReferenceEquals(this, other))
            {
                if (!ReferenceEquals(this, result))
                {
                    CopyTo(result);
                }

                Control.LinearAlgebraProvider.ScaleArray(2.0, _storage.Values, _storage.Values);
                return;
            }

            SparseMatrix left;

            if (ReferenceEquals(sparseOther, sparseResult))
            {
                left = this;
            }
            else if (ReferenceEquals(this, sparseResult))
            {
                left = sparseOther;
            }
            else
            {
                CopyTo(sparseResult);
                left = sparseOther;
            }

            var leftStorage = left._storage;
            for (var i = 0; i < leftStorage.RowCount; i++)
            {
                var endIndex = leftStorage.RowPointers[i + 1];
                for (var j = leftStorage.RowPointers[i]; j < endIndex; j++)
                {
                    var columnIndex = leftStorage.ColumnIndices[j];
                    var resVal = leftStorage.Values[j] + result.At(i, columnIndex);
                    result.At(i, columnIndex, resVal);
                }
            }
        }

        /// <summary>
        /// Multiplies this matrix with transpose of another matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        protected override void DoTransposeAndMultiply(Matrix<double> other, Matrix<double> result)
        {
            var otherSparse = other as SparseMatrix;
            var resultSparse = result as SparseMatrix;

            if (otherSparse == null || resultSparse == null)
            {
                base.DoTransposeAndMultiply(other, result);
                return;
            }

            resultSparse.Clear();

            var rowPointers = _storage.RowPointers;
            var values = _storage.Values;

            var otherStorage = otherSparse._storage;

            for (var j = 0; j < RowCount; j++)
            {
                var startIndexOther = otherStorage.RowPointers[j];
                var endIndexOther = otherStorage.RowPointers[j + 1];

                if (startIndexOther == endIndexOther)
                {
                    continue;
                }

                for (var i = 0; i < RowCount; i++)
                {
                    var startIndexThis = rowPointers[i];
                    var endIndexThis = rowPointers[i + 1];

                    if (startIndexThis == endIndexThis)
                    {
                        continue;
                    }

                    var sum = 0d;
                    for (var index = startIndexOther; index < endIndexOther; index++)
                    {
                        var ind = _storage.FindItem(i, otherStorage.ColumnIndices[index]);
                        if (ind >= 0)
                        {
                            sum += otherStorage.Values[index]*values[ind];
                        }
                    }

                    resultSparse._storage.At(i, j, sum + result.At(i, j));
                }
            }
        }

        /// <summary>
        /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to pointwise multiply with this one.</param>
        /// <param name="result">The matrix to store the result of the pointwise multiplication.</param>
        protected override void DoPointwiseMultiply(Matrix<double> other, Matrix<double> result)
        {
            result.Clear();

            var rowPointers = _storage.RowPointers;
            var columnIndices = _storage.ColumnIndices;
            var values = _storage.Values;

            for (var i = 0; i < RowCount; i++)
            {
                var endIndex = rowPointers[i + 1];
                for (var j = rowPointers[i]; j < endIndex; j++)
                {
                    var resVal = values[j]*other.At(i, columnIndices[j]);
                    if (resVal != 0d)
                    {
                        result.At(i, columnIndices[j], resVal);
                    }
                }
            }
        }

        /// <summary>
        /// Multiplies this matrix with another matrix and places the results into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to multiply with.</param>
        /// <param name="result">The result of the multiplication.</param>
        protected override void DoMultiply(Matrix<double> other, Matrix<double> result)
        {
            result.Clear();
            var columnVector = new DenseVector(other.RowCount);

            var rowPointers = _storage.RowPointers;
            var columnIndices = _storage.ColumnIndices;
            var values = _storage.Values;

            for (var row = 0; row < RowCount; row++)
            {
                var startIndex = rowPointers[row];
                var endIndex = rowPointers[row + 1];

                if (startIndex == endIndex)
                {
                    continue;
                }

                for (var column = 0; column < other.ColumnCount; column++)
                {
                    // Multiply row of matrix A on column of matrix B
                    other.Column(column, columnVector);

                    var sum = 0d;
                    for (var index = startIndex; index < endIndex; index++)
                    {
                        sum += values[index] * columnVector[columnIndices[index]];
                    }

                    result.At(row, column, sum);
                }
            }
        }

        /// <summary>
        /// Puts the strictly upper triangle of this matrix into the result matrix.
        /// </summary>
        /// <param name="result">Where to store the lower triangle.</param>
        private void StrictlyUpperTriangleImpl(Matrix<double> result)
        {
            var rowPointers = _storage.RowPointers;
            var columnIndices = _storage.ColumnIndices;
            var values = _storage.Values;
            var valueCount = _storage.ValueCount;

            for (var row = 0; row < result.RowCount; row++)
            {
                var startIndex = rowPointers[row];
                var endIndex = row < rowPointers.Length - 1 ? rowPointers[row + 1] : valueCount;
                for (var j = startIndex; j < endIndex; j++)
                {
                    if (row < columnIndices[j])
                    {
                        result.At(row, columnIndices[j], values[j]);
                    }
                }
            }
        }

        /// <summary>
        /// Pointwise divide this matrix by another matrix and stores the result into the result matrix.
        /// </summary>
        /// <param name="other">The matrix to pointwise divide this one by.</param>
        /// <param name="result">The matrix to store the result of the pointwise division.</param>
        protected override void DoPointwiseDivide(Matrix<double> other, Matrix<double> result)
        {
            result.Clear();

            var rowPointers = _storage.RowPointers;
            var columnIndices = _storage.ColumnIndices;
            var values = _storage.Values;
            var valueCount = _storage.ValueCount;

            for (var i = 0; i < RowCount; i++)
            {
                // Get the begin / end index for the current row
                var startIndex = rowPointers[i];
                var endIndex = i < rowPointers.Length - 1 ? rowPointers[i + 1] : valueCount;

                for (var j = startIndex; j < endIndex; j++)
                {
                    if (values[j] != 0d)
                    {
                        result.At(i, columnIndices[j], values[j]/other.At(i, columnIndices[j]));
                    }
                }
            }
        }