2013年10月8日 星期二

計算兩條直線的交點

開發工具: Microsoft Visual studio 2010 with OPENCV2.1
CPU: Intel Core Quad Q8300
直接看圖說故事:
 

以上為推導過程。
(2013/10/29更新) 上面的推導過程錯了一個地方,在於Vx是cos,Vy是sin,也就是上面的m(斜率)是錯誤的。
以下才是正確的推導過程















以下為程式碼:


#include 
#include 
#include 

#define SMALL 0.0000001        //定義一個很小的數字以防除以0

/*
   使用Matrix形式解,前提是兩條直線必有交點,line = { vx , vy , x0 , y0 }
*/
inline cv::Point2f twoLineIntersection1( const float line1[4] , const float line2[4] ) {
  //錯誤的
 /*float Vx1 = -line1[0] , Vx2 = -line2[0] ;
 float b1 = Vx1 * line1[2] + line1[1] * line1[3] , b2 = Vx2 * line2[2] + line2[1] * line2[3] ;
 float det = ( Vx1 * line2[1] + line1[1] * line2[0] ) ;
 return cv::Point2f( (line2[1] * b1 - line1[1] * b2) / det , ( line2[0] * b1 - line1[0] * b2 ) / det ) ;*/

 float b1 = line1[0] * line1[3] - line1[1] * line1[2] ;
 float b2 = line2[0] * line2[3] - line2[1] * line2[2] ;
 float det = 1.0f / (line1[0] * line2[1] - line1[1] * line2[0]) ;
 return cv::Point2f( (line2[0] * b1 - line1[0] * b2) * det , (line2[1] * b1 - line1[1] * b2) * det ) ;
}
/*
    使用y = mx + b形式解,前提以及直線的儲存方式與Matrix一樣
*/
inline cv::Point2f twoLineIntersection2( const float line1[4] , const float line2[4] ) {
 //錯誤的
 /*float m1  , m2 ;
 float b1 , b2 , x , y ;
 bool rot1 = false , rot2 = false ;
 if ( abs(line1[0]) <= SMALL )  //判斷vx有沒有小於設定的數值,如果有就標示那條線為垂直線
  rot1 = true ;
 if ( abs(line2[0]) <= SMALL )
  rot2 = true ;

 if ( rot1 | rot2 ) {        //如果兩條線之中有一條為垂直線,因為前提是兩條直線有交點,所以不會有兩條直線都是垂直線
  if ( rot1 ) {       //如果第一條為垂直線
   y = line1[3] ;
   m2 = line2[0] / line2[1] ;
   x = (y - (line2[3] - m2 * line2[2])) / m2 ;
  } else {            //如果第二條為垂直線
   y = line2[3] ;
   m1 = line1[0] / line1[1] ;
   x = (y - (line1[3] - m1 * line1[2])) / m1 ;
  }
 } else {
  m1 = line1[0] / line1[1] , m2  = line2[0] / line2[1] ;
  b1 = line1[3] - m1 * line1[2] , b2 = line2[3] - m2 * line2[2] ;
  x = ( b2 - b1 ) / ( m1 - m2 ) ;
  y = m1 * x + b1 ;
 }
 return cv::Point2f( x , y ) ;*/
 float m1  , m2 ;
 float b1 , b2 , x , y ;
 bool rot1 = false , rot2 = false ;
 if ( abs(line1[0]) <= SMALL )  //cos -> 0 => 垂直線,x = line1[2]
  rot1 = true ;
 if ( abs(line2[0]) <= SMALL )
  rot2 = true ;

 if ( rot1 | rot2 ) {
  if ( rot1 ) {
   x = line1[2] ;
   m2 = line2[1] / line2[0] ;
   y = m2 * x + (line2[3] - m2 * line2[2]) ;
  } else {
   x = line2[2] ;
   m1 = line1[1] / line1[0] ;
   y = m1 * x + (line1[3] - m1 * line1[2]) ;
  }
 } else {
  m1 = line1[1] / line1[0] , m2  = line2[1] / line2[0] ;
  b1 = line1[3] - m1 * line1[2] , b2 = line2[3] - m2 * line2[2] ;
  x = ( b2 - b1 ) / ( m1 - m2 ) ;
  y = m1 * x + b1 ;
 }
 return cv::Point2f( x , y ) ;
}

int main() {
 float line1[4] = { 0 , 1 , 1 , 1 } , line2[4] = { 1 , 0 , 3 , 3 };
 cv::Point2f p1 , p2 ;
 LARGE_INTEGER t1, t2, ts;
 QueryPerformanceFrequency(&ts);
 QueryPerformanceCounter(&t1);       //第一個方法的效能測試
 for ( int i = 0 ; i < 0xffffff ; ++i )
  p1 = twoLineIntersection1( line1 , line2 ) ;
 QueryPerformanceCounter(&t2);
 printf("(1)Lasting Time: %lf\n",(t2.QuadPart-t1.QuadPart)/(double)(ts.QuadPart));
 QueryPerformanceCounter(&t1);       //第二個方法的效能測試
 for ( int i = 0 ; i < 0xffffff ; ++i ) 
  p2 = twoLineIntersection2( line1 , line2 ) ;
 QueryPerformanceCounter(&t2);
 printf("(2)Lasting Time: %lf\n",(t2.QuadPart-t1.QuadPart)/(double)(ts.QuadPart));
 return 0 ;
}

經在下測試後,運算速度較快的是以Matrix的形式下去解,以程式上面設定的迴圈數執行時間大約是2.76s,而y = mx + b大約要5.2s。

造成這樣情況的原是,以y = mx + b的形式要判斷垂直線,所以這部分消耗太多CPU cycle,而導致運算速度變慢。

沒有留言: