(Offer to) move kabsch algorithm to CoordinateTransformations.jl?

[timholy]

I've found several cases where it would be nice to have the kabsch algorithm available for other applications of alignment. In JuliaGeometry/CoordinateTransformations.jl#87 I added an algorithm for constructing general affine transformation from point-pairs; any reservations about moving the kabsch implementation to CoordinateTransformations? I can git-credit you with the PR.

[tmcgrath325]

This seems like a good idea to me! I'd be happy to open a PR there, but before doing so, I have a couple of style/etiquette questions:

• Is there any particular way that the interface should look? It seems like CoordinateTransformations has either functions to generate AffineMap, LinearMap, etc, or functions that do operations on them. Would the Kabsch implementation be a replacement for the implementation you've already added when calling AffineMap, or should it exist separately under something called kabsch()?
• I'd been using GenericLinearAlgebra.svd() for compatibility with ForwardDiff rather than the implementation in LinearAlgebra. I'm not sure how common a use case this will be. Is it worth adding a new dependency to CoordinateTransformations for this?

If you'd prefer to add this yourself, there's no need to worry about git-credit. It's a pretty straightforward algorithm, and I didn't do anything particularly noteworthy code-wise.

[timholy]

I was hoping you'd see it this way! I'll confess that I do have a local branch of CoordinateTransformations that gets started on this. I could either push it or guide you through doing so. If you want to do it yourself, here are some highlights from what I've started:

• I called the algorithm kabsch just like it is here, and exported that name from CoordinateTransformations.
• I implemented the "point-pairs" interface of Construct affine transforms from point-pairs JuliaGeometry/CoordinateTransformations.jl#87 in kabsch
• I decoupled the package from an obligate dependency on GenericLinearAlgebra by allowing the user to supply their own svd algorithm:

function kabsch_centered(P, Q, w; svd::F = LinearAlgebra.svd) where F
    # function body
end

That way anyone who wants differentiability can get it, but it's not required for just computing the transformation.

Let me know if you want to do it or whether I should submit a PR with my branch.

I'm also struggling a bit with whether what I actually want is a generalization to allow isotropic scaling, i.e., a transformation AffineMap(s*R, t) where s is a scalar, R is a rotation matrix (possibly allowing reflection too??), and t a translation. That's not kabsch, of course, but it's still a constrained affine transformation. I'm not sure how useful that is, but from the standpoint of "do these two things line up visually?" I typically think of scaling & reflection as things I'd like to allow, but skewing is not.

[tmcgrath325]

That all sounds good to me. Since you've already got something that sounds ready to push, go for it!

Regarding isotropic scaling, check out this paper that extends the algorithm to consider scaling (to which I was directed by StackOverflow). I haven't followed along with all of the math in detail, but it looks like the result you're looking for is equation (61).

I don't have any ideas for allowing reflections, though, other than maybe trying separate runs of Kabsch after reflecting one of the sets of points.

[timholy]

JuliaGeometry/CoordinateTransformations.jl#97