topic Re: numpy linalg.lstsq - coordinate translations in Python Questions

numpy linalg.lstsq - coordinate translations

NeilAyres — Sun, 12 Dec 2021 07:56:43 GMT

I have some control points from a local grid to a known grid (a national grid system).

What I am trying to do is come up with the appropriate parameters to define a local prj for the input data.

This involves finding the origin (or rotation point), the amount of rotation, and any shifts or scaling in X & Y.

So, I have been messing with this code found here :

how-to-perform-coordinates-affine-transformation-using-python-part-2

So, what I have got is (almost the same as the referenced post 😞

import sys, os
import numpy as np
localData = [["AJP1", 5671.666, 19466.156], ["GCP CR", 9634.659, 13064.567], \
 ["MHP2", 10325.093, 7797.802], ["MILE 8", 6508.217, 7511.953]]
refData =[["AJP1", 157171.409, 75402.914], ["GCP CR", 158414.820, 67978.050], \
 ["MHP2", 157059.513, 62842.546], ["MILE 8", 153418.976, 64023.272]]

inArr = np.array([[l[1], l[2], 0.0] for l in localData])
refArr = np.array([[l[1], l[2], 0.0] for l in refData])
n = inArr.shape[0]
pad = lambda x: np.hstack([x, np.ones((x.shape[0], 1))])
unpad = lambda x: x[:,:-1]
X = pad(inArr)
Y = pad(refArr)
# Solve the least squares problem X * A = Y
# to find our transformation matrix A
A, res, rank, s = np.linalg.lstsq(X, Y)
transform = lambda x: unpad(np.dot(pad(x), A))
print "Target:"
print refArr
print "Result:"
print transform(inArr)
print "Max error:", np.abs(refArr - transform(inArr)).max()‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

which give me this :

>>> 
Target:
[[ 157171.409 75402.914 0. ]
 [ 158414.82 67978.05 0. ]
 [ 157059.513 62842.546 0. ]
 [ 153418.976 64023.272 0. ]]
Result:
[[ 157171.39578228 75402.91459879 0. ]
 [ 158414.84986487 67978.04864706 0. ]
 [ 157059.48564043 62842.54723945 0. ]
 [ 153418.98671243 64023.2715147 0. ]]
Max error: 0.0298648653843‍‍‍‍‍‍‍‍‍‍‍‍

So, I can see the fit is quite good.

But how do I use the solution from np.linalg.lstsq to derive the parameters I need for the projection definition of the localData?

In particular, the origin point 0,0 in the target coordinates, and the shifts and rotations that are going on here??

Tagging out very own numpy expert and all around math wiz Dan Patterson here.

Or is this not even possible with this approach.

My matrix math is extremely rusty to non-existent. But I am willing to try any methods offered

Thanks in advance.

Neil

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Sun, 12 Dec 2021 07:56:46 GMT

from a cursory look you are trying to solve the matrix, using homogenous coordinates to pick up the rotation terms. I can't deal with it right now, but I will do so at a later time. When WH helped me through some of this some years ago, a translation about the origin was made by subracting the mean X, and Y from the data set. This centers the point cloud about the origin. Solving for the least squares, essentially gives you the correlation coefficient of the points which is related to the axis of rotation of the points from which one gets the angle of rotation. Think of a perfectly straight line of points parallel to the x axis. Translated to the origin would center the points on the x axis with points to the left of the y-axis and some to the right. The deviation in the y direction would be zero, the correlation coefficient would be zero and therefore the axis of rotation would be zero.

I will leave the above paragraph for you to mull over while I attend to business then dig up my complete affine transformation code.

Others may wade in so Neil doesn't have to wait

here some bits

import numpy as np

def cent_pnt(pnts):
    """Return the centre of a point array"""
    pnts = np.asarray(pnts)
    return np.average(pnts, axis=0).squeeze()


def translate(pnts, cent=None):
    """Translate the points about the origin
    :Requires:
    :--------
    :  pnts - 2D sequence of point coordinates (N, 2)
    :  cent - a list, tuple, ndarray of xc, yc
    :
    :Returns:
    :-------
    :  p_trans - input points translated about the origin
    :
    :-------
    """
    pnts = np.asarray(pnts)
    if cent is None:
        p_trans = pnts - cent_pnt(pnts)
    else:
        p_trans = pnts - np.asarray(cent)
    return p_trans    
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Sun, 12 Dec 2021 07:56:49 GMT

Still busy, but this will delay you some more

Consider the following points:

>>> a = np.array([[1,1], [2,2], [3,3], [4,4]], dtype='float64')
>>> b = affine(a, angle=-22.5)[0]  # rotate to half-way to the x-axis
>>> c = b + [10,10]
>>>
>>> a
array([[ 1.00000000,  1.00000000],
       [ 2.00000000,  2.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.00000000,  4.00000000]])
>>> b
array([[-1.95984445, -0.81179415],
       [-0.65328148, -0.27059805],
       [ 0.65328148,  0.27059805],
       [ 1.95984445,  0.81179415]])
>>> c
array([[ 8.04015555,  9.18820585],
       [ 9.34671852,  9.72940195],
       [ 10.65328148,  10.27059805],
       [ 11.95984445,  10.81179415]])‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

This function (although only lightly tested) should translate and rotate the points from one form to another.

def warp(src=None, dest=None):
    """ move the src points based on the dest pnt affine parameters
    """
    a_slope, a_inter, a_cent = reg_params(src)
    b_slope, b_inter, b_cent = reg_params(dest)
    d_slope = b_slope - a_slope
    results = affine(src, d_slope)
    new_pnts = results[0] + b_cent
    return new_pnts, results‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

def reg_params(src=None):
    """ determine the correlation properties of a point set
    """
    a_cent = cent_pnt(src)
    a_t = translate(src, a_cent)
    xt = a_t[:,0]
    yt = a_t[:,1]
    #co_var = np.cov(m=src.T, y=None, rowvar=1, bias=0, ddof=None)
    #sigmX, sigmY = np.var(a, axis=0, dtype=None, out=None, ddof=0, keepdims=False)
    rho = np.corrcoef(x=xt, y=yt, rowvar=1, bias=0, ddof=None)  
    arr = np.vstack([xt, np.ones(len(xt))]).T
    #arr = xt.reshape(len(xt), 1)
    slope, inter = np.linalg.lstsq(arr, yt)[0]
    deg_slope = np.degrees(np.arctan2(slope,1.0))
    return deg_slope, inter, a_cent #, rho‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

I left in 3 commented lines for examining outputs of the covariance and x, y variance calculations if one wants to examine them. The other line, just shows an alternate method of producing the arr variable.

So if we put this together (with reg_params duplicated), then we have

def cent_pnt(pnts):
    """Return the centre of a point array"""
    pnts = np.asarray(pnts)
    return np.average(pnts, axis=0).squeeze()


def translate(pnts, cent=None):
    """Translate the points about the origin
    :Requires:
    :--------
    :  pnts - 2D sequence of point coordinates (N, 2)
    :  cent - a list, tuple, ndarray of xc, yc
    :
    :Returns:
    :-------
    :  p_trans - input points translated about the origin
    :
    :-------
    """
    pnts = np.asarray(pnts)
    if cent is None:
        p_trans = pnts - cent_pnt(pnts)
    else:
        p_trans = pnts - np.asarray(cent)
    return p_trans    


def rot_matrix(rotation=0):
    """Return the rotation matrix given points and rotation angle
    :Requires:
    :--------
    :  pnts - 2D sequence of point coordinates (N, 2)
    :  rotation - angle in degrees
    :
    :Returns:
    :-------
    :  rot_mat - rotation matrix for 2D transform
    :
    :Notes:
    :-----
    :  for homogenous coordinates use
    :    rot_mat = np.array([[c, -s, 0], [s, c, 0], [0, 0, 1]])
    :--------
    """
    rad = np.deg2rad(rotation)
    c = np.cos(rad)
    s = np.sin(rad)
    rot_mat = np.array([[c, -s], [s, c]])
    return rot_mat

       
def affine(pnts, angle=0, cent=None):
    """Produce an affine transformation of a point set, given an angle of
    :  rotation and a point to rotate about.
    :Requires:
    :--------
    :  pnts - a series of 2D points a list or ndarray (N,2)
    :  angle - rotation angle, in degrees, relative to the x-axis.
    :  cent - center of rotation, either provided or determined from the
    :         center of the point cloud. [x,y], np.array([x,y])
    :
    :Returns:
    :  p_rot - the rotated points.
    :
    :Notes:
    :-----
    :  The functions that are used are: cent_pnt, translate, rot_matrix.
    :  Once the points are translated, and the rotation matrix is obtained,
    :  the points are transposed to facilite broadcasting, the dot product 
    :  determined and the resultant transposed back to the input shape.
    :
    """
    p = np.asarray(pnts)
    if cent is None:
        cent = cent_pnt(pnts)
    p_trans = translate(p, cent=None)
    rot_mat = rot_matrix(rotation=angle)
    p_rot = np.dot(rot_mat, p_trans.T).T
    return p_rot, cent, rot_mat‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Notice I did a slice on the warp in this example, because I didn't need the rest of the returned values from the warp. Looks pretty good when you compare the from-to results with the inputs.

>>> warp(a, b)[0]
array([[-1.95984445, -0.81179415],
       [-0.65328148, -0.27059805],
       [ 0.65328148,  0.27059805],
       [ 1.95984445,  0.81179415]])
>>> warp(a, c)[0]
array([[ 8.04015555,  9.18820585],
       [ 9.34671852,  9.72940195],
       [ 10.65328148,  10.27059805],
       [ 11.95984445,  10.81179415]])
>>>
>>> # some more ...
>>>
>>> warp(b, c)[0]
array([[ 8.04015555,  9.18820585],
       [ 9.34671852,  9.72940195],
       [ 10.65328148,  10.27059805],
       [ 11.95984445,  10.81179415]])
>>> warp(c, a)[0]
array([[ 1.00000000,  1.00000000],
       [ 2.00000000,  2.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.00000000,  4.00000000]])
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

I hope this gets the idea across Neil. This is somewhat simplified and I will still look for the fuller code but the basics entail determining the nature of the parameters of the input points (center, slope etc), then translating, rotating and translating to new center.

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Tue, 27 Sep 2016 12:14:13 GMT

That's fantastic Dan. Plenty to chew over here, but once i have digested this and hopefully understand what's going on, I will attempt to push this through my problem data.

Will get back to you.

I presume the complicated stuff from WH is in your archive somewhere.

But, this has, I hope, wider application with the broader community here. At least anyone who is faced with the problem of trying to define an ArcGIS style local projection from a set of measured points.

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Tue, 27 Sep 2016 16:50:11 GMT

do note, that the approach that was used by the other script uses the transformation matrices of both point clouds since their correlations are not perfect. The examples that I gave have perfect correlations so the full transformation is not needed. My examples were illustrative of how the process works and not a full 'warp'... until I get a chance to get to my files.

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Wed, 28 Sep 2016 06:35:22 GMT

Thanks Dan.

Looking forward to your digging into your files.

Will have a play with these codes blocks in the next few days.

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Sun, 12 Dec 2021 07:56:52 GMT

Okay Dan,

Been messing with that code you posted.

It seems that it goes a bit pear shaped if the rotation angle is > 45deg.

Just added a few prints etc... (inside def warp to display the angles and b_cent)

Test 1 :

>>> 
Using Rotation 44.5, Shift [1000, 2000]
Array A 
[[ 1. 1.]
 [ 2. 2.]
 [ 3. 3.]
 [ 4. 4.]]
Array B 
[[-0.01851178 -2.12123957]
 [-0.00617059 -0.70707986]
 [ 0.00617059 0.70707986]
 [ 0.01851178 2.12123957]]
Array C 
[[ 999.98148822 1997.87876043]
 [ 999.99382941 1999.29292014]
 [ 1000.00617059 2000.70707986]
 [ 1000.01851178 2002.12123957]]
inWarp - Slope A 45.0, B 89.5, D 44.5
Shifts [ 0. 0.]
Warp A -> B
[[-0.01851178 -2.12123957]
 [-0.00617059 -0.70707986]
 [ 0.00617059 0.70707986]
 [ 0.01851178 2.12123957]]
inWarp - Slope A 45.0, B 89.5, D 44.5
Shifts [ 1000. 2000.]
Warp A -> C
[[ 999.98148822 1997.87876043]
 [ 999.99382941 1999.29292014]
 [ 1000.00617059 2000.70707986]
 [ 1000.01851178 2002.12123957]]
inWarp - Slope A 89.5, B 89.5, D 1.67688085639e-12
Shifts [ 1000. 2000.]
Warp B -> C
[[ 999.98148822 1997.87876043]
 [ 999.99382941 1999.29292014]
 [ 1000.00617059 2000.70707986]
 [ 1000.01851178 2002.12123957]]
inWarp - Slope A 89.5, B 45.0, D -44.5
Shifts [ 2.5 2.5]
Warp C -> A
[[ 1. 1.]
 [ 2. 2.]
 [ 3. 3.]
 [ 4. 4.]]
>>>‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Test 2 :

>>> 
Using Rotation 45.5, Shift [1000, 2000]
Array A 
[[ 1.  1.]
 [ 2.  2.]
 [ 3.  3.]
 [ 4.  4.]]
Array B 
[[ 0.01851178 -2.12123957]
 [ 0.00617059 -0.70707986]
 [-0.00617059  0.70707986]
 [-0.01851178  2.12123957]]
Array C 
[[ 1000.01851178  1997.87876043]
 [ 1000.00617059  1999.29292014]
 [  999.99382941  2000.70707986]
 [  999.98148822  2002.12123957]]
inWarp - Slope A 45.0, B -89.5, D -134.5
Shifts [ 0.  0.]
Warp A -> B
[[-0.01851178  2.12123957]
 [-0.00617059  0.70707986]
 [ 0.00617059 -0.70707986]
 [ 0.01851178 -2.12123957]]
inWarp - Slope A 45.0, B -89.5, D -134.5
Shifts [ 1000.  2000.]
Warp A -> C
[[  999.98148822  2002.12123957]
 [  999.99382941  2000.70707986]
 [ 1000.00617059  1999.29292014]
 [ 1000.01851178  1997.87876043]]
inWarp - Slope A -89.5, B -89.5, D -1.67688085639e-12
Shifts [ 1000.  2000.]
Warp B -> C
[[ 1000.01851178  1997.87876043]
 [ 1000.00617059  1999.29292014]
 [  999.99382941  2000.70707986]
 [  999.98148822  2002.12123957]]
inWarp - Slope A -89.5, B 45.0, D 134.5
Shifts [ 2.5  2.5]
Warp C -> A
[[ 4.  4.]
 [ 3.  3.]
 [ 2.  2.]
 [ 1.  1.]]
>>> ‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

You see everything gets neatly reversed. And the shifts from C -> A don't look right.

B -> C always works, because that's just a shift and no rotation.

Thanks for your cooperation and work on this so far.

Neil

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Wed, 28 Sep 2016 13:42:32 GMT

Neil... as I said... stripped down illustrative examples... not full all-encompassing working code (b-c has a rotation... it is just 0... of course there is 90-angle to look at.)

I have to finish my sabbatical report by Friday, ergo, some of the others that have viewed the thread are free to leap in... (we know who you are... we know you can do it)

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Wed, 28 Sep 2016 14:39:38 GMT

Sorry Dan, better let you get on with it...

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Thu, 06 Oct 2016 16:41:05 GMT

Dan_Patterson‌ could I try to revive your interest in this?

Hopefully you are now close to your "archive" and ready to share some of that mathematical wizardry from yourself or Bill Huber.

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Thu, 06 Oct 2016 18:41:07 GMT

Yes... I think I have one nostril above water...how about providing a coordinate list of 'from' and 'to' points so I can check with real world coordinates. If you want, post here and/or email me at my uni email on my profile page.

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Fri, 07 Oct 2016 07:08:07 GMT

Dan, thanks for your help on this...

Well a representative sample are those 4 points I posted originally. The local points are from a rotated project grid.

And the refPoints are as measured using the national grid system (which is normally oriented). Both are on the same classical datum for the area. The datum transformation dosn't come into this.

I have other examples here, but these are good enough I think.

The code I found give me a good idea of the overall fit of the points, and will enable me to identify any points which are way off because of measurement or transcription errors. I was just trying to come up with a method that would enable me to describe a local coordinate system - origin point, rotation and scale, to correctly define this system in an ArcGIS prj.

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Fri, 07 Oct 2016 07:11:10 GMT

And, in the original code I found, why do i need the extra 0.0 in the array before that pad lambda function?

Never use those, cos i can't get my brain around them....

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Sun, 12 Dec 2021 07:56:55 GMT

An affine transformation is a combination of a translation, scaling and rotation transformations, which can be expressed as:

    [1 0 tx]   [Sx 0  0]   [cos(a) -sin(a) 0]   [a b c]
A = [0 1 ty] . [0  Sy 0] . [sin(a)  cos(a) 0] = [d e f]
    [0 0 1 ]   [0  0  1]   [0       0      1]   [0 0 1]                 
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

the extra 0 is for use with homogenous coordinates. There are several ways of accomplishing the same thing

>>> a
array([[ 1.00000000,  1.00000000],
       [ 2.00000000,  2.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.00000000,  4.00000000],
       [ 5.00000000,  5.00000000]])
>>> b = affine(a, angle=-45)[0]
>>> c = affine(a, 45)[0]
>>> b
array([[ 0.17157288,  3.00000000],
       [ 1.58578644,  3.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.41421356,  3.00000000],
       [ 5.82842712,  3.00000000]])
>>> c
array([[ 3.00000000,  0.17157288],
       [ 3.00000000,  1.58578644],
       [ 3.00000000,  3.00000000],
       [ 3.00000000,  4.41421356],
       [ 3.00000000,  5.82842712]])‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

Where the simplied version drops the last row and right-most column as was done in the above example

'a' is the points at a 45 degree angle at a 45 degree (circles). 'b' is rotated -45 (parallel to x) and 'c' parallel to y.

All are centered about the point (3,3) which is the center of input points. All points are translated (moved to the

origin, then rotated, then translated back to their original center.

# array 'a' is pnts in the examples below
# rotation is the angle in degrees
#
#translate....
p_trans = pnts - cent_pnt(pnts) 
a_c = np.average(a, axis=0).squeeze() # cent_pnt is basically this
#
#rotation matrix....
rad = np.deg2rad(rotation)
c = np.cos(rad)
s = np.sin(rad)
rot_mat = np.array([[c, -s], [s, c]])
#
#rotate the translated points, get the dot product, take its translation and move back
p_rot = np.dot(rot_mat, p_trans.T).T + to_pnt‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

A rigid 2d Affine essentially.

>>> a  # input array
array([[ 1.00000000,  1.00000000],
       [ 2.00000000,  2.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.00000000,  4.00000000],
       [ 5.00000000,  5.00000000]])
>>> a.T  # translated .... from column to row for x, y values
array([[ 1.00000000,  2.00000000,  3.00000000,  4.00000000,  5.00000000],
       [ 1.00000000,  2.00000000,  3.00000000,  4.00000000,  5.00000000]])
>>> a_c  # the center of 'a'
array([ 3.00000000,  3.00000000])
>>> a_h  # the coordinates translated to the center
array([[-2.00000000, -2.00000000],
       [-1.00000000, -1.00000000],
       [ 0.00000000,  0.00000000],
       [ 1.00000000,  1.00000000],
       [ 2.00000000,  2.00000000]])
>>> a_h.T  # translate those to rows so that you can use the dot product
array([[-2.00000000, -1.00000000,  0.00000000,  1.00000000,  2.00000000],
       [-2.00000000, -1.00000000,  0.00000000,  1.00000000,  2.00000000]])
>>> 
>>> rm  # a rotation matrix for a 45 degree rotation
array([[ 0.70710678,  0.70710678],
       [-0.70710678,  0.70710678]])
>>>
>>> np.dot(rm, a_h.T).T   # from the calc. class we all skipped
array([[-2.82842712, -0.00000000],
       [-1.41421356, -0.00000000],
       [ 0.00000000,  0.00000000],
       [ 1.41421356,  0.00000000],
       [ 2.82842712,  0.00000000]])
>>>
>>> # as above, but with a translation back to the original center
>>> np.dot(rm, a_h.T).T + a_c  
array([[ 0.17157288,  3.00000000],
       [ 1.58578644,  3.00000000],
       [ 3.00000000,  3.00000000],
       [ 4.41421356,  3.00000000],
       [ 5.82842712,  3.00000000]])
‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍‍

You noticed that I skipped the fancy homogenous coordinate option since I am only doing a simple affine.

If I was going to do not only translation and rotation, but include scaling and shear, then I would have gone

the extra step of making the coordinates homogenous to facilitate the array multiplication since the final

array in the first example would be used once the coefficients of the affine array were determined.

MatPlotLib has a good example in their code in the Affine2D class which you can look at since it is installed with Arc*.

Anyway.... this might keep you entertained for a bit. So in short... if you are required to put an 'extra' number in the

coordinates chances are they are using homogenous coordinates. If you just want to work with plain old regular

coordinates, you simply do it piece-wise as I have shown... I will finish... but I have to get ready for classes

Re: numpy linalg.lstsq - coordinate translations

NeilAyres — Tue, 18 Oct 2016 13:37:04 GMT

Dan,

sorry, been a bit busy.

Just looked into matplotlib, but where do i find this Affine2D class?

Thanks,

Re: numpy linalg.lstsq - coordinate translations

DanPatterson_Retired — Tue, 18 Oct 2016 18:05:06 GMT

me too. but the quickest is to loo at the code itself in the matplotlib folder transforms.py to get help etc

>>> from matplotlib import transforms as tr
>>> dir(tr)
['Affine2D', 'Affine2DBase', 'AffineBase', 'Bbox', 'BboxBase', 'BboxTransform', 'BboxTransformFrom', 'BboxTransformTo', 'BboxTransformToMaxOnly', 'BlendedAffine2D', 'BlendedGenericTransform', 'CompositeAffine2D', 'CompositeGenericTransform', 'DEBUG', 'IdentityTransform', 'MaskedArray', 'Path', 'ScaledTranslation', 'Transform', 'TransformNode', 'TransformWrapper', 'TransformedBbox', 'TransformedPath', 'WeakValueDictionary', '__builtins__', '__cached__', '__doc__', '__file__', '__loader__', '__name__', '__package__', '__spec__', 'absolute_import', 'affine_transform', 'blended_transform_factory', 'composite_transform_factory', 'count_bboxes_overlapping_bbox', 'division', 'interval_contains', 'interval_contains_open', 'inv', 'ma', 'nonsingular', 'np', 'offset_copy', 'print_function', 'six', 'unicode_literals', 'update_path_extents', 'warnings']
>>> help(tr.affine_transform)
Help on built-in function affine_transform: