Coming out of Post 2 we have gained a decent understanding of how the shape of the PIXAR letters is constructed. But there’s still a lot of ground left to cover, and for some unfathomable reason I’m still sorta committed to keep this deep-dive going.

This is the third in a series of posts:

• The Postcard Pathtracer, Part 1: Deciphering, Continued
• The Postcard Pathtracer, Part 2: Of Letters, Coords & Squircles
• The Postcard Pathtracer, Part 3: Of Boxes, Lines & Rays (you are here)
• The Postcard Pathtracer, Part 4: Hemispheres on vector bases

— Overview —

In the previous two posts of this series I’ve been chipping away at a cleaned up version of the Vec class and the QueryDatabase function, trying to make sense of what even the smallest calculation/computation does.

— BoxTester class —

Before we wrap up what’s left of QueryDistance, let’s take a look at the BoxTest function again:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

 // Rectangle CSG equation. Returns minimum signed distance from
 // space carved by lowerLeft vertex and opposite corner vertex upperRight.
 float BoxTest(Vector3 position, Vector3 lowerLeft, Vector3 upperRight) {
     Vector3 toLowerLeft = position - lowerLeft;
     Vector3 toUpperRight = upperRight - position;

     return -std::min( // <- minimum distance out of all, i.e. to closest wall of the box
                 std::min( // <- minimum distance on x- OR y-axis
                     std::min(toLowerLeft.x_, toUpperRight.x_),  // <- minimum distance on x-axis
                     std::min(toLowerLeft.y_, toUpperRight.y_)   // <- minimum distance on y-axis
                 ),
                 std::min(toLowerLeft.z_, toUpperRight.z_)       // <- minimum distance on z-axis
             );
 }

The function takes three arguments: The position being tested and the lowerLeft & upperRight corners of the box we’re testing against. In a moment we’ll see multiple calls to BoxTest within the same line of code and it looks rather cluttered. So my idea here was to turn the BoxTest into a class that holds its defining vectors and does the same test:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21

 class BoxTester {
 public:
     Vector3 start;
     Vector3 end;

     BoxTester(Vector3 b, Vector3 e) : start(b), end(e) {}

     float Test(Vector3& position)
     {
         Vector3 toLowerLeft = position - start;
         Vector3 toUpperRight = end - position;

         return -std::min( // <- minimum distance out of all, i.e. to closest wall of the box
                     std::min( // <- minimum distance on x- OR y-axis
                         std::min(toLowerLeft.x_, toUpperRight.x_),  // <- minimum distance on x-axis
                         std::min(toLowerLeft.y_, toUpperRight.y_)   // <- minimum distance on y-axis
                     ),
                     std::min(toLowerLeft.z_, toUpperRight.z_)       // <- minimum distance on z-axis
                 );
     }
 };

In hindsight this might be one of those “You Ain’t Gonna Need It” things, since we’re currently not checking any box more than once… but it’s a fairly quick and painless change, so let’s roll with it.

Wrapping up QueryDistance

With the new class at our disposal we can return to the rest of QueryDistance. At the end of the last post we got all the distance checks against the letters covered.

What’s left after that are these checks:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21

float roomDist ;
roomDist = std::min(// min(A,B) = Union with Constructive solid geometry
                 //-min carves an empty space
                 -std::min(// Lower room
                      BoxTest(position, Vector3(-30, -.5, -30), Vector3(30, 18, 30)),
                      // Upper room
                      BoxTest(position, Vector3(-25, 17, -25), Vector3(25, 20, 25))
                 ),
                 BoxTest( // Ceiling "planks" spaced 8 units apart.
                   Vector3(fmodf(fabsf(position.x_), 8),
                       position.y_,
                       position.z_),
                   Vector3(1.5, 18.5, -25),
                   Vector3(6.5, 20, 25)
                 )
             );

if (roomDist < distance) distance = roomDist, hitType = HIT_WALL;

float sun = 19.9 - position.y_; // Everything above 19.9 is light source.
if (sun < distance)distance = sun, hitType = HIT_SUN;

