2

I have been scratching my head for days understanding how to fully describe a block to insert in DXF. For background I am interested in the rules and behaviour of how to do this, not simply trying to find an existing library.

I can create a block in dxf, lets say this is a simple cube. I can then in the Entities section insert this block. To do this I enter the coordinates relative to the object coordinate system codes:(10, 20, 30), and the normal vector of the object codes:(210,220,230). There is also a value to rotate about the normal vector code: 50.

So there are two things to do:

  1. Calculate the normal vector for my object
  2. Calculate the object co-ordinates, given the world coordinates for the object.

To calculate the normal vector, I use quaternions to calculate the rotation applied to the world z-axis (0,0,1) if a yaw, pitch, or roll angle is applied (if not, I simply use the world z-axis). I can use the same quaternions to calculate the arbitrary x and y axis using the same rotation for each of the real world x and y axes. For those looking for more info, this link has a very clear explanation of how to calculate them. https://danceswithcode.net/engineeringnotes/quaternions/quaternions.html

To help you, this calculator can help confirm if this has been implemented correctly or not by checking the results:

Rotation calculator

Then to calculate the object co-ordinates I simply cross multiply the world coordinate with each of the arbitrary axes/normal axes.

For most scenarios this appears to work. However I am left at a loss for additional rules and requirements which I cannot see documented anywhere:

  1. If only a yaw angle is applied (rotation about z), then the normal is still (0,0,1). We need to apply a rotation about z in the scenario that pitch and roll are zero.

  2. If I alter the pitch, the normal vector can change from a negative y-coordinate to a positive one. This change in slope appears to affect the direction the vector is headed and if positive, requires a rotation about the normal of 180 degrees to correct it.

  3. For some unknown reason when I apply a roll angle, my object is rotated around the normal by 90 degrees, and I need to apply a correction here.

I am struggling to find any more clear direction on this online. Does anyone have any thorough explanations that describe the behaviour above, or any pointer material?

samp17
  • 547
  • 1
  • 4
  • 16

1 Answers1

2

As suggested in another topic, a little matrix calculus can solve this issue.

Shortly, they are mainly two ways to describe the same Euler angles:

  • Tait-Bryan angles (yaw, pitch, roll) rotations about X, Y, Z axis;
  • Proper Euler angles (precession, nutation, intrisic rotation) rotations about Z, X, Z axis; AutoCAD uses the second one to describe a bloc orientation. The combination of precession and nutation provides the OCS transformation matrix and the intrinsic rotation is the rotation value of the block.

A little F# example with 3x3 matrices to describe 3d rotations.

type vector = float list

module Vector =
    // Computes the dot product of two vectors
    let dotProduct (v1: vector) (v2: vector) = List.map2 (*) v1 v2 |> List.sum

type matrix(rows: vector list) =
    // Gets the row vectors
    member _.Rows = rows

    // Transposes a matrix
    member _.Transpose = List.transpose rows |> matrix

    // Applies a matrix to a vector
    static member (*)(m: matrix, v: vector) =
        List.map (Vector.dotProduct v) m.Rows

    // Multipies two matrices
    static member (*)(m: matrix, q: matrix) =
        let trp = q.Transpose
        List.map (fun r -> trp * r) m.Rows |> matrix

// Describes a coordinate system data
type CoordinateSystem =
    { WcsToOcs: matrix
      Normal: vector
      Rotation: float }

// Matrix 3x3
module Matrix3x3 =
    // Gets the identity matrix
    let identity =
        matrix [ [ 1.0; 0.0; 0.0 ]
                 [ 0.0; 1.0; 0.0 ]
                 [ 0.0; 0.0; 1.0 ] ]

    // Gets the rotation matrix about X axis
    let xRotation a =
        matrix [ [ 1.0; 0.0; 0.0 ]
                 [ 0.0; cos a; -sin a ]
                 [ 0.0; sin a; cos a ] ]

    // Gets the rotation matrix about Y axis
    let yRotation a =
        matrix [ [ cos a; 0.0; sin a ]
                 [ 0.0; 1.0; 0.0 ]
                 [ -sin a; 0.0; cos a ] ]

    // Gets the rotation matrix about Z axis
    let zRotation a =
        matrix [ [ cos a; -sin a; 0.0 ]
                 [ sin a; cos a; 0.0 ]
                 [ 0.0; 0.0; 1.0 ] ]

    // Creates the matrix according to Yaw, Pitch and Roll values
    let createFromYawPitchRoll yaw pitch roll =
        zRotation yaw * yRotation pitch * xRotation roll

    // Gets the coordinate system data from a matrix 3x3
    let getCoordinateSystem (mat: matrix) =
        match mat.Rows with
        | [ [ m00; m01; m02 ]; [ m10; m11; m12 ]; [ m20; m21; m22 ] ] ->
            let nutation = acos m22

            if abs nutation < 1e-8 then
                { WcsToOcs = identity
                  Normal = [ 0.0; 0.0; 1.0 ]
                  Rotation = atan2 m10 m11 }
            else
                let precession = atan2 m02 -m12
                let spin = atan2 m20 m21

                let xform =
                    (zRotation precession * xRotation nutation)
                        .Transpose

                let normal = xform.Rows.Item 2

                { WcsToOcs = xform
                  Normal = normal
                  Rotation = spin }
        | _ -> invalidArg "mat" "Invalid 3x3 matrix"

