Как выполнить LU-декомпозицию с OpenCV? - PullRequest
Новая
       17

Как выполнить LU-декомпозицию с OpenCV?

1 голос
/ 22 октября 2011

Метод cvInvert () принимает флаг CV_LU, который выполняет факторизацию LU для инвертирования входной матрицы. Однако есть ли способ получить матрицы L и U, которые формируются во время этого вычисления? Кажется бессмысленным писать новую функцию для LU-декомпозиции, если OpenCV уже оптимизировал код для этого.

Ответы [ 2 ]

2 голосов
/ 11 июля 2013

Можно использовать функцию Cholesky(), предоставляемую OpenCV (используя 2.4.6), см. Исходный код modules\stitching\src\autocalib.cpp. Но для получения сопоставимого результата с Matlab требуется некоторое масштабирование: 

Mat chol = mat.clone();
if (Cholesky(chol.ptr<float>(), chol.step, chol.cols, 0, 0, 0))
{
    Mat diagElem = chol.diag();
    for (int e = 0; e < diagElem.rows; ++e)
    {
        float elem = diagElem.at<float>(e);
        chol.row(e) *= elem;
        chol.at<float>(e,e) = 1.0f / elem;
    }
}
2 голосов
/ 22 октября 2011

К сожалению, не похоже, что OpenCV предоставляет вам доступ к матрицам L и U. Здесь - как реализована функция. И, похоже, по причинам производительности LU-разложение выполняется на месте. Так что вам, вероятно, придется сделать это самостоятельно.

РЕДАКТИРОВАТЬ: Похоже, что после анализа того, как Matlab и Eigen выполняют LU-декомпозицию, вы можете получить их после вызова cvInvert. Матрица L - это строго матрица нижнего треугольника результата плюс матрица идентичности, а матрица U - матрица верхнего треугольника.

РЕДАКТИРОВАТЬ: Eigen на самом деле довольно хорошо интегрируется с OpenCV. И, похоже, у них реализован класс разложения LU здесь . Eigen уже является зависимостью для OpenCV, если вы создали его самостоятельно, он должен быть доступен для использования (он полностью реализован в заголовочных файлах, что делает его действительно простым в использовании). Существует также заголовок OpenCV, реализующий преобразование между собственными матрицами и матрицами OpenCV @ #include <opencv2/core/eigen.hpp>.

Однако, по крайней мере, в моей сборке SVN этот заголовок работал неправильно, поэтому я сделал свой собственный: 

#ifndef __OPENCV_CORE_EIGEN_HPP__
#define __OPENCV_CORE_EIGEN_HPP__

#ifdef __cplusplus

#include "opencv/cxcore.h"
#include <eigen3/Eigen/Dense>

namespace cv
{

template<typename _Tp, int _rows, int _cols, int _options, int _maxRows, int _maxCols>
void eigen2cv( const Eigen::Matrix<_Tp, _rows, _cols, _options, _maxRows, _maxCols>& src, Mat& dst )
{
    if( !(src.Flags & Eigen::RowMajorBit) )
    {
        Mat _src(src.cols(), src.rows(), DataType<_Tp>::type,
              (void*)src.data(), src.stride()*sizeof(_Tp));
        transpose(_src, dst);
    }
    else
    {
        Mat _src(src.rows(), src.cols(), DataType<_Tp>::type,
                 (void*)src.data(), src.stride()*sizeof(_Tp));
        _src.copyTo(dst);
    }
}

template<typename _Tp, int _rows, int _cols, int _options, int _maxRows, int _maxCols>
void cv2eigen( const Mat& src,
               Eigen::Matrix<_Tp, _rows, _cols, _options, _maxRows, _maxCols>& dst )
{
    CV_DbgAssert(src.rows == _rows && src.cols == _cols);
    if( !(dst.Flags & Eigen::RowMajorBit) )
    {
        Mat _dst(src.cols, src.rows, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        if( src.type() == _dst.type() )
            transpose(src, _dst);
        else if( src.cols == src.rows )
        {
            src.convertTo(_dst, _dst.type());
            transpose(_dst, _dst);
        }
        else
            Mat(src.t()).convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
    else
    {
        Mat _dst(src.rows, src.cols, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        src.convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
}

template<typename _Tp>
void cv2eigen( const Mat& src,
               Eigen::Matrix<_Tp, Eigen::Dynamic, Eigen::Dynamic>& dst )
{
    dst.resize(src.rows, src.cols);
    if( !(dst.Flags & Eigen::RowMajorBit) )
    {
        Mat _dst(src.cols, src.rows, DataType<_Tp>::type,
             dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        if( src.type() == _dst.type() )
            transpose(src, _dst);
        else if( src.cols == src.rows )
        {
            src.convertTo(_dst, _dst.type());
            transpose(_dst, _dst);
        }
        else
            Mat(src.t()).convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
    else
    {
        Mat _dst(src.rows, src.cols, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        src.convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
}


template<typename _Tp>
void cv2eigen( const Mat& src,
               Eigen::Matrix<_Tp, Eigen::Dynamic, 1>& dst )
{
    CV_Assert(src.cols == 1);
    dst.resize(src.rows);

    if( !(dst.Flags & Eigen::RowMajorBit) )
    {
        Mat _dst(src.cols, src.rows, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        if( src.type() == _dst.type() )
            transpose(src, _dst);
        else
            Mat(src.t()).convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
    else
    {
        Mat _dst(src.rows, src.cols, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        src.convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
}


template<typename _Tp>
void cv2eigen( const Mat& src,
               Eigen::Matrix<_Tp, 1, Eigen::Dynamic>& dst )
{
    CV_Assert(src.rows == 1);
    dst.resize(src.cols);
    if( !(dst.Flags & Eigen::RowMajorBit) )
    {
        Mat _dst(src.cols, src.rows, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        if( src.type() == _dst.type() )
            transpose(src, _dst);
        else
            Mat(src.t()).convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
    else
    {
        Mat _dst(src.rows, src.cols, DataType<_Tp>::type,
                 dst.data(), (size_t)(dst.stride()*sizeof(_Tp)));
        src.convertTo(_dst, _dst.type());
        CV_DbgAssert(_dst.data == (uchar*)dst.data());
    }
}

}

#endif

#endif

Надеюсь, это полезно для вас!

...