Physically based shading - ambient/indirect lighting

Real-time graphics deploys a variety of approximations to deal with the computational expense of simulating indirect lighting, trading off between runtime performance and lighting fidelity. This is an area of active research, with new techniques appearing every year.

Ambient lighting

At the very simplest end of the range, you can use ambient lighting: a global, omnidirectional light source that applies to every object in the scene, without regard to actual light sources or local visibility. This is not at all accurate, but is extremely cheap, easy for an artist to tweak, and can look okay depending on the scene and the desired visual style.

Common extensions to basic ambient lighting include:

Make the ambient color vary directionally, e.g. using spherical harmonics (SH) or a small cubemap, and looking up the color in a shader based on each vertex's or pixel's normal vector. This allows some visual differentiation between surfaces of different orientations, even where no direct light reaches them.
Apply ambient occlusion (AO) techniques including pre-computed vertex AO, AO texture maps, AO fields, and screen-space AO (SSAO). These all work by attempting to detect areas such as holes and crevices where indirect light is less likely to bounce into, and darkening the ambient light there.
Add an environment cubemap to provide ambient specular reflection. A cubemap with a decent resolution (128² or 256² per face) can be quite convincing for specular on curved, shiny surfaces.

Baked indirect lighting

The next "level", so to speak, of techniques involve baking (pre-computing offline) some representation of the indirect lighting in a scene. The advantage of baking is you can get pretty high-quality results for little real-time computational expense, since all the hard parts are done in the bake. The trade-offs are that the time needed for the bake process harms level designers' iteration rate; more memory and disk space are required to store the precomputed data; the ability to change the lighting in real-time is very limited; and the bake process can only use information from static level geometry, so indirect lighting effects from dynamic objects such as characters will be missed. Still, baked lighting is very widely used in AAA games today.

The bake step can use any desired rendering algorithm including path tracing, radiosity, or using the game engine itself to render out cubemaps (or hemicubes).

The results can be stored in textures (lightmaps) applied to static geometry in the level, and/or they can also be converted to SH and stored in volumetric data structures, such as irradiance volumes (volume textures where each texel stores an SH probe) or tetrahedral meshes. You can then use shaders to look up and interpolate colors from that data structure and apply them to your rendered geometry. The volumetric approach allows baked lighting to be applied to dynamic objects as well as static geometry.

The spatial resolution of the lightmaps etc. will be limited by memory and other practical constraints, so you might supplement the baked lighting with some AO techniques to add high-frequency detail that the baked lighting can't provide, and to respond to dynamic objects (such as darkening the indirect light under a moving character or vehicle).

There's also a technique called precomputed radiance transfer (PRT), which extends baking to handle more dynamic lighting conditions. In PRT, instead of baking the indirect lighting itself, you bake the transfer function from some source of light—usually the sky—to the resultant indirect lighting in the scene. The transfer function is represented as a matrix that transforms from source to destination SH coefficients at each bake sample point. This allows the lighting environment to be changed, and the indirect lighting in the scene will respond plausibly. Far Cry 3 and 4 used this technique to allow a continuous day-night cycle, with indirect lighting varying based on the sky colors at each time of day.

One other point about baking: it may be useful to have separate baked data for diffuse and specular indirect lighting. Cubemaps work much better than SH for specular (since cubemaps can have a lot more angular detail), but they also take up a lot more memory, so you can't afford to place them as densely as SH samples. Parallax correction can be used to somewhat make up for that, by heuristically warping the cubemap to make its reflections feel more grounded to the geometry around it.

Fully real-time techniques

Finally, it's possible to compute fully dynamic indirect lighting on the GPU. It can respond in real-time to arbitrary changes of lighting or geometry. However, again there is a tradeoff between runtime performance, lighting fidelity, and scene size. Some of these techniques need a beefy GPU to work at all, and may only be feasible for limited scene sizes. They also typically support only a single bounce of indirect light.

A dynamic environment cubemap, where the faces of the cubemap are re-rendered each frame using six cameras clustered around a chosen point, can provide decently good ambient reflections for a single object. This is often used for the player car in racing games, for instance.
Screen-space global illumination, an extension of SSAO that gathers bounce lighting from nearby pixels on the screen in a post-processing pass.
Screen-space raytraced reflection works by ray-marching through the depth buffer in a post-pass. It can provide quite high-quality reflections as long as the reflected objects are on-screen.
Instant radiosity works by tracing rays into the scene using the CPU, and placing a point light at each ray hit point, which approximately represents the outgoing reflected light in all directions from that ray. These many lights, known as virtual point lights (VPLs), are then rendered by the GPU in the usual way.
Reflective shadow maps (RSMs) are similar to instant radiosity, but the VPLs are generated by rendering the scene from the light's point of view (like a shadow map) and placing a VPL at each pixel of this map.
Light propagation volumes consist of 3D grids of SH probes placed throughout the scene. RSMs are rendered and used to "inject" bounce light into the SH probes nearest the reflecting surfaces. Then a flood-fill-like process propagates light from each SH probe to surrounding points in the grid, and the result of this is used to apply lighting to the scene. This technique has been extended to volumetric light scattering as well.
Voxel cone tracing works by voxelizing the scene geometry (likely using varying voxel resolutions, finer near the camera and coarser far away), then injecting light from RSMs into the voxel grid. When rendering the main scene, the pixel shader performs a "cone trace"—a ray-march with gradually increasing radius—through the voxel grid to gather incoming light for either diffuse or specular shading.

