Towards Mahony 2024 on-manifold Gaussian products

Project.toml

@@ -53,7 +53,7 @@ Rotations = "1"
 StaticArrays = "0.15, 1"
 Statistics = "1"
 TensorCast = "0.2, 0.3, 0.4"
-TransformUtils = "0.2.10"
+TransformUtils = "0.2.17"
 julia = "1.9"
 
 [extensions]

src/ApproxManifoldProducts.jl

@@ -1,39 +1,45 @@
 module ApproxManifoldProducts
 
-using Reexport
-@reexport using KernelDensityEstimate
+import Base: *, isapprox, convert, show, eltype, length
+
+using Logging
+using StaticArrays
+using LinearAlgebra
+import LinearAlgebra: rotate!, det
+
+using TensorCast
+using DocStringExtensions
+using Distributions
 using Random
+import Random: rand
+
+using Statistics
+import Statistics: mean, std, cov, var, entropy
 
-import TransformUtils as TUs
 import Rotations as _Rot
+using CoordinateTransformations
 
-import ManifoldsBase
-import ManifoldsBase: AbstractManifold, distance
 using RecursiveArrayTools: ArrayPartition
 export ArrayPartition
 
+import ManifoldsBase
+import ManifoldsBase: AbstractManifold, distance
 using Manifolds
 
-using DocStringExtensions
-
 using NLsolve
 import Optim
-using CoordinateTransformations
-using Requires
-using LinearAlgebra
-using TensorCast
-using StaticArrays
-using Logging
-using Statistics
-using Distributions
 
-import Random: rand
+import TransformUtils as TUs
+import TransformUtils: skew
 
-import Base: *, isapprox, convert, show, eltype, length
-import LinearAlgebra: rotate!, det
-import Statistics: mean, std, cov, var, entropy
+using Reexport
+@reexport using KernelDensityEstimate
 import KernelDensityEstimate: getPoints, getBW, evalAvgLogL, entropy, evaluate
 
+
+# FIXME ON FIRE OBSOLETE REMOVE
+using Requires
+
 const MB = ManifoldsBase
 const CTs = CoordinateTransformations
 const AMP = ApproxManifoldProducts
@@ -52,6 +58,7 @@ include("entities/KernelEval.jl")
 include("entities/ManellicTree.jl") # experimental
 include("entities/ManifoldKernelDensity.jl")
 
+include("services/ManifoldsOverloads.jl")
 include("CommonUtils.jl")
 include("services/ManellicTree.jl")
 

src/CommonUtils.jl

@@ -104,6 +104,8 @@ end
 # end
 
 
+
+
 """
     $SIGNATURES
 
@@ -132,8 +134,11 @@ function calcProductGaussians(
   #
   # calc sum of covariances  
   Λ = zeros(MMatrix{dim,dim})
+  # FIXME first transport (push forward) covariances to common coordinates
+  # see [Ge, van Goor, Mahony, 2024]
   Λ .= sum(Λ_)
 
+  # naive mean
   # Tangent space reference around the evenly weighted mean of incoming points
   u0 = mean(M, μ_)
 
@@ -144,7 +149,7 @@ function calcProductGaussians(
     Δuvee = vee(M, u0, log(M, u0, u))
     ΛΔμ += s*Δuvee
   end
-
+  Λ
   # calculate the delta mean
   Δμ = Λ \ ΛΔμ
 
@@ -179,8 +184,9 @@ calcProductGaussians(
 EXPERIMENTAL: On-manifold product of Gaussians.
 
 DevNotes
-- CHECK make sure this product is properly on manifold, Manifolds.jl likely already has solutions:
-  - https://juliamanifolds.github.io/Manifolds.jl/stable/features/distributions.html
+- FIXME do product of concentrated Gaussians on Lie group (approximation):
+  - See Section 3.2 and 4 of [Ge, van Goor, Mahony: A Geometric Perspective on using Gaussian Distributions on Lie Groups, 2024].
+  - Also see upstream utils, https://juliamanifolds.github.io/Manifolds.jl/stable/features/distributions.html
 - FIXME is parallel transport needed when multiplying with covariances from difffent tangent spaces?
 """
 function calcProductGaussians(