The last two lines are pretty straight-forward — anything above a certain position hits the sun or light source of the scene.

But before that there are the distance checks for the entire room boiled down into essentially a single line. To break that up I created a few BoxTesters:

1
2
3
4
5
6

// Lower room (the room itself)
BoxTester box_room(Vector3(-30, -.5, -30), Vector3(30, 18, 30));
// Upper room (The space where the ceiling planks are)
BoxTester box_roof(Vector3(-25, 17, -25), Vector3(25, 20, 25));
// Ceiling "planks", or one of them. Actual tested position will be modulo'd along an axis.
BoxTester box_planks(Vector3(1.5, 18.5, -25), Vector3(6.5, 20, 25));

The main room box_room has a footprint of 60 by 60 units centered around 0,0 and is roughly 18 units tall. Inset into the ceiling is another space box_roof that spans from height 17 to 20. As seen in the snippet before anything above y-position 19.9 is our light source already though, allowing light in through this space in the ceiling.

And lastly our third box box_planks sits at height 18.5 to 20 in the ceiling space. It spans the entire space on the z-Axis from -25 to 25, and sits at 1.5 to 6.5 on the x-Axis, near the center. The fact that it is on the positive side of the x-Axis is an intentional choice too, as we’re about to discover.

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11

Vector3 plank_test_pos(position);
plank_test_pos.x_ = fmodf(fabsf(plank_test_pos.x_), 8.0f);

float roomDist = std::min(// min(A,B) = Union with Constructive solid geometry
                     //-min carves an empty space
                     -std::min(
                          box_room.Test(position), // Lower room
                          box_roof.Test(position) // Upper room
                     ),
                     box_planks.Test(plank_test_pos) // Ceiling "planks" spaced 8 units apart.
                 );

Now, with all the Vector objects tucked away in the BoxTester instance, this all looks a lot cleaner. What’s going on here exactly?

Right in the middle there’s the Test() calls against the main two boxes, which as a reminder return a positive distance for anything outside of the box and a negative one for anything inside. As with most of our distance checks we’re only interested in finding the closest of these geometry objects, so we take the smaller out of the two. Then we negate the result, which reverses what’s considered solid/empty space: The space inside the boxes is empty now, returning positive distances towards the nearest wall. And lastly we take another minimum in comparison to the box_planks check, this one being a regular non-inverted box.

But hold on, we are using the Vector3 plank_test_pos on the last check! This is the same as our position-to-query that we’ve used throughout this function, but with a twist or two:

• On the x-Axis component we take the absolute value (fabsf), meaning anything between -6.5 < x < -1.5 within the right location will also be considered to be “inside” the plank!
• And secondly, we’re taking the modulo (fmodf) of the position with 8.0f — so after every “segment” of 8 units along the x-Axis we loop back to checking at 0.0f: Just as there’s a box at 1.5 - 6.5 there’ll be a box at 9.5 - 14.5 and 17.5 - 22.5.¹

So that use of fmodf and fabsf is a clever way of constructing some repeating geometry to make the scene more interesting.

All changes mentioned so far, which is mostly the introduction of the BoxTester class, are in this commit.

— New Class: LetterLine —

Next I decided that I still found the vector math we looked at in the last post to be rather unwieldy. Especially things like line_vector.Normalized() and having separate helper variables for line_begin, line_end & line_vector irked me.

So here I’m introducing a new class:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

class LetterLine {
public:
    Vector2 start;
    Vector2 end;
    Vector2 direction;
    Vector2 normal;
    float length;

    LetterLine(Vector2& b, Vector2& e) : start(b), end(e) {
        direction = e - b;
        normal = direction.Normalized();
        length = direction.Length();
    }
};

Looking back at it now, I realize that this isn’t exactly ideal. We have all those vectors in our line, but they’re public since I wanted to keep it simple… but they aren’t meant to be altered after creating an instance of this class, since that won’t update the direction, normal and length properties.

