Rigid Motions on 3D Digital Space

by Kacper Pluta, University Paris-Est Marne-la-Vallée, Gaspard Monge Laboratory of Computer Science

Image from page 252 of "Coloured illustrations of British birds, and their eggs" (1842) and Image from page 339 of "The natural history of British insects"

Technion 13/05/2018

Framework

[\mathbf{R}.\mathbf{x} + \mathbf{t}]

Input data

Output data

Rigid Motions on

\mathbb{Z}^2

rotation matrix

- \mathbf{R}

translation vector

- \mathbf{t}

\cos = \frac{a}{c}, \sin = \frac{b}{c},a,b,c \in \mathbb{Z}

\gcd(a,b,c) = 1.

Interesting properties when

and

integer point

- \mathbf{x}

What are Digitized Rigid Motions?

In general, such transformations are neither surjective nor injective.

Agenda

Introduction to the cubic grid and 3D rigid motions
Bijectivity of 3D digitized rotations
Algorithm for characterizing bijective 3D digitized rotations
Conclusion & Perspectives

Part I: 3D Bijective Digitized Rotations

Motivation & Framework
Alteration of digital images under digitized rigid motions
Problem as a quadrics arrangement
Algorithm for computing the arrangement
Conclusion & Perspectives

Part II: Classifying Rigid Motions of 3D Digital Images

In Collaboration With...

N. Passat of the University of Reims

G. Moroz of Inria Nancy

P. Romon of the University Paris-Est Marne-la-Vallée

Y. Kenmochi of the CNRS (University Paris-Est)

Part I: 3D Bijective Digitized Rotations

Introduction to 3D Digital Geometry

Image from page 27 of "A history of British birds. By the Rev. F.O. Morris .." (1862)

The Cubic Grid

The cubic lattice is then:

\mathbf{x}

Digitization Model

A digitization cell of

\mathbf{x}

denoted by

\mathcal{C}(\mathbf{x}).

Digitization Model

The digitization operator is defined as

\mathcal{D}: \mathbb{R}^3 \rightarrow \mathbb{Z}^3

such that

\forall \mathbf{x} \in \mathbb{R}^3, \exists! \mathcal{D}(\mathbf{x}) \in \mathbb{Z}^3

\mathbf{x} \in \mathcal{C}(\mathcal{D}(\mathbf{x})).

and

Spatial Rotations and Quaternions

\mathbb{H} = \{w + a i + b j + c k | a,b,c,d \in \mathbb{R} \}

with

i^2 = j^2 = k^2 = i j k = {-1}

Properties:

Hamilton product is not commutative
Any rotation can be written as

\mathcal{U} : \mathbf{x} \mapsto q \mathbf{x} q^{-1}

where

q \in \mathbb{H}, \mathbf{x} \in \mathbb{R}^3

We can derive a rotation matrix

where

K = w^2 + a^2 + b^2 + c^2.

The figure source: http://www.3dgep.com/understanding-quaternions/

U = \mathcal{D}\ \circ \mathcal{U}_{|\mathbb{Z}^3}

Properties

Non-injective

Non-surjective

Do not preserve distances

\mathcal{U}(\mathbb{Z}^3)

Digitized Rigid Motions

Remainders

3D digitized rotation is then
bijective when

\forall \mathbf{y} \in \mathbb{Z}^3\ \exists! \mathbf{x} \in \mathbb{Z}^3, S_q(\mathbf{x},\mathbf{y}) \in \mathcal{C}(\mathbf{0})

Bijectivity of 3D Digitized Rotations

Remainders

Equivalently and more computationally friendly (we will see later)

such that

S_{q}(\mathbb{Z}^3,\mathbb{Z}^3) \cap ((\mathcal{C}(\mathbf{0}) \cup q \mathcal{C}(\mathbf{0}) q^{-1}) \setminus (\mathcal{C}(\mathbf{0}) \cap q \mathcal{C}(\mathbf{0}) q^{-1})) = \varnothing

Group Spanned by Values of

S_q

Proposition

If all the generators of

\mathcal{G}_q

have only rational terms, then there exist vectors

which are the

minimal generators of

\mathcal{G}_q.

Group Spanned by Values of

S_q

Lemma

\mathcal{G}_q

is dense, then the

corresponding 3D digitized rotation is not bijective.

Group Spanned by Values of

S_q

Conjecture

\mathcal{G}_q

has a dense factor, the

corresponding digitized rotation is not bijective.

Lipschitz Quaternions and Rational Rotations

Proposition

There is a two-to-one correspondence between the set of Lipschitz quaternions

\mathbb{L} = \{ w + a i + b j + c k | w,a,b,c \in \mathbb{Z}\}

such that the

\gcd(w,a,b,c) = 1,

and the set of rational rotations.

where

K = w^2 + a^2 + b^2 + c^2.

\mathbf{R} =

Bijectivity Characterization Algorithm

Image from page 187 of "Coloured illustrations of British birds, and their eggs" (1842)

