OpenCV QR decomposition

openCV的QR分解

參考這篇的

http://suzuichibolgpg.blog.fc2.com/blog-entry-258.html

http://rosettacode.org/wiki/QR_decomposition

跟維基百科
https://en.wikipedia.org/wiki/QR_decomposition#Example_2

OpenCV裡面似乎沒有直接得到QR 分解的函式，輾轉找了一下最後在一個日文網誌上找到，試著照著上面的步驟寫了個自己的，有問題再留言跟我說。記得H.rows>H.cols

#include <iostream>

#include <fstream>

#include <sstream>

#include <opencv2/core/core.hpp>

#include <opencv2/opencv.hpp>

#include <opencv2/highgui/highgui.hpp>

using namespace cv;

using namespace std;

int main(){

Mat H = (Mat_<double>(3,3) << 12,-51,4,6,167,-68,-4,24,-41);

Mat Q;

Mat R;

QRDecomp( H, Q, R );

cout << Q << "\n";

cout << R << "\n";

cout << Q*R << "\n";

return 0;

}

void QRDecomp( Mat& m, Mat& Q, Mat& R)
{
 Mat e;
 Mat m2=m;  //current A matrix, it would be smaller and smaller for each loop
 Mat Qtemp=Mat::eye(m.rows,m.rows,CV_32F); //square identical matrix for temporarily saving the Q matrix in for loop
 Mat Q1;
 float sign;
 for(int i=0;i<m.rows-1;i++)
 {
  e=Mat::zeros(m2.rows,1,CV_32F);
  e.at<float>(0,0)=1;
  Mat a1=m2.col(0);
  float u_n=float(norm(a1,NORM_L2));
  sign=a1.at<float>(0.0)>0?1:-1;
  Mat u=m2.col(0)+sign*u_n*e;
  Mat v=2*u*u.t()/((float)(norm(u,NORM_L2))*(norm(u,NORM_L2)));
  Q1=Mat::eye(m2.rows,m2.rows,CV_32F)-v;
  cout<<Q1<<Q1.rows<<" Q1"<<endl;

  Mat Q1_2=Mat::eye(m.rows,m.rows,CV_32F);  //make a identical matrix and fill the QA matrix to the buttom-right
  Mat Q1_2_roi=Q1_2(Rect(i,i,m.rows-i,m.rows-i));  //cols=width row=height   To make another matrix for iterating
 
  addWeighted(Q1_2_roi,0,Q1,1,0,Q1_2_roi); 


  Mat QA=Q1*m2;

  {
   m2=QA(Rect(1,1,QA.cols-1,QA.rows-1)).clone(); //cols=width row=height   To make another matrix for iterating
  }

  Qtemp=Q1_2*Qtemp;

 }

 R=Qtemp*m;
 Q=Qtemp.t();
 cout<<R<<"R"<<endl;
 cout<<Q<<"Q"<<endl;

}