Potentially this would’ve been better:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14

class LetterLine {
public:
    const Vector2 start;
    const Vector2 end;
    const Vector2 direction;
    const Vector2 normal;
    const float length;

    LetterLine(Vector2& b, Vector2& e) : start(b), end(e),
                         direction(e - b),
                         normal(direction.Normalized()),
                         length(direction.Length())
    {}
};

I really just wanted objects with no logic to them, to just use them as containers for some useful properties of our lines. I also hoped calculating the extra two vectors and the length up-front would speed up the overall pathtracing, but I didn’t actually check if it makes any noticeable difference. /shrug

Here’s the new function to create LetterLine instances from our pixar_points array:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

std::vector<LetterLine> letterlines = {};

void SetupPathtracer() {

    for (int i = 0; i+3 < pixar_points.size(); i += 4) {
        Vector2 b = Vector2(pixar_points[i  ], pixar_points[i+1]);
        Vector2 e = Vector2(pixar_points[i+2], pixar_points[i+3]);

        LetterLine* l = new LetterLine(b, e);
        letterlines.push_back(*l);
    }
}

Not much to see here.

Now the fun part comes when we actually make use of these new “container objects”. For comparison, here’s the previous version of our code that we’ve been working with:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23

float QueryDatabase(Vector3 position, int &hitType) {
    //[...]

    for (int i = 0; i+3 < pixar_points.size(); i += 4) {
        Vector2 line_begin = Vector2(pixar_points[i], pixar_points[i + 1]);
        Vector2 line_end = Vector2(pixar_points[i + 2], pixar_points[i + 3]);
        Vector2 line_vector = line_end - line_begin;

        // factor to scale line_vector by to project position onto line
        // i.e. with this factor alone (line_begin + line_vector.Normalized() * scale_factor)
        // gives us the closest point on an infinite line
        float scale_factor = Vector2(pos_in_z0 - line_begin).DotProduct(line_vector.Normalized());

        // But now we clamp the scaling within 0.0f <= scale_factor <= line length
        scale_factor = std::min(std::max(scale_factor,  0.0f), line_vector.Length());

        Vector2 letter_distance = pos_in_z0 - (line_begin + line_vector.Normalized() * scale_factor);

        distance = std::min(distance, letter_distance.Length()); // compare distance.
    }

    //[...]
}

It’s already pretty clean, but we can go cleaner:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21

// Sample the world using Signed Distance Fields.
float QueryDatabase(Vector3 position, int &hitType) {
    //[...]

    for (LetterLine line : letterlines) {
        // factor to scale line.normal by to project position onto line
        // i.e. with this factor alone (line.start + line.normal * scale_factor)
        // gives us the closest point on an infinite line
        float scale_factor = Vector2(pos_in_z0 - line.start).DotProduct(line.normal);

        // But now we clamp the scaling within 0.0f <= scale_factor <= line length
        scale_factor = std::min(std::max(scale_factor,  0.0f), line.length);

        Vector2 letter_distance = pos_in_z0 - (line.start + line.normal * scale_factor);

        distance = std::min(distance, letter_distance.Length()); // compare distance.
    }

   //[...]
}

Just four lines of code within the loop, using two or three mathematical operators and variables each! Can’t get much simpler than this, if you ask me.

The corresponding commit is here, which includes the next section’s tiny set of changes.

— misc. changes —

std::vector of curves

Not far from our last change there’s also these:

1
2
3
4
5
6
7
8

         // Two curves (for P and R in PixaR) with hard-coded locations.
         Vector2 curves[] = {Vector2(-11, 6), Vector2(11, 6)};

         for (int i = 2; i--;) {
             Vector2 center_to_pos = pos_in_z0 - curves[i];

             //[...]
         }

While I was at it I decided to pull the curve vectors out of the function to make them easier to access/edit near the top. They got a new shiny std::vector as their home as well, which allows us to switch to a range-based for loop:

1
2
3
4
5
6