Turning into Integers

Since we have

q \cdot \mathbf{x} \cdot q^{-1}

where

q^{-1} = \frac{\bar{q}}{|q|^2}

and

\mathcal{G}_q

is rational, we can

multiply by

|q|^2

and obtain an integer lattice.

Reducing the Left Factor

We are allowed to divide on the left by q while keeping integer valued functions.

Characterization Algorithm Step 1

To check if

\mathbf{w} \in \mathbb{Z}^4

belongs to values

S'_{q},

first verify whether

w w_1 - a w_2 - b w_3 - c w_4 = 0

Characterization Algorithm Step 2

Then we solve the following Diophantine system:

\mathbf{Az} = \mathbf{w}

where

\mathbf{z}^t = (\mathbf{x}, \mathbf{y}) \in \mathbb{Z}^6

and

The complexity of reducing

\mathbf{A} \text{ is } \mathcal{O}(\log^2N(\mathbf{A})).

The final complexity is

\mathcal{O}(|q|^3 + \log^2N(\mathbf{A})).

Some Experimental Results

Image from page 45 of "The natural history of British birds, or, A selection of the most rare, beautiful and interesting birds which inhabit this country

The question is how to analytically characterize such quaternions?

Bijective Digitized Rotations

Square Lattice

A digitized rotation is bijective if and only if the corresponding Gaussian Integer

\alpha

is of the form

\alpha = (k+1) + k i, k \in \mathbb{Z}^+,

i.e. it corresponds to a twin primitive Pythagorean triple.

Hexagonal Lattice

A digitized rotation is bijective if and only if the corresponding Eisenstein Integer

\alpha

has generators of the form

s= 1, t > 0

s>0, t = s + 1.

Conclusion and Perspectives

Analytical characterization of bijective 3D digitized rotations
Proving our conjecture on semi-dense groups

Perspectives

Existence of bijective non-simple 3D digitized rotations
Algorithm for characterization of bijectivity of 3D digitized rotations
Implementation of the algorithm released under a BSD license

Conclusion

Part II: Quadrics Arrangement in Classifying Rigid Motions of a 3D Digital Image

Motivation

Spatial Rotations and Cayley Transform

\mathbf{R} = (\mathbf{I} - \mathbf{A})(\mathbf{I} - \mathbf{A})^{-1}

where

K = 1 + a^2 + b^2 + c^2

and

\mathbf{R} =

\mathbf{A} =

Image Patch

In general, we consider a finite set

\mathcal{N} \subset \mathbb{Z}^3,

called an image

patch whose center

\mathbf{c}

and radius

\mathcal{N}

are given by

\mathbf{c} = \frac{1}{|\mathcal{N}|} \sum_{\mathbf{v} \in \mathcal{N}}

and

r = \max_{\mathbf{v}\in \mathcal{N}} \|\mathbf{v} - \mathbf{c}\|,

respectively.

Definition

Image Patch

For example, let's consider the image patch

\mathcal{N} = \{ (1,0,0), (0,1,0), (0,0,1), (0,0,0), (-1,0,0), (0,-1,0), (0,0,-1) \}

where the center is

(0,0,0).

Alterations Step-by-Step

Critical Rigid Motions

Formulation

\mathbf{R}_i \mathbf{v} + t_i = k_i - \frac{1}{2}

where

\mathbf{v} \in \mathcal{N} \subset \mathbb{Z}^3, k_i \in \mathbb{Z} \cap [-r', r']

and

r + \sqrt{3}.

Critical Rigid Motions

Formulation

\mathbf{R}_i \mathbf{v} + t_i = k_i - \frac{1}{2}

where

\mathbf{v} \in \mathcal{N} \subset \mathbb{Z}^3, k_i \in \mathbb{Z} \cap [-r', r']

and

r + \sqrt{3}.

The parameter space

\Omega = \left\{(a,b,c,t_1,t_2,t_3) \in \mathbb{R}^6 | a,b,c \geq 0, -\frac{1}{2} < t_i < \frac{1}{2} \text{ for } i = 1,2,3\right\}

Problem as Arrangement of Hypersurfaces

The surfaces are given by

\mathbf{R}_i \mathbf{v} + t_i = k_i - \frac{1}{2}

Figure: Visualization 3D parameter space of 2D digitized rigid motions

Problem as Arrangement of Hypersurfaces

In particular, the surfaces are given by

Figure: Visualization 3D parameter space of 2D digitized rigid motions

\frac{1 - a^2}{1+a^2}v_1 - \frac{2a}{1+a^2}v_2 + t_1 = k_1 - \frac{1}{2}

\frac{2a}{1+a^2}v_1 + \frac{1-a^2}{1+a^2}v_2 + t_2 = k_2 - \frac{1}{2}

Problem as Arrangement of Hypersurfaces

Complexity of the arrangement

Overall number of hypersurfaces is

\mathcal{O}(r^4)

and the complexity of the

