Hello,

Like in title.

I know that you can use determinant for solving system of equations with Cramer's rule but I love to see any practical example of use in Houdini ?

Dot and Cross products have clear practical use in everyday work but determinant ? I don't see any need for use.


Thanks in advance.

As it happens I just browsed some literature and came across a use for it.

You can use it, for example, to calculate a covariance matrix to get better bounding boxes. A technique used in some games.

If I have time I'll look into it and see if I can replicate the algorithm in vopsop.

hey,

What literature ? Books ? Maybe something that I have and read but forgot about.


thanks for interest.

It's in here:

http://www.amazon.com/Game-Programming-Gems-v/dp/1584502959 [amazon.com]

Let us know if you implement it.

Soothsayer
It's in here:
http://www.amazon.com/Game-Programming-Gems-v/dp/1584502959 [amazon.com]
Let us know if you implement it.

Also here: http://www.amazon.com/Mathematics-Programming-Computer-Graphics-Development/dp/1584502770/ref=sr_1_1?ie=UTF8&s=books&qid=1268657048&sr=8-1 [amazon.com]

…interesting, but it seems they use PCA to compute best fit for an axis of bounding box. AFAIK computing pca in vex would be pretty difficult (whereas rather nicely in python with a help of a couple of modules).

SYmek
…interesting, but it seems they use PCA to compute best fit for an axis of bounding box. AFAIK computing pca in vex would be pretty difficult (whereas rather nicely in python with a help of a couple of modules).

Ok, I'm getting interested in this now.

I had a little go at it and calculated the covariance matrix with python but now I realize there's not much matrix functionality to do the rest of the stuff.

NumPy would be able to help me out there, right? I can't figure out how to import it into houdini python. Any ideas?

Actually, it's ok, I managed to get it working.

Below is a screenshot. The red green and yellow lines are the eigenvectors and as you can see they are orthogonal and align with the object. The only problem with this is that they can flip as the object changes and moves but I guess if you were to use them to orient a bounding box that wouldn't matter much (?)

Soothsayer
Actually, it's ok, I managed to get it working.)

Let us know if you implement it.

Soothsayer
NumPy would be able to help me out there, right? I can't figure out how to import it into houdini python. Any ideas?

I supose you are on Windows like me. I had the same problem.

Swann_

Soothsayer
I supose you are on Windows like me. I had the same problem.

Yes. windows. It's annoying isn't it? I had to put those numpy folders inside the houdini directory.

I think I may be able to get those bounding boxes, and will post the script then. I only interpreted the algorithm, not translated it literally, so it's a tiny bit different (maybe a bit slower, but still fast enough). I think I'll have to use those eigenvectors as my new coordinate system and then get the min and max positions of each direction and that should be it. I hope. I'll give it a go tomorrow.

Soothsayer
Actually, it's ok, I managed to get it working.

Good work Soothsayer!

Hmph, I can't really get rid of those axis flips. The procedure is referenced but not described in the book. It's something to do with the way the vectors are ordered, but I don't understand it well.

If anybody wants to have a go at it, below is the script. Just append it to some geometry.

It also creates 3 vector attributes that, if visualized in the display options, show the axis.

edit: cleaned up code


import numpy as np
# This code is called when instances of this SOP cook.
geo = hou.pwd().geometry()
# create lists containing x,y,z position
p=[ ,,]
points = geo.points()
for i in points:
p.append( i.position() )
p.append( i.position() )
p.append( i.position() )

# create a numpy array containing point positions
m = np.array([p,p,p])

#calculate covariance matrix and eigenvalues
matrix_covariant = np.cov(m)
eigen = np.linalg.eig(matrix_covariant)
eigenmatrix = eigen
eigenvalues = eigen

#eigenvectors are the basis vectors of the new coordinate system that is aligned with the object's “direction”
eigenvec = [eigenmatrix,eigenmatrix,eigenmatrix]

#I think we need to order those vectors (they swap) but this isn't having the desired effect
#eigenvec = [eigenvec for i,x in sorted(enumerate(eigenvalues), key=lambda pair: pair, reverse=True)]

#Calculate the center of the geometry
centroid = hou.Vector3([ sum(p)/len(p), sum(p)/len(p), sum(p)/len(p) ])

#Adds a box into the network
def makebox():
thisobj = hou.node(hou.pwd().path())
thispath = “/obj/”+thisobj.name()
thisnode = hou.node(thispath)
#Check if box already exists
if thisnode.glob(“myBoundingBox”) !):
pass
else:
col = hou.Color((0.6,0.1,0.1))
box = thisobj.createNode(“box”)
box.setName(“myBoundingBox”)
boxtransform = thisobj.createNode(“xform”)
boxtransform.setName(“myBox_xForm”)
boxtransform.setFirstInput( box )
boxtransform.setTemplateFlag(1)
boxtransform.setColor(col)
box.setColor(col)