// Two curves (for P and R in PixaR) with hard-coded locations.
// curve centers are between right ends of the bars.
std::vector<Vector2> curve_centers = {
  Vector2(-11, 6),
  Vector2( 11, 6)
};

1
2
3
4
5

         for (Vector2 c : curve_centers) {
             Vector2 center_to_pos = pos_in_z0 - c;

             //[...]
         }

centralize hard-coded vectors

And three more vectors to add to the list up top:

1
2
3

Vector3 camera_position = Vector3(-22, 5, 25);
Vector3 lookat_position = Vector3(-3, 4, 0);
Vector3 light_direction = Vector3(.6, .6, 1);

…

1
2

// OLD
Vector3 lightDirection(.6, .6, 1); // Directional light

1
2

// NEW
Vector3 lightDirection(light_direction); // Directional light

…

1
2
3

// OLD
Vector3 position(-22, 5, 25);
Vector3 goal = Vector3(-3, 4, 0) - position;

1
2
3

// NEW
Vector3 position(camera_position);
Vector3 goal = lookat_position - position;

Not much to add here, other than to wonder why I decided to use copy-constructors to make the local variables still. None of the three get modified anywhere other than the light direction getting normalized just to be safe.

Find the commit here.

— Function: RayMarching —

Finally it’s time to tackle the next function: RayMarching. Just to recall, this function performs the signed sphere marching or distance-aided ray marching. Again, other posts like the original one by Fabien explain the concept very well. In short: We’re casting a ray, and want to see what it hits and where. This is done by starting at the origin of the ray, following these steps:

1. Get distance to closest object via QueryDatabase (also updates hitType since it’s a reference parameter!)
2. Move forward by this distance, since the sphere of this radius is guaranteed to be empty
3. Repeat for a finite amount of steps and distance

There’s not much more to it than that and the implementation we’re working with is pretty compact:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17

// Perform signed sphere marching
// Returns hitType 0, 1, 2, or 3 and update hit position/normal
int RayMarching(Vector3 origin, Vector3 direction, Vector3 &hitPos, Vector3 &hitNorm) {
    int hitType = HIT_NONE;
    int noHitCount = 0;
    float d; // distance from closest object in world.

    for (float total_d=0; total_d < 100; total_d += d)
        if ((d = QueryDatabase(hitPos = origin + direction * total_d, hitType)) < .01
            || ++noHitCount > 99)
            return hitNorm =
                 Vector3(QueryDatabase(hitPos + Vector3(.01, 0, 0), noHitCount) - d,
                      QueryDatabase(hitPos + Vector3(0, .01, 0), noHitCount) - d,
                      QueryDatabase(hitPos + Vector3(0, 0, .01), noHitCount) - d).Normalized()
             , hitType; // Weird return statement where a variable is also updated.
    return 0;
}

Unfortunately (for readability at least!) there’s a lot going on with variables being modified in the if-condition and that big assignment to hitNorm within the return statement with the use of the comma operator (which we’ve come across in Part 2’s chapter on the curves)

These changes are not all that extensive however, so let’s look at the rewritten function:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25

// Perform signed sphere marching
// Returns hitType 0, 1, 2, or 3 and update hit position/normal
int RayMarching(Vector3 origin, Vector3 direction, Vector3 &hitPos, Vector3 &hitNorm) {
    int hitType = HIT_NONE;
    int noHitCount = 0;
    float d; // distance from closest object in world.

    // Signed distance marching
    for (float total_d=0; total_d < 100; total_d += d) {
        hitPos = origin + direction * total_d;
        d = QueryDatabase(hitPos, hitType);

        // negative distance will mean we're inside the object we hit, I guess
        // and we'll return right away so no worries there
        if (d < .01 || noHitCount > 99) {
            int unusedReturnRef = 0;
            hitNorm = Vector3(QueryDatabase(hitPos + Vector3(.01, 0, 0), unusedReturnRef) - d,
                              QueryDatabase(hitPos + Vector3(0, .01, 0), unusedReturnRef) - d,
                              QueryDatabase(hitPos + Vector3(0, 0, .01), unusedReturnRef) - d).Normalized();
            return hitType;
        }
        ++noHitCount;
    }
    return 0;
}

