gonum.org/v1/gonum@v0.14.0/spatial/r3/euler_example_test.go (about)

     1  // Copyright ©2021 The Gonum Authors. All rights reserved.
     2  // Use of this source code is governed by a BSD-style
     3  // license that can be found in the LICENSE file.
     4  
     5  package r3_test
     6  
     7  import (
     8  	"fmt"
     9  	"math"
    10  
    11  	"gonum.org/v1/gonum/num/quat"
    12  	"gonum.org/v1/gonum/spatial/r3"
    13  )
    14  
    15  // euler returns an r3.Rotation that corresponds to the Euler
    16  // angles alpha, beta and gamma which are rotations around the x,
    17  // y and z axes respectively. The order of rotations is x, y, z;
    18  // there are many conventions for this ordering.
    19  func euler(alpha, beta, gamma float64) r3.Rotation {
    20  	// Note that this function can be algebraically simplified
    21  	// to reduce floating point operations, but is left in this
    22  	// form for clarity.
    23  	var rot1, rot2, rot3 quat.Number
    24  	rot1.Imag, rot1.Real = math.Sincos(alpha / 2) // x-axis rotation
    25  	rot2.Jmag, rot2.Real = math.Sincos(beta / 2)  // y-axis rotation
    26  	rot3.Kmag, rot3.Real = math.Sincos(gamma / 2) // z-axis rotation
    27  
    28  	return r3.Rotation(quat.Mul(rot3, quat.Mul(rot2, rot1))) // order of rotations
    29  }
    30  
    31  func ExampleRotation_eulerAngles() {
    32  	// It is possible to interconvert between the quaternion representation
    33  	// of a rotation and Euler angles, but this leads to problems.
    34  	//
    35  	// The first of these is that there are a variety of conventions for
    36  	// application of the rotations.
    37  	//
    38  	// The more serious consequence of using Euler angles is that it is
    39  	// possible to put the rotation system into a singularity which results
    40  	// in loss of degrees of freedom and so causes gimbal lock. This happens
    41  	// when the second axis to be rotated around is rotated to 𝝿/2.
    42  	//
    43  	// See https://en.wikipedia.org/wiki/Euler_angles for more details.
    44  
    45  	pt := r3.Vec{1, 0, 0}
    46  
    47  	// For the Euler conversion function in this example, the second rotation
    48  	// is around the y-axis.
    49  	const singularY = math.Pi / 2
    50  
    51  	arb := math.Pi / 4
    52  
    53  	fmt.Printf("rotate around x-axis: %.2f\n", euler(arb, 0, 0).Rotate(pt))
    54  	fmt.Printf("rotate around y-axis: %.2f\n", euler(0, arb, 0).Rotate(pt))
    55  	fmt.Printf("rotate around z-axis: %.2f\n", euler(0, 0, arb).Rotate(pt))
    56  	fmt.Printf("rotate around x+y-axes: %.2f\n", euler(arb, arb, 0).Rotate(pt))
    57  	fmt.Printf("rotate around x+z-axes: %.2f\n", euler(arb, 0, arb).Rotate(pt))
    58  	fmt.Printf("rotate around y+z-axes: %.2f\n", euler(0, arb, arb).Rotate(pt))
    59  
    60  	fmt.Printf("rotate around y-axis to singularity: %.2f\n", euler(0, singularY, 0).Rotate(pt))
    61  	fmt.Printf("rotate around x+y-axes with singularity → gimbal lock: %.2f\n", euler(arb, singularY, 0).Rotate(pt))
    62  	fmt.Printf("rotate around z+y-axes with singularity → gimbal lock: %.2f\n", euler(0, singularY, arb).Rotate(pt))
    63  	fmt.Printf("rotate around all-axes with singularity → gimbal lock: %.2f\n", euler(arb, singularY, arb).Rotate(pt))
    64  
    65  	// Output:
    66  	//
    67  	// rotate around x-axis: {1.00 0.00 0.00}
    68  	// rotate around y-axis: {0.71 0.00 -0.71}
    69  	// rotate around z-axis: {0.71 0.71 0.00}
    70  	// rotate around x+y-axes: {0.71 0.00 -0.71}
    71  	// rotate around x+z-axes: {0.71 0.71 0.00}
    72  	// rotate around y+z-axes: {0.50 0.50 -0.71}
    73  	// rotate around y-axis to singularity: {0.00 0.00 -1.00}
    74  	// rotate around x+y-axes with singularity → gimbal lock: {0.00 0.00 -1.00}
    75  	// rotate around z+y-axes with singularity → gimbal lock: {0.00 0.00 -1.00}
    76  	// rotate around all-axes with singularity → gimbal lock: {0.00 0.00 -1.00}
    77  }