Most of these techniques are not widely used in games today due to problems scaling up to realistic scene sizes, or other limitations. The exception is screen-space reflection, which is very popular (though it's usually used with cubemaps as a fallback, for regions where the screen-space part fails).

As you can see, real-time indirect lighting is a huge topic and even this (rather long!) answer can only provide a 10,000-foot overview and context for further reading. Which approach is best for you will depend greatly on the details of your particular application, what constraints you're willing to accept, and how much time you have to put into it.


Test if polygon edge is facing another polygon

You can take the normal of the line segment you're testing and find the dot product of that normal and the normal of each segment of the polygon you're testing against. If it is positive, then they are facing the same way. If the dot product is negative, then the line segments are facing in opposite directions.

If they are facing in opposite directions, then they're either pointing at each other, or away from each other. You can figure out which by calculating a vector from a point on the first segment to a point on the other segment and seeing if it's facing the same direction as its normal (by again using the dot product). If so, they are facing each other, otherwise they're facing away.

You can skip all the tests above by first testing if the polygons are the same color and only doing the tests if they're not.


Depth of Field in Path Tracing: What do I do with the secondary ray?

The missing step

If you already understand how to generate a secondary ray, then you have already grasped the difficult part. All you need to do now is find the colour that this secondary ray results in. This is exactly the same process as using the primary ray to find a colour, in basic ray tracing.

After repeating this for a large number of secondary rays, the average of all those colours is the colour you assign to the pixel in question (you don't use the primary ray - its only purpose is generating the secondary rays). This is the only difference from basic ray tracing.

Why it works

My personal difficulty when trying to understand this was that I was trying to understand how the light gets from the different points on the lens back to the eye/camera. The thing to understand is that this is only an approximation of a lens, and the rays don't need to consider the eye/camera at all.

The lens is being approximated by a disc.

All these different secondary rays show you what you would see from different points on that disc, if you looked directly towards a certain point (the focal point).

For an object at the focal point, all these different viewpoints will see the same point on the object, and the object will look sharply defined.
For an object behind the focal point, all these different viewpoints will see slightly different positions on the object, and some of them will miss the object and see what is to the side of it. The edges of the object will be blended with the colour of whatever is to the side of it.
For an object in front of the focal point, again the different viewpoints will see slightly different positions on the object, giving blurred edges.

If a primary ray hits an object at the focal distance, but there is a closer object which is not quite in the way, but almost, then some of the secondary rays will hit the closer object before they reach the focal point. This will give a blurred edge to the closer object, partially obscuring the object at the focal point, even though that object itself is in sharp focus.

All of this is as desired. To this level of explanation, this mimics the effect of a real lens, without having to consider refraction at all.


Ray intersection with Sun positions throughout the year

The function that calculates sun position in the sky takes as input
local latitude, longitude and local date time. The output is a vector
which is sun position:

Inversing the above means using the ray direction(sun position in the sky) and retrieving back the local date time. If the date time is valid then the ray must hit the sun path throughout the year.

Here is the inverse function of the above:

public DateTime GetLocalDatetime(double haSolar, double vaSolar)
{
    double latitude = Helper.ToRadian(Latitude);
    double longitude = Helper.ToRadian(Longitude);

    double zenithTheta = Helper.PiOver2 - vaSolar;

    double sinSolarDeclination = (Math.Cos(Math.PI - haSolar) * Math.Cos(latitude) * Math.Sin(zenithTheta)) + Math.Sin(latitude) * Math.Cos(zenithTheta);
    double solarDeclination = Math.Asin(sinSolarDeclination);

    double cosH = (Math.Sin(vaSolar) - Math.Sin(latitude) * sinSolarDeclination) / (Math.Cos(latitude) * Math.Cos(solarDeclination));
    if (Math.Abs(cosH) > 1) return DateTime.MinValue;

    double H = Math.Acos(cosH);
    if (haSolar < 0) H *= -1;

    double declratio = solarDeclination / Helper.ToRadian(23.45);
    if (Math.Abs(declratio) > 1) return DateTime.MinValue;

    double beta = Math.Asin(declratio);
    int n = (int)Math.Round(beta / Helper.ToRadian(360d / 365d) + 81d, 1);

    DateTime localtime = new DateTime(DateTime.Now.Year, 1, 1).AddDays(n - 1);
    DateTime solartime = new DateTime(localtime.Year, localtime.Month, localtime.Day, 12, 0, 0);
    solartime = solartime.AddMinutes(4 * Helper.ToDegree(H));

    double eot = 9.87 * Math.Sin(2d * beta) - 7.53 * Math.Cos(beta) - 1.5 * Math.Sin(beta);
    double deltatimezone = 4d * (longitude - (Timezone * 15d));
    return solartime.AddMinutes(-(eot + deltatimezone));
}

And the image where each band represent each month in the year:


Alpha component in color buffer affects RGB components?

It looks like you're drawing with pre-multiplied alpha. This means that the red, green, and blue channels have been multiplied by the alpha channel. So your shader should look more like this:

precision mediump float;
uniform sampler2D uSampler;
varying vec2 vTextureCoord;
void main() {
    vec4 originalColor = texture2D(uSampler, vTextureCoord);
    vec4 tempColor = originalColor + vec4(0, 0, 0, -0.5);
    tempColor.rgb *=  tempColor.a;
    gl_FragColor = tempColor;
}

I believe you also need to set glBlendFunc(GL_ONE, GL_ONE_MINUS_SRC_ALPHA) to get the correct blending (and of course, glEnable(GL_BLEND)).


2D array texture with different filtering setting per texture?

Array textures do not contain textures; they are a single texture that contains a number of images, with each image being of the same size, format, mipmap depth, etc.

As such, the sampling parameters that are part of the texture apply to sample fetches for all images within the texture.

That having been said:

  In the shader I perform a manual GL_NEAREST lookup via textureFetch().

This means your FS knows when to look the texture up with linear filtering and when not to. As such, you could have two different array texture samplers in the shader, one of which you expect to use linear sampling and one which uses nearest.

Both of these samplers can be using the same texture object, just paired with a different sampler object (which contains all of the filtering parameters for the texture). One sampler would use nearest filtering, the other would use linear.


It's likely due to shadow rays. When the light source is exactly on the plane of the wall, a shadow ray traced to it may or may not be found to intersect the wall before it reaches the light. The result will depend on arithmetic roundoff error, and may appear as random noise, or stripes or some other artifact.

To fix it, you can either keep the light offset slightly from the plane, or set up the wall object not to cast shadows (assuming it doesn't need to in the scene).


It looks to me like you're resetting the distance calculation at each grid intersection? So every time you march another grid cell you're dividing the light intensity by $gridsize^2$ ? That would explain your observed results, since dividing by $1^2$ equals no change. It won't give you the correct result though, since to get the right attenuation you need to track the total distance back to the light source and square that.

Also note that the attenuation term used in rendering is usually a bit more complicated than just $1/D^2$ - generally it's a quadratic, $1/(aD^2 + bD + c)$. While not physically correct, it's more tunable and easier to get visually pleasing results. In your case, for example, using $1/(D^2 + 1)$ would prevent the edges immediately surrounding the light from being brighter than the light itself.


Without knowing the internals of Blender, I would say it uses shared normals for the quads, splitting the normals only on edges between the quads, not between the triangles. So your idea of uploading a normal per quad is probably close to the result you are looking for.

If you load the mesh as quads, you are set. If you are using an old version of OpenGL where quads still exist, you may even be able to use solution 3).

If you load the mesh as triangles, you will have to figure out which triangles originated in quads. You can likely do this simply by enforcing split normals beyond a fixed angle between adjacent triangles.

I am not sure there is a simple dynamic solution, where you just render any mesh. You may be able to use the geometry shader to look at the adjacent vertices, but the easiest version is to smooth the normals once before upload.


Hodge Decomposition in fluid simulations

Suppose I have a fluid which, like water, is incompressible; that is, its density is constant. I can find out how this fluid is flowing; let's say I inject dye at various points and follow how the dye particles move. I can estimate a velocity at each point, a velocity field $\mathbf{w}$. But water is incompressible; this means that the true velocity field $\mathbf{u}$ has to be mass-conserving. I want to "correct" my estimate $\mathbf{w}$ to find the closest mass-conserving velocity field.

The local accumulation or depletion of the fluid is computed with the divergence of the vector field, $\operatorname{div}\mathbf{w}$. (Stam writes the divergence of $\mathbf{w}$ as $\nabla\cdot\mathbf{w}$.) The divergence is a derivative, which is important here because that means it's linear; the divergence of the sum of two vector fields is the sum of their divergences:
$$ \operatorname{div}(\mathbf{u} + \mathbf{v}) = \operatorname{div}\mathbf{u} + \operatorname{div}\mathbf{v}. $$
We want to find a mass-conserving field; that is, we want to find a field with zero divergence. That means that we want to construct a field $\mathbf{v}$ with the same divergence as $\mathbf{w}$; then
$$ \operatorname{div}\mathbf{w} - \operatorname{div}\mathbf{v} = \operatorname{div}(\mathbf{w} - \mathbf{v}) = 0 $$
so we can get our mass-conserving field $\mathbf{u} = \mathbf{w} - \mathbf{v}$.

The only question is how we construct a field with a given divergence. For reasons I'll explain a few paragraphs down, we say that $\mathbf{v}$ is the gradient of some scalar field $q$, $\mathbf{v} = \operatorname{grad}q$. (Stam writes this as $\nabla q$.) Then we have
$$ \operatorname{div}\mathbf{w} = \operatorname{div}\mathbf{v} =
\operatorname{div}\operatorname{grad}q.$$
It's worth looking at what this says. The divergence of a vector field is a scalar field; at each point, it gives the local accumulation or depletion of the fluid. Both the divergence and the gradient are derivatives (which is why they can also both be written with $\nabla$). The equation is thus saying that the scalar field (the divergence of $\mathbf{w}$) is a particular second derivative of the scalar field $q$. This derivative is called the Laplacian (Stam writes it as $\nabla^2 q$). The differential equation relating scalar fields $q$ and $a$
$$ \operatorname{div}\mathrm{grad}\ q = a $$
is the Poisson equation. It's not so much that you use it to estimate $q$ as that it's the actual definition of $q$.

So, to "correct" our estimate $\mathbf{w}$, we can compute its divergence (a scalar field)
$$ f = \operatorname{div}\mathbf{w}, $$
solve the Poisson equation to find the scalar field $q$
$$ \operatorname{div}\operatorname{grad} q = f, $$
then subtract $\operatorname{grad}q$ from $\mathbf{w}$ to find our mass-conserving velocity field
$$ \mathbf{u} = \mathbf{w} - \operatorname{grad}q. $$

Above we just assumed that $\mathbf{v} = \operatorname{grad}q$ for some $q$. What we actually want is to find a field $\mathbf{v}$ that is in some sense the "smallest" field with the given divergence, so we can find the "closest" mass-conserving field $\mathbf{u} = \mathbf{w} - \mathbf{v}$; here I'll explain how the gradient of a scalar field satisfies this requirement.

We have a fixed amount of fluid moving through our space.
If our flow is mass-conserving, then we have a conserved amount of fluid, without local accumulation or depletion. If we have any motion at all, then, the fluid must move in cycles. In the same way as we measure the local accumulation with the divergence operator, we can measure the local cyclic behaviour of the flow by using another vector operator, the curl (which can be written $\nabla\times\mathbf{v}$). The curl at a point gives, broadly speaking, the axis about which the flow is circling near that point.

Now, consider the gradient of a scalar field, $\operatorname{grad} q$. (For instance, $q(x,y) = \sin x \cos y$.) The potential can be drawn as a height field:

Then the gradient of this field puts at each point a vector pointing in the steepest uphill direction:

But since you're always going "uphill" as quickly as you can, I think it's clear that there's no circulation in this vector field! In fact, for any scalar potential $q$,
$$ \operatorname{curl}\operatorname{grad}q = 0 ! $$

There's also a second similar identity: $\operatorname{div}\operatorname{curl}\mathbf{v} = 0$; unfortunately, I can't think of a good intuition for that one. (But look at this Math.SX answer for a more mathematical explanation.) We can use these two identities to explicitly write the Hodge decomposition: $\mathbf{w} = \mathbf{v} + \mathbf{u}$, where $\mathbf{v}$ has no curl, and $\mathbf{u}$ has no divergence, or
$$ \mathbf{w} = \operatorname{grad}q + \operatorname{curl}\mathbf{\nu},¹ $$
for some scalar field $q$ and vector field $\mathbf{\nu}$.

¹The Hodge decomposition is actually more complicated if your space isn't simple, i.e. finite and ball-like; on more complicated spaces, like say the surface of a torus, we have curl-free fields which aren't the gradient of anything, and divergence-free fields which aren't the curl of anything. For most fluid modelling in CG, though, this isn't usually an issue.


Calculating intersection of polyline and line

Given that you want to test for intersection against rays with many different starting points and directions, it's worth investigating raytracing-style acceleration structures, such as the bounding volume hierarchy (BVH). In 2D, this would look like a tree of axis-aligned bounding boxes that divide up the space. Each leaf node of this tree would store a reference to a small list of nearby segments of the polyline, as well as a bounding box that contains them. Higher nodes would store a bounding box that contains all their children.

There's a ton of material on BVH construction and traversal out there, so I'd recommend doing a few searches and reading up. It's not too difficult to implement a top-down construction method, by starting with the bounding box of the whole polyline and recursive splitting. Then, to use it for acceleration, when you do a ray query you can find which BVH nodes the ray touches and test only the segments in those nodes for intersections.

However, all that being said, I'd like to propose a different approach to your original problem of working with contour lines. It seems to me that this is a place where Delaunay triangulation can probably be very helpful. It basically constructs a triangle mesh that connects a set of input vertices and edges. For example, here is a contour map I found on the internet:

And now here it is with a constrained Delaunay triangulation, calculating using Jonathan Shewchuk's "Triangle" program:

As you can see, the triangulation gives you a lot of information about neighboring regions in the contour map, so I suspect that you could find a way to get what you want from it somehow. For instance, you could find all the triangles that touch the two contour lines you're interested in, and look at their edges, midpoints, and so on.


Programmatically generating vertex normals

Computing the normal from vertex positions is quite simple using the vector cross product.

The cross product of two vectors $u$ and $v$ (noted $u \times v$, or sometimes $u \wedge v$) is a vector perpendicular to $u$ and $v$, of length $||u \times v|| = ||u|| \cdot ||v|| sin(\theta)$, with $\theta$ the angle between $u$ and $v$. The direction of the vector will depend on the order of the multiplication: $u \times v$ is the opposite of $v \times u$ (the two directions perpendicular to the plane).

If you are not familiar with the cross product, I invite you to read about it and get comfortable with it. Normals will then seem simple.

Flat shading normals

If you have a triangle $ABC$, $AB \times AC$ is a vector perpendicular to the triangle and with a length proportional to its area. Since the normal is the unit vector perpendicular to the triangle's plane, you can get the normal with:

$N = \dfrac{AB \times AC}{||AB \times AC||}$

In code, this would look like n = normalize(cross(b-a, c-a)) for example. Just apply this over all your faces and you will have your normals per face.

For each triangle ABC
    n := normalize(cross(B-A, C-A))
    A.n := n
    B.n := n
    C.n := n

Note that this assumes vertices are not shared between triangles. I am not familiar with the Kinect API; it's quite possible that they are shared, in which case you would have to duplicate them, or move on to the next solution:

Smooth shading normals

After lighting with normals computed as above, you will notice that triangle edges are apparent. If this in not desirable, you can compute smooth normals instead, by taking into account all the faces that share a same vertex.

The idea is that if a same vertex is shared by three triangles $T1$, $T2$ and $T3$ for example, the normal $N$ will be the average of $N1$, $N2$ and $N3$. Moreover, if $T1$ is a big triangle and $T2$ is a tiny one, you probably want $N$ to be more influenced by $N1$ than by $N2$.

Remember how the cross product is proportional to the area? If you add up the cross products then normalize the sum, it will do exactly the weighted sum we want. So the algorithm becomes:

For each vertex
    vertex.n := (0, 0, 0)

For each triangle ABC
    // compute the cross product and add it to each vertex
    p := cross(B-A, C-A)
    A.n += p
    B.n += p
    C.n += p

For each vertex
    vertex.n := normalize(vertex.n)

This technique is explained in longer detail this article by Iñigo Quilez: clever normalization of a mesh.

For more on normals, see also:

How to Calculate Surface Normals for Generated Geomtery


What's the point of "work groups" for OpenGL compute shaders?

Work groups share resources. The most notable resources they share are barriers and LDS (Local Data Storage aka shared memory in GL lingo, aka Thread Group Shared Memory).

See Efficient Compute Shader Programming for an example where LDS/TGSM is used to accelerate a simple gaussian blur.

It is an interesting thought exercise of how a simple algorithm (weight average N values and output the result) can be rewritten to better leverage the hardware. Basically with a kernel radius of 5, Pixel A will load pixels [A; E] to average them, pixel B loads pixels [A; F] to average, pixel C loads pixels [A; G], and so on until we hit pixel F which will load pixels [B; J] to average them. So the paper asks the question: why do we load pixels A through F from main RAM so many times? (note: GPUs DO have caches, but if multiple threads are requesting all the same regions is likely many of these loads will miss the cache or have bank conflicts) and instead tries to load A through F ONCE (each thread loads only one pixel) and store it into TGSM memory which is much faster (and can also be accessed directly by the HW). It's basically very similar to manual cache management.

Also if the workgroup happens to be <= than the GPU's native workgroup size, barriers are much cheaper.

Note not all GPUs are created equal. For example current Mali's GPU do not have LDS, and shared memory is "emulated" via RAM fetches, thus trying to use LDS like a cache is only going to hurt performance.


Generating the end caps of a cylinder using triangle strips for use in a geometry shader

If the end of your cylinder is convex, with vertices $V_0, V_1, V_2,...V_{N-1}$ then just replace your fan order (i.e. as above) with this strip $V_0, V_1, V_{N-1}, V_2, V_{N-2}, V_3...$, i.e. just zig-zag across the polygon.


http://www.cs.sandia.gov/~samitch/papers/vor_final.pdf There isn't the computer algorithm in any programming languague, but u should be able to replicate it easily with some Plane Reflections and the algorithms from the above link.


Is denoising ALWAYS about doing a low pass filter / blur?

No, but this is the most obvious technique. A good denoiser isn't just a filter that runs on the image, but actually performs the reconstruction; i.e. it's a function from random samples to an image, not a function from an image to an image.

  Or are there other ideas and techniques for removing noise?

Yes. For years, the main family of denoise techniques worked in two parts. First, a line detection algorithm would detect the true edges in the image. Then, the output in this would set weights in a filter kernel to perform some kind of guided blur, so that the blur would be strongest (i.e. remove more frequencies) where there are no edges, and weakest (i.e. allow high frequencies) where there were edges. This would remove more noise from gradients while preserving edges.

Nowadays, deep learning has provided a new family of techniques. Using the same technique that's producing the noise (i.e. path tracing if you want to denoise path-tracer renders), produce lots of image pairs where one has high noise (few samples) and one has low noise (high sample rate). Feed these pairs into a convolutional network, to train it to produce the low-noise images from the high-noise images. Then feed your high-noise image into the trained network, and enjoy the low-noise output. Variations on this technique include having the network produce a filter kernel for each image patch, instead of learning the output image directly; or classifying an image patch according to its structure (what edges are present), so similar patches can be handled similarly. Some networks can even take the raw samples as input and directly produce a denoised image; others can guide your adaptive sampling algorithm.


Is there some kind of Bresenham algorithm or equivalent for scanline rendering a rotated ellipse?

The classic book Computer Graphics: Principles and Practice (second edition) by Foley, van Dam, et al. describes such an algorithm in section 19.2.

The explanation in the book seems to come from an MSc thesis, Raster Algorithms for 2D Primitives by Dilip Da Silva.

See also these papers:

Curve-drawing algorithms for raster displays by Van Aken and Novak (ACM TOG, 1985).
A tracking algorithm for implicitly defined curves by Chandler (IEEE Computer Graphics and Applications, 1988).


A few more options:

don't use a scenegraph. They don't map well to the data a gpu needs for rendering. A flat array of object tends to work much better for the cache as well
separate the scenegraph into render-related data and non-render-related data. Then only double buffer the render-related data.

The render related data would end up in a flat array containing per rendered object: reference to the mesh, the final model matrix, a reference to the texture

Flat arrays are much easier to double buffer than node-based data structures.


Short answer: Yes, It can be done. But no one does so.

Long answer:
Scene graphs can be stored and processed on a GPU using OpenCL/WebCL. But it is not practical to do so. Updating scene graphs (a tree not in flat memory) on a GPU is slow, and processing it on a GPU is also slow (again, the tree is not in flat memory), while computing transformation matrices is cheap on a CPU. Thus no one is doing this.


Dealing with OpenGL versions across multiple computers

After running across an aside in a tutorial and doing some further digging, I think I understand this now.

According to the OpenGL ABI specification, only the functions defined in OpenGL 1.2 need to be present in the GL/gl.h header file. (At least for Linux - other platforms may be different. This answer indicates that Windows may expose only OpenGL 1.1 functions by default.) Functionality from more recent OpenGL versions (even commonly used functions) may not even be defined in the header files. -- This seems to be the case for glDrawElementsInstanced. In the headers of one Linux version I looked at, by default glDrawElementsInstanced isn't exposed by the the headers, even though OpenGL 4.5 is supported.

Instead, the functionality from newer OpenGL versions are apparently handled through the OpenGL extension system. (Even if it's not strictly an "extension", but part of standard functionality.) That is, each individual user is nominally responsible for loading the function address into the appropriate function variable using the appropriate (platform-specific) *GetProcAddress() function from the GL/glext.h header. The usage is to pass the appropriate string to the function, which will return a function pointer to the appropriate function (if present), or a null pointer if the functionality can't be found. (See the tutorial section for more detail on usage.)

This is, obviously, somewhat less than satisfactory if you're using modern OpenGL features. There is a preprocessor define GL_GLEXT_PROTOTYPES which, when defined prior to including GL/glext.h, will enable the definitions of a large number of standard functions. However, this may not work robustly if you're running on a computer which doesn't implement all of those functions which you're using. (You're likely to get missing symbol error during shared library loading.)

The better way to handle things is to use the OpenGL Extension Wrangler (GLEW). This hides the heavy lifting of the whole *GetProcAddress() process, defining symbols like glDrawElementsInstanced directly in the glew.h header. So long as you call glewInit() after initializing OpenGL but before using any version dependent function, these functions should be initialized appropriately if available, or to a null pointer if not.

You can either test for the function pointer being non-null, or use the GLEW-provided version test variables like GLEW_VERSION_3_1, which are booleans (not compile-time defines) indicating if that OpenGL version is supported on the computer the program is running on. Alternatively, you can test on a feature-by-feature level by using glewIsSupported()/glewGetExtension(). You can look up the respective strings for individual extensions at the OpenGL Registry. It doesn't lend itself to forward searching, but may permit use of features prior to the nominal version of OpenGL in use.

So why wasn't this clearer to me before? I think it's because most of the tutorials and documentation about OpenGL simply assume that you're using GLEW. Despite its name, GLEW isn't just about extensions, it also abstracts much of the version-to-version compatibility issues. As long as you put use behind version checks and don't actually attempt to call a non-supported function, GLEW hides all the pesky details about functions which aren't actually being implemented.


Where should I project a polygon corner when it is behind me?

Trying to explain a problem seems to help the thought process. This is what eventually worked for me:

I realized that I can easily find the point (in 3D space) where a line between two of the corners brakes through one of the walls in my pyramid field of view, given the coordinates of the out-of-view corner and one of the visible ones. It is only a matter of finding the intersection between two lines. The position of the intersection point gives away the height on the frame wall.

Self answer, yes, but nothing came to mind before posting at SE.


What can be the cause of some strange ray shadows artefacts on a perfectly flat surface?

It appeared that it was exactly the problem explained in the link of my original post (Cause of shadow acne).
The shadow function offset has been increased and the problem is solved!


After a deformation operation on polygons, how can I check for and fix inverted polys?

If you have an oriented triangle mesh (a non oriented surface would be the Moebius Strip for instance), you can check the triangle inversion by looking at the signed area of triangles. The signed area is positive when the vertices are in CCW order and negative if they are in CW order. If your triangle vertices are all sorted in CCW then all triangles with negative area are inverted. Fixing that depends on the deformation technique, I don't know your case.

Checking self-intersection is more involved. You should actually do Collision Detection of the mesh with itself. Collision detection is fast but need to implement proper data structures such as AABB Tree for narrow the amount of triangle-triangle intersection tests. Collision Detection is usually performed using very low-poly version of the meshes. In some computer games, the bounding volume hierarchy (spheres, AABBs, etc) that approximate the mesh is the only thing considered for collision, no the mesh itself.

Once a collision is detected, there are several ways to response, applying a  bouncing force for instance, is common.


Ray tracing scene with 2 lights and sphere

I think what you're describing is what happens when you mix colours in gamma space. You should be using linear colour space when mixing.

When reading colours of lights and materials, convert them from gamma->linear using the inverse gamma function.

colourLin = pow(colourLight, 2.2);
... do lambert calculation here and accumulate in 'colour' variable ..

// On output (linear->gamma conversion)
gl_FragColor = pow(colour, 1.0/2.2);

There is an in depth chapter in GPU Gems on this topic: The Importance of Being Linear.


Because the geometry will be subdivided into quads. Anyway the primary reason for this is if you have a function with 2 directions it is easer for humans to think it out andcanticipate results if its a quad.

Software renderers dont have to rely on triangles. And did not for a long time. But then if you microdice your surface who cares its not like you can see the micropolygon anyway.

For more on this read this post
You might want to control the triangulation yourself. Interpolation does not happen the same way regardless of how the triangles are laid out. So being able to triangulate and then flip offending triangles is useful.


Anti-aliasing - Controlling colour density of pixel that comes within the "rectangle" surrounding my line

You have a background color and the color of your line, an antialiased line has additional lines drawn on either side of the first that are simply part way between the color and brightness of your origonal line and the background.
If the line you are drawing is horizontal or vertical the additional antialiasing lines (however many, depending upon the intended resulting line thickness) should be closer in color and brightness to the intended color - there is absolutely minimal antialiasing on such lines, that there is a tiny amount is due to offset subpixel rendering on some but not all matrixes.
This shows the difficulty:

When the lines are not horizontal or vertical it gets more complicated ...
If the additional lines drawn on each side of your line are drawn:

with color and brightness more closely matching the background's color then your line remains sharper and the background color near your line is fuzzier.
or

with color and brightness more closely matching your line's color then your line becomes fuzzier (wider) and the background color near your line remains sharper.

That means that if you intend to draw a thin line, a single pixel wide line, you would use the first option. (Make the result skinny).

If you want a two pixel wide line you would choose the second option but weigh the color of the background more towards the background's color and brightness. (Make the result not as skinny as a single line, with each of the two additional lines forming almost ~50% of the total width).

If you want a three pixel wide line you still use the second option but would instead weigh the color of the background more towards the line's color and brightness. (Make the result clearly much fatter than a single line).

For a four pixel wide line go back to the first method where you draw two lines of the intended color and brightness surrounded by two lines using the second method, making the two lines wider. (Sort of a blending of the first and second methods)

For five pixel wide lines the center is exactly your intended color, the two nearest surrounding line are very close to the color and brightness of the center line, with the outer two fairly close to the color of the center line (but less than the ones immediately next to the center line)

For a six pixel wide line there are two central lines of the intended drawing color and brightness, surrounded by two very similar lines, surrounded by another two less similar to the central two's color and brightness; but still very different from (lighter or darker) than the background color.

A seven pixel line or greater follows number 5 or 6 above.

The weighing of the additional lines more towards the intended color gives a sharper line, if you thought of it from a sideview analogy it would be top hat shaped. If you weighed additional lines more towards the background color the  (impossible) sideview analogy is that the profile of the 2D line would be rounded or triangular. If you are drawing 2D shapes it is in your interest to weigh towards the intended color. If you are rendering a 3D shape on a 2D screen then favoring the background color blends towards subsequent lines of a different Z depth.
A simplistic but effective algorithm with good speed is the Xiaolin Wu improvement to Bresenham's line drawing algorithm to implement antialiasing.
// Pascal
function plot(x, y, c) is
    plot the pixel at (x, y) with brightness c (where 0 ≤ c ≤ 1)

// integer part of x
function ipart(x) is
    return floor(x)

function round(x) is
    return ipart(x + 0.5)

// fractional part of x
function fpart(x) is
    return x - floor(x)

function rfpart(x) is
    return 1 - fpart(x)

function drawLine(x0,y0,x1,y1) is
    boolean steep := abs(y1 - y0) > abs(x1 - x0)

    if steep then
        swap(x0, y0)
        swap(x1, y1)
    end if
    if x0 > x1 then
        swap(x0, x1)
        swap(y0, y1)
    end if

    dx := x1 - x0
    dy := y1 - y0
    gradient := dy / dx
    if dx == 0.0 then
        gradient := 1.0
    end if

    // handle first endpoint
    xend := round(x0)
    yend := y0 + gradient * (xend - x0)
    xgap := rfpart(x0 + 0.5)
    xpxl1 := xend // this will be used in the main loop
    ypxl1 := ipart(yend)
    if steep then
        plot(ypxl1,   xpxl1, rfpart(yend) * xgap)
        plot(ypxl1+1, xpxl1,  fpart(yend) * xgap)
    else
        plot(xpxl1, ypxl1  , rfpart(yend) * xgap)
        plot(xpxl1, ypxl1+1,  fpart(yend) * xgap)
    end if
    intery := yend + gradient // first y-intersection for the main loop

    // handle second endpoint
    xend := round(x1)
    yend := y1 + gradient * (xend - x1)
    xgap := fpart(x1 + 0.5)
    xpxl2 := xend //this will be used in the main loop
    ypxl2 := ipart(yend)
    if steep then
        plot(ypxl2  , xpxl2, rfpart(yend) * xgap)
        plot(ypxl2+1, xpxl2,  fpart(yend) * xgap)
    else
        plot(xpxl2, ypxl2,  rfpart(yend) * xgap)
        plot(xpxl2, ypxl2+1, fpart(yend) * xgap)
    end if

    // main loop
    if steep then
        for x from xpxl1 + 1 to xpxl2 - 1 do
           begin
                plot(ipart(intery)  , x, rfpart(intery))
                plot(ipart(intery)+1, x,  fpart(intery))
                intery := intery + gradient
           end
    else
        for x from xpxl1 + 1 to xpxl2 - 1 do
           begin
                plot(x, ipart(intery),  rfpart(intery))
                plot(x, ipart(intery)+1, fpart(intery))
                intery := intery + gradient
           end
    end if
end function


There are no easy answers here, if you are serious then it is possible to replace cubemaps with hemispherical maps, which result in just 2 offscreen renders each instead of 6. They are not popular because the math and setup are nontrivial. RTR4 has references to papers. If you get that working do a writeup here please! Look for Dual-paraboloid shadow mapping.
Another solution is to only update your maps every other frame, or every third frame. So if you have say 12 maps, then update 4 of them each frame round robin style. Also frugal culling can greatly reduce the number of faces that have to be rendered but quality is reduced. I see so many implementations that render all the geometry in the scene to every face when much of that geometry is outside the frustum of that face. This alone can increase performance but scene complexity will play a major role in overall performance.
As far as infinite reflections is concerned, think about that for a second, infinite rays is going to take infinity to render no matter what your technique is. If you manage to trap photon's in a perfect reflector in real life then you can't observe it because observing it destroys the perfect reflector. I'm not trying to be sarcastic, I'm just giving a reality check.
But simple techniques can get good reflections in real time. Screen Space reflections are tricky little buggers to get right but when done well can yield very nice results. You don't see many implementations around because there are a lot of corner cases to consider. Again if you get a good implementation of this working, please share it here!
What you are implying in your post sounds like you would, at the very least, have to do reflections against all the cube faces, this too is not realistic, even on the best hardware, your computes are better spent someplace else, like ray tracing...


Why if we rotate an object by R and then translate it with T, from the object perspective, we're actually applying $-T$ and then $-R$

If you have object$\rightarrow$world space transformation matrix:
$$M=T*R$$
then inverse (i.e. world$\rightarrow$object) of this transformation matrix is:
$$M^{-1}=(T*R)^{-1}=R^{-1}*T^{-1}$$
So the order of multiplication of matrix inverses is reversed per matrix inversion rules.


Euler rotations work around the axes. This is represented in 3 angles (rot_x, rot_y, rot_z). The order of the rotation needs to be determined as the order affects the result.

For example, if we have a point (1, 0, 0) and we rotate 90 degrees around the Y and then 90 degrees around the Z, then the result will be

v1 = (1,0,0)
v2 = rotate_y(v1, 90) # v2 = (0,0,1)
v3 = rotate_z(c2, 90) # v3 = (0,0,1)

When doing the other order first Z and then Y, the result is different:

v1 = (1,0,0)
v2 = rotate_z(v1, 90) # v2 = (0,1,0)
v3 = rotate_y(c2, 90) # v3 = (0,1,0)

In computer graphics systems, quaternion rotations are often used, which allows rotation around any axis. This then solves many of the problems with Euler systems.


I did some research and found the answer I was looking for.  The three most common ways to interpolate vectors are:

Slerp - short for "spherical interpolation", this is the most correct way, but is also the costliest.  In practice you likely do not need the precision.
lerp - short for "linear interpolation", you just do a regular linear interpolation between the vectors and use that as a result.
nlerp - short for "normalized linear interpolation" you just normalize the result of a lerp.  Useful if you need your interpolated vector to be a normalized vector.

In practice, lerp/nlerp are pretty good at getting a pretty close interpolated direction so long as the angle they are interpolating between is sufficiently small (say, 90 degrees, or 120 like Simon mentions in his comment), and nlerp is of course good at keeping the right length, if you need a normalized vector.  If you want to preserve the length while interpolating between non normalized vectors, you could always interpolate the length and direction separately.

Here is an example of the three interpolations on a large angle.  Dark grey = start vector, light grey = end vector.  Green = slerp, blue = lerp, orange = nlerp.

Here is an example of a medium sized angle (~90 degrees) interpolating the same time t between the angles:

Lastly, here's a smaller angle (~35 degrees).  You can see that the results of lerp / nlerp are more accurate as the angle between the interpolated vectors gets smaller.

If you do lerp or nlerp, you can definitely do both bilinear as well as bicubic interpolation.

Using slerp, you can do bilinear interpolation, but I'm not sure how bicubic would translate.

I generated these images by taking screenshots from an interactive shadertoy demo I made to demonstrate these differences.  You can see that shadertoy here:
https://www.shadertoy.com/view/4sV3zt


Translate an arcball camera upwards or left-right

Here's a basic camera with panning (some parts removed for brevity). It's not arcball but should get the idea across. Basically it just transforms some directions of a plane in view space to one that is in world space.

rotate and pan are called on mouse move with the change in mouse position.

// Updates the view matrix based on parameters
void updateCam()
{
    // Calc eye position on sphere
    glm::vec3 eyePos = glm::vec3(
        rho * sin(phi) * cos(theta),
        rho * cos(phi),
        rho * sin(phi) * sin(theta));

    setEyePos(eyePos + focalPt);
    setFocalPt(focalPt);
    updateLookAt(); // Uses eyePos and focalPt to calculate a view matrix
}

// Input is a change in position
void rotate(glm::vec2 diff)
{
    theta += diff.x * rotateSpeed;
    phi -= diff.y * rotateSpeed;

    // Clamp to prevent user flipping camera at top of sphere
    phi = clamp(phi, 0.01f, 3.14f);
    updateCam();
}
void pan(glm::vec2 diff)
{
    // Up and right are the directions to pan when camera is facing forward
    glm::mat3 invView = glm::inverse(view);
    glm::vec3 dx = invView * glm::vec3(1.0f, 0.0f, 0.0f);
    glm::vec3 dy = invView * glm::vec3(0.0f, 1.0f, 0.0f);

    // Scale shift with scale so shift is relative to how far you are zoomed out
    focalPt += (dy * diff.y - dx * diff.x) * shiftSpeed;
    updateCam();
}


What algorithm is used in the animation of the uTorrent 's "About window"?

As you discovered and mentioned in your self-answer, the pattern in the background appears to be a sum of sinusoidal gradients.

However, the example linked to in your answer is more complicated than that used by µTorrent. The background of the About window appears to be a static pattern, rather than the animated sinusoidal pattern used in the plasma post.

Several sinusoidal gradients have been summed to give a single image, and the illusion of movement is given by simply cycling the colours in that one image, rather than generating a number of different images. This is most noticeable if you focus on the centre of one of the rings of colour. In the µTorrent pattern you will notice that each ring stays in one place, and has colour flowing either into it or out of it. In contrast, the rings of colour in the fully animated pattern move around, occasionally dividing or merging.

The simplified approach used by µTorrent is reminiscent of animations used in the past when recalculating the sinusoidal patterns each frame was not realistic.


Trouble making path tracer multithreaded, read access violation

I found out what was wrong. It wasn't related to multithreading at all.
I was taking in a vector in a constructor hitable_list(std::vector<hitable*> listVector) and using a pointer to the first element of the vector elsewhere. So as soon as the vector went out of scope the pointer became invalid.

To solve this issue I just removed this constructor and whenever I need to pass in a vector I make sure that the vector stays in scope always and directly pass in the pointer to the first element of the vector, and it's size.


Phong: Real-time rendering of Metal and Plastic Materials

Diffuse colours on materials typically come from within the material, while the specular colour is from the very surface. Coloured plastic materials are made by embedding particles of dye inside a colourless medium, so the diffuse colour is the colour of the dye, while the specular colour is white from the colourless surface. With metals, all of the reflection is at the surface, so the diffuse and specular colours are the same.

There are many different kinds of paint, but most commonly they're made like plastics: putting dye particles inside a white or colourless medium. Emulsion, acrylic, and enamel paints are all made this way. For this reason, painted or varnished metal surfaces have white specular as if they were plastic. There are some kinds of paint (metallic paints being the most obvious example) that aren't made this way, so it's not always true.

Ambient doesn't correspond to a physical property of the material: it's a cheap substitute for global illumination. It represents all of the indirect light that falls on the shading point and is reflected in the viewing direction. There's no physically correct way to set it: you should just do what looks best in your scene. That said, I'd offer a rule of thumb that shinier objects probably want less ambient, because less of the indirect light is going to be reflected diffusely to where you can see it. However, you might consider not using ambient at all, and pointing dim fill lights at places that will come out too dark without indirect light.


How to crop with edge-detection using imagemagick

Try -trim instead of -shave:

convert original.png -trim trimmed.png

See the docs for this and other options.


Why do adjacent triangles never overlap when rasterized?

This is the same problem as discussed in What are sample gaps during scan conversion?

Briefly, rasterisation - at least with the majority of hardware systems - tests at a single 'infinitesimal' point for each pixel to determine if that pixel is inside a given triangle.

If the sample point lies exactly on an edge, then additional tie-breaking rules, often referred to as "fill rules", come in to play. Typically these will be something equivalent to, say, only consider the edge points on the top/left** edges of a triangle to be "inside".

**Other rules are available ;-)


PBRT path tracing: negation of ray direction before next event estimation

The incoming direction $\omega_o$ (coming from eye) is flipped in pbrt:

"Another convention we will follow is that the incident light direction $\omega_i$ and the outgoing viewing direction $\omega_o$ will both be normalized and outward facing after being transformed into the local coordinate system at the surface. By convention, the surface normal always points to the "outside" of the object...".


If I understand you correctly, you've added a new compute shader. In DX12, you'll need at least one separate PSO for each shader, since the shader code is part of the PSO.


Mix shader looks wrong on my path tracer

As suggested in the comments, I created an exact scene in Blender using the same .obj. After adding a mix shader node and setting one side to diffuse and the other to a perfect mirror, you can see that the results are almost identical. This would indicate that my implementation is correct.

Here is a side by side comparison. Mine is on the left, Blender's is on the right:


What algorithm to reassemble a broken image?

If the pieces are fairly close to where they were, ICP could be of help (under the assumption that the pieces are planar polygons).


Metal and Vulkan are low-level 3D graphics APIs, in much the same way that OpenGL and Direct3D are 3D graphics APIs. Metal is Apple's API, usable on iOS and MacOS. Vulkan was developed by the Khronos Group, who also oversees the various forms of OpenGL. Vulkan is directly usable on all non-Apple platforms, and there is the MoltenVK project which provides a translation layer from Vulkan to Metal (since the APIs are not that different), thus allowing Vulkan applications to run even on Apple platforms.

Metal, Vulkan, and Direct3D v12 represent a dramatic departure API-wise from prior low-level graphics APIs. They use a command-buffer-based interface, where you record rendering commands into a command buffer object, then execute those commands at some later period. The more immediate style of APIs like OpenGL and D3D-pre-12 have no intermediary command buffer. When you call a rendering function, it is given to the GPU "right then" and that's the end of it.

These styles of API are important because they allow threaded creation of rendering commands. In the immediate APIs, you can only effectively have one thread issuing rendering commands. In command-buffer APIs, multiple threads can build multiple CBs simultaneously, and then one thread submits them all as a single bunch of work. This allows better utilization of CPU resources for CPUs with lots of cores.


R16_FLOAT, SampleCmp and CheckFormatSupport

Because I was using the CheckFormatSupport wrong. The second parameter is supposed to be an input.

Here is what it should be like.

bool Renderer::InitRenderer(HWND hwnd)
{
    if (!InitD3D11App(hwnd)) return false;
    InitViewport();

    UINT formatSupport;
    HRESULT hr = d3d11Device->CheckFormatSupport(DXGI_FORMAT_R16_FLOAT, &formatSupport);
    if (!SUCCEEDED(hr) || !(formatSupport & D3D11_FORMAT_SUPPORT_SHADER_SAMPLE_COMPARISON))
    {
        printf("Not supported!\n");
        return false;
    }

    return true;
}

https://www.gamedev.net/forums/topic/700419-r16_float-samplecmp-and-checkformatsupport/?tab=comments#comment-5397459


How to do texturing with OpenGL direct state access?

You still have to bind the desired texture to the texture unit to use it for rendering. In your current code, you're not specifying which texture to use for rendering, so the GL driver doesn't know which one you want to use and is defaulting to "no texture".

There are a few ways to bind a texture. You can keep using glActiveTexture+glBindTexture as in your non-DSA code, or you can use either glBindTextures or glBindTextureUnit.

For example, with glBindTextureUnit, do this:

glBindTextureUnit(0, textureHandle);

The first parameter is the texture unit to bind to.

By the way, instead of using glTextureImage2DEXT (which is from the old and deprecated EXT_direct_state_access extension), it's probably better to stick with core OpenGL calls. You can initialize the texture by using glTextureStorage2D to allocate the memory and set the format, then glTextureSubImage2D to upload the data.

(Here's an article with more on the differences between EXT_direct_state_access and the newer core version of the functionality, if you're interested.)


How to avoid slowdown with 25-30 students running simple GPU kernels on 4 GeForce GTX 650 Ti s?

You need nvidia-docker. If your ssh server is in a Docker instance, you can assign each one to a different GPU. Then just set things up so only one team is sshing into each Docker instance. You don't need to do anything fancy with cudaSetDevice(), or risk students using the wrong device and messing things up for others, because each instance looks like it only has one GPU.


This is pretty straightforward matrix math and can be found on Wikipedia. It's quite elaborate, though, if you want to learn the math behind it and the other types of projection that exist. I will encourage you to checkout either Real-Time Rendering or Computer Graphics, Principles and Practice.


How to generate binary gltf array buffers?

Technically any file format where you can dump the entire file into a VBO and then render from that will work for the .bin files. Unfortunately those formats are less well known than they should be.

Then you only need to adjust the bufferview elements to reference the proper subset of the file


Here's a more complete answer from the extent of my knowledge. A sort of brief overview you can google more about. I'm sure I didn't cover all the different techniques but most branch off the mentioned ones. Note: if it's a liquid we are more concerned with getting an accurate surface representation.

Surface Extraction Methods: The most common way would be to rasterize the particles into a 3d image/grid using a basis function, unless they're already in a grid. Then extract the polygonal surface from the image using marching cubes or some variant of it. Then render that polygon. There are some methods to extract polygonal surfaces directly from particles without rasterizing.
Ray Marching with an Implicit Surface: Again we rasterize the particles into an image but quickly convert into an implicit surface in our 3d image then ray march render the implicit surface, very fast.
Splat Volume Rendering (doing it in screen space): Splatting is good for visualizing extremely large datasets and used to be used a lot for volume rendering because of how fast it is. Basic idea is to rasterize all or most points as discs from the volume. You can do many effects, lighting, etc. Nvidia's screen space fluids is a form of this (SSF). Simple, fast, quality is acceptable for many cases.
Slicing/Texture Based Volume Rendering: Here we overlay 2d semi transparent images. You can just overlay some precomputed textures on various particles or compute a series of sliced 2d images a set distance apart in the 3d image/volume given the camera angle. More expensive. Gives better results than splatting.
Ray cast Volume Rendering: Similar to ray marching, we march along the ray from every pixel but instead of looking for an implicitly defined surface we sample the volume directly and map it too color and transparency then composite it back along the ray. For lighting we can do it similar to surface based rendering but use the gradient of the image as the normal. There are many other tricks.
Ray traced Volume Rendering: Basically the same as the ray casting but we bounce the rays around inside of the volume. They're distributed according to a phase function. I'm sure all the clouds in your favorite animated movies were rendered this way. Disney did this for their MPM simulated Snow in Frozen. The ray tracing approach is the most expensive of all of them, more expensive then surface based ray tracing. Not for real time as of writing this.

Note: The first two methods aren't for gasses. The last 3 are best for gasses due to their lack of surface representation.

I think it's also worth noting the value of having a stored representation for post processing effects. For manual deformations or effects by an artist. Procedural effects like spray, foam, bubbles, etc. Or just some extra geometric processing (smoothing or whatever).


Ray-tracing the Cornell Box results in really inconsistent image

Two symptoms

There appear to be two problems with the image.

The background is showing through along the line between adjacent triangles.
The colour displayed is not always from the closest intersection.

Note that the background colour is white, rather than black, due to the line:

                PutPixelSDL(screen, x, y, vec3(1, 1, 1));

Why this background shows through at the triangle edges I do not know. However, the second problem is identified below, so I recommend you fix that first and then the first problem can be analysed in isolation without the additional distractions.

Cause of the problem with ClosestIntersection

The ClosestIntesection function returns true when an intersection is found, and updates an Intersection object called closestIntersection. However, the intersection stored here is not always the closest to the camera.

The reason for this is here:

                closestIntersection.distance = glm::distance(start, x);

This appears in two places so it will need to be fixed in both.

start is cameraPos, presumably in world coordinates. x appears to be in the local coordinates of the triangle being intersected (one of the triangle vertices being the origin). The distance between two points in different coordinate systems is meaningless, so that in different regions of a triangle it will appear to be behind or in front of another triangle for no immediately apparent reason.

If x is first converted to world coordinates then the distance will be correct and the intersection chosen will be the closest to the camera.


Rendering 3-D point clouds obtained from depth maps as a 2.5D video

Depending on your application, you could use image-based solutions that are much cheaper than meshing a uniformly sampled set of points. For example, you could use the depth map as a height map.

There are many techniques that allow you to render it in an efficient way, the more realistic the more difficult to code as a shader (but not too much). From easiest to harder:

height maps
normal maps
displacement mapping
parallax mapping
geometry shaders to recreate the actual mesh from a set of primitives and the depth map
and so on.


I think the most efficient way might be either

return float(shadowValue > shadowTexZ);

or

return step(shadowTexZ, shadowValue);


Looking for a cylinder / axis-aligned box volume intersection test

As mentioned the problem of cylinder box intersection reduces to the problem of line-parallelogram intersection, and line segment-cylinder intersection.

I will assume that you know how to intersect a line and a parallelogram. Here is a sketch for solving the line segment-cylinder intersection:

The canonical equation of an elliptic cylinder aligned with the $Z$ axis and having center $(c_x, c_y, c_z)$ is given as:
$$\frac{(x-c_x)^2}{a^2} + \frac{(y-c_y)^2}{b^2} = 1$$
The parametric form of a ray is given as:
$$\vec{r}(t) = \vec{o} + t\vec{d}$$
Then by plugging r into the canonical equation one gets:
$$\frac{(o_x + td_x -c_x)^2}{a^2} + \frac{(o_y + td_y -c_y)^2}{b^2} = 1$$
Solve the quadratic equation for $t$ and you have your intersection. If the discriminant is negative then there is no intersection. To get a rotated cylinder, the easiest way is to apply the inverse rotation to the ray (the equation may be modified however too). Note that for a line segment, you simply need to check whether $t \in [t_0, t_1]$, where $\vec{a} = \vec{o} + t_0\vec{d}$ is the first point of the line segment, and $\vec{b} = \vec{o} + t_1\vec{d}$ is the second one.

EDIT: To clarify you further questions in the comments:

1) a,b in the denominator are the radii along the $X$ and $Y$ axes respectively. In the trivial case where $a=b=r$ we have the canonical equation for a circular cylinder:
$$(x-c_x)^2 + (y-c_y)^2 = r^2$$

2) For the intersection between the cylinder axis and the box a rotation is only required if your box is OOBB and you want to use your AABB intersection code.

3) The solution of the quadratic equation is trivial (note, I renamed $a,b$ to $r_1,r_2$ to further emphasize that they are the radii):