Assigning to d now happens in line 11 above, before the if-statement. I added in a comment to remind myself why we’re only seemingly checking against the small positive distance of 0.1, but if d is a negative number that means we’re inside an object and should definitely return from the function as well.

The int noHitCount is used to count how many iterations we’ve gone through, and I moved the increment operator used on that back towards the end of the loop, instead of doing it right within the comparison. In principle we could also use noHitCount as our for-loop iteration variable, and then the total_d += d would be inside the loop and total_d >= 100 would be within the if statement.

Technically we do something slightly different depending on whether we step into the if statement and return, or finish the for loop and return.

• within the if statement we calculate hitNorm (more on that coming up) and return hitType
• after the for loop we simply return which happens to be HIT_NONE. The reference hitNorm is still unchanged from what was passed in Does this difference matter? …It might. But these are edge cases: We hit total_d >= 100 when our ray has traveled a long way without hitting anything — and this is a ray that is not being reflected anywhere yet, its direction never changes. Which leads me to believe this might be impossible within our current scene. If we did have a ray this long, it’d return HIT_NONE, which as we’ll see ends the traced path in Trace().

In the other case of noHitCount > 99 we actually return as if we’d hit something, and hitType will invariably have been set by QueryDistance — looking back through that it cannot possibly be HIT_NONE anymore! We will be returning the normal vector towards, and type of, the closest thing in our scene. Which will continue the traced ray. This case of noHitCount reaching 100 is actually fairly common, and is an area where a lot of computation is spent: Imagine our ray is passing by a surface really closely but still outside of the d < .01. Since we’re already close to the surface, our sphere is as small as d is and we’ll be inching closer and closer to the surface at an ever decreasing step size.

I’ve tried visualizing this at some point, although I have to admit I don’t fully remember if this is the original code or if I already attempted to optimize the edge cases in any way.

If I recall correctly the middle image is a mapping of the biggest number of iterations (=noHitCount) that occured on that pixel, even across reflections, while the bottom image is a map of the max num of iterations on the first ray sent into the scene?…

But in any case it makes sense to return something useable if we exceed our maximum iterations before exceeding maximum distance, since this most likely means we’ve been approaching a surface. If we’ve traveled a long total distance and “not much happened” we can’t really be sure there’s anything around anywhere.

A cursory check with a log statement before the return 0 just now contributed no output, so that last return statement really doesn’t get hit in our current setup. If I reduce the maximum ray distance from 100 down to 40 I get something like this:

…and a 136 MB log file because of course I forgot to take out the debug logging statement… Whoops!

I’ve gone off on a bit of a tangent there, but there’s one thing I still wanted to highlight within the RayMarching function: Calculating the normal vector to the hit.

1
2
3
4

int unusedReturnRef = 0;
hitNorm = Vector3(QueryDatabase(hitPos + Vector3(.01, 0, 0), unusedReturnRef) - d,
                  QueryDatabase(hitPos + Vector3(0, .01, 0), unusedReturnRef) - d,
                  QueryDatabase(hitPos + Vector3(0, 0, .01), unusedReturnRef) - d).Normalized();

First off, QueryDatabase still takes that second parameter to hand back a hitType — which is why I added the variable unusedReturnRef just to signal that we don’t care about that value here at all. In the original code noHitCount was passed in by reference here, simply because that variable was available and not needed anymore either.

Second of all, QueryDatabase only returns a distance value for a coordinate, not a direction. That’s why we need to query three times with an offset per axis to calculate a normal vector. The reason why this works has kinda thrown my brain for a loop for a while, but the easiest analogy I can think of is a finite distance approximation of the slope of a 2d graph/function. And what we have in our case is a multi-variable function. We could also think of it as getting the gradient of a vector field, or the gradient of an isosurface within that… but right now I’m not even sure I’m using these terms correctly anymore.

