दो 3 डी वेक्टर

Question

द्वारा परिभाषित घूर्णन से 3x3 रोटेशन मैट्रिक्स की गणना करने के लिए कुशल तरीका, जबकि मुझे इसके 2 समाधान मिल गए हैं, लेकिन अगर मैं इस ऑपरेशन को करने के लिए अच्छी तरह से ज्ञात विधि है तो यह उत्सुक था क्योंकि यह काफी आम काम लगता है।

यहाँ 2 स्पष्ट तरीकों psudo-कोड रहे हैं ...

एक्सिस कोण

यह काफी तार्किक है, लेकिन (रूपांतरण मैट्रिक्स कोण गणना और अक्ष कोण में) sin कॉल दो बार और एक बार cos ।

Matrix3x3 rotation_between_vectors_to_matrix(const Vector v1, const Vector v2) 
{ 
    angle = v1.angle(v2); 
    axis = v1.cross(v2); 

    /* maths for this is well known */ 
    Matrix3x3 matrix = axis_angle_to_matrix(axis, angle); 

    return matrix; 
} 

संपादित: सबसे सरल कार्य है, एक काफी धीमी है लेकिन जैसा कि यहाँ उत्तर में बताया गया है: कोण की गणना angle_sin और angle_cos रही, अक्ष लंबाई और v1,v2 डॉट से से बचा जा सकता क्रमशः उत्पाद।

दो मैट्रिक्स

के बीच

अंतर यहाँ एक और तरीका मैंने पाया जो वैक्टर से दो 3x3 मैट्रिक्स निर्माण करती है और अंतर देता है।

हालांकि यह धीमा है तो अक्ष/कोण गणना जिसे अनुकूलित किया जा सकता है (ऊपर उल्लिखित)।

नोट। यह मानता है कि दोनों वैक्टर सामान्यीकृत हैं, मैट्रिक्स कॉलम-प्रमुख (ओपनजीएल) है।

Matrix3x3 rotation_between_vectors_to_matrix(const Vector v1, const Vector v2) 
{ 
    Matrix3x3 m1, m2; 

    axis = v1.cross(v2); 
    axis.normalize(); 

    /* construct 2 matrices */ 
    m1[0] = v1; 
    m2[0] = v2; 

    m1[1] = axis; 
    m2[1] = axis; 

    m1[2] = m1[1].cross(m1[0]); 
    m2[2] = m2[1].cross(m2[0]); 

    /* calculate the difference between m1 and m2 */ 
    m1.transpose(); 

    Matrix3x3 matrix = m2 * m1; 

    return matrix; 
} 

क्या इस गणना को करने के बेहतर तरीके हैं?

संपादित करें: इस प्रश्न का उद्देश्य माइक्रो-ऑप्टिमाइज़ और प्रत्येक विधि को बेंचमार्क नहीं करना है। इसके बजाए - अगर मैं कुछ पूरी तरह से अलग और बेहतर तरीका है जो मुझे नहीं पता था तो मैं उत्सुक था।

---

नोट: मैं उद्देश्यपूर्ण सह-रैखिक वैक्टर के लिए पतित मामले (जहां अक्ष शून्य लंबाई है) के लिए चेक बाहर छोड़ दिया, उदाहरण सरल रखने के लिए।

Answer 1

आपके द्वारा पोस्ट की गई दोनों विधियों को अनुकूलित किया जा सकता है।

विधि 1

इसके बजाय acos का उपयोग कर दो वैक्टर के बीच के कोण को खोजने के लिए की

, करने के लिए एक बेहतर बात बिल्कुल कोण खोजने से बचना है। कैसे? रॉड्रिग्स द्वारा axis-angle formula केवल sin θ, cos θ और 1 - cos θ की आवश्यकता है, इसलिए वास्तविक कोण ढूंढना अनावश्यक है।

हम जानते हैं कि v1 और v2 इकाई वैक्टर हैं; v1 · v2 = |v1| |v2| cos θ|v1| = |v2| = 1, v1 · v2 से हमें सीधे cos θ देता है, 1 - cos θ महंगा नहीं है। v1 × v2 = |v1| |v2| sin θ n = sin θ n, जहां nv1 और v2 पर एक यूनिट वेक्टर लंबवत है, |v1 × v2| क्रॉस उत्पाद की परिमाण सीधे sin θ देगी।

