# Coordinate transformation from three known points
**URL:** https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153
**Category:** General Usage
**Tags:** question
**Created:** [April 7, 2020, 8:59am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153 "2020-04-07T08:59:34Z")
**Posts on this page:** 19
**Page:** 1
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 7, 2020, 8:59am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/1 "2020-04-07T08:59:34Z")
I have 3 coordinates, `moving_points`. I know where each of these is supposed to really be, `fixed_points`. I want to transform them (and many other coordinates) to these known locations (involving translation, rotation, but also a bit of stretching).
How can I create a transformation so that when it’s applied to these 3 coordinates they become a set of known other coordinates?
MWE:
```julia
fixed_points = [[0,0], [0,50], [130,50]]
moving_points = [rand(2) for _ in 1:3]
trans = # Probably something involving CoordinateTransformations, the fixed_points, and moving_points
@assert fixed_points ≈ trans(moving_points)
```
---
### Author:  [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [April 7, 2020, 9:11am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/2 "2020-04-07T09:11:38Z")
If your “fixed points” are y\_1, y\_2, y\_3, and your “moving points” are x\_1, x\_2, x\_3, you are probably looking for A,b such that
Ax\_i + b = y\_i \quad i=1,2,3
which is a linear system, so just formulate it as one and use `\`.
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 7, 2020, 9:22am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/3 "2020-04-07T09:22:25Z")
Do you mean this?
```julia
fixed_points = [[0,0], [0,50], [130,50]]
X = hcat(vcat(fixed_points...), ones(2*3))
moving_points = [rand(2) for _ in 1:3]
y = vcat(moving_points...)
p = X\y
```
which wouldn’t work for any of the points. Or this:
```julia
Ab = [hcat(m, [1,1])\f for (m, f) in zip(moving_points, fixed_points)]
```
but then, how do I apply that…? As you can see, I’m a wee bit lost.
I was under the impression that `trans` would be a 3×3 transformation matrix…
---
### Author:  [@mschauer](https://discourse.julialang.org/u/mschauer)
#### Post date: [April 7, 2020, 10:48am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/4 "2020-04-07T10:48:44Z")
Trans would be 3x3 and solve
```julia
M*[xs ys ones(3)] = [xsnew ysnew ones(3)] ?
```
Edit: Yes
```julia
julia> M = [xsnew ysnew ones(3)]/[xs ys ones(3)]
3×3 Array{Float64,2}:
1.02264 -0.615986 0.593346
0.140419 0.425081 0.434501
2.92372 -1.79821 -0.125504
julia> [xsnew ysnew ones(3)]
3×3 Array{Float64,2}:
0.485407 0.575846 1.0
0.232837 0.587207 1.0
0.905755 0.955727 1.0
julia> M*[xs ys ones(3)]
3×3 Array{Float64,2}:
0.485407 0.575846 1.0
0.232837 0.587207 1.0
0.905755 0.955727 1.0
```
---
### Author:  [@Skoffer](https://discourse.julialang.org/u/Skoffer)
#### Post date: [April 7, 2020, 10:54am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/5 "2020-04-07T10:54:40Z")
This one works (of course, it can be written in more efficient way, it’s only proof of concept)
```julia
X = permutedims(hcat(map(x -> [x..., 1.0], moving_points)...))
Y = hcat(map(x -> x[1], fixed_points), map(x -> x[2], fixed_points))
c = X \ Y
c = c'
A = c[1:2, 1:2]
b = c[:, 3]
```
And transformation can be written like this
```julia
mov_points = hcat(moving_points...)
A * mov_points .+ b
```
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 7, 2020, 10:58am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/6 "2020-04-07T10:58:13Z")
That works! Awesome! I can’t help imagine that there is (or should be) a mechanism for this somewhere in [GitHub - JuliaGeometry/CoordinateTransformations.jl: A fresh approach to coordinate transformations...](https://github.com/JuliaGeometry/CoordinateTransformations.jl)…
Here is something:
```julia
using CoordinateTransformations
f = Translation(b)∘LinearMap(A)
@assert f.(moving_points) ≈ fixed_points
```
---
### Author:  [@Skoffer](https://discourse.julialang.org/u/Skoffer)
#### Post date: [April 7, 2020, 11:06am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/7 "2020-04-07T11:06:03Z")
You can use directly `f = AffineMap(A, b)`, there is no need to compose two transformations.
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 7, 2020, 11:30am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/8 "2020-04-07T11:30:44Z")
Here’s an OK implementation that works:
```julia
using CoordinateTransformations
function CoordinateTransformations.AffineMap(ps::Vector{Pair{T,R}}) where {T, R}
X = vcat(hcat(last.(ps)...), ones(1, length(ps)))'
Y = hcat(first.(ps)...)'
c = (X \ Y)'
A = c[:, 1:2]
b = c[:, 3]
AffineMap(A, b)
end
fixed_points = [[0,0], [0,50], [130,50]]
moving_points = [rand(2) for _ in 1:3]
f = AffineMap(Pair.(fixed_points, moving_points))
@assert fixed_points ≈ f.(moving_points)
```
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 7, 2020, 1:45pm UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/9 "2020-04-07T13:45:12Z")
OK… I’m stumped again.
So while this works great for **three** coordinates, I’m trying to make it work for only **two** coordinates. But the problem is that this implementation “flattens” all the coordinates I apply it to:
```julia
ps = [[0,0] => [1,1], [10,0] => [9,4]] # using just two coordinates
f = AffineMap(ps)
@assert first.(ps) ≈ f.(last.(ps)) # it "works"
n = 1000
for (x,y) in zip(rand(-20:20, n), rand(-20:20, n))
_, y2 = f([x,y]) # but no matter what the moving point is
@assert y2 == 0 # its `y` value is always zero
end
```
I’m sure this makes sense because there is no information about the y axis in those two coordinates (both have zero y-values). But I’m now wondering how to make this work. Now the transformation is a translation followed by a rotation and finally stretching just one of the axis (in this case, the x axis).
---
### Author:  [@mschauer](https://discourse.julialang.org/u/mschauer)
#### Post date: [April 7, 2020, 3:14pm UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/10 "2020-04-07T15:14:04Z")
With two coordinates, you need an additional constraint, for example that the map is only rotation and scaling, but no shearing. Scaling can be computed from the change in the distance of the points, rotation from change in the angle of the line connecting the two points.
---
### Author:  [@Skoffer](https://discourse.julialang.org/u/Skoffer)
#### Post date: [April 7, 2020, 4:17pm UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/11 "2020-04-07T16:17:28Z")
Wow, I have never worked with `CoordinateTransformations` before, but it is really very powerful.
Solution to this problem is in the sequence of linear transformations. We should start from moving points to the origin, than rotate them to align along `x` axis, than scale, than rotate backward in such a way so points will be collinear and finally shift them to `fixed_point` (it’s not very clear from this description, try to read the code and draw necessary transformations). And great thing: these are all linear transformations and `CoordinateTransformations` internally concatenate them together into single `AffineMap` so no Linear Algebra is needed, it’s pure geometry!
```julia
using CoordinateTransformations
using Rotations
fixed_points = [[1, 1], [9, 4]]
moving_points = [rand(2) for _ in 1:2]
z1 = moving_points[2] .- moving_points[1]
z2 = fixed_points[2] .- fixed_points[1]
α = norm(z2)/norm(z1)
θ1 = PolarFromCartesian()(z1).θ
θ2 = PolarFromCartesian()(z2).θ
t1 = Translation(-moving_points[1]...)
lm1 = LinearMap(RotZ(-θ1)[1:2, 1:2])
lm2 = LinearMap([α 0; 0 1])
lm3 = LinearMap(RotZ(θ2)[1:2, 1:2])
t2 = Translation(fixed_points[1]...)
am = t2 ∘ lm3 ∘ lm2 ∘ lm1 ∘ t1 # This is a single affine map!
@assert fixed_points ≈ am.(moving_points)
```
---
### Author:  [@mschauer](https://discourse.julialang.org/u/mschauer)
#### Post date: [April 8, 2020, 7:55am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/12 "2020-04-08T07:55:25Z")
> [@Skoffer](#):
>
> lm2 = LinearMap([α 0; 0 1])
Perhaps lm2 = LinearMap([α 0; 0 α])?
---
### Author:  [@Skoffer](https://discourse.julialang.org/u/Skoffer)
#### Post date: [April 8, 2020, 8:08am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/13 "2020-04-08T08:08:42Z")
> [@yakir12](#):
>
> finally stretching just one of the axis (in this case, the x axis).
It can be, but OOP explicitly said “finally stretching just one of the axis (in this case, the x axis).”
As a side note, one can remove `Rotations` dependence, and do instead
```julia
nz1 = z1 ./norm(z1)
lm1 = LinearMap([nz1[1] -nz[2]; nz[2] nz[1])
```
and the same for `lm3`
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 8, 2020, 8:49am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/14 "2020-04-08T08:49:20Z")
Thanks a lot!!! I’ve managed to avoid this kind of transformation altogether.
But! the bigger problem of creating a transformation from a set of fixed and moving coordinates is still open (see issues [30](https://github.com/JuliaGeometry/CoordinateTransformations.jl/issues/30) and [33](https://github.com/JuliaGeometry/CoordinateTransformations.jl/issues/33) in `CoordinateTransformations.jl`). Matlab solves it nicely, see “More About” and then “Transformation Types” [here](https://se.mathworks.com/help/images/ref/fitgeotrans.html#More%20About). Perhaps an additional package that does exactly that would be beneficial. It would be extra nice if said package would solve it for when there are more than the minimum number of coordinate pairs and the resulting transformation object would be an optimum considering all of the input pairs.
---
### Author:  [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)
#### Post date: [April 8, 2020, 9:23am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/15 "2020-04-08T09:23:13Z")
> [@yakir12](#):
>
> the bigger problem of creating a transformation from a set of fixed and moving coordinates is still open (see issues [30](https://github.com/JuliaGeometry/CoordinateTransformations.jl/issues/30)
I think you were suggested a perfectly fine solution in that issue:
>
>
> I might have missed how to do this (and am glad to add this to the documentation…), but how would one get a transformation from a (large) set of locations, say on an image, and the corresponding real-world coordinates?
For exactly determined linear problems (like the original one here), just use the linear solver. For overdetermined linear problems you can use eg least squares (again `\`), and for nonlinear transformations just minimize a cost function.
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [April 8, 2020, 9:31am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/16 "2020-04-08T09:31:25Z")
Yea, for sure. I need to get on that.
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [April 8, 2020, 9:42am UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/17 "2020-04-08T09:42:30Z")
I’m surprised no one has mentioned [Orthogonal Procrustes problem - Wikipedia](https://en.wikipedia.org/wiki/Orthogonal_Procrustes_problem)
with a julia blog post here [Page Redirection](https://simonensemble.github.io/2018-10-27-orthogonal-procrustes/)
Another related post on the topic, although quite old at this point
[Calculating rigid transformation between two point set](https://discourse.julialang.org/t/calculating-rigid-transformation-between-two-point-set/6199)
---
### Author:  [@yakir12](https://discourse.julialang.org/u/yakir12)
#### Post date: [September 15, 2020, 2:54pm UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/18 "2020-09-15T14:54:08Z")
Does this orthogonal Procrustes problem include all the possible transformations, namely translation, scale, shear/stretching or “just” rotation? Cause otherwise that blog post would basically solve camera calibration for me…
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [September 15, 2020, 4:07pm UTC](https://discourse.julialang.org/t/coordinate-transformation-from-three-known-points/37153/19 "2020-09-15T16:07:37Z")
It’s no longer an orthogonal transformation if it includes shearing?
In general, if the correspondence between the points in the two sets to be aligned is known, the problem is relatively easy for linear operators. If the correspondence is not known it becomes a bit messier.