$$\frac{(o_x + td_x -c_x)^2}{r_1^2} + \frac{(o_y + td_y -c_y)^2}{r_2^2} = 1$$
$$r_2^2(o_x + td_x - c_x)^2 + r_1^2(o_y + td_y -c_y)^2 = r_1^2r_2^2$$
$$(r_2^2d_x^2+r_1^2d_y^2)t^2 - 2(r_2^2d_x(c_x-o_x) + r_1^2d_y(c_y-o_y))t + r_2^2(c_x-o_x)^2 + r_1^2(c_y-o_y)^2 - r_1^2r_2^2 = 0$$

Let us substitute:
$$A = (r_2^2d_x^2+r_1^2d_y^2)$$
$$B = (r_2^2d_x(c_x-o_x) + r_1^2d_y(c_y-o_y))$$
$$C = r_2^2(c_x-o_x)^2 + r_1^2(c_y-o_y)^2 - r_1^2r_2^2$$

Then you have the quadratic equation:

$$At^2 -2Bt + C = 0$$
$$D = B^2-AC$$

If $D<0$ there is no intersection between the ray $\vec{o}+t\vec{d}$ and the cylinder, otherwise:

$$t_1 = \frac{B-\sqrt{D}}{A}, t_2 = \frac{B+\sqrt{D}}{A},$$