अब हमारे पास sin θ और cos θ है, हम रॉड्रिग्स फ़ोरमला का उपयोग करके सीधे रोटेशन मैट्रिक्स बना सकते हैं; here's a simple implementation।

विधि 2

के बाद आप मैट्रिक्स के रूप में दो orthonormal फ्रेम का निर्माण किया है, तो आप दूसरे पक्षांतरित आप कर बच सकते हैं। सबूत यहाँ है।

कहो A और B दो मैट्रिक्स हो, जब से तुम A से B हम कुछ मैट्रिक्स X जो जब A साथ गुणा B दे देंगे की जरूरत करने के लिए बारी बारी से करना चाहते हैं:

XA = बी

एक्स = BA⁻¹

यह सब आपको चाहिए; जब आप X से A को पूर्व-गुणा करते हैं तो आपको B मिल जाएगा। हालांकि, क्या आप पाते हैं Y Y

वाई = AB⁻¹

वाई बी = एक

तो फिर तुम Y स्थानांतरित Y⁻¹ यानी

Y⁻¹YB पाने के लिए है = ⁻¹ ए

बी = Y⁻¹A

दो इनवर्क्स (यहां स्थानांतरित) करने के बजाय, आप केवल उपर्युक्त विधि कर सकते हैं जिसमें केवल एक ही ट्रांज़ेक्शन शामिल है।

मैं अभी भी कहूंगा कि उनके अनुकूलित रूपों में विधियों को बेंचमार्क किए बिना, हम यह नहीं कह सकते कि विधि 2 विधि 1 से तेज है। इसलिए मैं वास्तव में आपको दो तरीकों के बीच बेंचमार्क करने के लिए आग्रह करता हूं (कुछ गैर-तुच्छ भार के साथ) और फिर एक निष्कर्ष निकालें।

Answer 2

अक्ष कोण विधि, सबसे तेज तरीका है heres सी कोड मैं कुशल अक्ष के लिए/कोण 3x3 मैट्रिक्स रूपांतरण के लिए के साथ आया था।

यह सह-रैखिक मामलों के लिए भी जांच करता है।

नोट: यदि आपके पास अपनी गणित लाइब्रेरी है, तो आप शायद rotation_between_vecs_to_mat3 पूर्णता के लिए शामिल किसी भी संबंधित कार्यों के बिना काम कर सकते हैं।

#include <math.h> 
#include <float.h> 


/* -------------------------------------------------------------------- */ 
/* Math Lib declarations */ 

static void unit_m3(float m[3][3]); 
static float dot_v3v3(const float a[3], const float b[3]); 
static float normalize_v3(float n[3]); 
static void cross_v3_v3v3(float r[3], const float a[3], const float b[3]); 
static void mul_v3_v3fl(float r[3], const float a[3], float f); 
static void ortho_v3_v3(float p[3], const float v[3]); 
static void axis_angle_normalized_to_mat3_ex(
     float mat[3][3], const float axis[3], 
     const float angle_sin, const float angle_cos); 


/* -------------------------------------------------------------------- */ 
/* Main function */ 
void rotation_between_vecs_to_mat3(float m[3][3], const float v1[3], const float v2[3]); 

/** 
* Calculate a rotation matrix from 2 normalized vectors. 
* 
* v1 and v2 must be unit length. 
*/ 
void rotation_between_vecs_to_mat3(float m[3][3], const float v1[3], const float v2[3]) 
{ 
    float axis[3]; 
    /* avoid calculating the angle */ 
    float angle_sin; 
    float angle_cos; 

    cross_v3_v3v3(axis, v1, v2); 

    angle_sin = normalize_v3(axis); 
    angle_cos = dot_v3v3(v1, v2); 

    if (angle_sin > FLT_EPSILON) { 
axis_calc: 
     axis_angle_normalized_to_mat3_ex(m, axis, angle_sin, angle_cos); 
    } 
    else { 
     /* Degenerate (co-linear) vectors */ 
     if (angle_cos > 0.0f) { 
      /* Same vectors, zero rotation... */ 
      unit_m3(m); 
     } 
     else { 
      /* Colinear but opposed vectors, 180 rotation... */ 
      ortho_v3_v3(axis, v1); 
      normalize_v3(axis); 
      angle_sin = 0.0f; /* sin(M_PI) */ 
      angle_cos = -1.0f; /* cos(M_PI) */ 
      goto axis_calc; 
     } 
    } 
} 