arrangement is bounded by number of hypersurfaces to the power of dimensionality of the space

\mathcal{O}(r^{24}).

Dimension Reduction

Image taken from page 371 of 'Lahore to Yārkand. Incidents of the route and natural history of the countries traversed by the expedition of 1870, under T. D. Forsyth, etc'

Uncoupling the Parameters

k_i - \frac{1}{2} < \mathbf{R}_i\mathbf{v} + t_i < k_i + \frac{1}{2}

\mathcal{U}(\mathbf{v})

Uncoupling the Parameters

k_i - \frac{1}{2} -\mathbf{R}_i\mathbf{v} < t_i < k_i + \frac{1}{2} -\mathbf{R}_i\mathbf{v}

\mathcal{U}(\mathbf{v})

Uncoupling the Parameters

\mathbf{k}'_i - \frac{1}{2} - \mathbf{R}_i \mathbf{v}' < \mathbf{k}_i + \frac{1}{2} - \mathbf{R}_i \mathbf{v}.

\exists \mathbf{k}, \mathbf{k}'

a,b,c

A configuration exists if there are

\forall \mathbf{v}, \mathbf{v}' \in \mathcal{N}_r

such that

Uncoupling the Parameters

Thanks to the rational expressions in Cayley Transform we obtain
polynomials of degree 2:

q_i[\mathbf{v},k_i](a,b,c) = (1+a^2+b^2+c^2)(2k_i-1-2\mathbf{R}_i\mathbf{v})

Uncoupling the Parameters

We then study arrangement of quadrics defined as

Q_i[\mathbf{v},\mathbf{v}',k_i,k'_i](a,b,c) = q_i[\mathbf{v},k_i](a,b,c) + 2(1+a^2+b^2+c^2) - q_i[\mathbf{v}',k'_i](a,b,c)

Using Sylvester's Law of Inertia we can prove that the only real Euclidean types are hyperboloid of one sheet and hyperbolic paraboloid.

Reducing Set of Quadrics

Discarding quadrics

For example, if we consider the image patch

\mathcal{N} = \{ (1,0,0), (0,1,0), (0,0,1), (0,0,0), (-1,0,0), (0,-1,0), (0,0,-1) \},

we directly obtain 441 quadrics while only 81 are necessary to consider.

The discarded quadrics are:

of non-real Euclidean type i.e., empty or complex
quadrics equivalent up to a multiplication by a constant

What do we want to compute?

At least one sample point for each 3-dimensional connected component of the set

Q_i[\mathbf{v},\mathbf{v}',k_i,k'_i](a,b,c) \neq 0.

Computing Arrangement of Quadrics

Coloured illustrations of British birds, and their eggs (1842)

Relations with Previous Studies

B. Mourrain, J. P. Tecourt, and M. Teillaud: On the Computation of an
Arrangement of Quadrics in 3D, Computational Geometry 30(2), Elsevier (2005)

What is different?

Use of non-generic directions
Support for asymptotic critical values
We only store 3D sample points

Global Idea of the proposed Algorithm

Along non-generic direction, detect and sort all the events in which the topology of the arrangement changes

Between two consecutive events, place a plane, intersect it with the quadrics, and compute one point in each connected component bounded by conics

a_1

a_2

Detection of Events

First intersections for an image patch of

r = 1,2,3,

respectively.

Detection of Events

Q_i(s,b,c) = \partial_b Q_i(s,b,c) = \partial_c Q_i(s,b,c) = 0

Detection of Events

Q_i(s,b,c) = Q_j(s,b,c) = (\nabla Q_i \times \nabla Q_j)_1(s,b,c) = 0

Detection of Events

Q_i(s,b,c) = Q_j(s,b,c) = Q_k(s,b,c) = 0

Detection of Events

An asymptote lives in a quadric.

Detection of Events

A curve given by an intersection of two quadrics is an asymptote.

a = \frac{1}{2}

Recovering the Translational Part

Image from page 114 of "Coloured illustrations of British birds, and their eggs" (1842)

Recovering the Translational Part

Sample points of the translational part can be computed from the
sample points of the previous step from

k_i - \frac{1}{2} -\mathbf{R}_i\mathbf{v} < t_i < k_i + \frac{1}{2} -\mathbf{R}_i\mathbf{v}

Implementation

Maple 2015 and later
Support for Maple Grid framework
Computation can be restarted e.g., after a crash
Partial support for clusters compatible with POSIX
Support for FGb library written by Jean-Charles Faugère
We found and reported a few bugs in Maple and a bug in Linux
Released under Revisited BSD License

https://github.com/copyme/RigidMotionsMapleTools

Conclusion and Perspectives

Optimization of the implementation by exploring symmetry of the variables
Identification of image patches which break connectivity under 3D digitized rigid motions

Perspectives

A new algorithm for computing arrangement of quadrics in 3D
Use of non-generic directions
Taking in to account asymptotic critical values

Conclusion

Thank you for your attention!