2013-04-21 112 views
0

我有一个椭圆形的粒子系统,我试图计算每个粒子的长轴和距下一个粒子中心距离的矢量之间的夹角。我希望这张照片应该很好地说明。 Elliptical particles vector angle计算中间矢量角

在这个例子中,我在(0,0)和(10,-10)初始化两个粒子。每个粒子相对于笛卡尔系统的角度是45度。所以这两个是相互平行的。我想计算我在图片中标记为90度角的角度。

对于这一点,我做了这个功能,使用点积:

double temp = 0; 
    double result = 0; 
    double cosPhi = 0; 
    double dotProd = 0; 
    double cosTheta = cos(current->getTheta()); 
    // If Theta = PI/2, there is a tiny difference between this rounded value 
    // and the exact value of PI/2 due to double precision representation. This 
    // gives cos(PI/2) = 6*10^-17 instead of 0. So, manually set cos(PI/2) = 0. 
    if(almostEqual(cosTheta,0)) 
     cosTheta = 0; 
    double cellLengthX = current->getLength()*cosTheta; 
    double cellLengthY = current->getLength()*sin(current->getTheta()); 
    double dist = Distance(current, next); 
    if(dist == 0) 
     return 0.0; 
    //cout << "Dist " << dist << endl; 
    // calculate the vector of distance 
    // next - current, because the end of the vector is the next cell, 
    // while the start of the vector is the current cell. 
    double distX = next->getCurrX() - current->getCurrX(); 
    double distY = next->getCurrY() - current->getCurrY(); 
    //cout <<"DistX " << distX <<" DistY " << distY << endl; 
    // Dot product 
    dotProd = cellLengthX*distX + cellLengthY*distY; 
    //cout << "DotP " << dotProd << endl; 
    // angle from formula : vector_a*vector_b = |a||b|cos(phi) 
    cosPhi = dotProd/(current->getLength()*dist); 
    cout << "aCos " << acos(cosPhi) << endl; 
    if(cosPhi == 1) // 0 degrees 
    { 
     return 0; 
    } 
    else if(cosPhi == -1) // 180 degrees 
    { 
     return M_PI; 
    } 
    else if(cosPhi == 0) // 90 or 270 degrees 
    { 
     if(distX>0 || distY>0) 
      return M_PI_2; // 90 
     else 
      return 3*M_PI_2; // 270 
    } 
    else if(-distX<=0 && -distY <0) // 1st quadrant -distX<=0 
    { 
     //cout << "Angle 1st :" << acos(cosPhi) << endl; 
     return acos(cosPhi); 
    } 
    else if(-distX>=0 && -distY >0) // 3rd quadrant -distX>=0 
    { 
     temp = acos(cosPhi); 
     result = M_PI_2 + temp; 
     //cout << "Angle 3rd :" << result << endl; 
     return result; 
    } 
    else if(-distX >0 && -distY <=0) // 2nd quadrant -distY <=0 
    { 
     //cout << "Angle 2nd :" << acos(cosPhi) << endl; 
     return acos(cosPhi); 
    } 
    else if(-distX <0 && -distY >=0) // 4th quadrant -distY >=0 
    { 
     temp = acos(cosPhi); 
     result = 3*M_PI_2 + temp; 
     //cout << "Angle 4th :" << result << endl; 
     return result; 
    } 

这种方法原则上工作,但也有一些是错误的,当我在哪个象限判断分配的acos(cosPhi)回报。正如我从其他人那里读到的(发现我自己),acos存在问题,因为您不知道在哪个象限中分配最终角度。

所以我用ATAN2做了一个方法(Y,X):

double x1, y1, x2, y2, phi1, phi2, theta, omega; 
    omega = 0; 
    theta = current->getTheta(); 
    phi1 = 0; 
    phi2 = 0; 
    x1 = current->getLength()*cos(theta); 
    y1 = current->getLength()*sin(theta); 
    x2 = next->getCurrX() - current->getCurrX(); 
    y2 = next->getCurrY() - current->getCurrY(); 
    //phi1 = atan2(y1,x1); 
    //phi2 = atan2(y2,x2); 
    omega = atan2(y2-y1,x2-x1);//phi2 - phi1; 
    cout << phi1 << " " << phi2 << " " << omega << endl; 
    return omega; 

但我得到的回复是第一个粒子对第二-90度(我想这是确定)和-180度不好的粒子。任何帮助将非常感激。

+0

你考虑载体的产品,而不是双重交叉? $ | axb | = | a || b | sin(theta)$ – Bernhard 2013-04-21 13:53:46

+0

你说你想要“每个粒子的长轴与距下一个粒子中心的距离的向量之间的夹角”,但在你的例子中,第二个粒子没有下一个粒子。 – 2013-04-21 14:01:53

+0

我正在计算点积:'vector_a * vector_b = | a || b | cos(phi)',而不是叉积。但是,也许我可以尝试交叉产品。 – Dima1982 2013-04-21 14:03:57

回答

1

这是很难确定,根据您的描述,但是...试试这个:

double x2, y2, theta, omega; 

theta = current->getTheta(); 

x2 = next->getCurrX() - current->getCurrX(); 
y2 = next->getCurrY() - current->getCurrY(); 

omega = atan2(y2,x2) - theta; 
cout << omega << endl; 
return omega; 
+0

你是对的!所以,这基本上做的是找到向量d(在图中)的绝对角度(根据Cartesian(0,0)),然后减去当前粒子的角度?如果我没有说清楚,我很抱歉。你无法想象今天你成为同胞的人有多快乐! :) 谢谢! – Dima1982 2013-04-21 19:51:56