# Tag: math

Entries for tag "math", ordered from most recent. Entry count: 58.

Uwaga! Informacje na tej stronie mają ponad 5 lat. Nadal je udostępniam, ale prawdopodobnie nie odzwierciedlają one mojej aktualnej wiedzy ani przekonań.

Pages: 1 2 3 4 5 6 ... 8

# Properties of Pascal Triangle

Thu
12
Nov 2009

Pascal Triangle is a mathematical object that looks like triangle with numbers arranged the way like bricks in the wall. Each next row has one more number, ones on both sides and every inner number is the sum of two numbers above it. It can span infinitely.

Despite simple algorithm this triangle has some interesting properties. First, if you draw only odd numbers, you get a fractal - Sierpinski Triangle.

Second, graph of numbers in each row resembles Gaussian distribution function. The lower row you take from the triangle (containing more numbers), the more precise graph you get.

I've made these images with Octave - a free Matlab alternative. I like it and I think concepts behind this Matlab programming language can be quite useful for doing computations, but doing graphics is painfully slow. Whole seconds to draw several dozens of texts or rectangles? WTF?!

# Calculating Linear and Quadratic Equation Coefficients

Sun
30
Aug 2009

Some time ago I've written about formulas to calculate coefficients of linear and quadratic equation [pl] when having 2 or 3 given points (x,y). Yesterday I suddenly noticed that my code needs to do it every frame so I need functions to calculate these coefficients. Below you can find the code of my functions. They are templates, so they work with many types including float, double, as well as vectors of any dimmension and other types which have addition, subtraction and assignment operators, as well as multiplication and divistion by float scalar.

# Reflect and Refract Functions

Mon
24
Aug 2009

Every game developer must have some knowledge about maths and especially geometry. Math library for game developers always contains structures like vector, matrix and functions like vector normalization or matrix inversion. But high level shader languages, designed to do geometrical computations, have many more interesting intrinsic functions. Some of more complex ones are reflect and refract, which calculate vector reflection and refraction in relation to the normal of a surface. How to implement them in C++?

Reflection works according the simple rule that "angle of incidence equals angle of reflection". Although there are no angles in the formula - all the necessary work is done by dot product.

```void Reflect(VEC3 &out, const VEC3 &incidentVec, const VEC3 &normal)
{
out = incidentVec - 2.f * Dot(incidentVec, normal) * normal;
}```

For reflect function, the formula is given in th HLSL documentation, but for refract function it is not. It took me some time to find the implementation, but finally I've found it in the GLSL Specification (thanks KriS!). So here is the code converted to C++. Eta is the refraction index (e.g. 1.5f).

```inline void Refract(
VEC3 &out, const VEC3 &incidentVec, const VEC3 &normal, float eta)
{
float N_dot_I = Dot(normal, incidentVec);
float k = 1.f - eta * eta * (1.f - N_dot_I * N_dot_I);
if (k < 0.f)
out = VEC3(0.f, 0.f, 0.f);
else
out = eta * incidentVec - (eta * N_dot_I + sqrtf(k)) * N;
}```

These functions are part of my library: CommonLib.

# Three Ways to Calculate Mean and Variance

Mon
20
Jul 2009

There is a free e-book on DSP (Digital Signal Processing) called The Scientist and Engineer's Guide to Digital Signal Processing (by Steven W. Smith, Ph.D.). As I have a holiday now, I've started reading it and I've found some interesting information even in introductory chapters. For example, there are three algorithms to calculate mean and variance. Let's say we have some data in a vector of bytes and we want to calculate its statistics.

```std::vector<unsigned char> Bytes;
// Fill Bytes...
uint N = Bytes.size();```

```float Mean = 0.f, Variance = 0.f;
for (uint i = 0; i < N; i++)
Mean += (float)Bytes[i];
Mean /= (float)N;
for (uint i = 0; i < N; i++)
Variance += ((float)Bytes[i] - Mean) * ((float)Bytes[i] - Mean);
Variance /= (float)(N-1);```

Now the incremental one, which can return current mean and variance at any time while you can also post next piece of data. All we need to implement it is to keep track of current sample number, sum of samples and sum of squared samples. I've created template class for that, which can be parametrized by float, double, int or other numeric types.