give you the two intersections (it is one intersection if $D=0$).

Now let your segment vertices be $\vec{v}_a = \vec{o}+t_a\vec{d}$ and $\vec{v}_b = \vec{o} + t_b\vec{d}$. If $t_a \leq t_1 \leq t_b$ then you have an intersection at $t_1$, if $t_a \leq t_2 \leq t_b$ then you have an intersection at $t_2$, otherwise there is not intersection.


Maximal and minimal no. of angles obtained from clipping a convex polygon with n angles

Since this is a homework question, I'll give hints rather than a numerical answer.

Clipping a convex polygon

Think about how many times the polygon can cross each edge of the rectangle. In general, a polygon can cross one of the edges of the rectangle an arbitrarily large number of times. How is this number reduced if the polygon is convex?

This will give you the maximal number of new angles that can be added per rectangle edge, which should lead you to the answer.

Clipping a non-convex polygon

As for your guess of 2n, even for a non-convex polygon it may not be able to reach this high. If one vertex is outside the rectangle, then in order to create 2 new vertices both of its neighbouring vertices must be inside the rectangle. This means you cannot create 2n new vertices this way.

The exception is when an edge of the polygon crosses 2 edges of the rectangle, for example at a corner of the rectangle. This allows creating two new vertices from a single edge. So for small enough n, it is possible to clip to a polygon with 2n vertices, but for larger n this is not possible.


