3D L system orientation problem

Hi,

I’m reading a book called “Algorithmic beauty of plants” (you can read the book online at this link) where author explains how to define tree structures using L system. To understand this better i’m trying to implement this in open frame works. I even got the 2D trees working.
see this tree i successfully generated :

2D example tree using L system(image)

To proceed into 3D there are simple rules

F : move forward and draw line forward
f : move forward without drawing line forward
+ : turn left by angle, using rotation matrix Ru(-angle)
- : turn right by angle, using rotation matrix Ru(angle)
& : pitch down by angle, using rotation matrix Rl(-angle)
^ : pitch up by angle, using rotation matrix Rl(angle)
\ : roll left by angle, using rotation matrix Rh(-angle)
/ : roll right by angle, using rotation matrx rh(angle)
| : turn around using rotation matrix Ru(-180)
[ : Push the current state of the turtle onto a pushdown
    operations stack.The information saved on the stack contains the
    turtle's position and orientation, and possibly other attributes
    such as the color and width of lines being drawn.
] : Pop a state from the stack and make it the current state
    of the turtle.No line is drawn, although in general the
    position of the turtle changes.
Ru(angle) : 
            cos(angle)  sin(angle)      0
           -sin(angle)  cos(angle)      0
               0            0           1
Rl(angle):
           cos(angle)      0      -sin(angle)
               0           1           0
           sin(angle)      0      cos(angle)
Rh(angle)
               1           0           0
               0      cos(angle)   -sin(angle)
               0      sin(angle)    cos(angle)

from the book he gave a example tree and its rules for construction

3D tree and its rules for construction(image)

But when i try to construct the tree the top portion of the tree is pointing towards down instead of pointing up like in the above image. I haven’t implemented the coloring yet so my output is just white lines but still the structure should look same at least. But its not and i cant figure out where i did wrong.

My output for the same formulae(image)

Here is the code i wrote for creating lstring

in ofApp.cpp

string str = lString.GetLString(5, 8, 
		string("plant "), 
		string("plant = internode + [ plant + flower ] - - // [ - - leaf ] internode [ + + leaf ] - [ plant flower ] + + plant flower "), 
		string("internode = F seg [// & & leaf ] [// ^ ^ leaf ] F seg "), 
		string("seg = seg F seg "), 
		string("leaf = [' { +f-ff-f+ | +f-ff-f } ]"),
		string("flower = [ & & & pedicel ' / wedge //// wedge //// wedge //// wedge //// wedge ]"),
		string("pedicel = FF"),
		string("wedge = [' ^ F ] [ { & & & & -f+f | -f+f } ]"));

	n1 = lString.GenerateQuad( str, 18.0f );

and the code inside LString class

string LString::GetLString(int pIteration, int pRulesCount, ...){
	std::vector<string> rules;

	va_list ap;
	va_start(ap, pRulesCount);

	// first argument will be axiom
	string axiom = va_arg(ap, string);

	// after axiom remaining all arguments will be rules
	for (int i = 0; i < pRulesCount - 1; i++) {
		string rule = va_arg(ap, string);
		rules.push_back(rule);
	}
	va_end(ap);

	string result = axiom;
	
	// looping for number of iterations
	for (int itr = 0; itr < pIteration; itr++){

		// replacing the variables with equations from the string
		for (int ruleCount = 0; ruleCount < pRulesCount - 1; ruleCount++){
			
			string rule = rules[ruleCount];
			string variable = "";
			string equation = "";
			bool isFindingVariable = true;

			// seperating the variable and equation from
			// the rule
			for (int i = 0; i < rule.size(); i++){
				if (rule[i] != '='){
					if (isFindingVariable){
						variable += rule[i];
					}
					else{
						equation += rule[i];
					}
				}
				else{
					isFindingVariable = false;
				}
			}

			// now replacing the variable with equation
			result = replace(result, variable, equation);
		}
	}

	return result;
}

string LString::replace(string pSource, string pDelimiter, string pReplaceString){
	string result="";
	for (int i = 0; i < pSource.size(); i++){
		string tmp = pSource.substr( i, pDelimiter.size() );
		if (tmp.compare(pDelimiter) == 0){
			result += pReplaceString;
			i = i + pDelimiter.size() - 1;
		}
		else{
			result += pSource.substr(i, 1);
		}
	}

	return result;
}