```template <typename T>
class IncrementalMeanAndVarianceCalc
{
public:
IncrementalMeanAndVarianceCalc() : m_Sum(T()), m_SqSum(T()), m_Count(0) { }
void Clear() { m_Sum = m_SqSum = T(); m_Count = 0; }

{
m_Sum += v;
m_SqSum += v*v;
m_Count++;
}
void Add(const T *values, uint valCount)
{
for (uint i = 0; i < valCount; i++)
{
m_Sum += values[i];
m_SqSum += values[i]*values[i];
}
m_Count += valCount;
}

bool IsEmpty() { return m_Count == 0; }
uint GetCount() { return m_Count; }
T GetSum() { return m_Sum; }
T GetSqSum() { return m_SqSum; }

T GetMean()
{
assert(!IsEmpty());
return m_Sum / m_Count;
}
T GetVariance(bool varianceBiased = true)
{
if (varianceBiased)
return (m_SqSum - m_Sum*m_Sum/m_Count) / (m_Count-1);
else
return (m_SqSum - m_Sum*m_Sum/m_Count) / m_Count;
}
void GetMeanAndVariance(T &outMean, T &outVariance, bool varianceBiased = true)
{
assert(!IsEmpty());
outMean = m_Sum / m_Count;
if (varianceBiased)
outVariance = (m_SqSum - m_Sum*m_Sum/m_Count) / (m_Count - 1);
else
outVariance = (m_SqSum - m_Sum*m_Sum/m_Count) / m_Count;
}

private:
T m_Sum, m_SqSum;
uint m_Count;
};```

Finally, there is an algorithm to calculate mean and variance from histogram. It is very efficient method. If we have only 256 possible values for each sample, we don't have to do binning and calculating histogram is very simple:

```uint Histogram[256] = { 0 };
for (uint i = 0; i < N; i++)
Histogram[Bytes[i]]++;```

Now, the mean and variance can be calculated this way:

```float Mean = 0.f, Variance = 0.f;
for (uint i = 0; i < 256; i++)
Mean += (float)i * Histogram[i];
Mean /= (float)N;
for (uint i = 0; i < 256; i++)
Variance += (i - Mean)*(i - Mean) * Histogram[i];
Variance /= N - 1;```

# KD-Tree

Thu
09
Jul 2009

People at forums are usually advised to learn BSP, quadtree or octree as space partitioning data structure, but there are many more interesting structures than these three. This time my choice for home project is KD-tree. It's actually something between BSP and octree. Just like BSP it's a binary tree (each non-leaf node has two child nodes) and optimal splitting plane is estimated each time by special algorithm, but splitting planes are always aligned to one of three main axes and thus each node can be described by an AABB (axis-aligned bounding box), just like in octree.

As my tree is designed to manage objects and not geometry, nothing can be split and some of objects may intersect splitting planes. How to deal with them? I simply assign them to the parent node, so not only leaves are allowed to contain list of objects. To avoid too many small objects intersecting splitting planes to degrade performance by falling into top level nodes (it's called "sticky planes" or something like that :) I adopted "loose octree" idea to my KD-tree. It simply means I extend each node’s bounding box so that each node's children slightly overlap each other and small objects intersecting splitting plane fall into one of the children.

My KD-tree is also dynamic, which means it reorganizes itself as objects get added, removed and moved in the tree. It's actually quite simple. Each time an object is added to a node, that node can be split into two children if it's object number exceeds constant limit. Similarly node can be merged by deleting it's children each time an object is removed from one of its children, if the number of objects in that node and its children drops under constant minimum.

For additional performance, I allocate tree nodes from my own "Free List" memory pool and I keep objects connected to each node as doubly-linked list.

I also came up with an idea how to easily visualize quality of my space partitioning technique. I keep track of the number of tree nodes and objects on each depth level. This way I can tell from these several numbers whether tree is more "tall" or "wide" and whether most of objects stay in leaves instead of some top-level nodes.

Here are some screenshots and a video from my recent code:

Comments | #rendering #video #gallery #algorithms #math

# Krzywe i gówienko na kole

Thu
23
Apr 2009

Matematycy mają mądre nazwy na różnego rodzaju krzywe. Czytałem sobie o nich ostatnio na Wikipedii. Nie chodzi mi tu o dobrze znane programistom grafiki krzywe Beziera czy B-spline.

Na przykład taka krzywa Lissajous powstaje przez wykreślenie w 2D albo 3D śladu ruchu punktu wg funkcji sinus, z określoną amplitudą, okresem i fazą w każdej osi. Dawno temu napisałem wygaszacz ekranu rysujący coś takiego - ProgrameX Linear Screensaver - nie wiedząc nawet, że to się tak nazywa.

Z kolei Roulette to ogólna klasa krzywych wykreślanych przez punkt przyczepiony do jakiejś figury, która "toczy się" wzdłuż innej figury. Na przykład jeśli okrąg toczy się po prostej, to jakiś punkt skojarzony z tym okręgiem tworzy Trochoid. Jeśli ten punkt leży na tym okręgu, to jest przypadek szczególny - cykloida (Cycloid). Skojarzenie nazwy z cyklistami zapewne nie jest przypadkowe - krzywa przypomina kształt, jaki zakreśla w powietrzu gówienko przylepione do koła jadącego roweru :)