Why did AMD's putting an SSD on a GPU board reduce latency so much?

Answering from a somewhat layperson's perspective (I'm not a hardware expert, but I have built my own computer), I'll analogize, using the situation of getting some ingredients for cooking food.

Using "main system" RAM/storage is like you looking for an ingredient, not finding it in your cupboards, and, say, walking over to your neighbor and asking if they have any.  It turns out they do; they bring you it, you thank them and return to your house and use it.

Using "on board" RAM/storage is like you looking for an ingredient, finding it in your cupboard, and using it.

With on board RAM/storage, the GPU just has to do less work to obtain the information it needs.  The information is also simply closer, both physically and when considering the amount of hardware and software the information needs to go through.  With the on board RAM/storage, the GPU needs only look up where the information is, and then retrieve that information from RAM/storage, which then travels through maybe a couple inches of wire to the GPU's processing units.

With main system RAM/storage, there is a larger distance of wire for the signal to travel through which, although it may be only a minuscule fraction of a second, can still add up over thousands of calls to retrieve information.  Also, there are more potential bottlenecks, with all the different pieces that are interconnected.

Lastly, with on board RAM/storage, AMD has the opportunity to optimize the system, since they're controlling precisely hardware is being used.  This is just like how Apple tailors Mac OS to its Mac computers.