#supposed to transform the bounding box accordingly. Not working well yet.
def movebox():
makebox()
transform = hou.node(“/obj/grid_object1/myBox_xForm”)
pt =
pr =

pt.set( centroid )
pt.set( centroid )
pt.set( centroid )

# 3 vector attributes to visualize the new coordinate system. Not necessary but good to see what's going on.
# visualize axis_x, axis_y, axis_z in display properties to see them
axisx = geo.addAttrib(hou.attribType.Point,“axis_x”,[eigenvec, eigenvec, eigenvec])
axisy = geo.addAttrib(hou.attribType.Point,“axis_y”,[eigenvec, eigenvec, eigenvec])
axisz = geo.addAttrib(hou.attribType.Point,“axis_z”,[eigenvec, eigenvec, eigenvec])

###———————————


movebox()

Soothsayer
The procedure is referenced but not described in the book. It's something to do with the way the vectors are ordered, but I don't understand it well.

It's in detailed explained here:

http://books.google.pl/books?id=bfcLeqRUsm8C&pg=PA218&dq=covariance+matrix+in+games&cd=1#v=onepage&q=&f=false [books.google.pl]

I don't have time to look at it, but it seems your code follows the recipe (some optimisation on hom side and array initialisation are possible though). Maybe bounding box construction should be approached differently? Apparently they look for points with min/max position along eigenvector and from that they build bounding box.

Yeah….hmmm….if I'm not mistaken, I tried it in this part of my code (you need to un-comment it)

#I think we need to order those vectors (they swap) but this isn't having the desired effect
#eigenvec = [eigenvec for i,x in sorted(enumerate(eigenvalues), key=lambda pair: pair, reverse=True)]

But it isn't working. Hmmm, hmmm, hmmm.


Edit: I rewrote the covariance matrix calculations and changed the way it calculates the eigenvectors. Now my eigenvalues and covariance matrix are in line with the book-examples. The vectors spanning the coordinate system have stopped rotating, but they still mirror reverse. If that is a bug, or expected behavior I don't know yet. At least it is more stable, and if not animated should give correct results.

I'm closer and closer to getting this working.

The bounding box now translates properly. It rotates as well. The bounding isn't perfect yet but it kind of works.

edit: deleted. problem solved.

Here is what it does:

http://www.vimeo.com/10253538 [vimeo.com]

And the whole thing, if you want to try it, is here:


import numpy as np

# This code is called when instances of this SOP cook.
geo = hou.pwd().geometry()

# create lists containing x,y,z position
p=[ ,,]
points = geo.points()
for i in points:
p.append( i.position() )
p.append( i.position() )
p.append( i.position() )

# create a numpy array containing point positions
m = np.array([p,p,p])

#Calculate the center of the geometry
length = len(p)
centroid = hou.Vector3([ sum(p)/length, sum(p)/length, sum(p)/length ])

#calculate covariance matrix
var_x = np.var(p,ddof = 0)
var_y = np.var(p,ddof = 0)
var_z = np.var(p,ddof = 0)
cov_xy = np.cov(p,p,bias = 1)
cov_xz = np.cov(p,p,bias = 1)
cov_yz = np.cov(p,p,bias = 1)

matrix_covariant = np.array([
,
,
])

#Calculate the eigenvalues and eigenvectors
eigen = np.linalg.eigh(matrix_covariant)
eigenmatrix = eigen
eigenvalues = eigen

#eigenvectors are the basis vectors of the new coordinate system that is aligned with the object's “direction”
eigenvec = [eigenmatrix,eigenmatrix,eigenmatrix]

#Order the eigenvectors so that if eigenvalues correspond to eigenvector then
#eigenvalues>eigenvalues>eigenvalues
eigenvec = [eigenvec for i,x in sorted(enumerate(eigenvalues), key=lambda pair: pair, reverse=True)]


#Adds a box into the network
def makebox():
obj = hou.node(“/obj”)
#check for existence
if hou.node(“/obj”).glob(“Bounding*”) !):
pass
else:
bounding_obj = obj.createNode(“geo”,“Bounding”, run_init_scripts=False)
bounding_obj.createNode(“box”)


#Transform the bounding box to the calculated eigenvectors. Not working well yet.
def movebox():
makebox()
obj = hou.node(“/obj”)
bbox_object = obj.node(“Bounding”)
scale =

