0

I want to compute the gradient (direction and magnitude) of an overdefined plane (> 3 points), such as e.g. four x, y, z coordinates:

[0, 0, 1], [1, 0, 0], [0, 1, 1], [1, 1, 0]

enter image description here

My code for computing the gradient looks like this (using singular value decomposition from this post, modified):

import numpy as np
def regr(X):
    y = np.average(X, axis=0)
    Xm = X - y
    cov = (1./X.shape[0])*np.matmul(Xm.T,Xm)   # Covariance
    u, s, v = np.linalg.svd(cov)               # Singular Value Decomposition
    return u[:,1]

However as a result I get:

[0, 1, 0]

which does not represent the gradient nor the normal vector. The result should look like this:

[sqrt(2)/2, 0, -sqrt(2)/2]

Any ideas why I am getting a wrong vector?

blackst0ne
  • 3,134
  • 4
  • 16
  • 28
  • 1
    Are you sure a `dtype` is being used that can hold floating point numbers? It's not just `int` values, right? – Random Davis Nov 10 '20 at 21:02
  • @RandomDavis I have just checked it, it can hold floating point numbers. It doesn't look like a numerical but a conceptual problem to me... – blackst0ne Nov 10 '20 at 21:14
  • @blackst0ne your output result is normal to that surface. However, I am not sure what you mean by slope on a 2-D surface. Could you please elaborate on that so we can help better? Thank you. – Ehsan Nov 11 '20 at 00:12
  • @Ehsan thanks for your comment. Sorry for the spongy expression, I meant the gradient. And the result [0, 1, 0] is not the normal vector to the plane specified by the four points. However Joe could provide me some useful links. – blackst0ne Nov 12 '20 at 14:03

1 Answers1

1

You can use the function numpy.gradient for that. Here is some math background how it works.

In short:

The directional derivative will be the dot product of the gradient with the (unit normalized) vector of the direction.

An easier solution is to use numdifftools, especially the convenience function directionaldiff. An example can be found here. And if you look at the source code you can see the relation to the math mentioned above.

Joe
  • 6,758
  • 2
  • 26
  • 47
  • I just read your links in more depth, and unfortunately they do not solve my problem. I have an overdefined plane (> 3 points), and that's why I tried the singular value decomposition. Your solution works fine for well defined planes, however not for overdefined ones (e.g. 4 points or 9 points). – blackst0ne Nov 12 '20 at 14:33
  • You could fit your points to a plane equation and feed that into the function. – Joe Nov 12 '20 at 15:44