Honestly I was digging through some interesting math-related Wikipedia entries because I also wondered how inaccurate this method is, how taking more samples or sampling at different distances would affect it, and so on. But that could use another post to dig into at some later date, since I got these current ones to finish.

All that aside, with the final commit for today we have today’s final code base ironed out.

— Interim result —

And once again for completeness' sake, here’s the entire body of code we got so far, minus the Urho3D scaffolding:

 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241

float randomVal() { return (float) rand() / RAND_MAX; }

// Box CSG equation. Returns minimum signed distance from space
// carved by start vector and opposite corner vector end.
class BoxTester {
public:
    Vector3 start;
    Vector3 end;

    BoxTester(Vector3 b, Vector3 e) : start(b), end(e) {}

    float Test(Vector3& position)
    {
        Vector3 toLowerLeft = position - start;
        Vector3 toUpperRight = end - position;

        return -std::min( // <- minimum distance out of all, i.e. to closest wall of the box
                    std::min( // <- minimum distance on x- OR y-axis
                        std::min(toLowerLeft.x_, toUpperRight.x_),  // <- minimum distance on x-axis
                        std::min(toLowerLeft.y_, toUpperRight.y_)   // <- minimum distance on y-axis
                    ),
                    std::min(toLowerLeft.z_, toUpperRight.z_)       // <- minimum distance on z-axis
                );
    }
};

// Lower room (the room itself)
BoxTester box_room(Vector3(-30, -.5, -30), Vector3(30, 18, 30));
// Upper room (The space where the ceiling planks are)
BoxTester box_roof(Vector3(-25, 17, -25), Vector3(25, 20, 25));
// Ceiling "planks", or one of them. Actual tested position will be modulo'd along an axis.
BoxTester box_planks(Vector3(1.5, 18.5, -25), Vector3(6.5, 20, 25));

#define HIT_NONE 0
#define HIT_LETTER 1
#define HIT_WALL 2
#define HIT_SUN 3

class LetterLine {
public:
    Vector2 start;
    Vector2 end;
    Vector2 direction;
    Vector2 normal;
    float length;

    LetterLine(Vector2& b, Vector2& e) : start(b), end(e) {
        direction = e - b;
        normal = direction.Normalized();
        length = direction.Length();
    }
};


std::vector<int> pixar_points = {
         -13, 0,  -13, 8, // P stem
         -13, 4,  -11, 4, // P top bar
         -13, 8,  -11, 8, // P mid bar

          -7, 0,   -5, 0, // I bottom bar
          -6, 0,   -6, 8, // I stem
          -7, 8,   -5, 8, // I top bar

          -3, 0,    1, 8, // X slash
          -3, 8,    1, 0, // X backslash

           3, 0,    5, 8, // A slash
           5, 8,    7, 0, // A backslash
           4, 4,    6, 4, // A bar

           9, 0,    9, 8, // R stem
           9, 4,   11, 4, // R mid bar
           9, 8,   11, 8, // R top bar
          10, 4,   13, 0  // R leg
        };

std::vector<LetterLine> letterlines = {};

// Two curves (for P and R in PixaR) with hard-coded locations.
// curve centers are between right ends of the bars.
std::vector<Vector2> curve_centers = {
  Vector2(-11, 6),
  Vector2( 11, 6)
};

Vector3 camera_position = Vector3(-22, 5, 25);
Vector3 lookat_position = Vector3(-3, 4, 0);
Vector3 light_direction = Vector3(.6, .6, 1);

void SetupPathtracer() {

    for (int i = 0; i+3 < pixar_points.size(); i += 4) { // letters.Length()
        //Vector3 begin = Vector3(letters.At(i) - 79, letters.At(i+1) - 79);
        //Vector3 end = Vector3(letters.At(i+2) - 79, letters.At(i+3) - 79);
        Vector2 b = Vector2(pixar_points[i  ], pixar_points[i+1]);
        Vector2 e = Vector2(pixar_points[i+2], pixar_points[i+3]);

        LetterLine* l = new LetterLine(b, e);
        letterlines.push_back(*l);
    }
}