OpenGL with SFML, create an n-pointed star?

Try starting at the origin, then dividing 360 degrees by n to find some useful angles for the points of your star.
I hope this drawing will help you to find a solution.

Now that you have the angle theta, and the values of alpha and x are chosen by you, you have the tools to calculate all required vertices to make the star using trigonometry.


implementing right handed coord system APIs on top of a left-handed rendering engine

What you probably need is this 4 x 4 matrix:

$$\mathtt{T} = \begin{bmatrix}1 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 \\ 0 & 0 & -1 & 0 \\ 0 & 0 & 0 & 1 \\  \end{bmatrix}$$
which is flipping the z-axis. The inverse of this particular matrix is the matrix itself:
$\mathtt{T}^{-1} = \mathtt{T}$

For 3D points $\mathbf{x}$ multiply it from the left:
$$\mathbf{x}' = \mathtt{T}\,\, \mathbf{x}$$
and for 3D transformations multiply it from the right:
$$\mathtt{P}' = \mathtt{P} \,\,\mathtt{T}$$

Thus, you have not changed anything in the projection because:
$$\mathtt{P}'\,\,\mathbf{x}' = \mathtt{P} \,\,\mathtt{T} \, \,\mathtt{T}\,\,  \mathbf{x} =  \mathtt{P} \,\,(\mathtt{T}^{-1} \, \,\mathtt{T})\,\,  \mathbf{x} = \mathtt{P}\,\,\mathbf{x}$$

Quaternions $q=(s,x,y,z)^\top$ could be convert into
a 3D transformation matrix with:

$$\mathtt{P} = \begin{bmatrix} 1-2y^2-2z^2 & 2xy-2sz& 2xz+2ys & 0\\
              2xy+2sz & 1-2x^2-2z^2& 2yz-2sx & 0\\
              2xz-2sy & 2yz+2sx& 1-2x^2-2y^2 & 0\\  0 & 0 & 0 & 1\end{bmatrix}$$
and could then be treated like other 3D transformations.

However, 3D engines typically apply transformation in a transformation hierarchy. If this is the case, make sure you are only modifying the 3D transformation "closest" to the 3D point in the hierarchy.

Sources:
Graphics Programming I - 3D Transformations | Thorsten Thormählen


This is a very common question that most people in your situation (interested or fascinated shall I say in Computer Graphics) have and a problem they want to figure out. The solution is pretty simple, though not necessarily that straightforward to implement.

I can't make a specific answer because it depends on the OS you are using, and as I said this is not something you can do in 2 lines of code.

The answer is to use some sort of library that:

allows you to create a window on your system
allows you to display some content (an image) in the window you created and displayed to the screen

The first task is generally quite complex regardless of the OS you use (MacOS, Windows, Linux). For the second step, you need to generally go through some sort of graphics API (like OpenGL).