Quad LString::GenerateQuad(string pLString, float pAngle){

	Quad quad;

	// this unit vector is only for direction
	ofVec3f unitDirectionVector(0.0f, -1.0f, 0.0f);
	LLine prevLine(ofVec3f(0.0f, -1.0f, 0.0f), ofVec3f(0.0f, 0.0f, 0.0f));

	float m_CurrentLineLength = 10.0f;
	float rotationAngle = pAngle * (M_TWO_PI / 360.0f);

	LStringStack.clear();

	for (int i = 0; i < pLString.size(); i++){
		char operation = pLString[i];

		if (operation == 'F'){
			// this is for batching the 'F'. if there are four 'F'
			// instead of drawing 4 small lines drawing one big line is
			// better for performance
			/*int numOfF = 0;
			for (int j = i; j < pLString.size(); j++){
				if (pLString[j] == 'F')	numOfF++;
				else break;
			}
			i = i + numOfF - 1;*/

			ofVec3f start = prevLine.end;
			ofVec3f end = start + m_CurrentLineLength  * unitDirectionVector;

			LLine tmpLine(start, end);

			prevLine = tmpLine;

			quad.AddLine(tmpLine);
		}
		else if (operation == 'f'){
			ofVec3f start = prevLine.end;
			ofVec3f end = start + m_CurrentLineLength * unitDirectionVector;

			LLine tmpLine(start, end);
			prevLine = tmpLine;
		}
		else if (operation == '+'){
			//  cos(angle)  sin(angle)      0
			// -sin(angle)  cos(angle)      0
			//     0            0           1
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = cosf( -rotationAngle );     rotationMatrix.b = sinf( -rotationAngle );     rotationMatrix.c = 0;
			rotationMatrix.d = -sinf( -rotationAngle );    rotationMatrix.e = cosf( -rotationAngle );     rotationMatrix.f = 0;
			rotationMatrix.g = 0;                          rotationMatrix.h = 0;                          rotationMatrix.i = 1;

			//unitDirectionVector = unitDirectionVector.rotateRad(-rotationAngle, ofVec3f(0.0f, 0.0f, 1.0f));
			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
		}
		else if (operation == '-'){
			//  cos(angle)  sin(angle)      0
			// -sin(angle)  cos(angle)      0
			//     0            0           1
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = cosf(rotationAngle);       rotationMatrix.b = sinf(rotationAngle);     rotationMatrix.c = 0;
			rotationMatrix.d = -sinf(rotationAngle);      rotationMatrix.e = cosf(rotationAngle);     rotationMatrix.f = 0;
			rotationMatrix.g = 0;                         rotationMatrix.h = 0;                       rotationMatrix.i = 1;

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
			//unitDirectionVector = unitDirectionVector.rotateRad(rotationAngle, ofVec3f(0.0f, 0.0f, 1.0f));
		}
		else if (operation == '&'){
			//  cos(angle)      0      -sin(angle)
			//      0           1           0
			//  sin(angle)      0      cos(angle)
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = cosf(-rotationAngle);     rotationMatrix.b = 0;      rotationMatrix.c = -sinf(-rotationAngle);
			rotationMatrix.d = 0;                        rotationMatrix.e = 1;      rotationMatrix.f = 0;
			rotationMatrix.g = sinf(-rotationAngle);     rotationMatrix.h = 0;      rotationMatrix.i = cosf(-rotationAngle);

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
			//unitDirectionVector = unitDirectionVector.rotateRad(PI-rotationAngle, ofVec3f(0.0f, 1.0f, 0.0f));
		}
		else if (operation == '^'){
			//  cos(angle)      0      -sin(angle)
			//      0           1           0
			//  sin(angle)      0      cos(angle)
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = cosf(rotationAngle);     rotationMatrix.b = 0;     rotationMatrix.c = -sinf(rotationAngle);
			rotationMatrix.d = 0;                       rotationMatrix.e = 1;     rotationMatrix.f = 0;
			rotationMatrix.g = sinf(rotationAngle);     rotationMatrix.h = 0;     rotationMatrix.i = cosf(rotationAngle);

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
			//unitDirectionVector = unitDirectionVector.rotateRad(rotationAngle, ofVec3f(0.0f, 1.0f, 0.0f));
		}
		else if (operation == '\\'){
			//      1           0           0
			//      0      cos(angle)   -sin(angle)
			//      0      sin(angle)    cos(angle)
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = 1;      rotationMatrix.b = 0;                        rotationMatrix.c = 0;
			rotationMatrix.d = 0;      rotationMatrix.e = cosf(-rotationAngle);     rotationMatrix.f = -sinf(-rotationAngle);
			rotationMatrix.g = 0;      rotationMatrix.h = sinf(-rotationAngle);     rotationMatrix.i = cosf(-rotationAngle);

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
		}
		else if (operation == '/'){
			//      1           0           0
			//      0      cos(angle)   -sin(angle)
			//      0      sin(angle)    cos(angle)
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = 1;    rotationMatrix.b = 0;                      rotationMatrix.c = 0;
			rotationMatrix.d = 0;    rotationMatrix.e = cosf(rotationAngle);    rotationMatrix.f = -sinf(rotationAngle);
			rotationMatrix.g = 0;    rotationMatrix.h = sinf(rotationAngle);    rotationMatrix.i = cosf(rotationAngle);

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
		}
		else if (operation == '|'){
			//  cos(angle)  sin(angle)      0
			// -sin(angle)  cos(angle)      0
			//     0            0           1
			ofMatrix3x3 rotationMatrix;
			rotationMatrix.a = cosf(-PI);     rotationMatrix.b = sinf(-PI);    rotationMatrix.c = 0;
			rotationMatrix.d = -sinf(-PI);    rotationMatrix.e = cosf(-PI);    rotationMatrix.f = 0;
			rotationMatrix.g = 0;             rotationMatrix.h = 0;            rotationMatrix.i = 1;

			unitDirectionVector = RotateVector(rotationMatrix, unitDirectionVector);
		}
		else if (operation == '['){
			LStringStack.push_back( LLineStack(prevLine, unitDirectionVector) );
		}
		else if (operation == ']'){
			LLineStack tmpLine = LStringStack.back();
			LStringStack.pop_back();
			prevLine = tmpLine.line;

			unitDirectionVector = tmpLine.direction;
		}
	}

	return quad;
}
// multiply the vector with rotation matrix and return the resultant vector
ofVec3f LString::RotateVector(ofMatrix3x3 pRotationMatrix, ofVec3f pVector){
	ofVec3f result;
	result.x = pRotationMatrix.a * pVector.x + pRotationMatrix.b * pVector.y + pRotationMatrix.c *  pVector.z;
	result.y = pRotationMatrix.d * pVector.x + pRotationMatrix.e * pVector.y + pRotationMatrix.f *  pVector.z;
	result.z = pRotationMatrix.g * pVector.x + pRotationMatrix.h * pVector.y + pRotationMatrix.i *  pVector.z;

	return result;
}