// Sample the world using Signed Distance Fields.
float QueryDatabase(Vector3 position, int &hitType) {
    float distance = 1e9;
    Vector2 pos_in_z0(position.x_, position.y_); // Flattened position (z=0)

    for (LetterLine line : letterlines) {
        // factor to scale line.normal by to project position onto line
        // i.e. with this factor alone (line.start + line.normal * scale_factor)
        // gives us the closest point on an infinite line
        float scale_factor = Vector2(pos_in_z0 - line.start).DotProduct(line.normal);

        // But now we clamp the scaling within 0.0f <= scale_factor <= line length
        scale_factor = std::min(std::max(scale_factor,  0.0f), line.length);

        Vector2 letter_distance = pos_in_z0 - (line.start + line.normal * scale_factor);

        distance = std::min(distance, letter_distance.Length()); // compare distance.
    }

    for (Vector2 c : curve_centers) {
        Vector2 center_to_pos = pos_in_z0 - c;

        if(center_to_pos.x_ > 0) {
            float comparison_dist = std::abs(center_to_pos.Length() - 2.0f);
            distance = std::min(distance, comparison_dist);
        }
    }

    distance = powf(powf(distance, 8) + powf(position.z_, 8), .125) - .5;
    hitType = HIT_LETTER;

    Vector3 plank_test_pos(position);
    plank_test_pos.x_ = fmodf(fabsf(plank_test_pos.x_), 8.0f);

    float roomDist = std::min(// min(A,B) = Union with Constructive solid geometry
                        //-min carves an empty space
                        -std::min(
                             box_room.Test(position), // Lower room
                             box_roof.Test(position) // Upper room
                        ),
                        box_planks.Test(plank_test_pos) // Ceiling "planks" spaced 8 units apart.
                    );
    if (roomDist < distance) distance = roomDist, hitType = HIT_WALL;

    float sun = 19.9 - position.y_ ; // Everything above 19.9 is light source.
    if (sun < distance)distance = sun, hitType = HIT_SUN;

    return distance;
}

// Perform signed sphere marching
// Returns hitType 0, 1, 2, or 3 and update hit position/normal
int RayMarching(Vector3 origin, Vector3 direction, Vector3 &hitPos, Vector3 &hitNorm) {
    int hitType = HIT_NONE;
    int noHitCount = 0;
    float d; // distance from closest object in world.

    // Signed distance marching
    for (float total_d=0; total_d < 100; total_d += d) {
        hitPos = origin + direction * total_d;
        d = QueryDatabase(hitPos, hitType);

        // negative distance will mean we're inside the object we hit, I guess
        // and we'll return right away so no worries there
        if (d < .01 || noHitCount > 99) {
            int unusedReturnRef = 0;
            hitNorm = Vector3(QueryDatabase(hitPos + Vector3(.01, 0, 0), unusedReturnRef) - d,
                              QueryDatabase(hitPos + Vector3(0, .01, 0), unusedReturnRef) - d,
                              QueryDatabase(hitPos + Vector3(0, 0, .01), unusedReturnRef) - d).Normalized();
            return hitType;
        }
        ++noHitCount;
    }
    return 0;
}