So there is hope for you. One thing you can do is use some other library that works as a wrapper around 1) and 2) In other words they help you especially with the creation of the window, and then if you know OpenGL you can do the rest of the work. You will need to create some texture that you will fill with the content of your image as it gets rendered and display this texture to the screen (that's pretty basic in OpenGL though you need to know this API in order to do so).

GLFW is a great library that will allow you to do that. It's simple to use and compile and install, etc. You can alternatively use Qt but I don't recommend it unless you want to become a professional programmer.

I had some documentation about the whole process but it's offline right now. I might put it back in a few weeks but you will need to wait I am afraid.

EDIT

I am editing the answer to respond to your comment. The algorithm would look like this:

Color *buffer = Color  char [w * h];
// create a texture that you can map onto the quad stretched over the area
// of the window. You will need to create a quad rendered in the space of the
// window, aka you in 2D not 3D.
...
// now render your image using you ray-tracer
for (j = 0; j < h; ++h) {
    for (i =0; i < w; ++i) {
        buffer[j * w + i] =  trace(orig, dir, ...);
    }
    // it is not very efficient to update your window each time you have
    // a new pixel so do this every row or every 10 rows ... you get the idea
    // copy content of your buffer into the GL texture
    // I am not using the right calls here, too lazy to find the right ones
    glCopyTexture(myGLTextureId, buffer, (sizeof(char) * 3 * w * h);
    // do the texture binding, etc... all GL or DX stuff
    ...
    // this will update your window with the content of the buffer
    glSwapBuffer();
}
// now save content of buffer to image file
saveToDisk(buffer);

Windows? Arg(. You will have to use DirectX calls indeed but they 2 APIs are very close.


OpenGL Unable to retrieve the values of data sent to shader using textures

What GL_LUMINANCE as a texture format actually does is store one component of data per pixel that is used for the red, green and blue color components. So, a classic greyscale image (It doesn't use it for the alpha component, though, that's why there was also GL_LUMINANCE_ALPHA for two compnent storage and GL_INTENSITY that broadcasts it to all components).

When using shaders you don't actually need this distinction anymore, though, since you can do the swizzling yourself in the shader. So you just use a proper typed 1-component internal format, most probably GL_R8 (i.e. 8-bit normalized fixed point red) in your case, but depending on your exact data type requirements really. And for setting it you would (in that case) use GL_RED as the format of the data in your RAM storage. So something like:

glTexImage2D(GL_TEXTURE_2D, 0, GL_R8, <width>, <height>, 0, GL_RED, <type>, <data>);

And in the shader you can do with that red color component whatever you want. If you just need a single value per pixel, you use that. If you need a grey color (i.e. GL_LUMINANCE behaviour), you swizzle it to the RGB components, and if you need it as color and opacity, you swizzle it into all four components. Or whatever else you want to do.

float data = texture(...).r;
vec4 luminance = texture(...).rrra;
vec4 intensity = texture(...).rrrr;

Though, this requires your shader to be somewhat aware of what texture you're working with. But if you want to use a more generic shader working with a 4-component color from an arbitrary texture, this might be too much of a hassle. In this case what you can also do is have the texture itself control its own swizzling behaviour. This is an OpenGL 3.3+ feature and can be set with glTexParameter and for good old GL_LUMIMANCE behaviour this would be:

glTexParameter(..., GL_TEXTURE_SWIZZLE_R, GL_RED);
glTexParameter(..., GL_TEXTURE_SWIZZLE_G, GL_RED);
glTexParameter(..., GL_TEXTURE_SWIZZLE_B, GL_RED);
glTexParameter(..., GL_TEXTURE_SWIZZLE_A, GL_ONE);

For a broadcast into all components, i.e. GL_INTENSITY behaviour, you can use the shortcut

glTexParameter(..., GL_TEXTURE_SWIZZLE_RGBA, GL_RED);

So then you can retrieve a normal 4-component color from your texture in the shader and get the 1-component it's storing conveniently broadcast into the RGBA components of the result.


Refraction in a Ray Tracer: What do with an intersection within the medium?

this is an interesting question (and I am actually an author on Scratchapixel so I can maybe help on that one)).

Things go as follows:

you cast the primary ray into the scene
the ray hits the glass which is a refractive-reflective/transparent material
you compute and cast two rays from the point of intersection: a reflective ray and a refractive ray
if the refractive ray hits an object contained in the transparent object (the glass) then you need to shade that object too. So if it's a plastic material, compute the result of that shader (including the lighting, so loop over all lights and add their contribution, etc.)
the color at the intersection point can finally be computed as a mix between the reflection and the refraction color (where as the mix is driven by the fresnel effect).

Now that works but this is essentially wrong.

Why? Because in fact light rays too are bent by say glass and water. So in fact when you cast a shadow ray from a point on the pen in the direction of the light, this ray direction is in fact wrong, since you should account for refraction and reflection (the ray should be bent when it goes from one medium to another in this case from water to air) but this is impossible to do with this algorithm/approach.

Though in general this gives acceptable results so no one cares much but if your goal is really a physically plausible image then this approach is wrong.

A solution to this maybe is to use photon mapping, where photons of light are cast from the light into the scene in a pre-pass. Thus when a light-photon hits the glass/water surface it will bent due to the law of refraction and eventually it will be "deposited" onto the pen. This is a bit hard to explain here but hopefully there will be a lesson on Scratchapixel about this in the future. Hope it helps.

EDIT

> Thanks for clarifying, I had one other issue though. In the source code provided on the website for that lesson there doesn't seem to be code present to spawn a ray when the light ray leaves the glass cylinder. Is this unnescessary or is it actually already present? If it is could you point out where in the code?

Hum it seems like my explanation wasn't clear if you ask this question. When the refractive ray hits the surface of the pen you need to compute the color of the pen at that intersection point (where the refractive ray intersected the pen's geometry). Now the pen is made out of wood for example which for the sake of simplicity is a diffuse material only (no specular). So in order to compute the color of the pen at the intersection point, you simply execute the traditional illumination loop:

Vec3f computePenColor(const Vec3f& objectColor, const Vec3f& intersectionPoint, normalAtIntersection) {
    Vec3f shadedPointColor = 0;
    // loop over all lights in the scene to cadd their contribution
    for (i = 0; i < scene.lights.size(); ++i) {
        Light currentLight = scene.lights[i];
        // let's say this light source is a point light source
        Vec3f lightDirection = currentLight.position - intersectionPoint;
        // this cast a shadow ray from the shaded point to the light position
        // if this ray intersects an object along the way then this point is in the shadow of this light and the function returns true
        bool shadow = castShadowRay(lightDirection.length(), lightDirection.normalize());
        // this point is not in the shadow of the light so add its contribution
        if (!shadow) {
            // the object is diffuse so apply cosine law only
            shadedPointColor += max(0, normalAtIntersection.dot(lightDirection)) * currentLight->intensity * currentLight->color;
        }
    }
    return shadedPointColor * objectColor;
}

So maybe as you can see more clearly from this example, is that when I say "cast a ray from the point of intersection to the light" I mean cast a shadow ray. This shadow ray helps to determine if the shaded point is in the shadow of that light. If not, then you can add the contribution of that light to the shaded point.

Now you can also see that the direction of this shadow ray is computed as the position of the point on the pen (the intersection point or shaded point) to the light position (we assume a point light source here for simplicity). So this is a straight line between the point on the pen and the light.

But this where things are wrong. This shoudn't be a straight line, since when the "shadow ray" (which is also a light ray but we go from the point on the pen to the light rather than from the light to the pen, but in essence they are the same thing) exits the water volume it should be bent due to refraction (and part of it should be reflected to). You can compute this new direction if you wish, but then if the shadow ray is bent when it leaves the volume of water then it won't travel in the direction of the light anymore (as showed in the pic below), so you can't use it as a shadow or light ray any longer.

So this algorithm that consists of looping over the lights and casting shadow ray can't be used if you really want to account for light ray bending due to reflection and refraction. Though as I said few people care about that because in reality, ignoring this fact, doesn't make much of a visual difference in 99% of cases. And as I pointed out, the solution to get it right if really needed is something like photon mapping where photon of lights are cast from the light into the scene into a pre-pass, which allows to follow the paths of light rays into the scene as they are being bent by refractive/reflective surfaces. But photon mapping is a rather complex algorithm to implement.


Screen space coordinates to Eye space conversion

Object Space → World Space → Eye Space → Clip Space → Normalized Device Space → Window Space

You get from Object Space to World Space by multiplying by the "World
Matrix" (This may actually be multiple matrices depending
on the parent-child object hierarchy of the scene.  Compose them all together to get a single "World Matrix" for a given object.)
You get from World Space to Eye Space by multiplying by the "View Matrix".
You get from Eye Space to Clip Space by multiplying by the "Projection Matrix".
You get from Clip Space to Normalized Device Space by performing the "perspective divide" (divide by w in clip space coordinates).  This is a non-linear transformation requiring a divide, and is not represented by a matrix.  Information can be destroyed in this step.  If w is divided by w you just get 1, and there's no way to recover the original w.

You get from Normalized Device Space to Window Space by performing the Viewport and Depth Range transform. (Simple linear relationship)

The World Matrix, View Matrix, and Projection Matrix are all just 4x4 matrices and can be inverted to get the reverse transformation.

So, to go backwards...

You get from Window Space to Normalized Device Space by performing the inverse of the Viewport and Depth Range transform.
You can get from Normalized Device Space to Clip Space by performing the inverse of the perspective divide (multiply by w).  This is a non-linear transformation requiring a multiply, and is not represented by a matrix.  In general this operation doesn't magically recover 3D data from 2D data.  You must know something about the depth externally to be able to recover the 3D position of a point on the 2D screen.
You can get from Clip Space to Eye Space by multiplying by the inverse of the "Projection Matrix"
You can get from Eye Space to World Space by multiplying by the inverse of the "View Matrix".
You can get from World Space to Object Space by multiplying by the inverse of the "World Matrix".

May I suggest reading this explanation of the OpenGL Transformation.


How to get assembly code from GLSL Shader

In unextended OpenGL that's indeed a little difficult, if possible at all. There is really no actual notion of "GLSL assembly code", or not one that's remotely vendor-independent, let alone a way to retrieve that.

There are a few extensions that might help you get that information, though, albeit in a not entirely intuitive or even platform-dependent way (but well, you're talking about micro-optimizations, so we're firmly in platform-dependent territory anyway).

