Adds Hartley's data normalization method

Includes methods for applying Hartley's data normalization on a set of  N dimensional (N > 1) points. Also includes docstrings and basic tests.

README.md

@@ -1 +1,3 @@
 # MultipleViewGeometry
+
+Under construction..

src/MultipleViewGeometry.jl

@@ -1,5 +1,24 @@
+__precompile__()
+
 module MultipleViewGeometry
 
+using Compat
+
+include("types.jl")
+include("math_aliases.jl")
+
+# Types exported from `types.jl`
+export HomogeneousPoint
+
+# Functions exported from `operators.jl`.
+export 𝑛
+
+# Functions exported from `hartley_transformation.jl`.
+export hartley_normalization, hartley_transformation
+
+include("operators.jl")
+include("data_normalization/hartley_transformation.jl")
+
 # package code goes here
 
 end # module

src/data_normalization/hartley_transformation.jl

@@ -0,0 +1,125 @@
+"""
+    hartley_transformation(pts::AbstractArray{T}) where T<:HomogeneousPoint
+
+Returns a matrix which can be used to map a set of ``d``-dimensional Cartesian
+points  which are represented by ``\\text{length-}(d+1)`` homogeneous coordinates into a
+data-dependent coordinate system. In the data-dependent coordinate system the
+origin is the center of mass (centroid) of the points  and the root-mean-square
+distance of the points to the origin is equal to ``\\sqrt{d}``.
+
+# Details
+
+A point in ``\\mathbb{R}^d`` with
+Cartesian coordinates  ``\\left(m_1, m_2, \\ldots, m_d \\right)`` can also be
+expressed in homogeneous coordinates with the vector ``\\mathbf{m} =
+\\left[m_1, m_2, \\ldots, m_d , 1 \\right]^\\top``.
+
+Suppose one has a set ``\\left \\{ \\mathbf{m}_n \\right \\}_{n = 1}^{N} `` of
+Cartesian points which are represented by homogeneous coordinates.
+Let
+```math
+\\overline{\\mathbf{m}} = \\frac{1}{N} \\sum_{n = 1}^{N} \\mathbf{m}_n
+\\quad  \\text{and} \\quad
+\\sigma = \\left( \\frac{1}{d \\times n} \\sum_{n = 1}^{N}  \\left \\| \\mathbf{m}_n -
+\\overline{\\mathbf{m}} \\right \\|^{2} \\right)^{1/2}
+```
+represent the centroid of the points and the root-mean-square distance of the
+points to the centroid, respectively.
+
+This function returns the matrix
+```math
+\\mathbf{T} =
+\\begin{bmatrix}
+\\sigma^{-1} & 0              &   0     & \\ldots        & -\\sigma^{-1} \\overline{m}_1 \\\\
+           0 & \\sigma^{-1}   &   0     & \\ldots        & -\\sigma^{-1} \\overline{m}_2 \\\\
+           0 & 0              & \\ddots &  0             &  \\vdots                      \\\\
+           \\vdots & \\vdots  &   0     &  \\sigma^{-1}  & -\\sigma^{-1} \\overline{m}_d \\\\
+           0 & 0              &   0     &           0    &                             1
+\\end{bmatrix}
+```
+such that a transformed point ``\\tilde{\\mathbf{m}}_n = \\mathbf{T} \\mathbf{m}_n``
+has a root-mean-square distance to the origin of a new coordinate system
+equal to ``\\sqrt{d}``.
+
+
+"""
+function hartley_transformation(pts::AbstractArray{T}) where T<:HomogeneousPoint
+    if isempty(pts)
+        throw(ArgumentError("Array cannot be empty."))
+    end
+    # Convert list of homogeneous coordinates into an npts x ndim matrix.
+    npts = length(pts);
+    ndim = length(pts[1].coords);
+    array = reinterpret(Float64,pts,(npts,ndim))
+    𝐌 = transpose(reshape(array,(ndim,npts)))
+    _hartley_transformation(view(𝐌,:,1:ndim-1))
+
+end
+
+function _hartley_transformation(𝐌::AbstractArray{T}) where T<:Number
+    if isempty(𝐌)
+        throw(ArgumentError("Array cannot be empty."))
+    end
+    𝐜 = mean(𝐌,1)
+    ndim = length(𝐜)
+    # Compute root mean square distance of each point to the centroid.
+    σ = √((1/length(𝐌)) * ∑((𝐌 .- 𝐜).^2))
+    σ⁻¹ = 1./σ
+    𝐓 = [σ⁻¹*eye(ndim) -σ⁻¹*transpose(𝐜);
+         zeros(1,ndim)                1]
+
+end
+
+"""
+    hartley_normalization(pts::AbstractArray{T}) where T<:HomogeneousPoint
+
+Maps a set of ``d``-dimensional Cartesian points  which are represented by
+``\\text{length-}(d+1)`` homogeneous coordinates into a data-dependent coordinate system.
+In the data-dependent coordinate system the origin is the center of mass
+(centroid) of the points  and the root-mean-square distance of the points to the
+origin is equal to ``\\sqrt{d}``.
+
+# Details
+
+A point in ``\\mathbb{R}^d`` with
+Cartesian coordinates  ``\\left(m_1, m_2, \\ldots, m_d \\right)`` can also be
+expressed in homogeneous coordinates with the vector ``\\mathbf{m} =
+\\left[m_1, m_2, \\ldots, m_d , 1 \\right]^\\top``.
+
+Suppose one has a set ``\\left \\{ \\mathbf{m}_n \\right \\}_{n = 1}^{N} `` of
+Cartesian points which are represented by homogeneous coordinates.
+Let
+```math
+\\overline{\\mathbf{m}} = \\frac{1}{N} \\sum_{n = 1}^{N} \\mathbf{m}_n
+\\quad  \\text{and} \\quad
+\\sigma = \\left( \\frac{1}{d \\times n} \\sum_{n = 1}^{N}  \\left \\| \\mathbf{m}_n -
+\\overline{\\mathbf{m}} \\right \\|^{2} \\right)^{1/2}
+```
+represent the centroid of the points and the root-mean-square distance of the
+points to the centroid, respectively.
+
+This function returns a new set of points ``\\left \\{ \\tilde{\\mathbf{m}}_n \\right \\}_{n = 1}^{N} ``
+where ``\\tilde{\\mathbf{m}}_n = \\mathbf{T} \\mathbf{m}_n`` for each ``n``, and
+```math
+\\mathbf{T} =
+\\begin{bmatrix}
+\\sigma^{-1} & 0              &   0     & \\ldots        & -\\sigma^{-1} \\overline{m}_1 \\\\
+           0 & \\sigma^{-1}   &   0     & \\ldots        & -\\sigma^{-1} \\overline{m}_2 \\\\
+           0 & 0              & \\ddots &  0             &  \\vdots                      \\\\
+           \\vdots & \\vdots  &   0     &  \\sigma^{-1}  & -\\sigma^{-1} \\overline{m}_d \\\\
+           0 & 0              &   0     &           0    &                             1
+\\end{bmatrix}.
+```
+These new points have the property that their root-mean-square distance to
+the origin of the coordinate system is equal to ``\\sqrt{d}``.
+
+
+"""
+function hartley_normalization(pts::AbstractArray{T}) where T<:HomogeneousPoint
+    𝐓 = hartley_transformation(pts)
+    map!(pts , pts) do p
+        𝐦 = collect(p.coords)
+        𝐦 = 𝑛(𝐓 * 𝐦)
+        HomogeneousPoint(tuple(𝐦...))
+    end
+end