Źródło: Wikipedia

Jeśli okrąg porusza się po wewnętrznej stronie większego okręgu, powstaje hipotrochoida (Hypotrochoid). Jej szczególnym przypadkiem (punkt leży na okręgu) jest hipocykloida (Hypocycloid), a z kolei jej przykładem jest czteroramienna gwiazdka - Astroid. Podobnie, okrąg poruszający się po zewnętrznej stronie większego okręgu tworzy epitrochoidę (Epitrochoid), której szczególnym przypadkiem (punkt na okręgu) jest epicykloida (Epicycloid), której przykładem z kolei może być przypominająca kształtem serce kardioida (Cardioid).

Źródło: Wikipedia

# Programy matematyczne

Sun
01
Mar 2009

Czasem trzeba coś policzyć. Do prostych obliczeń wystarczy systemowy kalkulator. Pewne obliczenia na wektorach i kolorach daje się zrobić za pomocą mojego GameDev Calc. Czasami potrzebne są jednak bardziej zaawansowane funkcje. Jaki program matematyczny jest dobry? Niedościgniony jest podobno Matlab, ale on niestety nie należy do darmowych. Na szczęście są darmowe programy, które do wielu rzeczy z powodzeniem wystarczą.

Pierwszy z nich to Scilab. Używa składni podobnej do Matlaba i potrafi robić dużo rzeczy. Na przykład aby znormalizować wektor i pomnożyć go przez macierz:

```v=[1 2 3]
vn=v/norm(v)
M=[1 0 0; 0 0 1; 0 1 0]
v2=vn*M```

Wed
14
Jan 2009

W programowaniu bardzo często stosuje się funkcję liniową lub kwadratową. Przykładowo, jeśli mgła ma się zaczynać w głębokości Min i kończyć w głębokości Max, to jej intensywość od głębokości można wyrazić prostym wzorem:

`FogIntensity = saturate(Depth * FogScale + FogBias);`

Problem w tym, żeby znaleźć współczynniki tej funkcji. Do tego przydają się wzory, które wyliczają współczynniki dla funkcji przechodzącej przez dane punkty. Potrafi to robić mój GameDev Calc, ale żeby policzyć je w swoim programie albo na kartce, warto mieć pod ręką te wzory.

Funkcja liniowa przechodząca przez dwa punkty (x1, y1), (x2, y2) ma wzór:

Co w przełożeniu na kod daje:

```float W = p2.x - p1.x;
if (W == 0.f) Error();
float a = (p2.y - p1.y) / W;
float b = (p2.x * p1.y - p2.y * p1.x) / W;```

Z kolei funkcja kwadratowa przechodząca przez trzy punkty (x1, y1), (x2, y2), (x3, y3) ma wzór:

Co daje trochę dłuższy kod:

```float x1 = p1.x, x2 = p2.x, x3 = p3.x;
float y1 = p1.y, y2 = p2.y, y3 = p3.y;
float W =
x1 * x1 * x2 + x3 * x3 * x1 +
x2 * x2 * x3 - x1 * x1 * x3 -
x2 * x2 * x1 - x3 * x3 * x2;
if (W == 0.f) Error();
float a =
y1 * x2 + y3 * x1 + y2 * x3 -
y1 * x3 - y2 * x1 - y3 * x2;
float b =
x1 * x1 * y2 + x3 * x3 * y1 +
x2 * x2 * y3 - x1 * x1 * y3 -
x2 * x2 * y1 - x3 * x3 * y2;
float c =
x1 * x1 * x2 * y3 + x3 * x3 * x1 * y2 +
x2 * x2 * x3 * y1 - x1 * x1 * x3 * y2 -
x2 * x2 * x1 * y3 - x3 * x3 * x2 * y1;
a /= W;
b /= W;
c /= W;```