First of all, there's GL_ARB_get_program_binary (core since 4.1), which let's you retrieve a platform-dependent binary version of the program (after linking). At least on NVidia hardware (maybe others, too, but I wouldn't know) this binary blob also contains ASCII assembly listings of all the contained shaders in NVidia-extended ARB syntax (continuously refined in the GL_NV_gpu_programX extensions). Those listings are quite easy to find in the binary when opening it with a normal text editor and the syntax is quite straight forward if you know your way around shaders and general assembly programming principles. As said, other vendors might offer similar listings inside their binaries.

Then there's also external programs that might compile your shader into some intermediate form, the most platform-independent of those being SPIR-V, for which Khronos-approved external compilers exist. Being an intermediate language this is somewhat of a "high-level assembly language" with an "assembly-like" structure but with many high-level abstractions that builds a bit of a compromise between the easy interpretability of assembly and the easy programability and platform-independence of a high-level language (think of Java bytecode).

With the GL_ARB_spirv extension (core since 4.6) you can directly use (externally precompiled) SPIR-V programs in OpenGL in place of GLSL shaders. I don't know if there is currently a way to let your driver compile GLSL into SPIR-V and then retrieve that. Otherwise some external SPIR-V module created from a reference compiler would not necessarily be equal to the machine code that your driver's actual GLSL compiler produces. It might also be that such an intermediate representation might already be too high-level for analysing micro-optimizations anyway.

All in all, and especially if you just want to use it for debugging/analyzing purposes on a dedicated machine, trying to look for a whatever-natured assembly listing in your driver's binary blob retrieved with glGetProgramBinary might be the easiest and most "close-to-the metal" way to get some insight into the compiled product. I use it occasionally just out of interest.

But of course no answer to this question would be complete without a warning that you shouldn't make too much out of it either. As on the CPU, generally trusting your compiler is a good default approach and you shouldn't let micro-optimizations like this guide your practice too much, let alone rely on any of that guidance being remotely portable (let alone assume that the assembly listed in the binary is actually what's uses as machine code). But if you want to know what your driver might make of it, it can be quite interesting and helpful, especially when deciding between two equally meaningful ways to write some function that might have unexpected differences in the "assembly".


How to calculate matching roundness of two offset rectangles?

What you (probably) want to achieve is something like this:

When having a closer look at one of the corners and add a few lines, we see this:

The black lines indicate that the center points of the circles along the borders of the red and blue boxes is the same. If the outer radius of the red box, for example, is $50px$, and the distance between the outer borders of the red and blue boxes is $25px$, this leads to a border radius of $50px - 25px = 25px$ for the blue box. The distance, obviously, is dependent on the stroke width of your red stroke.

Note: When creating the above images, I simply drew two filled boxes with rounded corners. So first a filled red box and then a smaller filled blue box above it. This also avoids some pixels being left white due to rounding errors (no pun intended).


My issue was in the call to glTextureStorage2D(texture, 1, GL_RGB, width, height);

From the man page

  internalformat must be one of the sized internal formats given in Table 1 below, one of the sized depth-component formats GL_DEPTH_COMPONENT32F, GL_DEPTH_COMPONENT24, or GL_DEPTH_COMPONENT16, one of the combined depth-stencil formats, GL_DEPTH32F_STENCIL8, or GL_DEPTH24_STENCIL8, or the stencil-only format, GL_STENCIL_INDEX8

Changing GL_RGB to GL_RGB8 fixed my issue.


Indeed there are/were 'just in time' renderers.

For example, Dreamcast (PowerVR CLX2) had a mode where the 'frame buffer' only had to be a few rows of tiles in size. The system would be rendering into one row while the display was reading out of some subset of the others. It would then cycle through the buffers.

In order to do this, however, you do need a display-list rendering system. A standard 'immediate mode' renderer is not appropriate.

Of course, scene complexity needs to be considered carefully for this sort of scheme to work. In particular, if one tiny section of the screen contains particularly complex geometry, the renderer might not finish it in time before the display needs to read it. Allocating more memory for the partial screen buffer allows for more variability in scene complexity.

I also have a feeling/seem to recall some old (e.g 60s/70s?) flight simulators might have rendered 'on the fly' because memory (for a frame buffer) was once fiendishly expensive.

Alternatively, the PowerVR PCX1/PXC2 ports of Tomb Raider had an optional single buffered rendering mode. (This was primarily because one of the guys doing the port (cough, me, cough) only had a 2MB 2D graphics card which could only allocate a single 1024x768 framebuffer).  Again, since the system was tile-based, it didn't suffer from the nasty artefacts you'd get with single buffering on a traditional Z-buffer system.

UPDATE:
FWIW there's another related question on "Racing the Beam" here


What are the benefits of downscaling on modern games?

What you are describing is supersampling, a technique for combating aliasing (a.k.a. jaggies), which is a term from signal theory that also applies to graphics. It's a spatial domain approach. Basically it means that the sampling resolution (render target pixels in this case) makes signals with different frequencies be reconstructed as the same frequency. Or, in other words, the scene is sampled too infrequently for sub-pixel detail to be properly represented.

It is especially prevalent in modern games since the recent advances in GPU computing power have made that quality/performance tradeoff (that you deem unjustified, but unless we're talking VR, everything above 60 FPS is just buttery-smooth and undistinguishable) even possible, and because the (fairly) recently introduced physically-based shading models are much more susceptible to aliasing from sources other than just undersampled geometry: a group of artifacts collectively called shading aliasing, with specular aliasing (a.k.a. fireflies) being especially visible. Epic Games' Brian Karis explains some of that in the opening of his presentation on temporal AA in UE4 (temporal AA is just another approach to anti-aliasing, one that distributes samples in time instead of space).

To sum up: it's the brute force approach to fixing aliasing. An anti-aliasing technique that simply captures more samples of the scene, and thus captures more sub-pixel detail that would otherwise be lost.

It is also worth mentioning that while it appears superficially similar to Multi-Sampled Anti-Aliasing, it differs in the sampling pattern (SSAA simply increases the pixel grid density, whereas MSAA has a sort of sheared grid pattern) and the sampled signal resolution (SSAA is a full colour capture, MSAA cheats slightly by sampling triangle coverage more finely than colour).


Why the BRDF of specular reflection is infinite in the reflection direction?

In the BRDF curve, the total area under the curve is the albedo: what fraction of incident light is reflected in total (as opposed to being absorbed or transmitted). For a perfect reflector, the area under the curve sums to 1, because it reflects all of the incident light. This area limits how high the curve can be. For example, a perfect diffuse reflector has a BRDF that's a flat line at $\frac{1}{2\pi}$, because it reflects the light evenly across the whole hemisphere, the solid angle of a hemisphere is $2\pi$, and the area has to be 1.

In a perfect specular reflector, the area under the curve is still 1: it still reflects all of the incident light. But this time, it's reflecting all of that light in one direction. The BRDF curve is a single spike. If you draw a graph with an infinitely thin spike, and the area under the graph has to be 1, then the spike must also be infinitely high. It's the only way to make the area correct.

Obviously a real material can't be a perfect specular reflector, so it will have two differences. First, the area under the curve will be less than 1,  because some light is absorbed. This doesn't matter much to us. More importantly, the reflection peak can't be infinitely thin: the reflection will be blurred ever so slightly. This means the peak won't be infinitely high: the wider the peak, the less high it has to be to still have an area of 1.

The remark about "$\delta$-functions" is because of a function (properly speaking, it's not a function but a distribution) called the Dirac $\delta$ (delta). This distribution is an infinitely thin, infinitely tall spike at $x = 0$, whose area is 1, just like our perfect specular BRDF.


Calculating the gradient of a triangular mesh

A function $f$ on your triangulated surface is an assignment of real numbers on the vertices of the triangulation and therefore $f \in \mathbb{R}^{|\mathcal{V}|}$. According to the continuous theory, the gradient should be a linear (differential) operator that takes a function and spits out a vector field. Hence, in the discrete theory, the gradient should be a linear operator that takes a function $f \in \mathbb{R}^{|\mathcal{V}|}$ and spit out some kind of a discrete version of a vector field. In your case, the vector field is given as an assignment of a vector to each triangle of the triangulation. That is why the discrete gradient operator should be a linear operator, i.e. a linear transformation
$$\text{grad} : \mathbb{R}^{\mathcal{|V|}} \to \mathbb{R}^{3|\mathcal{F}|}$$
which is represented by a matrix $\text{grad} \in \mathbb{R}^{3|\mathcal{F}| \times |\mathcal{V}|}$. In the continuous theory, the gradient is determined only by local data, i.e. the behaviou of the function near the point at which the gradient is being calculated and by the geometry near that point. Therefore, in the discrete theory, the gradient is determined only by the geometry of one triangle and the values of the function at the triangle's vertices. Let the triangle $j \in  \mathcal{F}$ have vertices denoted (according to the instructions' convention) by $1, 2, 3$. The labels are set up so that the vertices $1,2,3$ define the labels of the opposing edges $\vec{e}_{j_1}, \vec{e}_{j_2}, \vec{e}_{j_3}$, i.e. for the vertex  $i$ the edge opposite to the vertex $i$ is denoted by $\vec{e}_{j_i}$, where the latter vector is oriented so that the vectors form a closed chain of vectors summing to zero in a counterclockwise direction.
Hence, the gradient can be defined as
$$(\text{grad}\,f)(j) = \frac{1}{2\text{Area}(j)}\, \sum_{i=1}^{3} f_i\, J\vec{e}_{j_i}$$ Or to write it in coordinates,
$$\begin{bmatrix} (\text{grad}\,f)(j)_1\\ (\text{grad}\,f)(j)_2\\(\text{grad}\,f)(j)_3\end{bmatrix} = \frac{1}{2\text{Area}(j)}\left( f_1\, \begin{bmatrix} (Je_{j_1})_1\\ (Je_{j_1})_2\\(Je_{j_1})_3\end{bmatrix} +
f_2\, \begin{bmatrix} (Je_{j_2})_1\\ (Je_{j_2})_2\\(Je_{j_2})_3\end{bmatrix} + f_3\, \begin{bmatrix} (Je_{j_3})_1\\ (Je_{j_3})_2\\(Je_{j_3})_3\end{bmatrix}\right)$$ and in matrix form
$$\begin{bmatrix} (\text{grad}\,f)(j)_1\\ (\text{grad}\,f)(j)_2\\(\text{grad}\,f)(j)_3\end{bmatrix} = \frac{1}{2\text{Area}(j)}\, \begin{bmatrix} (Je_{j_1})_1 &  (Je_{j_2})_1 & (Je_{j_3})_1 \\ (Je_{j_1})_2 & (Je_{j_2})_2 & (Je_{j_3})_2\\(Je_{j_1})_3 & (Je_{j_2})_3 & (Je_{j_3})_3 \end{bmatrix} \begin{bmatrix} f_1\\ f_2\\f_3\end{bmatrix} $$
So, for the matrix $E \in \mathbb{R}^{3|\mathcal{F}|\times|\mathcal{V}|}$ is defined as follows. Let $j$ be a triangle and $i$ one of its vertices. Take the edge $e_{j_i}$ of triangle $j$ opposite to the vertex $i$, and consider the vector $\vec{e}_{j_i}$ along the edge, starting from one vertex, different from $i$, to the other vertex, different from $i$, oriented so that the the vector follows counterclockwise orientation of the boundary of the triangle $j$. Then rotate the vector $\vec{e}_{j_i}$ by $90^{\circ}$ so that it points inside the triangle $j$. Denote this rotated vector by  $$J\vec{e}_{j_i} = \begin{bmatrix} (Je_{j_i})_1\\ (Je_{j_i})_2\\(Je_{j_i})_3\end{bmatrix}$$ Then fill up the following entries of the matrix $E \in \mathbb{R}^{3|\mathcal{F}|\times|\mathcal{V}|}$:
\begin{align}
&E(j,i) = (Je_{j_i})_1\\
&E(j + |\mathcal{F}|, i) = (Je_{j_i})_2\\
&E(j + 2|\mathcal{F}|, i) = (Je_{j_i})_3
\end{align}
You can rotate the vector $\vec{e}_{j_i}$ by using the cross-product operation twice. Take two of the edge-vectors of the triangle $j$, say $\vec{e}_{j_1}$ and $\vec{e}_{j_2}$, and calculate the unit normal vector to the triangle $j$:
$$\vec{n}_j = \frac{1}{\|\vec{e}_{j_1} \times \vec{e}_{j_2}\|}\,\Big(\vec{e}_{j_1} \times \vec{e}_{j_2}\Big)$$ It doesn't matter whcih two out of the three edge-vector you would choose, the result is always the same unit normal $\vec{n}_j$. Then
$$J\vec{e}_{j_i} = \vec{n}_j \times \vec{e}_{j_i}$$