src/math_aliases.jl

@@ -0,0 +1,2 @@
+const  √ = sqrt
+const  ∑ = sum

src/operators.jl

@@ -0,0 +1,37 @@
+"""
+    𝑛(v::Vector{T}) where T<:Number
+
+Scales a length-``n`` vector ``v``  such that the last component
+of the vector is one, provided that the last component is not zero. If the last
+component is zero then the vector is left unchanged.
+
+# Details
+Suppose the length-``n`` vector ``v`` represents the homogeneous coordinates  of
+a point in a projective space. The corresponding Cartesian coordinates  usually
+just the first ``n-1`` numbers of homogeneous coordinates divided by  the last
+component. So if the last component is one, then the first  ``n-1`` homogeneous
+coordinates can be interpreted as Cartesian.  The exceptional case is when the
+last component of the homogenenous coordinates is zero. These homogeneous
+coordinates are associated with so-called *points at infinity* and have no
+Cartesian counterparts.
+
+# Example
+```julia
+h = [4, 4 , 2]
+c = 𝑛(h)
+
+3-element Array{Float64,1}:
+ 2.0
+ 2.0
+ 1.0
+```
+
+
+"""
+function 𝑛(v::Vector{T}) where T<:Real
+    if v[end] != 0 && v[end] != 1
+        v = v ./ v[end]
+    else
+        v
+    end
+end