/* -------------------------------------------------------------------- */ 
/* Math Lib */ 

static void unit_m3(float m[3][3]) 
{ 
    m[0][0] = m[1][1] = m[2][2] = 1.0; 
    m[0][1] = m[0][2] = 0.0; 
    m[1][0] = m[1][2] = 0.0; 
    m[2][0] = m[2][1] = 0.0; 
} 

static float dot_v3v3(const float a[3], const float b[3]) 
{ 
    return a[0] * b[0] + a[1] * b[1] + a[2] * b[2]; 
} 

static void cross_v3_v3v3(float r[3], const float a[3], const float b[3]) 
{ 
    r[0] = a[1] * b[2] - a[2] * b[1]; 
    r[1] = a[2] * b[0] - a[0] * b[2]; 
    r[2] = a[0] * b[1] - a[1] * b[0]; 
} 

static void mul_v3_v3fl(float r[3], const float a[3], float f) 
{ 
    r[0] = a[0] * f; 
    r[1] = a[1] * f; 
    r[2] = a[2] * f; 
} 

static float normalize_v3_v3(float r[3], const float a[3]) 
{ 
    float d = dot_v3v3(a, a); 

    if (d > 1.0e-35f) { 
     d = sqrtf(d); 
     mul_v3_v3fl(r, a, 1.0f/d); 
    } 
    else { 
     d = r[0] = r[1] = r[2] = 0.0f; 
    } 

    return d; 
} 

static float normalize_v3(float n[3]) 
{ 
    return normalize_v3_v3(n, n); 
} 

static int axis_dominant_v3_single(const float vec[3]) 
{ 
    const float x = fabsf(vec[0]); 
    const float y = fabsf(vec[1]); 
    const float z = fabsf(vec[2]); 
    return ((x > y) ? 
      ((x > z) ? 0 : 2) : 
      ((y > z) ? 1 : 2)); 
} 

static void ortho_v3_v3(float p[3], const float v[3]) 
{ 
    const int axis = axis_dominant_v3_single(v); 

    switch (axis) { 
     case 0: 
      p[0] = -v[1] - v[2]; 
      p[1] = v[0]; 
      p[2] = v[0]; 
      break; 
     case 1: 
      p[0] = v[1]; 
      p[1] = -v[0] - v[2]; 
      p[2] = v[1]; 
      break; 
     case 2: 
      p[0] = v[2]; 
      p[1] = v[2]; 
      p[2] = -v[0] - v[1]; 
      break; 
    } 
} 

/* axis must be unit length */ 
static void axis_angle_normalized_to_mat3_ex(
     float mat[3][3], const float axis[3], 
     const float angle_sin, const float angle_cos) 
{ 
    float nsi[3], ico; 
    float n_00, n_01, n_11, n_02, n_12, n_22; 

    ico = (1.0f - angle_cos); 
    nsi[0] = axis[0] * angle_sin; 
    nsi[1] = axis[1] * angle_sin; 
    nsi[2] = axis[2] * angle_sin; 

    n_00 = (axis[0] * axis[0]) * ico; 
    n_01 = (axis[0] * axis[1]) * ico; 
    n_11 = (axis[1] * axis[1]) * ico; 
    n_02 = (axis[0] * axis[2]) * ico; 
    n_12 = (axis[1] * axis[2]) * ico; 
    n_22 = (axis[2] * axis[2]) * ico; 

    mat[0][0] = n_00 + angle_cos; 
    mat[0][1] = n_01 + nsi[2]; 
    mat[0][2] = n_02 - nsi[1]; 
    mat[1][0] = n_01 - nsi[2]; 
    mat[1][1] = n_11 + angle_cos; 
    mat[1][2] = n_12 + nsi[0]; 
    mat[2][0] = n_02 + nsi[1]; 
    mat[2][1] = n_12 - nsi[0]; 
    mat[2][2] = n_22 + angle_cos; 
}