#define vectors and create target transformation matrix
v0 = eigenvec #p is the x axis ?
v1 = eigenvec
v2 = eigenvec
v3 = centroid

Mp = hou.Matrix4((
(v0,v0,v0,0),
(v1,v1,v1,0),
(v2,v2,v2,0),
(v3,v3,v3,1))) #translation (centroid)

#set the bounding box transformation matrix to the target transformation matrix
bbox_object.setWorldTransform(Mp)

#find the bounding box size. Still a bit rough
p_trans=[ , , ]
for p in points:
pos = p.position()
pos = pos*Mp # transform to new space
p_trans.append( pos )
p_trans.append( pos )
p_trans.append( pos )


diff = [max(p_trans) - min(p_trans), max(p_trans) - min(p_trans), max(p_trans) - min(p_trans) ]


scale.set( diff )
scale.set( diff )
scale.set( diff )



# 3 vector attributes to visualize the new coordinate system.
axisx = geo.addAttrib(hou.attribType.Point,“axis_x”,[eigenvec, eigenvec, eigenvec])
axisy = geo.addAttrib(hou.attribType.Point,“axis_y”,[eigenvec, eigenvec, eigenvec])
axisz = geo.addAttrib(hou.attribType.Point,“axis_z”,[eigenvec, eigenvec, eigenvec])

###———————————

Below is the near final version. The bounding box is now calculated correctly as far as I can see. It's also faster than before, I think.

The points are projected onto the eigenvectors and the min max of the length of these is used to determine the bounding boxes.

video:

http://www.vimeo.com/10274660 [vimeo.com]


Edit: cleaned up code


import numpy as np

# This code is called when instances of this SOP cook.
geo = hou.pwd().geometry()

# create lists containing x,y,z position
p=[ ,,]
points = geo.points()
for i in points:
p.append( i.position() )
p.append( i.position() )
p.append( i.position() )

# create a numpy array containing point positions
m = np.array([p,p,p])

#Calculate the center of the geometry
length = len(p)
centroid = hou.Vector3([ sum(p)/length, sum(p)/length, sum(p)/length ])

#calculate covariance matrix
var_x = np.var(p,ddof = 0)
var_y = np.var(p,ddof = 0)
var_z = np.var(p,ddof = 0)
cov_xy = np.cov(p,p,bias = 1)
cov_xz = np.cov(p,p,bias = 1)
cov_yz = np.cov(p,p,bias = 1)

matrix_covariant = np.array([
,
,
])

#Calculate the eigenvalues and eigenvectors
eigen = np.linalg.eigh(matrix_covariant)
eigenmatrix = eigen
eigenvalues = eigen

#eigenvectors are the basis vectors of the new coordinate system that is aligned with the object's “direction”
eigenvec = [eigenmatrix,eigenmatrix,eigenmatrix]

#Order the eigenvectors so that if eigenvalues correspond to eigenvector then
#eigenvalues>eigenvalues>eigenvalues
eigenvec = [eigenvec for i,x in sorted(enumerate(eigenvalues), key=lambda pair: pair, reverse=True)]


#Adds a box into the network
def makebox():
obj = hou.node(“/obj”)
#check for existence
if hou.node(“/obj”).glob(“Bounding*”) !):
pass
else:
bounding_obj = obj.createNode(“geo”,“Bounding”, run_init_scripts=False)
bounding_obj.createNode(“box”)

#Transform the bounding box to the calculated eigenvectors.
def movebox():
makebox()
obj = hou.node(“/obj”)
bbox_object = obj.node(“Bounding”)

#calculate scale vector (by projecting points onto eigenvectors and taking the max and min)
basis = [hou.Vector3( eigenvec ),hou.Vector3( eigenvec ),hou.Vector3( eigenvec )]
tlength=[ , , ]
for i in points:
p_vector = i.position()
difference_vector = p_vector - centroid
for j in basis:
c = basis.index(j)
tlength.append( difference_vector.dot(j) )

s = [ np.ptp( tlength), np.ptp( tlength), np.ptp( tlength) ] #ptp gets the range of an axis

# set the scale matrix
Ms = hou.Matrix4((
(s,0,0,0),
(0,s,0,0),
(0,0,s,0),
(0, 0, 0, 1)))

#define vectors and create target transformation matrix
v0 = eigenvec
v1 = eigenvec
v2 = eigenvec
v3 = centroid

Mp = hou.Matrix4((
(v0,v0,v0,0),
(v1,v1,v1,0),
(v2,v2,v2,0),
(v3,v3,v3,1))) #translation (centroid)


#set the bounding box transformation matrix to the target transformation matrix
Mp = Ms*Mp
bbox_object.setWorldTransform(Mp)

###———————————


movebox()