stuff growing and crawling - Hyphae - OTL version

I've been playing with path finding algorithms recently and I came up with an OTL that creates a venation structure given root points, seed points and an optional collision surface.

This can be used to create venations on, inside surfaces, and I guess some more cool stuff that at the moment I cannot even think of !

This is a quick demonstration video of the otl (I speeded it up 2x so you don't fall asleep while you watch)

A few days ago I was reading this article cause I was interested in reproducing Jack Frost effect from Rise of the Guardians.
This line of the article really caught my attention:

Our developer on that would take a model, whatever the frost needed to grow on, and he would run a cellular automata simulation across it

Cellular automata ? what the hell is that ? So I started researching a bit and as it always happens...the by-products of the research are usually even more interesting than the destination.

So I came across this video:

Hyphae - growth process diagram in 2D from Nervous System on Vimeo.

Check out this blog out when you have time, there is a lot of really interesting stuff, presented in a really nice way. 

Researching a bit more I found the details of the algorithm in this amazing paper:

http://algorithmicbotany.org/papers/venation.sig2005.pdf

(thank you Adam Runions, Martin Fuhrer, Brendan Lane, Pavol Federl, Anne−Gaëlle Rolland−Lagan, and Przemyslaw Prusinkiewicz for this amazing paper !!!)

This paper blew my mind. Bye Jack Frost (for now), welcome Hyphae ! 

So I quickly implemented a prototype in Houdini.

The Algorithm

Given this initial data:

• hormone points (H)
• roots points (R)

The target of the algorithm is to add points to R till the branches reach all the points in H. When one of the R branches reaches one of the H points, that H point will be removed by the group H and not used anymore for the next iteration.

So the way I implemented it in Houdini is the following:

1. for each point h in H
   1. find the closest point cp in R and store (on h) its index (integer CLIDX) , location (vector CLPOS) and normalized vector pointing to it from h (vector CLDIR)
   2. if the cp is too close to h flag it for deletion (int IGNORE=1)
2. delete all points in H where IGNORE==1
3. in H create groups GR_* based on CLIDX (using partition sop)
4. for each group gr in GR_*
   1. find the average CLDIR and normalize it --> AVGCLDIR
   2. create a new point r=CLPOS-AVGCLDIR
   3. merge the new r point into R
5. goto 1

Test 1

"let's see if it works" 

H (hormones) = 100 points scattered by a sphere volume centered at the origin 

R (roots) = 1 point at the origin

 Test 2 

"hell yeah it does ! ...how about more roots ??" 

H = same as before 

R = 3 points (manually placed randomly in the space not too far from the sphere I used to generate H)

Test 3 

"this is amazing. How about a lot more roots ?"

H = same as before 

R = points scattered from the surface of a sphere larger than the one I used to generate H

Uhmm...that is not so cool. Why ?

Problem 1: 
Well, to start we notice that not all of the points R actually grow. This is because of the step 1.1 in the algorithm ! Only one point between all R is chosen to be the closest one to a specific point in H. So some point of R might not be the closest to ANY H since the beginning.

Solution: 

1. we could decrease the number of R but that's cheating , isn't it.

2. we can increase the radius and number of the H distribution so that the outer H points are closer to the roots.

Problem 2: 
The veins grow radially towards the center, but as soon as they find a hormone, they stop without splitting any new branch. So we don't get that awesome intricate complexity we like so much and that doesn't let us sleep on night. The reason for this is point 2 in the algorithm. What happens is that we have so many veins now, that they take care of 'killing' the few H very quickly (and more or less at the same time, cause the distribution of the chasers (R) is similar to the distribution of the chased (H), since both are spherical distributions).

Solution: 

We need to add more H. The best option is to add more density the more we move towards the center of the spherical distribution of H. A quick way to do that is use Poly to VDB sop, check "fill interior" and scatter from this.

Problem 3: 
Lines are too straight ! We like chaos, variation !

Solution:
How about adding in the algorithm some noise ?
Let's modify the step 4.B in the algorithm adding some noise fuction of the position.
Something like this:
4.2 --> r=CLPOS-AVGCLDIR+noise(CLPOS)

Test 4

H = 5000 points scattered from VDB volume with interior filled.
R = same as test 3
I've added a noise function as described in the solution of problem 3 above.

Collision Objects

This algorithm allows creating any sort of volumetric venation , which is really cool.

So I wanted to test it on a real complex model, and see if I could create the effect of venations growing and crawling over the surface of an object.

So I took the model of an arm from internet , I scattered points over the surface of the arm (H points), and I placed one single root point on the tip of the index finger. When I ran the simulation I noticed that, even if the result was mind blowing, there was a big problem : the venation was growing even "inside" the arm not only on the surface. I had to expect this of course, since the each loop iteration of the algorithm is performed on ALL the H points. So it doesn't really matter if the root is initially on the index finger tip : it will be affected even from the H points on the fore arm. So the venation would grow regardless.

You can see it even from another point of view : the algorithm is not aware about the surface of any model ! All we feed into him are point clouds.

So I decided to include in the algorithm the concept of Collision Objects: when a root point enters in a collition object (SDF sample value < 0) then it will be pushed out of the collision object along the normal of the collision object in that point (SDF gradient value).

So all we need is modify once more the step 4.2 of the algorithm:

4.2 --> r=CLPOS - AVGCLDIR + noise(CLPOS) + CVEC

where

CVEC = if ( sample (CLPOS,cobj) < 0 , -sample(CLPOS, cobj)*gradient(CLPOS,cobj) , 0 )

cobj is the collision object

Now CVEC component will make sure that the roots will not grow into the collision object.

What if we want to make sure the roots never leave the surface ?

For instance in this case, we don't want to see roots growing in the air. We want them to stick on the surface.

That's easy to fix : we just have to remove the condition to the CVEC component.
In other words, now the roots will be pushed outside of the collision object if they're inside, and towards the collision object if they're outside.

CVEC = -sample(CLPOS, cobj)*gradient(CLPOS,cobj)

This will make sure that the roots stick to the collision object.

In the following example, I duplicated the arm , poly-extruded it out so that the copy is slightly bigger than the arm. Then I selected the hand and a few more polygons on the arm and poly-extruded it again so that the selection is inside the original geometry. Then used it (the copy) as a collision object.
This way I created very quickly a small path on the arm to allow the root to reach the hand which is free to be filled with veins.
This is the collision object (in red) and the arm (in yellow):

This is a quick test


IMPORTANT NOTE:
The Hyphae algorithm is great to create vines / venations structures using STATIC point clouds. What I mean is that, moving one single H or R point in the whole system can radically change the outcoming shape.

So, if you plan to have H, or R or both, animated and expect to have a smooth animation of the resulting venation system, this is not going to happen. The vines will be radically different between two adjacent frames. (this is great if you have to create some electricity effect, or lightnings, or voltaic arcs a la Frankenstein etc..).

The workaround to smoothly animate the resulting venation is to generate it on static point clouds, then deform it using some other deformer (Lattice for instance).

stuff growing and crawling - L-Systems

In this nth experiment I wanted to create a way to wrap a growing geometry on the surface of another geometry without relying upon UVs on the target geometry.

The OTL I created takes 3 inputs:

1. lsystem (this will be wrapped on top of 2)
2. poly geometry
3. one point

L-System (points and edges)

The L-System provided as first input can be any L-System where the root point of a new branch overlaps the point it generates from. In other words if you 'window' select the first point of a branch in your L-System you should end up selecting 2 points. Inside the OTL the L-System is pre-processed adding point attributes (like an attr that specifies if a certain point is a branch root, or the direction of each segment) that will make traversing it's hierarchy later on way faster. Something like this:

Wrapped Geometry (any geometry)

The OTL wraps the original L-System around this geometry. It can be any geometry. This geometry will be converted into a VDB SDF in the OTL, sample and gradient volume vops are used to quickly find the closest point on this surface. Because of this approach, the OTL doesn't require any UV coordinate to locate positions on this geometry and it's super fast since the SDF is calculate only once at the beginning of the whole simulation and contains already all the info to calculate the closest point on the surface with just a couple of multiplications and sums.

Start Point (one point)

The closest point of the Start Point on the Wrapped Geometry is where the root of the 'wrapped' L-System will sit. This point must be as close as possible to the the Wrapped Geometry surface (in order to have it within the Wrapped Geometry SDF VDB voxels).

How it works

simple otl interface

The original L-System is flattened on the XY plane before being processed. The root (the first point calculated) of the wrapped L-System (WL) is the closest point of the Start Point (SP) on the Wrapped Geometry (WG).The initial direction is given by the vector specified in the OTL UI (see image above).Every 'next' point position is calculated starting from the position of the previous point. For this reason I used a foreach loop making sure to uncheck the parameter "merge", since I want the same data to be processed and provided to the next cycle after adding a new point.To calculate the closest point on the surface from a certain point in the world space, I've used VDB SDF. SDF are already a blessing provided by Houdini but they've always had a certain level of imprecision. VDB SDF are way faster to calculate, and much more precise (thank you DreamWorks Animation for this gift !). The following is the VOP SOP structure to calculate the closest point on surface, given an SDF volume and a sample location in world space:

In this example I've created the wrapped L-System, and then used a Ray Sop to make sure each point really was on the surface, once per frame (Ray Sop is very fast).

I made sure to preserve the Width attribute calculated by default in the L-System Sop and I used it in the PolyWire Sop to drive the width of the Tubes as explained in the PolyWire documentation:

....This node works like the Wireframe node, but this node creates more complex tube geometry from curves, with smoother bends and intersections than Wireframe, especially for L-systems. The four numerical parameters support all the local variables of the Point operation, plus the LSYSTEM specific variables of $WIDTH, $SEGS, $DIV, $LAGE, $GEN, and $ARC....

WL has been calculated on a static WG (I choose the first frame of the hand animation). Then I used a Lattice Sop to wrap WL on the animated version of WG. Before applying the Lattice Sop I've calculated and stored the density of the points on a point attribute on the animated WG via VOP SOP where I used Point Cloud vex nodes to access the same geometry and return the number of points within a certain radius. Then I've used the density attribute to drive the size of the metaballs internally used by Lattice Sop.

Note: L-Systems are perfect to create intricate veins systems. Because L-Systems are a parallel rewriting systems, their size might get huge very fast. For instance in the previous example the original L-System had 95 generations and ~66k points. Nevertheless the generation of the L-System itself was very fast (less than 1s per frame for 95 generations).

OTL processing time to calculate the Wrapped L-System in this example:

Frames 1-84 (83 frames) : less than 2 minutes
Frames 85-91 (7 frames) : ~1 minute
Frames 92-96 (4 frames) : ~1 minute
Frames 97-99 (2 frames) : ~1 minute
Frames 100-102 (2 frames) : ~1 minute
Frames 103-104 (1 frame): ~1 minute
Frame 116 : ~2 minutes
Frame 122 : ~3 minutes
Frame 131 : ~5 minutes
Frame 138: ~10 minutes
Frame 145: ~21 minutes

Time to calculate 145 frames ~ 3h 40m.

on an Intel Core i7 3.0Ghz - 4 cores - 8 logical All VOP SOPs are set to 1 thread per proc.

As you can see , processing times increase exponentially when the number of L-System generation starts getting crazy so if the idea is, for instance, to cover a full animated character with veins, the best approach is divide and conquer, instead of using one single massive L-System.

Now, I guess this is just a way to accomplish this effect, I am sure there are many others I've not thought about. So, feel free to comment and pitch ideas !