I uploaded the whole project in github. You can download the project here.
There must be something very small i’m missing. And I looked at this code for past one week and couldn’t find what i did wrong. :confounded:
Thanks in advance.

one thing that can trip folks up is that our top left corner is 0,0, which means angular rotation (unless you alter it) instead of going counter clockwise goes clockwise. I haven’t looked at your code in much detail, but wonder if some general rotational inversion might be causing this?

Hi Plato,

It so happens that I’ve also been working on an L-systems engine based in part on the ABOP book. I plan to release it as an ofx addon within the next month or two, once I’ve finished adding support for parametric L-systems.

I haven’t studied your code too closely but I suspect your problem stems from the fact that you’re only maintaining a single direction vector for your 3D turtle. You need to keep track of three vectors for its orientation: a forward, an up, and a left vector, and your rotation operations need to adjust all three of them. (Each operation will adjust two of the three vectors.)

@zach
Yeh i considered that. that is the reason i inverse the rotations.

@Dewb
hanks for reply, there is something new i didn’t know. Since rotation always happens on one axis remaining other two axis will be effected by that rotation. I thought having a vector and rotating it will give the same effect. If I have 3 unit vectors in all three directions then should the next point be calculated like this

ofVec3 end;
end.x = start.x + (unitVectorAlongX * m_CurrentLineLength).x
end.y = start.y + (unitVectorAlongY * m_CurrentLineLength).y
end.z = start.z + (unitVectorAlongZ * m_CurrentLineLength).z

and with every operation rotate all three unit vectors instead of one vector?

The three vectors aren’t always along X, Y, and Z (although they might start out that way.) They’re forward, up, and left to specify the current orientation of the turtle. Each operation will rotate two of the vectors around the third. Turning left rotates forward and left around up, pitching down rotates forward and up around left, etc.

The up and left vectors are to make sure the notion of “turning left” is relative to the turtle, not to the world. When you move forward to calculate the next position, you’ll still just use the forward unit vector.

Great. I worked a lot on unity so I’m very familiar with concepts of local and global coordinates. But I never thought that all these operations should be done on local coordinates. I changed my code so that the operations should be done on local coordinates (forward, left and up). Output looks a lot like how it is supposed to look.

Expect for the leaves and flowers. Flowers looks like this.

In the book author says

Production p4 specifies the leaf as a filled polygon with six edges. Its boundary is formed from the edges f enclosed between the braces { and } (see Chapter 5 for further discussion).

I guess i have to read till chapter 5 to get this correct. Because I still haven’t implemented the { and } operations yet. or is there still something wrong I’m doing here?

Yep, that looks correct; you’re rendering the branches and stamen filaments correctly, you’re just missing the leaves and petals. You will indeed need { } polygons implemented to see those in this ruleset.