src/types.jl

@@ -0,0 +1,3 @@
+struct HomogeneousPoint{T <: AbstractFloat,N}
+    coords::NTuple{N, T}
+end

test/data_normalization_tests.jl

@@ -0,0 +1,43 @@
+using MultipleViewGeometry, Base.Test
+
+
+# Tests for a set of two-dimensional Cartesian points represented by homogeneous
+# coordinates.
+pts2D = map(HomogeneousPoint,
+        [(-10.0, -10.0, 1.0),
+         (-10.0,  10.0, 1.0),
+         ( 10.0, -10.0, 1.0),
+         ( 10.0,  10.0, 1.0)])
+
+@test hartley_normalization(pts2D) == map(HomogeneousPoint,
+                                        [(-1.0,-1.0, 1.0),
+                                         (-1.0, 1.0, 1.0),
+                                         (1.0, -1.0, 1.0),
+                                         (1.0,  1.0, 1.0)])
+
+ @test hartley_transformation(pts2D) == eye(3)
+
+
+# Tests for a set of three-dimensional Cartesian points represented by homogeneous
+# coordinates.
+pts3D = map(HomogeneousPoint,
+           [(-10.0, -10.0, -10.0, 1.0),
+            (-10.0, -10.0,  10.0, 1.0),
+            (-10.0,  10.0, -10.0, 1.0),
+            (-10.0,  10.0,  10.0, 1.0),
+            ( 10.0, -10.0, -10.0, 1.0),
+            ( 10.0, -10.0,  10.0, 1.0),
+            ( 10.0,  10.0, -10.0, 1.0),
+            ( 10.0,  10.0,  10.0, 1.0)])
+
+@test hartley_normalization(pts3D) == map(HomogeneousPoint,
+                                         [(-1.0,-1.0, -1.0, 1.0),
+                                          (-1.0,-1.0,  1.0, 1.0),
+                                          (-1.0, 1.0, -1.0, 1.0),
+                                          (-1.0, 1.0,  1.0, 1.0),
+                                          (1.0, -1.0, -1.0, 1.0),
+                                          (1.0, -1.0,  1.0, 1.0),
+                                          (1.0,  1.0, -1.0, 1.0),
+                                          (1.0,  1.0,  1.0, 1.0)])
+
+@test hartley_transformation(pts3D) == eye(4)

test/operators_tests.jl

@@ -0,0 +1,8 @@
+using MultipleViewGeometry, Base.Test
+
+# Vectors are scaled so that the last component is unity.
+v = [2,2,2]
+@test 𝑛(v) == [1.0, 1.0, 1.0]
+# Vectors which represent points at infinity are unchanged.
+v = [2,2,0]
+@test 𝑛(v) == [2,2,0]

test/runtests.jl

@@ -1,5 +1,6 @@
 using MultipleViewGeometry
 using Base.Test
 
-# write your own tests here
-@test 1 == 2
+
+@testset "Operators Tests" begin include("operators_tests.jl") end
+@testset "Data Normalization Tests" begin include("data_normalization_tests.jl") end