[CoordinateTransformations.jl + Meshes.jl + Rotations.jl] Design for mixed 3d homogeneous and cartesian meshes

Specific Domains Visualization
1

TL;DR

Cartesian meshes seems to be the standard in Julia (Meshes.jl, Rotations.jl…). I read files using 4x4 matrices, should I:

  1. define my own transformation in CoordinateTransformations.jl that takes a Point3 in and a 4x4 and return a Point3 (don’t know how to do it),
  2. stick to homogeneous coordinates everywhere,
  3. find a way to decompose the 4x4 matrix into a composition of 3x3 matrices (don’t know how to do it for rotation*scaling) ?

Longer version with details

I am building a visualization library for 3d plants called PlantGeom.jl (not registered yet). The design is mainly based on AmapStudio and OpenAlea if anyone is interested.

Now we have a (quite) standard file format called OPF (Open Plant Format) for exchange of plants between platforms and models. The file is divided into mainly two parts: reference meshes, and the topology description. The reference meshes are just a list of meshes with an ID (+materials). The topology is a tree-like structure where each node is typically an organ (a leaf, a phytomer…). The node holds the information about which reference mesh is used to represent the organ in 3d (its ID) + an homogeneous (4x4) transformation matrix that is applied to the mesh to get the effective mesh of the organ in the plant world. This design was adopted because plants usually have a limited number of organs, and only their scale and position changes (so having reference mesh * homogeneous matrix saves a lot of space).

Now PlantGeom is not meant for reading OPF files only, I plan to develop functionalities for building plants from scratch too, using reference meshes, and letting the user use transformations to place the organs in the plant.

I am far from a specialist on 3D, so I try to use existing packages as much as possible. I found out that Meshes.jl is becoming the new standard for meshes, and that CoordinateTransformations.jl and Rotations.jl for transformations. So i try using them.

Now here is my problem: the OPF files uses homogeneous transformation matrices in the 4x3 form (I transform them into 4x4 when parsing). But I feel that the design around the other packages is more oriented toward Cartesian matrices (3x3) (or at least Rotation.jl do). (Any idea why?)

I would prefer to stick to whatever is the standard in the Julia community to avoid spending a lot of time on the 3D aspects now and mostly to benefit from the future improvements in those packages “freely” + other packages eventually.

At the moment each node in my tree-like data structure (see MultiScaleTreeGraph.jl for more info) holds a geometry attribute that is a mutable struct with the reference mesh (an immutable struct), the transformation, and optionally the resulting mesh (computed only when needed) + some other unimportant stuff (see geometry struct).

So here is what came to my mind as possible solutions:

  1. use cartesian points and cartesian transformations, except when we have a 4x4 matrix (from reading an OPF file). In this case put the 4x4 matrix in a LinearMap().
  2. use homogenous matrices and points everywhere. Meaning when we apply a new transformation, we just update the homogenous matrix direclty.
  3. decompose the 4x4 matrix and use a composition of 3x3 transformations. This could be elegant but I couldn’t find out how to properly extract the rotation/scaling matrices from the 4x4 matrix.

Point 1 is appealing because then I can save memory space with all my points (only 3 values instead of 4) + I can use Rotations.jl and CoordinateTransformations.jl (and bonus: my transformations are reversible!). The downside is I would have to transform my cartesian 3d points into homogeneous points before applying the LinearMap(4x4) transformation, and I don’t know how to do that.

Point 2 is in my opinion the regular way, but maybe not the preferred way ?

Point 3 is appealing too but I would have to define my own transformations, and it would use more memory for my meshes (one value more per point).

Sorry for the long post! Any opinion that could help me find the right way ?

2

I believe version 2 is standard for a good reason, the cpu will prefer operating on 4x4 matrices since this size aligns very well with the typical SIMD lane width. The memory savings you’d get from storing a rotation matrix and vector separately is typically not beneficial from a speed perspective. In the end, it comes down to how your typical usage looks though.

1 Like
3

OK thanks @baggepinnen ! I’ll use version 2. The thing is I don’t think I can use CoordinateTransformations.jl and Rotations.jl anymore if I do that right ?

4

I’m not familiar with CoordinateTransformations, but you always have the option to extract the rotation matrix should you need to use Rotations.jl for something. If there is a lot of functionality you find useful in those two packages, 2.) might not be the better choice for you?

5

Well the thing is I like the way CoordinateTransformations works :slight_smile: : you can inverse a transformation using inv(), and also (from the docs) " Composition may be intelligent, for instance by precomputing a new Translation by summing the elements of two existing Translation". Which is nice.

When you say “you always have the option to extract the rotation matrix”, you mean use Rotation.jl and put the result in a 4x4 matrix ? Or decompose the 4x4 matrix into 3x3 matrices and add the new rotation? Because I couldn’t find the proper way to decompose the 4x4 matrix. I tried this but it doesn’t work:

# Use the 3x3 matrix inside the 4x4 matrix as the rotation x scaling matrix:
cartesian_rotation = homogeneous_matrix[1:3, 1:3]

# The translation is just the 3 first values from the last column:
cartesian_translation =
    SMatrix{3,3}(
        hcat(
            one(cartesian_rotation)[1:3, 1:2],
            homogeneous_matrix[1:3, 4]
        )
    )

# Build the transformation using CoordinateTransformations's LinearMap and composition (∘):
transformation_3x3 =  LinearMap(cartesian_translation) ∘ LinearMap(cartesian_rotation)

But the result gives me a wrong transformation (my mesh is looking all weird). I also tried changing the order of my transformations but it’s not any better.

6

The problem is with the way you’ve constructed your cartesian_translation — it’s not a LinearMap at all, it’s a Translation. Both linear maps and translations are types of affine maps, so you can use the AffineMap type to combine them, or you can combine them with composition.

Here’s an example of how to deconstruct your 4x4 matrix into the necessary parts:

julia> m44 = [3 2 1 4;
              1 3 1 2;
              1 2 3 1;
              0 0 0 1]
4×4 Matrix{Int64}:
 3  2  1  4
 1  3  1  2
 1  2  3  1
 0  0  0  1

The 3x3 linear map and translation vector:

julia> m33 = m44[1:3,1:3]
3×3 Matrix{Int64}:
 3  2  1
 1  3  1
 1  2  3

julia> t = m44[1:3,4]
3-element Vector{Int64}:
 4
 2
 1

Putting these together, you can see you get the same result using AffineMap or 4x4 matrix multiplication:

julia> AffineMap(m33, t)([1,2,3])
3-element Vector{Int64}:
 14
 12
 15

julia> m44 * [1,2,3,1]
4-element Vector{Int64}:
 14
 12
 15
  1

Alternatively, you can construct the AffineMap from the translation and linear map parts. The library will eagerly compose these together into a single AffineMap data structure (the same one we manually constructed above):

julia> Translation(t) ∘ LinearMap(m33)
AffineMap([3 2 1; 1 3 1; 1 2 3], [4, 2, 1])
2 Likes
7

Haaa yes that’s it ! What was I thinking ? Thank you so much for your answer it works perfectly now !

1 Like