Vector3 Trace(Vector3 origin, Vector3 direction) {
    Vector3 sampledPosition(1, 1, 1);
    Vector3 normal(1, 1, 1);
    Vector3 color(0, 0, 0);
    Vector3 attenuation(1, 1, 1);
    Vector3 lightDirection(light_direction); // Directional light
    lightDirection.Normalize();

    for (int bounceCount = 3; bounceCount--;) {
        int hitType = RayMarching(origin, direction, sampledPosition, normal);
        if (hitType == HIT_NONE) break; // No hit. This is over, return color.
        if (hitType == HIT_LETTER) { // Specular bounce on a letter. No color acc.
            direction = direction + normal * ( normal.DotProduct(direction) * -2);
            origin = sampledPosition + direction * 0.1;
            attenuation = attenuation * 0.2; // Attenuation via distance traveled.
        }
        if (hitType == HIT_WALL) { // Wall hit uses color yellow?
            float incidence = normal.DotProduct(lightDirection);
            float p = 6.283185 * randomVal();
            float c = randomVal();
            float s = sqrtf(1 - c);
            float g = normal.z_ < 0 ? -1 : 1;
            float u = -1 / (g + normal.z_);
            float v = normal.x_ * normal.y_ * u;
            direction = Vector3(v,
                          g + normal.y_ * normal.y_ * u,
                          -normal.y_) * (cosf(p) * s)
                      +
                      Vector3(1 + g * normal.x_ * normal.x_ * u,
                          g * v,
                          -g * normal.x_) * (sinf(p) * s) + normal * sqrtf(c);
            origin = sampledPosition + direction * .1;
            attenuation = attenuation * 0.2;
            if (incidence > 0 &&
                 RayMarching(sampledPosition + normal * .1,
                             lightDirection,
                             sampledPosition,
                             normal) == HIT_SUN)
                color = color + attenuation * Vector3(500, 400, 100) * incidence;
        }
        if (hitType == HIT_SUN) { //
          color = color + attenuation * Vector3(50, 80, 100); break; // Sun Color
        }
    }
    return color;
}

477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512

   //  int w = 960, h = 540, samplesCount = 8;
   Vector3 position(camera_position);
   Vector3 goal = lookat_position - position;
   goal.Normalize();

   // this vector is actually to the right of goal... flipped signs of x and z for left.
   // Up is down. Hence the subtractions instead of additions for trace directions
   // Note to self: Yes, it was rendering upside-down, but only because originally it
   // wrote values for y = h .. 0 and x = w .. 0, so from bottom right to top left, sequentially
   Vector3 left = Vector3(-goal.z_, 0, goal.x_).Normalized() / w_;

   // Cross-product to get the up vector
   Vector3 up = goal.CrossProduct(left);

   // this used to be the argument to Trace
   Vector3 trace_dir{};

   for (int x = 0; x < w_; x++) {
       Vector3 color;
       for (int p = samplesCount_; p--;) {
           // this was the second Trace argument as both lines in one:
           // !(goal + left * (x - w_ / 2 + randomVal()) + up * (y_ - h_ / 2 + randomVal()))
           trace_dir = goal + left * (x - w_ / 2 + randomVal()) + up * (y_ - h_ / 2 + randomVal());
           trace_dir.Normalize();

           color = color + Trace(position, trace_dir);
       }

       // Reinhard tone mapping
       color = color / samplesCount_ + Vector3(1.0f, 1.0f, 1.0f) * 14. / 241.0;
       Vector3 o = color + Vector3(1.0, 1.0, 1.0);
       color = color / o;

       //[...]

   }

(full main.cpp here)

As always, thank you for reading along so far. I started digging into this at the start of 2020, and now the year is almost over. I’ve had episodes where I’d work on these posts, but a lot of the time I’ve been putting it of… even though my lowest bar was to at least get all four planned posts out before the year ends.

In the meantime I also got sidetracked by implementing the pathtracer in CUDA and building supporting features around that to tinker with it further. That’s been great fun to see the whole thing running in real time and with camera controls and more, although the code is far messier than what I went for in this series!

Here’s a quick attempt at embedding a webm-video, hope that works out. It’s not very smooth especially while capturing with Peek. The GUI that didn’t want to cooperate with the low refresh rate is Nuklear which is a single header immediate-mode GUI that is otherwise pleasantly easy to deal with.

Maybe at some point I’ll figure out a more sustainable way of blogging about these things, where it comes more naturally alongside the tinkering. I might make that my goal for 2021.

If you have any corrections, constructive criticism or comments, please post them in the comments section below, or write to me via e-mail or Twitter! Feedback highly welcome. Stay safe out there.

~rf