Prime Factorization

This script uses a brute force approach.

Number	Primes

Matrix to Euler Angles

Every orientation in 3D space can be described using the product \(T\) of three rotation matrices, usually around the main axes \(R_X\), \(R_Y\) and \(R_Z\). Given such a matrix \(T\), it is often necessary to decompose it into three seperate angles - the Euler Angles. These can be used as input arguments for a sequence of basic rotation operations, so that the following equation holds:

\begin{align} T = R_i(\alpha_0) R_j(\alpha_1) R_k(\alpha_2) \end{align}

The choice of \(i, j, k\) determines the rotation convention. One common possibility would be Yaw-Pitch-Roll (XYZ, depending on the local coordinate system) but the following implementation works for arbitrary conventions.

Mathematics

The decomposition of \(T\) depends highly on the definition (internal implementation) of \(R_i\) so make sure you know what your matrix class is doing. The standard definition of rotation matrices around X, Y and Z looks like this:

\begin{align} R_X=\begin{pmatrix} 1 & 0 & 0 \\ 0 & A & -B \\ 0 & B & A \end{pmatrix} \end{align}

\begin{align} R_Y=\begin{pmatrix} C & 0 & D \\ 0 & 1 & 0 \\ -D & 0 & C \end{pmatrix} \end{align}

\begin{align} R_Z=\begin{pmatrix} E & -F & 0 \\ F & E & 0 \\ 0 & 0 & 1 \end{pmatrix} \end{align}

Here \(A, C, E\) denote the cosine and \(B, D, F\) the sine of the unknown Euler Angles. Let's take a look at \(T\) for the convention YXY.

\begin{align} T=R_{YXY}= \begin{pmatrix} CC-DDA & DB & CD+DAC \\ DB & A & -CB \\ -CD-DAC & CB & -DD+CCA \end{pmatrix} \end{align}

Apparently the element (indexing starting with 0) at \(T(1,1)\) is the cosine of the Euler Angle \(\alpha_1\) for the middle rotation matrix \(R_X\). This is true for all of the 12 possible conventions. However the location of this pivot element is not always the same. The coordinates are determined by the first and third rotation component by using the simple conversion \(X \rightarrow 0\), \(Y \rightarrow 1\) and \(Z \rightarrow 2\). The next step is to compute \(\alpha_0\) and \(\alpha_2\) using the elements horizontal and vertical of the pivot. Knowing \(A\) one can easily compute \(\alpha_1\) and therefore \(B\). You can now use the horizontal neighbors of the pivot element to get \(\alpha_2\) and the vertical neighbors to get \(\alpha_0\). Since the applied convention uses the same axis for the first and third rotation there are no \(E, F\) components in this matrix.

Gimbal Lock

Now what happens if \(B\) is close to zero in the above case? This means a so called Gimbal Lock has occured because the divisions in the row and column of the pivot element are undefined. This can happen if the orientation of \(T\) is so that one angle becomes ambiguous which means it can be set to any value, e.g., \(\alpha_0 = 0\). Since one Euler Angle is now zero the corresponding rotation matrix yields \(R_i=I\) and the new reduced transformation matrix looks like this:

\begin{align} T=R_{XY}= \begin{pmatrix} C & 0 & D \\ DB & A & -CB \\ -DA & B & CA \end{pmatrix} \end{align}

Interpreting the axes as coordinates again a new pivot element is found at \(T(0,1)\), which is always zero. Once again it's horizontal neighbors can be used to determine the remaining Euler Angle \(\alpha_2\). Keep in mind that the order of cosine and sine elements may vary depending on the convention applied. As previously mentioned there are 12 possibilities.

:: CODE BLOCK ::

\(R_X\) - pivot	\(R_Y\) pivot	\(R_Z\) pivot
Y-X-Y	X-Y-X	X-Z-X
Y-X-Z	X-Y-Z	X-Z-Y
Z-X-Y	Z-Y-X	Y-Z-X
Z-X-Z	Z-Y-Z	Y-Z-Y

Listing: A - "Conventions"

If you want to construct an Euler Angle decomposition for a special convention just follow these instructions to extract the necessary equations. I recommend a symbolic calculation, e.g., in MatLab to make this task easier.

City Gen

City Gen was a group project of 2008 based on OpenGL using advanced techniques, mainly Shader (implemented Phong shading) and L-Systems. It demonstrates the procedural generation of a virtual city with regard to several environmental parameters.

The scene consists of an initially planar mesh which can be modeled using the elevation tool by simply brushing over the target area. Brush diameter and pressure can be adapted. Figure A shows the subsequent steps in particular the application of the population map, the automatic generation of streets and finally the extraction of the blue blocks representing buildings and other structures.

Examples 1 to 3 show possible output after the street generation. The street generator uses a two phase L-System. In a preprocessing step up to 3 population centers are calculated. These centers are then connected using a "highway"-style L-System. Finally each street node spawns new street segments based on the population density, terrtain height and conflicting street segments. New street segments vary in length and angle and are placed so that they fit best into the existing street network for example by avoiding dead ends.

A more detailed description may follow in the future. For now, take a look at the source which will be made accessible if requested.

Binary Filter

Have you ever wondered how the sequence of natural numbers looks like if you filter out all the numbers that don't match a given number of ones or zeros in its binary representation? Well, i did, a while ago and that's why i came up with this little piece of code, the Binary Filter. Although it may not seem like it, the solution for this particular problem was somewhat tricky to obtain, hence i decided to share it with the world.

The key idea to my solution is the observation that the problem "Give me all numbers with \(k\) ones in its binary representation" is very similar to calculating the binomial coefficient. Consider the following table containing all 5 digit numbers with 3 ones.

:: TABLE BLOCK ::

\(i\)	\(N_{10}\)	\(N_2\)
0	28	11100
1	26	11010
2	25	11001
3	22	10110
4	21	10101
5	19	10011
6	14	01110
7	13	01101
8	11	01011
9	7	00111

Listing: A - "Binary Filter in action"

Naturally, the number of items equals the binomial coefficient 5 over 3.
In order to calculate the correct number at an index \(i \in [0,9]\), we need to know how many ones are at the top and how many zeros are at the bottom for the first column. This happens to be the quotient of \(3 \over 5\) times the binomial coefficient 5 over 3. This is true for any combination of \(n\) and \(k\). Finally, after knowing whether we're in the upper or lower half of the table, we can dismiss the first column and analyze the remaining \(n-1\) columns.

The Code

As usual, here the necessary code snippets for your own project. You're going to need an implementation of the binomial coefficient. Instead of using the standard version, you might want to use this optimized version without explicitly calling the factorial function.

:: CODE BLOCK ::

int binomial(int _n, int _k)
{
    if(_k < _n-_k) _k = _n-_k;
    int c = 1;

    for(int i=0; i<_k; ++i){
        c = c*(_n-i);
        c = c/(i+1);
    }

    return c;
}

Listing: B - "Binomial Coefficient"

And finally the actual binary filter function in a recursive implementation. The function will return the \(i\)-th, \(n\) digit number with \(k\) ones in its binary representation. The order of values is descending.

:: CODE BLOCK ::

int bfilter(int _n, int _k, int _i)
{
    if(_n < =0) return 0;

    //Number of "1"s
    int n_x = (_k*binomial(_n, _k))/_n;

    //_i is in "1" block
    if(_i < n_x){
        return pow(2, (_n-1)) + bfilter(_n-1, _k-1, _i);
    }
    //_i is in "0" block
    else{
        return bfilter(_n-1, _k, _i-n_x);
    }
}

Listing: C - "Binary Filter"

Check out the ▷docs for further information.