// Testing
module test =
    let radians x = x * System.Math.PI / 180.0

    // Input yaw, pitch, roll angles and WCS point
    let yaw = radians 10.0
    let pitch = radians 20.0
    let roll = radians 30.0
    let wcsPoint = [ 8.0; 5.0; 3.0 ]

    // Computation of the coordinate system
    let ocs =
        Matrix3x3.createFromYawPitchRoll yaw pitch roll
        |> Matrix3x3.getCoordinateSystem

    let ocsPoint = ocs.WcsToOcs * wcsPoint

    // Print results
    printfn "Normal X (210): %f" (ocs.Normal.Item 0)
    printfn "Normal Y (220): %f" (ocs.Normal.Item 1)
    printfn "Normal Z (230): %f" (ocs.Normal.Item 2)
    printfn "Rotation (50): %f" ocs.Rotation
    printfn "OCS point X (10): %f" (ocsPoint.Item 0)
    printfn "OCS point Y (20): %f" (ocsPoint.Item 1)
    printfn "OCS point Z (30): %f" (ocsPoint.Item 2)

A C# implementation:

using static System.Math;
using static System.Console;

namespace CsharpMatrix3d
{
    class Program
    {
        /// <summary>
        /// Console testing
        /// </summary>
        static void Main()
        {
            double yaw = D2R(10.0);
            double pitch = D2R(20.0);
            double roll = D2R(30.0);
            var wcsPoint = new Vector(8.0, 5.0, 3.0);

            var ocs = ObjectCoordinateSystem.FromYawPitchRoll(yaw, pitch, roll);
            var ocsPoint = ocs.WorldToPlane * wcsPoint;

            WriteLine($"Normal X (210): {ocs.Normal.X}");
            WriteLine($"Normal Y (220): {ocs.Normal.Y}");
            WriteLine($"Normal Z (230): {ocs.Normal.Z}");
            WriteLine($"Rotation (50): {ocs.Rotation}");
            WriteLine($"OCS point X (10): {ocsPoint.X}");
            WriteLine($"OCS point Y (20): {ocsPoint.Y}");
            WriteLine($"OCS point Z (30): {ocsPoint.Z}");
        }
        static double D2R(double x) => x * PI / 180.0;
    }

    /// <summary>
    /// Provides properties and method for an object coordinate system
    /// </summary>
    struct ObjectCoordinateSystem
    {
        public ObjectCoordinateSystem(Matrix3x3 m)
        {
            double nutation = Acos(m.Row2.Z);
            if (Abs(nutation) < 1e-8)
            {
                Normal = new Vector(0.0, 0.0, 1.0);
                Rotation = Atan2(m.Row1.X, m.Row1.Y);
                WorldToPlane = Matrix3x3.Identity;
            }
            else
            {
                var precession = Atan2(m.Row0.Z, -m.Row1.Z);
                WorldToPlane = (Matrix3x3.ZRotate(precession) * Matrix3x3.XRotate(nutation)).Transpose();
                Normal = WorldToPlane.Row2;
                Rotation = Atan2(m.Row2.X, m.Row2.Y);
            }
        }

        public Vector Normal { get; }
        public Matrix3x3 WorldToPlane { get; }
        public double Rotation { get; }

        public static ObjectCoordinateSystem FromYawPitchRoll(double yaw, double pitch, double roll) =>
            new ObjectCoordinateSystem(
                Matrix3x3.ZRotate(yaw) * 
                Matrix3x3.YRotate(pitch) * 
                Matrix3x3.XRotate(roll));
    }

    /// <summary>
    /// Provides methods for vector calculus
    /// </summary>
    struct Vector
    {
        public Vector(double x, double y, double z)
        { X = x; Y = y; Z = z; }

        public double X { get; }
        public double Y { get; }
        public double Z { get; }

        public double DotProduct(Vector v) =>
            X * v.X + Y * v.Y + Z * v.Z;

        public static Vector operator *(Matrix3x3 m, Vector v) =>
            new Vector(m.Row0.DotProduct(v), m.Row1.DotProduct(v), m.Row2.DotProduct(v));
    }

    /// <summary>
    /// Provides methods for matrix calculus
    /// </summary>
    struct Matrix3x3
    {
        public Matrix3x3(Vector row0, Vector row1, Vector row2)
        { Row0 = row0; Row1 = row1; Row2 = row2; }

        public Vector Row0 { get; }
        public Vector Row1 { get; }
        public Vector Row2 { get; }

        public static Matrix3x3 Identity =>
            new Matrix3x3(
                new Vector(1.0, 0.0, 0.0),
                new Vector(0.0, 1.0, 0.0),
                new Vector(0.0, 0.0, 1.0));

        public Matrix3x3 Transpose() =>
            new Matrix3x3(
                new Vector(Row0.X, Row1.X, Row2.X),
                new Vector(Row0.Y, Row1.Y, Row2.Y),
                new Vector(Row0.Z, Row1.Z, Row2.Z));

        public static Matrix3x3 XRotate(double a) =>
            new Matrix3x3(
                new Vector(1.0, 0.0, 0.0),
                new Vector(0.0, Cos(a), -Sin(a)),
                new Vector(0.0, Sin(a), Cos(a)));

        public static Matrix3x3 YRotate(double a) =>
            new Matrix3x3(
                new Vector(Cos(a), 0.0, Sin(a)),
                new Vector(0.0, 1.0, 0.0),
                new Vector(-Sin(a), 0.0, Cos(a)));

        public static Matrix3x3 ZRotate(double a) =>
            new Matrix3x3(
                new Vector(Cos(a), -Sin(a), 0.0),
                new Vector(Sin(a), Cos(a), 0.0),
                new Vector(0.0, 0.0, 1.0));

        public static Matrix3x3 operator *(Matrix3x3 a, Matrix3x3 b)
        {
            var trp = b.Transpose();
            return new Matrix3x3(trp * a.Row0, trp * a.Row1, trp * a.Row2);
        }
    }
}
gileCAD
  • 2,295
  • 1
  • 10
  • 10