Feeding acceleration structure to the GPU then raytracing it in a compute shader?

If using the new real-time ray tracing functionality is an option for you, it would be worth a try! It might be that getting it to run is less work than porting your CPU data structure into a GPU-compatible format. Furthermore, using the hardware functionality would result in much better performance since BVH traversal (and triangle intersection, if that is relevant to you) are accelerated in hardware.

Official real-time ray tracing support has been announced by Khronos (in the form of an KHR extension) lately, see Ray Tracing In Vulkan. With the Vulkan API version 1.2, the option for "Ray Queries" has been added. It could be relevant for your use case since it allows to trace rays against the BVH from arbitrary shaders (i.e., e.g., also from compute shaders or fragment shaders).

Using VK_KHR_ray_tracing (or NVIDIA's precursor VK_NV_ray_tracing), you would have to update the acceleration structure on the GPU every time after you have updated your CPU data structure. The API supports updating acceleration structures every frame. In NVIDIA's extension this is done via the VkBool32 update flag in vkCmdBuildAccelerationStructureNV. In the KHR-extension it is handled a bit differently, but in a very similar fashion: The update-flags are contained in the VkAccelerationStructureBuildGeometryInfoKHR structures.


How to unproject cursor with orthographic projection

I am still not 100% sure, if I understood your question, because of this sentence:

Now however with orthographic projection the far and near plane are of the same size, so we can't calculate the direction of the cursor this way. The direction is going to be equal the world's z axis. (I haven't had much sleep so I hope this makes sense).

If I misunderstood you, let me know in the comments and I will adjust or remove my answer.
However, if I understood your intention correctly and you want to cast a ray through your frustum (for example to pick objects), then your statement is wrong. The direction will be equal to the view spaces' negative z-direction, not the world spaces'. So all you need to do is to transform your direction vector or the near and far plane points to world space. To proof, that this works, I have implemented everything in a Python script, that you will find at the end of this answer. If you have a Python interpreter with MatPlotLib and NumPy installed, you can modify the setup parameters and experiment a little bit yourself.
So let's have a look at the relevant implementation. First, we calculate the mouse position in clip space and the corresponding 2 points on the near and far plane.
mouse_pos_x_clip = mouse_pos_x_screen / screen_width * 2 - 1
mouse_pos_y_clip = mouse_pos_y_screen / screen_height * 2 - 1

mouse_pos_near_clip = np.array([mouse_pos_x_clip, mouse_pos_y_clip, -1, 1], dtype=float)
mouse_pos_far_clip = np.array([mouse_pos_x_clip, mouse_pos_y_clip, 1, 1], dtype=float)

Now we get the involved matrices. My notation here is as follows: I use two characters after M_ that are abbreviations of the involved spaces. The first character is the source and the second the target space. The characters are c for clip space, v for view space, and w for world space. So M_vc is the view space to clip space transformation a.k.a the projection matrix.
M_wv = get_world_to_view_matrix(camera_pitch, camera_yaw, camera_position)
if perspective:
    M_vc = get_perspective_mat(field_of_view, z_near_plane, z_far_plane, aspect_ratio)
else:
    M_vc = get_orthogonal_mat(frustum_width, frustum_height, z_near_plane, z_far_plane)

M_vw = np.linalg.inv(M_wv)
M_cv = np.linalg.inv(M_vc)

Now I simply use the correct transformation matrices to transform from clip to world space. Note that the perspective projection needs a division by w after the transformation to view space. This is not necessary for the orthographic projection, but performing it does not affect the result.
mouse_pos_near_view = np.matmul(M_cv, mouse_pos_near_clip)
mouse_pos_far_view = np.matmul(M_cv, mouse_pos_far_clip)

if perspective:
    mouse_pos_near_view= mouse_pos_near_view / mouse_pos_near_view[3]
    mouse_pos_far_view= mouse_pos_far_view / mouse_pos_far_view[3]

mouse_pos_near_world = np.matmul(M_vw, mouse_pos_near_view)
mouse_pos_far_world = np.matmul(M_vw, mouse_pos_far_view)

This is, as far as I can see, identical to your first code section. Now let's have a look a the result for perspective and orthographic projection with the following setup parameters:
screen_height = 1080
screen_width = 1980

mouse_pos_x_screen = 500
mouse_pos_y_screen = 300

camera_position = [3, 0, 1]
camera_yaw = 20
camera_pitch = 30

z_near_plane = 0.5
z_far_plane = 3

# only orthogonal
frustum_width = 3
frustum_height = 2

# only perspective
field_of_view = 70
aspect_ratio = screen_width / screen_height

The screen space and clip space values are identical for both projections:

The red line connects the two points on the near and far plane. The red dot is the point on the near plane, which is your "screen" in a 3d space. The green lines mark the borders of the frustum. In clip space, it is obviously just a cube. An important thing to realize is, that clip space is defined in a left-handed coordinate system while the other coordinate systems are usually right-handed (have a look at the images in this link). I mention it since I had some problems with the plots until I realized that.
Now for perspective projection, I get the following plots:

The blue dot is the camera position.
If I just exchange the perspective matrix with an orthographic projection matrix, the results look like this:

As you can see, the approach you used in your first code section works independently of the chosen projection. I don't know why you thought it wouldn't. An assumption of mine is, that you might have made a small mistake during the implementation of the orthographic projection matrix. For example, if you accidentally flipped rows and columns (transposed) of the orthographic projection matrix, you get total crap like this:

I know this looks like a wrong implementation of the perspective projection, but this is what I get when I transpose the orthographic projection matrix before the multiplication.
So make sure you use the correct orthographic projection matrix (source):
$$
\begin{bmatrix}
\frac{2}{w}&0&0&0\\
0&\frac{2}{h}&0&0\\
0&0&\frac{-2}{f-n}&-\frac{f+n}{f-n}\\
0&0&0&1
\end{bmatrix}
$$
Here $w$ is the frustum width, $h$ the frustum height, $f$ the far plane z-value and $n$ the near plane z-value. This is the representation if you use column vectors and left multiplied matrices. For row vectors and right multiplied matrices, you need to transpose it.
Your second approach:
Vector4 cScreen = Vector4(cursorNormX, cursorNormY, 0, 0);
Vector4 cView = Inverse(projection)*cScreen;
cView = Vector4(cView.x, cView.y, 0, 0);
Vector4 cWorld = Inverse(View) * cView;
cWorld = Vector4(cWorld.x, cWorld.y, 0, 0);

Vector3 cursorPos = cWorld.xyz();

has multiple issues and all are related to the z- and w-components of your vectors. Basically, you need to do the same transformations as in your first approach. So use Vector4 cScreen = Vector4(cursorNormX, cursorNormY, -1, 1); as initial vector.
One problem of the line cView = Vector4(cView.x, cView.y, 0, 0); is, that your z-component should be identical to your near plane value and not zero. You might get away with that since it would just shift your point a little bit in the camera viewing direction in world space, but more problematic is that you set w to 0. This makes it impossible to apply any translation to the vector by $4 \times 4$ matrix multiplication. So when you transform to world space, you will always end up with a point that treats the camera to be located at the coordinate system origin, regardless of its true position. So you need to set the w-component to 1. However, if the previous lines are correct, you should automatically get the correct z- and w-values which makes this line obsolete.
Lastly, the line cWorld = Vector4(cWorld.x, cWorld.y, 0, 0); doesn't make much sense either to me. Your camera is somewhere in 3d world space. Why do you remove the z-component you previously calculated? With this, you move the point into the XY-plane for no reason. Just remove this line.
To get the cameras viewing direction without using the far plane point, just multiply the vector [0, 0, -1, 0] with the view-to-world matrix (M_vw). In this case, the w-component really has to be 0, since you do not want to apply translations to it because it is a direction vector and no point. The z-component needs to be -1 because the camera looks into the negative direction by definition. Notice, that the transformed vector is usually not of unit length anymore. So you might want to normalize it.
Additional Note
For an orthographic projection, there is no need to calculate the inverse projection matrix. You can simply calculate the x and y values directly by something like this (untested pseudo-code):
x_view = (x_screen / screen_width - 0.5) * frustum_width
y_view = (y_screen / screen_height - 0.5) * frustum_height

Then you get the screen space point in view space by setting (untested pseudo-code):
point_view_near = [x_view, y_view, -z_near, 1]

Watch out to use the negative near plane z-value! For the far plane, you can do the same.
Full Python script
import numpy as np
from mpl_toolkits.mplot3d import Axes3D  # noqa: F401 unused import
import matplotlib.pyplot as plt

# setup --------------------------------------------------------------------------------

screen_height = 1080
screen_width = 1980

mouse_pos_x_screen = 500
mouse_pos_y_screen = 300

camera_position = [3, 0, 1]
camera_yaw = 20
camera_pitch = 30

# ----------------
# projection setup
# ----------------
perspective = False # set 'False' for orthogonal and 'True' for perspective projection

z_near_plane = 0.5
z_far_plane = 3

# only orthogonal
frustum_width = 3
frustum_height = 2

# only perspective
field_of_view = 70
aspect_ratio = screen_width / screen_height

# functions ----------------------------------------------------------------------------

def render_frustum(points, camera_pos, ax, right_handed=True):
    line_indices = [
        [0, 1],
        [0, 2],
        [0, 4],
        [1, 3],
        [1, 5],
        [2, 3],
        [2, 6],
        [3, 7],
        [4, 5],
        [4, 6],
        [5, 7],
        [6, 7],
    ]
    for idx_pair in line_indices:
        line = np.transpose([points[idx_pair[0]], points[idx_pair[1]]])
        ax.plot(line[2], line[0], line[1], "g")
    if right_handed:
        ax.set_xlim([-5, 5])
    else:
        ax.set_xlim([5, -5])
    ax.set_ylim([-5, 5])
    ax.set_zlim([-5, 5])
    ax.set_xlabel("z")
    ax.set_ylabel("x")
    ax.set_zlabel("y")
    ax.plot([-5, 5], [0, 0], [0, 0], "k")
    ax.plot([0, 0], [-5, 5], [0, 0], "k")
    ax.plot([0, 0], [0, 0], [-5, 5], "k")
    if camera_pos is not None:
        ax.scatter(
            camera_pos[2], camera_pos[0], camera_pos[1], marker="o", color="b", s=30
        )

def render_ray(p0,p1,ax):
    ax.plot([p0[2], p1[2]], [p0[0], p1[0]], [p0[1], p1[1]], color="r")
    ax.scatter(p0[2], p0[0], p0[1], marker="o", color="r")

def get_perspective_mat(fov_deg, z_near, z_far, aspect_ratio):
    fov_rad = fov_deg * np.pi / 180
    f = 1 / np.tan(fov_rad / 2)

    return np.array(
        [
            [f / aspect_ratio, 0, 0, 0],
            [0, f, 0, 0],
            [
                0,
                0,
                (z_far + z_near) / (z_near - z_far),
                2 * z_far * z_near / (z_near - z_far),
            ],
            [0, 0, -1, 0],
        ]
    )

def get_orthogonal_mat(width, height, z_near, z_far):
    r = width / 2
    t = height / 2

    return np.array(
        [
            [1 / r, 0, 0, 0],
            [0, 1 / t, 0, 0],
            [
                0,
                0,
                -2 / (z_far - z_near),
                -(z_far + z_near) / (z_far - z_near),
            ],
            [0, 0, 0, 1],
        ]
    )

def get_rotation_mat_x(angle_rad):
    s = np.sin(angle_rad)
    c = np.cos(angle_rad)
    return np.array(
        [[1, 0, 0, 0], [0, c, -s, 0], [0, s, c, 0], [0, 0, 0, 1]], dtype=float
    )

def get_rotation_mat_y(angle_rad):
    s = np.sin(angle_rad)
    c = np.cos(angle_rad)
    return np.array(
        [[c, 0, s, 0], [0, 1, 0, 0], [-s, 0, c, 0], [0, 0, 0, 1]], dtype=float
    )

def get_translation_mat(position):
    return np.array(
        [
            [1, 0, 0, position[0]],
            [0, 1, 0, position[1]],
            [0, 0, 1, position[2]],
            [0, 0, 0, 1],
        ],
        dtype=float,
    )

def get_world_to_view_matrix(pitch_deg, yaw_deg, position):
    pitch_rad = np.pi / 180 * pitch_deg
    yaw_rad = np.pi / 180 * yaw_deg

    orientation_mat = np.matmul(
        get_rotation_mat_x(-pitch_rad), get_rotation_mat_y(-yaw_rad)
    )
    translation_mat = get_translation_mat(-1 * np.array(position, dtype=float))
    return np.matmul(orientation_mat, translation_mat)

# script -------------------------------------------------------------------------------

mouse_pos_x_clip = mouse_pos_x_screen / screen_width * 2 - 1
mouse_pos_y_clip = mouse_pos_y_screen / screen_height * 2 - 1

mouse_pos_near_clip = np.array([mouse_pos_x_clip, mouse_pos_y_clip, -1, 1], dtype=float)
mouse_pos_far_clip = np.array([mouse_pos_x_clip, mouse_pos_y_clip, 1, 1], dtype=float)

M_wv = get_world_to_view_matrix(camera_pitch, camera_yaw, camera_position)
if perspective:
    M_vc = get_perspective_mat(field_of_view, z_near_plane, z_far_plane, aspect_ratio)
else:
    M_vc = get_orthogonal_mat(frustum_width, frustum_height, z_near_plane, z_far_plane)

M_vw = np.linalg.inv(M_wv)
M_cv = np.linalg.inv(M_vc)

mouse_pos_near_view = np.matmul(M_cv,mouse_pos_near_clip)
mouse_pos_far_view = np.matmul(M_cv,mouse_pos_far_clip)

if perspective:
    mouse_pos_near_view= mouse_pos_near_view / mouse_pos_near_view[3]
    mouse_pos_far_view= mouse_pos_far_view / mouse_pos_far_view[3]

mouse_pos_near_world = np.matmul(M_vw, mouse_pos_near_view)
mouse_pos_far_world = np.matmul(M_vw, mouse_pos_far_view)

# calculate view frustum ---------------------------------------------------------------

points_clip = np.array(
    [
        [-1, -1, -1, 1],
        [ 1, -1, -1, 1],
        [-1,  1, -1, 1],
        [ 1,  1, -1, 1],
        [-1, -1,  1, 1],
        [ 1, -1,  1, 1],
        [-1,  1,  1, 1],
        [ 1,  1,  1, 1],
    ],
    dtype=float,
)

points_view = []
points_world = []
for i in range(8):
    points_view.append(np.matmul(M_cv, points_clip[i]))
    points_view[i] = points_view[i] / points_view[i][3]
    points_world.append(np.matmul(M_vw, points_view[i]))

# plot everything ----------------------------------------------------------------------

plt.figure()
plt.plot(mouse_pos_x_screen,mouse_pos_y_screen, marker="o", color="r")
plt.xlim([0, screen_width])
plt.ylim([0, screen_height])
plt.xlabel("x")
plt.ylabel("y")
plt.title("screen space")

plt.figure()
ax_clip_space = plt.gca(projection="3d")
render_ray(mouse_pos_near_clip, mouse_pos_far_clip, ax_clip_space)
render_frustum(points=points_clip, camera_pos=None, ax=ax_clip_space, right_handed=False)
ax_clip_space.set_title("clip space")

plt.figure()
ax_view = plt.gca(projection="3d")
render_ray(mouse_pos_near_view, mouse_pos_far_view, ax_view)
render_frustum(points=points_view, camera_pos=[0, 0, 0], ax=ax_view)
ax_view.set_title("view space")

plt.figure()
ax_world = plt.gca(projection="3d")
render_ray(mouse_pos_near_world, mouse_pos_far_world, ax_world)
render_frustum(points=points_world, camera_pos=camera_position, ax=ax_world)
ax_world.set_title("world space")

plt.show()


ScratchAPixel : Understanding how to use Perlin's permutations hash table

First, (1/2) I didn't understand how the hash table is actually used. The workflow I think I understood is the following :

  For an input (coordinate), we find the both nearest pre-defined coordinates (i.e. : those present as indices in A). It's a multiple of
  255.
  The lower one is used to lookup in B
  We get a number between [0;255] that we use to lookup in A
  We get a number between [0;1] : the noise.

  But where is the interpolation between the both pre-defined
  coordinates' noise values ?

There are 4 of most of what you wrote here, then the interpolation. For example, if your 2D input is (2.7, 6.4):

You would truncate that to (2,6).
Then perform the permutation and noise lookup steps for (2,6), (3,6), (2,7), (3,7).
And finally perform a 2D interpolation (lerp) of the 4 lookups with the remainder of your input: (0.7, 0.4).

As for the second part of your question, the table is just a creative way to get a "(int, int) -> int" hashing function to go from input coordinates to an index to pick a value in the random data table. You could use any generic hash function instead. In fact, the method used is fairly obsolete for CPUs as it does not lend itself to SIMD implementation and has horrible memory dependencies. It might still have value for GPUs as they have different behavior. But it will have a visible period of 256 which can cause problems in images. So it is really a poor hash function but for the CPUs back then, it was fast and used little memory which is likely why it was picked.


The two are different because they do different things. Clipping and culling aren't the same thing at all. If an edge crosses a clip plane, the primitive has to be clipped to be within that plane. This means generating new primitives from the old one.
If any vertex of a primitive is outside of a cull plane, then the entire primitive is culled.
You use the cull planes when you want the VS to be able to cull out primitives based on some calculation. You use clip primitives when you want primitive rendering to be bounded to some non-clip-space region (though line and point width can violate those borders).


Perspective projection near plane not in world coordinates

The perspective projection matrix transforms from camera space to clip space. The near-plane value is the position of the near plane in camera space. So in case that your world space to camera space transformation does not use any scaling, only rotations and translations, then if you place the camera $z$ units away from the billboard surface and look directly onto it, a near plane value of $z$ will cause the billboard plane to lie on the near plane. Due to rounding errors, it might get clipped away, so your idea of setting it to $z-offset$ should solve the problem.
Why it doesn't work is hard to tell without any code. Maybe your implementation of the perspective matrix is wrong. Maybe you have a hidden scaling when you transform between world and camera space. To debug this you can try to multiply $z$ with values between 0 and 1 until you get the desired result. The resulting factor might be a clue. For example if you get $0.5$ we know that there is either a 2 missing or too much in your transformation matrices. Also, make sure that your far plane value is larger than your near plane value (even though I think from your description that this isn't the problem - but who knows)
In case you are using OpenGL, here is a link to how the perspective projection matrix has to look like. Maybe you have to transpose the matrix. That depends on your vector layout - column vectors or row vectors.
Here is a link to an answer of mine to another question. It has a Python script at the end that you can use to visualize how the different values affect the view frustum in the different spaces. Some pictures are also included in the answer. Maybe that is already enough to help you out.


What is the algorithm used for phong shading normal interpolation?

EDIT: I misunderstood your question, so please disregard this answer.

The Phong shading model assumes a surface normal $\mathbf{n}$ at the shading point $\mathbf{x}$ but it is independent of its type, so either the shading normal or the geometric normal will work. Moreover, this reflectance model requires an incoming and outgoing directions. In that sense, you can't just "interpolate the normals" to get the desired effect since you need to take into account the positions of the camera and the light, as well as the roughness parameters.

The (normalized) Phong BRDF is a sum of a diffuse and specular term:

$$
f_r(\mathbf{x},\omega_i,\omega_o) = \frac{\rho_d}{\pi} + \rho_s \frac{n+2}{2\pi}\max(0,\cos^n\alpha),
$$
where

$\rho_d$ is the diffuse albedo, i.e. the fraction of the incoming energy that is reflected diffusely,
$\rho_s$ is the specular albedo, i.e. the fraction of the incoming energy that is reflected specularly,
$n$ is the Phong exponent; higher values yield more mirror-like specular reflection, and
$\alpha$ is the angle between the perfect specular reflection direction $\omega_r$ and the outgoing direction $\omega_o$.

You can compute $\omega_r$ from the surface normal $\mathbf{n}$ and the lighting direction $\omega_i$, which is just

$$
\omega_r = 2\mathbf{n}(\mathbf{n}\cdot\omega_i) - \omega_i.
$$

Typically you pick the albedos such that $\rho_d + \rho_s \leqslant 1$ to conserve energy.


Your first image employs an orthographics projection, while the second uses a perspective projection. You can look up the perspective matrix derivation in: http://www.songho.ca/opengl/gl_projectionmatrix.html Note that there's a division at the end to account for perspective projection being a non-linear transformation.


Two questions on Gamma correction and Tonemapping

the colors are distorted and much brighter corresponding to (186,186,186) color value. If I don't apply gamma encoding then the image file displays the expected color that is (127,127,127).

No, the gamma-correct color is not distorted; your expectations are. If your lighting, meshes, textures, and other scene elements were not built to be properly linear and gamma correct to the screen, then every decision you have made has been through a distorted lens. If you don't like the effect you're getting with a properly gamma-correct lighting system, that's the fault of all of those decisions you've been making, not with gamma correction.

Half of the full lighting intensity in a linear colorspace does correspond to greater than half in a gamma-correct colorspace. And those gamma-correct values are what you ought to be passing to the screen. If that's not what you're shoving into the display, then you're doing lighting wrong.


Is there a way to resize an image without creating new colors in the output image?

A nearest-neighbor resize will do what you ask but I am not certain it will give you much smaller tiles.

It is quite likely that the tiles made from the original image are smaller simply because they have less high frequency content and thus compress far better.

You can look at it the other way around: there are more terrain features in the tiles of the resized image so it takes more space to store them. That's a very gross simplification of the truth, of course. And it assumes a reasonably efficient compression algorithm.


Does Vulkan have the equivalent of OpenGL's separable shader program and pipeline?

Does Vulkan have the equivalent of the OpenGL separable shader objects, or are there any plans for it? They seem just so much more flexible.

It would be "so much more flexible". But Vulkan is a low-level rendering API. It's job is to give you access to the hardware as close to the metal as possible; user-convenience isn't even in the top-5 goals for the API.
Putting all of the context state into one gigantic, immutable (more or less) object gives the implementation all the power it needs to make your code fast. Forcing implementations to deal with stitching together shaders would make implementations more complicated (and another goal of Vulkan is to be relatively easy to implement).
The best way to solve this (where possible) is to just have fewer shader variations. Use uniforms/push-constants to detect what variations you want to use, then have conditionals based on them. Obviously this can't be done in every situation, but the uber-shader approach is the way to handle it to the degree possible/reasonable.


What are Spherical Harmonics & Light Probes?

Spherical harmonics

If you know what a Fourier transform is, you already almost know what spherical harmonics are: they're just a Fourier transform but on a spherical instead of a linear basis. That is, while a Fourier transform is a different way of representing a function $f(x)$, spherical harmonics are the analogous thing for polar functions $f(\theta, \phi)$.

If you don't know what a Fourier transform is, you kinda need to know before you can understand spherical harmonics. The Fourier transform lets you represent a signal as a series of sine and cosine waves, each with twice the frequency of the last. That is, you can represent the signal as its average, plus a sine wave whose wavelength is the same as the length of the signal, plus a sine wave twice that wavelength, and so on. Because the Fourier transform fixes you to these particular wavelengths, you only need to record the amplitude of each one.

We commonly use Fourier transforms to represent images, which are just 2D digital signals. It's useful because you can throw away some of the sine waves (or reduce the precision with which you store their amplitude) without significantly changing what the image looks like to human eyes. OTOH, throwing away pixels changes the look of the image a lot.

In a sampled signal like an image, if you use the same number of sine waves as there were samples (pixels) in the original image, you can reconstruct the image exactly, so once you start throwing away any frequencies, you're making the image take less storage.

Spherical harmonics are just like Fourier transforms, but instead of sine waves, they use a spherical function, so instead of linear functions (such as images), they can represent functions defined on the sphere (such as environment maps).

Light probes

Just like how a standard image records all the light reaching a certain point through the image plane, a light probe records all the light reaching a certain point from all directions. They first came out of movie effects. If you want to add a computer-generated object to a real-world scene, you need to be able to light the synthetic object with the real-world lighting. To do that, you need to know what light reaches the point in the scene where the synthetic object will be. (N.B. Although I say "lighting", you're recording an image of all the light, so it can be used for reflections as well.)

Because you can't have a camera with a spherical lens that records all the light reaching a single point from all directions, you record this by taking normal photographs of a spherical mirror, and then reprojecting the images onto a sphere.

Outside of movie effects, it's more common to use light probes generated from an artificial scene. Imagine you have some expensive algorithm to compute global illumination (GI) in a scene, and you also have some smaller objects moving around in this scene (such as a game level with players in it). You can't run the whole GI algorithm every time any object moves, so you run it once with the static scene, and save light probes taken at various points in the level. Then you can get a good approximation to the GI by lighting the player with whichever light probe they're closest to.

Using them together

Generally you want to filter out sharp edges in your global illumination anyway, so you want a way to represent them that's compact and easily lets you throw away high frequencies. That's what spherical harmonics are really good at! That's why you'll hear these two terms used together a lot.

You compute light probes with your expensive GI algorithm - typically in the level-design tool, or maybe once per second (instead of once per frame) if you want to include your dynamic objects to. You store those cheaply with spherical harmonics - 16 floats is enough for pretty high quality lighting, but not reflections. Then for each dynamic object you want to light, you pick the nearest light probe (or linearly interpolate several together) and use it as a uniform or constant input to your shader. It's also commonplace to use spherical harmonics to represent ambient occlusion data, and it's very cheap to convolve that with the light probe, though there's some complexity around rotating spherical harmonic functions.


Maybe. HDR monitors are typically going to be some approximation of the DCI-P3 color space, while Adobe's wide-color rendering is the Adobe RGB color space. The two are similar but not quite the same; in particular, Adobe RGB extends a bit further into the green-cyan range than DCI-P3 does:

That said, these are only idealized color spaces, and the gamut of a given actual monitor is going to differ somewhat from these. That's why you'll see monitors advertised as having "99% Adobe RGB coverage" or "99% DCI-P3" or some such.
So, a typical HDR monitor can probably reproduce most of the Adobe space but maybe not quite all of it, but in practice it might be just about as good as an actual Adobe-supporting monitor.
That said, I don't know whether Photoshop and the drivers/OS/etc will actually cooperate to let you turn on Adobe RGB rendering on a DCI-P3 HDR monitor. If Photoshop is outputting 10-bit color in Adobe RGB space, but the monitor doesn't actually take that format, then it would need to be converted somewhere along the line to a color space and encoding that the monitor can take as input. That can certainly be done by the graphics card/drivers/OS, but I don't know whether the support for doing such conversions is actually there.


Rendering in radiometric units or photometric?

It's fine to use photometric units as an overall scale for setting light brightnesses. However, there's a technical subtlety you should be aware of. I'll quote from a blog post I wrote on the subject last year:

With RGB images, it’s important to recognize that our display devices behave more radiometrically than photometrically. A red pixel value of 255, and a green pixel value of 255, both result in about equal amounts of radiant flux (watts) being generated by a pixel on a screen—not equal amounts of luminous flux. By the same token, digital cameras capture pixel values that correspond to radiant flux, not luminous flux.
That’s why we need to use luma coefficients when converting images to grayscale, or calculating the brightness of a pixel, to get a perceptually accurate result; and it also means that rendering RGB images proceeds more naturally in radiometric units than in photometric ones.

In other words, the wavelength dependence of photometric units is different to what you might be expecting. In ordinary RGB color spaces, white is (1, 1, 1) and has a roughly flat radiometric spectrum; but in a supposed "photometric RGB", (1, 1, 1) would not be white; it would be a purpleish color, with less energy in the green range and more in the red and blue ranges. A similar problem would afflict spectral renderers trying to measure all their bins in wavelength-dependent photometric units, but even worse, as the radiance needed to generate a given luminance diverges toward either end of the visible spectrum, where the human luminous efficiency function goes to zero.
So, if you wish to use photometric units, IMO it's better to "cheat" a bit and not use true wavelength-dependent photometric units, but just use some fixed wavelength (such as 555 nm green, which is the peak of the human luminous efficiency function), or perhaps an average over the spectrum, as a reference unit, and apply that single unit to measure all wavelengths. This will get you in less trouble when importing RGB colors and spectra from other sources, and when generating them as outputs.


The block should total up to 8320 bytes.

No, it shouldn't. You use a uint[16] array. Well, std140 layout makes it clear that the array stride for any array is always rounded up to the stride of a vec4. So that array takes up 256 bytes.

You should instead use a uvec4[4] instead, which you would index from like a multidimensional array.


Practically speaking, probably not - try printing the size of each type using sizeof. Probably GLuint is 32 bits and size_t is 64 bits - so 32 bits are thrown out; but float and GLfloat are probably both 32 bit IEEE 754 numbers. "Probably" sprinkled above, because technically speaking the sizes of C ints, floats and the float format are compiler dependent, but GL type sizes (and formats) should be standard (see https://www.opengl.org/wiki/OpenGL_Type). Chances are, GLfloat is just typedef'd as a float.

You may be able to turn off this warning in your compiler - I choose to explicitly cast ("GLuint size = GLuint(vec.size())" works) so I am reminded of what's happening.


Flat shading is the simplest shading model.  Each rendered polygon has a single normal vector; shading for the entire polygon is constant across the surface of the polygon.  With a small polygon count, this gives curved surfaces a faceted look.

Phong shading is the most sophisticated of the three methods you list.  Each rendered polygon has one normal vector per vertex; shading is performed by interpolating the vectors across the surface and computing the color for each point of interest.  Interpolating the normal vectors gives a reasonable approximation to a smoothly-curved surface while using a limited number of polygons.

Gourard shading is in between the two: like Phong shading, each polygon has one normal vector per vertex, but instead of interpolating the vectors, the color of each vertex is computed and then interpolated across the surface of the polygon.

On modern hardware, you should either use flat shading (if speed is everything) or Phong shading (if quality is important).  Or you can use a programmable-pipeline shader and avoid the whole question.


How to assign/calculate triangle texture coordinates

UV or barycentric coordinates represent a point on a triangle as a weighted average of the three vertexes. You can use those coordinates to smoothly interpolate any value over the triangle, like texture coordinates or the normal vector. So think of it like a two-step process. We map the world coordinates of the intersection point to uv coordinates in "triangle space", and then use those to map to another space, in this case "texture space". Your ray hits triangle ABC and you determine the uv coordinates of the intersection point, let's say (0.3, 0.5, 0.2). Intuitively, we're sort of saying the intersection point is 30% A, 50% B, and 20% C. Each vertex would have a texture coordinate, i.e. a pixel in the texture image, associated with it (that's part of modeling / setting up your scene, not something calculated at runtime). Mix the texture coordinates at A, B and C in that ratio and you'll get the point on the texture image associated with that point on the triangle.

This scratchapixel page has some good diagrams and clear explanation.

edit: I'm sorry, I missed the part about how to calculate the coordinates in the first place. The page I linked explains it better than I could anyway.


Calculate Transformation between two point sets - but with constrained Degrees of Freedom

Every 3D rotation has an associated space that can be generated by 3 basis vectors, which can be seen as the axis of the coordenate system of this space. These vectors are also the 3 columns of the transformation matrix. In order to decrease the dof of the rotation, you should provide one or more of these basis vectors and solve for the rest of the matrix. The same principle applies to the translation vector, you should provide constants for one or more distances of the translation. Of course, if you restrict the transformation too much you will end with an unsolvable exact system or a transformation that take the first triangle too far away from the second one, if you are using a non-linear solver.


How should I generate the view position for a light at infinity when creating depth maps?

Instead of using a perspective projection you would use a orthographic projection.

Then the trick is to position the bounding box to in front of the normal camera.

An additional option is to add a skew/shear operation so it maps the light direction to the vertical and keeps the horizontal ground plane horizontal.

This avoids leaving a large dead zone underground near the light source.


Let's go through the paper you mentioned and try to list the factors which will affect your decision on whether to use GPU or CPU for simulation.
I will also add some comments about my opinion for each factor.

Computation bound

A GPU program runs in a SIMT(Single instruction, multiple threads) model.
In this model, if some threads stall due to memory access or divergence,
thread scheduler will schedule another set of threads to run
to hide the memory latency.

As we can see it requires adequate computation load to make this mechanism work.
If one implementation is memory bound, we hardly benefit from the SIMT programming model.

Fortunately, fluid simulation is a computation bound task for most of the cases.

Locality

In the paper you mentioned we can observe perfect locality:
for each grid point we only have to examine its neighbors.

The actual memory access pattern depends on the finite difference scheme you use.
But for all grid based simulation methods(Eulerian or Lattice-Boltzmann method), I think the locality always presents.

Divergence

A warp (a group of 32 threads) in a GPU program can be executed in parallel if they follow the same code path. Otherwise the program diverges and the execution of each branch will be serialized. Divergence will hurt the performance.

In particle based simulations (such as SPH method) we usually have to query the nearest N particles. A tree-like data structure will be used.

In GPU a tree data structure may introduce divergence because particles may fall into different branches.

The paper you mentioned uses grid based methods, where all threads execute along the same code path (Jacobi or Gauss-Seidel iteration), so divergence is not a problem here.

Interactivity

In the paper it uses static boundary conditions.
For this case you just need to upload the boundary conditions to GPU memory only once.

If you have dynamic boundary conditions, for example, which are specified by the user at each frame, the overhead of memory transfer between GPU and CPU should be taken into consideration.

Spatial complexity

Uniform grid based simulations cannot handle very large computation volume, because the spatial complexity is cubic.

There are adaptive grid methods, which use relatively higher resolution for the computation region we are most interested in. But that would also bring in divergence.

Numerical issues

For fluid simulation you cannot decide the number of steps solely based on the performance/frame rate.

In grid based methods, the Courant–Friedrichs–Lewy (CFL) condition puts a constraint on how large simulation step you can have.
Given a fixed time step, too large or too small steps would cause instability and blow up the simulation.

In particle based methods and some stable grid based methods (such as the semi-Lagrangian method used in the paper), we can use relatively larger simulation steps.

But we still cannot increase the step size with no bound because of the accuracy.
For example in grid based methods a large step size would introduce too much numerical viscosity.

Summary

As a summary, if your simulation has strong convergence and locality, prefer a GPU implementation; otherwise, prefer a CPU one.

In practice, I will always implement a CPU version because it is easier to debug and it will be useful as a baseline.
Then I will profile it to locate the bottleneck and decide whether to have a GPU version.


Updating the position/rotation/scale in the vertex buffer will use less memory bandwidth—that’s 9 values, where a transform matrix would be either 12 or 16—and the GPU can parallelize the work of turning the position/rotation/scale into a transform for the vertices. I’d go with 1.


I apologise for only skimming through your question/answer and the links you gave but I thought it would be a good idea to post what needs to be done to correctly do any blending and/or filtering on images with an alpha channel.
The first thing to note is that all calculations should be done in a linear colour space, but the image data in PNGs and JPEGs is most likely to be in sRGB which is non-linear. You should convert this data to a linear form before doing any calculations. For example, the average of black (0,0,0) and white (255,255,255) in sRGB is not (127,127,127) but (off the top of my head) something like (186,186,186). The filtering of images, especially those with neighbouring, complementary, saturated colours (e.g. Yellow&Blue), can look quite odd if you don't take gamma into account.
The second part of this theory dates back to Porter and Duff's seminal 1984 paper, "Compositing Digital Images", which introduces the concept of pre-multiplied alpha. Pre-multiplied alpha basically means that you multiply the RGBs by the alpha value before doing calculations.
The Porter-Duff paper states...

There are a few reasons to dismiss this proposal,
the most severe being that all compositing operations
will involve multiplying the 1 in the red channel by
the .5 in the alpha channel to compute the red contribution
of this object at this pixel. The desire to avoid this
multiplication points up a better solution, storing the
pre-multiplied value in the color component, so that
(.5,0,0,.5) will indicate a full red object half covering a
pixel.

... which suggests this was initially done just to reduce computation costs, but in Jim Blinn's article "Compositing-Theory", he demonstrates that the correct way to do any filtering or compositing operations is to first convert everything to premultipled.
Once you have linearised and then pre-multiplied your pixel data, do the filtering/blending and, finally, undo the pre-multiply and reapply the gamma function to map back to sRGB.


How to prevent moire artifacts in this light casting algorithm?

In principle you avoid using scatter (casting) behavior with GPU. They have offered random output coordinate write out since only shader model 5 as a need for extreme situations. But you should as general rule write your GPU code in a "gather" fashion.

The difference: the hardware threads are logically soft-locked to one output position in the render target. The scheduler decides to what rectangle (or cube) in the target buffers, the kicked thread group will output results.

So you should work around the designated destination, and figure out the start; instead of working from some start and computing a dynamic destination.

This way not only will you please the hardware by avoiding contention, and race conditions completely; but also you'll avoid holes.


Wrong tangential directions for cubic bezier curve

Using what Rahul kindly suggested in the comments, I realised that my tangent can be rotated in the opposite direction for two cases:

if dy/dt and dx/dt are both negative, yet dy/dx ends up being positive
if dy/dt is negative and dx/dt is positive, so the overall dy/dx is negative

This is because my angle calculation, which is done by float angle = atan(dydx) * 180/M_PI; will only calculate acute angles from the horizontal, and hence miss out the opposite cases.

Thanks a bunch!


Since you're asking in this Stack Exchange, I assume you intend to build a VR app yourself that utilizes this scenario.

While it is perfectly possible to utilize PC for rendering and stream the output to a mobile VR headset, it will come at the relatively high cost of increased latency (see e.g. here for why this is important; the tl;dr is that high latency → nausea). This is because wireless networks have intrinsically higher latency than wired ones, and you need to do the following in a loop:

Collect head tracking (gyroscope/accelerometer) data from the VR headset.
Send the head tracking data over the network to the PC.
Render a frame on the PC.
Stream the frame from the PC to the VR headset. This involves the following substeps:

Encode video and audio on the PC.
Transmit the encoded media frame to the VR headset over the network.
Decode video and audio on the VR headset.

As I mentioned previously, steps 2 and 4 introduce additional latency (probably from several to a dozen milliseconds, from stream encoding/decoding and network transmission), which – depending on your specific application and the individual tolerance of the user – may or may not be a deal breaker.


In order to dissociate the rotation from the CPU speed, you need to use some sort of timing in your function. I recommend passing in a time value of some sort. The easiest thing is to pass in the number of seconds since the app started. How you do that depends on your particular OS, which you didn't specify.

Once you have a timing constant, you can generate a triangle wave to control the rotation. This formula will oscillate between 0 and 1 over time:

double rotation = 2.0 * fabs(fmod(time, 1.0) - 0.5);

You can make it go from 0 to 10 by multiplying the result by 10:

rotation *= 10.0;

Then you plug that into your glRotated() call:

glRotated(rotation, 0, 1, 0);

If you want to make the motion a little smoother, you could replace the triangle wave with a sine wave, like this:

double rotation = 0.5 * (sin(time * 2.0 * M_PI + M_PI / 2.0) + 1.0);


Point light illuminating only half of a scene

The normals were indeed wrong. Which was crazy, because I had even checked in Blender that they were correct (and they were). Turns out, when loading the mesh, I was passing the normal as color and vice versa. Valletin's normal visualization shader was great (https://vallentin.dev/2015/02/16/visualizing-normals).

Thanks for your help!


In forward rendering what's the go to technique to smooth out rigid shadow due to low resolution shadowmap texture?

Percentage Closer Filtering with some amount of blur (ie. a minimum filter size in texture space) is the most basic method I know of. For example, see GPU Gems. My experience is that it is also quite robust.


What you are looking for is a way to get access to individual pixels in an image, in a way that you can modify those pixels with CPU code.

Once you have that you'll want to find software rendering / software rasterization methods.

A great place to start would be reading up in the various bresenham algorithms.

There is a lot of info out there on that. Here is one article for instance:
http://blog.demofox.org/2015/01/17/bresenhams-drawing-algorithms/


Trying to implement Microfacet BRDF but my result images are wrong

Problem fixed by RichieSams, trichoplax and xpicox. Thanks all of you for the answers.

I lower the roughness, change the color of material and reversed the ViewDirection then finally I start to see proper specular :).

Fixed Code:

#version 330 core

in vec3 Position;
in vec2 TexCoord0;
in vec3 Normal;
in vec3 Tangent;
out vec4 FinalColor;

uniform vec3 uCameraPosition;

#define PI 3.1415926f
#define EPSILON 10e-5f
#define saturate(value) clamp(value, 0.0f, 1.0f);

float BRDF_Lambert(float NdotL)
{
    return NdotL;
}

// Got these BRDF formulas Moving Frostbite to PBR slide by Sebastien Lagarde & Charles de Rousiers
// http://www.frostbite.com/wp-content/uploads/2014/11/course_notes_moving_frostbite_to_pbr_v2.pdf
float BRDF_D_GGX(float NdotH, float Roughness)
{
    float Roughness2 = Roughness * Roughness;
    float f = (NdotH * Roughness2 - NdotH) * NdotH + 1.0f;
    return Roughness2 / (f * f + EPSILON);
}

float BRDF_F_FresnelSchlick(float LdotH, float F0)
{
    float f = F0 + (1.0f - F0) * pow((1.0f - LdotH), 5);
    return f;
}

float BRDF_G_SmithGGXCorrelated(float NdotL, float NdotV, float Roughness)
{
    float Roughness2 = Roughness * Roughness;
    float GV = NdotL * sqrt((-NdotV * Roughness2 + NdotV) * NdotV + Roughness2);
    float GL = NdotV * sqrt((-NdotL * Roughness2 + NdotL) * NdotL + Roughness2);

    return 0.5f / (GV + GL + EPSILON);
}

float BRDF_Specular(float NdotV, float NdotL, float NdotH, float LdotH, float Roughness, float F0)
{
    float D = BRDF_D_GGX(NdotH, Roughness);
    float F = BRDF_F_FresnelSchlick(LdotH, F0);
    float G = BRDF_G_SmithGGXCorrelated(NdotL, NdotV, Roughness);
    return (D * F * G) / PI;
}

void main()
{
    vec3 normal = normalize(Normal);

    vec4 BaseColor = vec4(1.0f, 0.0f, 0.0f, 1.0f);
    vec4 SpecularColor = vec4(1.0f, 1.0f, 1.0f, 1.0f);

    vec3 LightDirection = normalize(vec3(0, 4, 4) - Position);
    vec3 ViewDirection = normalize(uCameraPosition - Position);
    vec3 HalfVector = normalize(ViewDirection + LightDirection);
    float Roughness = 0.04f;

    float RefractiveIndex = 0.24f; // RI for Gold materials. I got this from http://refractiveindex.info/
    float F0 = pow(((1.0f - RefractiveIndex) / (1.0f + RefractiveIndex)), 2);

    float NdotL = saturate(dot(LightDirection, normal));
    float NdotV = abs(dot(ViewDirection, normal)) + EPSILON; // Avoid artifact - Ref: SIGGRAPH14 - Moving Frosbite to PBR
    float LdotH = saturate(dot(LightDirection, HalfVector));
    float NdotH = saturate(dot(normal, HalfVector));

    float DiffuseFactor = BRDF_Lambert(NdotL);
    float SpecularFactor = 0.0f;
    if(DiffuseFactor > 0.0f)
    {
        SpecularFactor = BRDF_Specular(NdotV, NdotL, NdotH, LdotH, Roughness, F0);
    }
    FinalColor = BaseColor * DiffuseFactor + SpecularColor * SpecularFactor;
}

Final Image:


Polygons versus curve primitives in software rendering

Most software rendering engines dice the parametric primitives to micropolygons, usually on the fly as needed. In essence this reduces the needed complexity to determine intersections. The surface will still look smooth, since each polygon is smaller than a pixel.

This allows for

Data caching.
Discrete geometric derivates.
Displacement is easy to implement.
Same rendering code to be used on most primitive types.

Image 1: Micropolygon render of trimmed NURBS surfaces using Renderman from 2003. No polygon is visible as each sample gets more or less a separate polygon.

But yes it is possible to use iteration method to solve surfaces from implicit formulations. It is typically used in raymarching contexts. The entire render engine could work by doing this its essentially a tradeof between memory and computation.


Silhouette curve for isometric surface of revolution

The silhouette curve of a quadric surface is a conic (ellipse, parabola, hyperbola, etc.). Suppose you have the quadric $\mathbf{X}^T \mathbf{M} \mathbf{X} = 0$ and your eye is at the point $\mathbf{Y}$. Then the silhouette curve is the intersection of the quadric and the so-called polar plane of $\mathbf{Y}$, which is $\mathbf{Y}^T \mathbf{M} \mathbf{X} = 0$.

For anything more complex than quadrics, things get nasty in a hurry, and you won't find any analytic solutions; you'll need numerical methods to trace the curve. Even the silhouette of a torus is pretty complex. So, in essence, the best you can do is compute discrete points along the curve, and use these to construct either a polyline or a spline of some sort. But, of course, you can always get a polyline just by tesselating the surface and finding the silhouette of the tesselation.

One approach would be to approximate your surface of revolution by quadrics, but that may not be any better than tesselating. You'd get a piecewise conic curve, instead of a polyline, which might be an improvement.

Note that if your surface is not G2, then it's silhouette will not be G1 (i.e. it may have sharp corners). Tracing a curve with corners is tricky.

There are silhouette generation functions in geometry kernels like Parasolid and ACIS. Maybe in OpenCascade, too.

To illustrate the difficulty, suppose we have a surface with parametric equations $(u,v) \mapsto \mathbf{S}(u,v)$, and we are viewing it from an eye- point $\mathbf{E}$. Let $\mathbf{N}(u,v)$ be the surface normal at parameter values $(u,v)$, which we can calculate as the cross product of partial derivatives: $\mathbf{N}(u,v) = \mathbf{S}^u(u,v) \times \mathbf{S}^v(u,v)$.
If a point $\mathbf{P}$ is on the silhouette curve, then the surface normal at this point is perpendicular to the vector $\mathbf{P} - \mathbf{E}$. So, the silhouette curve is the locus of points where
$$
(\mathbf{S}(u,v) - \mathbf{E}) \cdot \mathbf{N}(u,v) = 0
$$
Substituting for $\mathbf{N}(u,v)$, this becomes
$$
(\mathbf{S}(u,v) - \mathbf{E}) \cdot (\mathbf{S}^u(u,v) \times \mathbf{S}^v(u,v) ) = 0
$$
If $\mathbf{S}(u,v)$ is even a very simple surface, like a torus or a biquadratic or bicubic patch, the silhouette equation will not have any analytic solution.


What's the current state-of-the-art algorithm for ray-tracing height-fields?

For the current state-of-the-art, look for this paper:
"Maximum Mipmaps for Fast, Accurate, and Scalable Dynamic Height Field
Rendering", Tevs et al. 2008

The basic idea is to skip a lot of space by having knowledge of the maximum value over large areas of terrain. If the ray stays above that, skip to the next large area.

If you look at Figure 8, you'll see a comparison of basic linear stepping vs. maximum mipmaps. Linear stepping results in 200 steps, which can be done real-time on modern gpus but is still actually slow. Max mipmaps do the same in about 10 steps, all in shader.


How does everyone count "Rays per Second"?

I've seen both, unfortunately. I'm not fond of rays per second as meaning exclusively primary rays and I'd suggest "paths per second" or better yet "samples per second" instead. "Complete ray" is not a common term: a ray is a (potentially unbounded) line segment and a sequence of rays is a path.

Rays per second in your second sense of total ray casts is not well-specified for a path tracer: do shadow rays count, for instance? It's a useful metric for evaluating an acceleration structure or intersection testing framework (i.e. Embree or OptiX) but I'd avoid it for a renderer.

Finally, be aware that samples per second still isn't a great metric of actual performance since sample quality will vary wildly depending on implementation details. It's probably the best thing you can do starting out though, as the better solutions involve fairly complex variance estimates.


Brightness and contrast VS Gamma, is it possible to identify the gamma "correction amount"?

If I get correctly what you are asking you basically just need to find the G in this equation:

$$Image_{out} = Image_{in}^G$$

This could be easily solved as

$$G = \frac{\log{Image_{out}}}{\log{Image_{in}}}$$

Because usually gamma is applied in a uniform fashion on the image, you can just pick any two non zero pixel values (one for source and one for destination) to find out the gamma value applied.

EDIT: As @ChristianRau pointed out, if you don't actually know if the image has been modified with gamma correction, you have to take an higher amount of samples and trying to fit a gamma function on those. If a fit to a gamma function results in too many significant outliers, then probably gamma correction wasn't the function applied.


How to render and store a color gradient in a ramp texture for 16bit dataset?

It may be easier to compute the gradient on the CPU, even with 16bit precision the CPU can grind out a texture quickly. Then just upload that to the GPU as an image. Also, it doesn't sound like users will be able to pick a new color set every frame so its sounds a lot more natural to have a work flow of: Pick colors -> Generate Image -> upload to the GPU -> use new palette to render.
Are the picked colors RGB choices? If yes then it does make sense to spread the color palette out as wide as possible (maybe even more then 16bits) 65536 pixels is just a drop in the bucket. And a nice wide color range could lead to less banding in the final image.
At some  level the code will have to narrow the range down to whatever the current precision is. If that is the table, and the table has 65536 positions, then the table might contain a lot of repeating colors, but it is easy to compute the table with a new range.


Fast clipping without clearing stencil buffer

It sounds like what you need is a scissor test. It's specifically designed for clipping against an unrotated rectangle, and should be faster than messing about with the stencil buffer.


Why do 3d modeling packages use a single normal per polygon and how is this viable for smooth shading

I think the problem comes from this starting assumption:

If we have a quad that is non-planar.

A non-planar quad does not have a normal. It's not a flat surface, so you can't talk about what its normal would be. You can't talk about the normal of a sphere, a cone, or any other non-flat surface. You can talk about what the normal might be at any particular point on that surface, but there is no normal for the surface in its entirety.
You are asking a question for which there is no correct answer. But there are many incorrect answers, and they're all equally wrong. So there's no reason to pick one wrong answer over another.
If you're talking about computing vertex normals, a non-planar quad is still not a valid surface. The two triangulations for a non-planar quad represent very different surfaces, which would have very different vertex normals. Neither answer can claim legitimacy over the other, so both answers are equally wrong.
So why choose that particular method? Because it's no more wrong than any other method you might pick.
It does not make sense to talk about the vertex or face normals for non-planar polygons. So if you want an answer that makes sense, fix that problem first.


A half-edge is an edge split along its length, and having a directional component, that is, a beginning vertex and an end vertex. Where two polygons share an edge, each polygon gets a single half-edge between the same two vertices, which will have opposite directions if the winding order is consistent. These half-edges will have references to one another as two halves of a pair.

The full half edge data structure stores for each half-edge:

The polygon it belongs to
Previous and next half edges within the polygon
Its pair half-edge in a neighbouring polygon
Its origin vertex

This structure allows you to extract all sorts of connectivity information from a mesh, such as which edges or polygons lie around a particular vertex, simply by traversing the half-edges. There's a good article explaining it here.


Comprehensive explanation of conservative rasterization algorithm

I'm not sure what exactly you mean by

  how to visualize a plane that use w value as one of this coordinates

but here's a sketch that will, hopefully, clarify this sentence:

  a line through each edge of the input triangle as a plane in homogeneous (x c , y c , w c )-space

I've used uppercase letters for clip-space, and lowercase letters for window space.

Let's say AB is an edge of your original triangle.
It projects to an edge ab in window coordinates.
Point A is defined using [XA, YA, ZA, WA] coordinates. After projection it will have coordinates [xA, yA] = [XA/WA, YA/WA] and the value ZA/WA will be used to determine which object is closer to the camera at this point.
So your [X, Y, Z, W] always turns into a [X/W, Y/W] in the end. You can see that Z doesn't play any role in where the point is going to end up on the screen. So we just throw it away for now and interpret the remaining [X, Y, W] as 3D coordinates of the points.

Now we have point A = [XA, YA, WA] and B = [XB, YB, WB].

The projected point a will lie on the line that connects [0, 0, 0] and A.
Point b will lie on the line [0, 0, 0] and B. Three points define a plane, and a normal to the plane defined by 3 points can be found as a cross product of two vectors in this plane (OA and OB, for example).
"For a convex polygon (such as a triangle), a surface normal can be calculated as the vector cross product of two (non-parallel) edges of the polygon."
More on why in Wikipedia here and here.

Generally this moves your problem from 2D domain in the projected coordinates into a 3D domain of not-yet-projected points. The 2D and 3D are "dual" to each other. For example, if a point lies in the plane of OAB, it for sure lies on the line through ab (and vice versa).

The thing you need to remember is that W is your 3'rd coordinate, not Z! I know, it's very counter-intuitive, but it is what it is, Z only carries depth information, it does not affect position of objects on the screen.

I hope I cleared it up for you a bit.


What algorithm(s) are behind Google Street View's 3d mapping mode?

Google said how in 2012, avoiding explaining too many details

The never-ending quest for the perfect map

  ...new imagery rendering techniques and computer vision that let us
  automatically create 3D cityscapes, complete with buildings, terrain
  and even landscaping, from 45-degree aerial imagery

and here it's a video of the feature.

So we know that "Computer Vision" is used. I guess the edges are detected, along with the shadows, to build the elevation of the object. Then the segmented image is used as a texture on the built mesh.

However, i saw this paper linked, and i think a similar technique is used in Google to extract the 3d mesh from different 45° images.

There are a lot of papers (a book too) about "building extraction/detection".

Detection of buildings in urban and suburban areas using very high resolution satellite
Google Scholar: Building Detection
Google Scholar: Building Detection from Satellite Images

src: How does the new google maps make buildings and cityscapes 3D?


get partial view points from 3d surface mesh

In your case you would test every point to see if it lies inside your camera frustum, this is basically a set of point-vs-plane tests. This is not a fast way to do do things but the simplest to understand and implement.


What is the meaning of words 'stochastic', 'temporal' in computer graphics?

A stochastic process is one where you gather your samples randomly. So stochastic AA is antialiasing where you gather multiple samples for the same pixel with small random changes.
Temporal refers to a process over time. Temporal coherence is for instance mentioned when you have no flickering artifacts, or aliasing that you see during movement. It means that frames are coherent in a certain time window.


So my understanding of VBO they're essentially Buffer objects (bunch of bytes to be sent to the GPU) with the difference that VBO are specifically dedicated to Vertex Data. Is this correct?

A "VBO" is not really a thing, and I really wish tutorials would stop pretending that they are.
You are correct in that there exist buffer objects. These are linear arrays of unformatted data that the GPU stores and can access for various processes. One of those "various processes" that the GPU can use them for is providing array data for vertex rendering operations.
But the buffer itself isn't special. It doesn't become a different thing just because you use it for vertex data, and there is no rule that says that a buffer you're currently using for vertex data cannot be used for other buffer usages too.
Note that this is different from how textures work. When you create a texture object and bind it as a 2D texture, that texture object forever becomes a 2D texture. You cannot later bind it as a 1D or a cubemap or whatever. Buffer objects aren't like that; where you bind them defines only how you're currently using them.

What does it mean "any subsequent vertex attribute calls from that point on will be stored inside the VAO"

This is yet another annoying tutorial trope: that objects "store calls" or are macros or other such nonsense.
Objects in OpenGL don't store calls; they store state. They're like structs; they have data members. Binding an object to the OpenGL context allows OpenGL operations that access certain state to find that state within the bound object. So if you call a function that modifies certain state that is from a bound object, that function will modify the state within the object.
When you bind a VAO, you are binding an object that stores the state used for vertex arrays. When you call a function that modifies the state for vertex arrays, that function will modify the state stored in the VAO currently bound. When you render, the rendering system will use the vertex array state found in the currently bound VAO.


Here is an interesting answer from Styves on gamedev.net.

  Calculating bent normals is just an extension of AO calculation where
  the direction of each sample is also averaged along with the occlusion
  amount.

  The bent normals from substance painter seem correct. In any region
  that there is no AO, your normal won't be bent in any direction,
  therefore in tangent space it's pointing straight up along the surface
  normal. This is why you only see details in the ears.

  I have no idea how Unity handles the bent normals but the usual way is
  to, instead of applying AO as a standard multiplier onto the lighting
  result (which gives this weird muddy look), use the dot product of the
  bent normal and the light instead. This is how we handle SSDO in
  CryENGINE for example.

  I've only ever used a texture-based version of this technique once,
  however in this version the bent normals also contained information
  from the original tangent space normal map, and was used directly
  during lighting in place of the original.


Some questions about working in 0-255 integers instead of single precision floating point colors

First, a quick correction:

  storing color values as one byte integers I can have a combination of 255^3

You mean $256^3$. There are 256 eight-bit integers (including zero), not 255. But here comes the more fundamental thing you need to understand.

  If I store the red color as 128, this ends up being a floating point value of about 0.502

It doesn't end up being a floating-point value at all: it ends up as an 8 bit (or possibly 10 bit) integer on your display. Whatever extra precision you had in the floats to start with is thrown away when the final colour goes into the framebuffer. So the rest of what I'm going to say shouldn't influence your choice of what to store, because everything will end up as 8 bit colour channels.

Now that doesn't mean there's never any value in using a higher-precision number earlier on. More precision in the intermediate data means you'll get the full precision in your output format. e.g. if you start out with these normalized integers, and you multiply them all by 2 in the shader, you'll only end up with the even-numbered integers in your output: you lost half of your precision.

As for the question you really asked, about visible differences, it's easy to test that in live conditions. Just make a grey colour gradient from $(0, 0, 0)$ to $(1, 1, 1)$ in whatever representation you like (floats or normalized integers) and show it in place of your real computation. If you see Mach bands on your gradient, then yes, you can tell the difference between adjacent intensities. In general, for adults with normal sight and 8 bit-per-channel RGB in common lighting conditions, there is a visible difference in a grey gradient. That doesn't mean there's a visible difference between every adjacent pair: the eye is less sensitive to contrast in red and blue than in green, and less sensitive to changes in hue, so it strongly depends on the colours involved.

However, there are techniques to make the differences invisible most of the time. Check out my answer explaining what sRGB is for on SO for some examples.


Subdivision scheme where the faces and edges have weights (not necessarily scalar weights)

What you are looking for is semi-sharp creases. You can find it in section 3 of this paper: https://graphics.pixar.com/library/Geri/paper.pdf
Basically, each edge is given a sharpness value $s$. This can either be a integer or a floating point value. This value signifies how many subdivision steps this edge will be subdivided using sharp-subdivision rules. When subdividing an edge with sharpness $s$, it subdivide into two edges with sharpness level $s-1$. Until $s$ is 0 and the edge will be subdivided as normal.For fractional levels the last subdivision with $s \in [0,1]$ the vertices will obtain a blended position between the subdivided $p$ and the sharp position $p_s$, i.e. $(1-s)p + s p_s$.
The sharp subdivision rules boil down to the following. Vertices at the end of an sharpness tagged edge will remain in the same position and vertices on the sharp edge are refined using the stencil [1, 6, 1]/8 applied to previous, current and next points on the tagged curve respectively. New points are inserted at midpoints of edges.
This scheme can be generalized to faces by propagating the sharpness value of a face to its edges.


The texture is probably generated by picking a random angle per pixel, and populating the image with its sine and cosine, remapped into [0, 1]:

$$\theta \sim [0, 2\pi] \quad \to \quad \begin{bmatrix} \tfrac{1}{2} \cos \theta + \tfrac{1}{2} \\ \tfrac{1}{2} \sin \theta + \tfrac{1}{2} \\ 0 \end{bmatrix}$$

I did a quick test and generated something that looks pretty similar (enlarged 4x):

If you want to be extra fancy, then instead of picking independent random angles for every pixel, you might try stratified sampling or a low-discrepancy sequence. That may give better final results for whatever effect you're using the jitter texture for, due to better local sample distribution.


Drawing a plane super-imposing over another plane

The problem that you encounter is called "z-fighting" and it occurs when two objects with the same depth (z value) are rendered.

In this example, we take two objects called ObjA and ObjB. They are the same, but have different colors. Due to the very small error that you get with fixed precision floating point numbers, ObjA may have a lower fragment depth than ObjB for a given pixel and thus ObjA's color is saved to that pixel, but ObjB could have a lower fragment depth for another given pixel and ObjB's color is then saved to that pixel. Which object will have a lower depth value is quite random and causes that flickering.

There are two ways to fix this. Which one works best for you, depends on the program.

Option #1

Disabling depth testing on the second surface that is being drawn. The first surface obviously draws well, and the problem comes when drawing that second surface. By disabling the depth test, OpenGL just draws every fragment to the framebuffer and you do not have that problem with "which fragment has the lowest depth".

If you only disable the depth test for the second surface and enable it again for the other surfaces (and render the second surface after the first surface), then it should work as the first surface still writes the depth values.

This approach starts to break, when the second surface either does not fully align with the first surface or the second surface is shaped in a way so that it overlaps itself. (When rendering a sphere, for each pixel you always render two fragments, the front side and the back side of the sphere. The back side would then not be culled away.)

Option #2

Adding a bias to the depth values of the first surface. The problem is that both surfaces have the same depth value for a given pixel, but if you offset the depth values for one of them, then you do not have that problem.

Just like gl_FragColor which allows you to set the color in the framebuffer, gl_FragDepth allows you to set the depth in the depth buffer. You could get the normal depth value from gl_FragCoord (as the documentation says) and add a small amount to it. Doing this, does not actually change anything of the mesh. The mesh does not move around and it does not go further into z-space. The shader simply writes a different depth value to the depth buffer. This depth buffer is only used for depth testing and can thus easily be changed to get the effect we want. This only needs to be done for the first surface and not the second surface or any other surface.

When the bias is too large, artifacts can visibly show up around intersections with the surface, but when the bias is small it should not be noticeable.


Zero-copy texture reading from CPU memory

Is such a method of zero-copy texture loading possible/optimal on modern discrete graphics hardware?

You say that as though it were optimal back then. The Intel 740 was not exactly warmly embraced by the market, as noted by the Wikipedia article in question:

  The AGP Texture concept soon proved to be a tremendous error in design, because the card had to constantly access the textures over a channel that was upwards of eight times slower than RAM placed on the graphics card itself.

The disparity between available memory bandwidth/performance and available compute resources has only increased in the last 20 years. So, will it be "optimal" to have a discrete GPU read from main memory?

No.

Could it be "possible"? Not in terms of marketability. Even $50 and below discrete GPUs carry their own on-board RAM. Anything much slower/smaller than that isn't worthwhile, because pretty much everyone already has Intel or AMD's on-die GPUs for when performance isn't of importance. If you have a discrete board, it will always be better to put some RAM on the board with its own bus than to have the board try to read most of its data directly from main memory.

  Furthermore, could zero-copy texture access functionality be exposed by using mapped buffers in OpenGL?

No. Mapped buffers would be buffers, not textures. They're different kinds of storage objects with (hypothetically) distinct allocations behind them. PBO is a means of asynchronously transferring data from a buffer to a texture. That is a transfer of data, not an adoption of a memory allocation.

Intel has an extension for mapping a texture's storage, but that requires that the texture be in a linear format. It's not widely supported, and likely never will be.

Vulkan provides for the possibility of something similar, but using linear texture as sources in sampler read operations from a shader is not required to be supported. You'd have to ask each implementation whether it allows it and in which formats, and for which usages are allowed.


It looks to me like you want a buffer containing a single array of integers (not a buffer containing multiple arrays, whatever that would mean). So, you should be able to just do this:

RWBuffer<int> indices;

then access it like this:

indices[int(id.y)]


Confusion around Lambert's Cosine Law in Ray Tracing in One Weekend

It is cancelled out by the probability density function in the estimator. The pdf in their case is exactly: $\frac{\cos\theta}{\pi}$, which is in the denominator: att = albedo * cos_theta / pdf = albedo * pi. They have absorbed pi in the albedo. Note that there's an update on github to the code, since random_in_unit_sphere() actually generates a $\cos^3\theta$ distribution, they fixed it to random_on_unit_sphere() which generates a $\cos\theta$ distribution. So technically you are working with a $\cos^2\theta$ brdf here and not constant.


Should animations be locked to FPS or Fixed time

You should use a fixed time for this. Consider the following: Your computer is super fast, you have high end hardware. Now you animate via FPS and make your animation look great. Let's say you had 120 fps.

Since you did great work, you want to show me. I have bad hardware, everything takes ages. I have to work with just 30 fps - sufficient to call it real-time, but not too nice. It will take four times to completely play your animation on my computer, because you made it follow the fps. Therefore it doesn't convince me at all.

This is even more problematic, if you do this in a multiplayer scenario. You could already start a new action after the animation, whereas I have to wait for the animation to finish. Therefore, there is an inconsistency in game state as well as an unfair game, because you can do things much faster.


OpenGL layered rendering with multiple additional clip planes per view

The GS is responsible for assigning a particular primitive to a particular layer. The GS is also responsible for writing values for gl_ClipDistance for a particular primitive. So the GS is perfectly capable of doing both of those things: assigning a primitive to a layer and doing the clip-plane computations for that layer.
And if you need to "deactivate" a clip plane for a particular primitive assigned to a particular layer, just assign the gl_ClipDistance[i] for that to 1.0 for all of its vertices.


From a quick gander through the examples it looks like QGLViewer is a thin wrapper around a opengl window using the legacy api.

So what you really need is a way to load a .off file to render in legacy opengl. it looks like a very simple format, even simpler than .obj. So simple that IMO it is worth it for you to figure out on your own.


How to open and render in a second window with OpenGL / JOGL with NEWT

So after letting this lie for quite a while, I just came back to it. Turns out, the error was in my code with a faulty false return value which would lead to a System.exit. I have no idea how I didn't catch this the last time round, but it now works with Try 2. That is:

public class SubWindow implements GLEventListener
{
    GLWindow glWindow;
    Program m_program;

    public void createSubWindow(int windowWidth, int windowHeight, boolean isFullscreen, String windowTitle
            , Program program)
    {
        GLCapabilities caps = new GLCapabilities(GLProfile.get(GLProfile.GL3));
        caps.setBackgroundOpaque(false);
        glWindow = GLWindow.create(caps);

        glWindow.setTitle(windowTitle);
        glWindow.setSize(windowWidth, windowHeight);
        glWindow.setUndecorated(false);
        glWindow.setPointerVisible(true);
        glWindow.setVisible(true);

        m_program = program;
        glWindow.addGLEventListener(this);
        glWindow.addKeyListener(m_program);
        glWindow.addMouseListener(m_program);

        Animator animator = new Animator();
        animator.add(glWindow);
        animator.start();
    }

Where the Program is one of my own classes handling the logic and the SubWindow class additionally implements the required init, display, etc. methods.


The camera does not need to move for this problem to exist. You can see the mixed polygons as in your linked image even with a static camera.

Things are worse with a moving camera because it makes the position of the polygons change, which can lead to different rounding. This means it may not be the same polygon which comes up in front at a given location for every frame. Even a very small change in the camera could swap which polygon is visible. So you get flickering which is visually a lot worse because it draws attention.

One possible solution is to progressively blend surfaces together when their Z values become too similar. However, I'm not certain how easily that can be implemented in traditional Z-buffer rendering.


Modern Screen Space Ambient Occlusion Techniques

This bachelor thesis briefly reviews six SSAO techniques. It could be a good start.

CryEngine 2 AO
StarCraft II AO
HBAO
Volumetric Obscurance
Alchemy AO
Unreal Engine 4 AO

http://frederikaalund.com/a-comparative-study-of-screen-space-ambient-occlusion-methods/


Which technique should be used for rendering under water images, Ray tracing or Radiosity?

Radiosity is only efficient for diffuse to diffuse rendering. Path tracing is more versatile on the type of light transport that could be computed. The main issue with underwater rendering is the caustics generated from specular from the water-air interface. These caustics are very important to generate a realistic underwater image (except for deepwater rendering here the sunlight is completely attenuated).

The problem of rendering caustics from path tracing is the lack of non-efficient sampling strategies to generate a valid path. Density estimation [1] (a.k.a. photon mapping) could be a good candidate for rendering these difficult caustics, but a large scene will require light subpaths guiding to get a reasonable rendering time. Metropolis light transport (MLT) [2,3] provide this guiding automatically and could be used in a density estimation framework [4].

Recently, path guiding techniques [5,6] was proposed to render these difficult paths efficiently with unidirectional path tracing. The simplest approach to use/implement is the technique proposed by Muller et al. [5]. A Mitsuba code plugin is available on GitHub: https://github.com/Tom94/practical-path-guiding

Note that, inside water, there might are volume scattering events [7] to take into account, depending on the amount of particle floating.

[1] "Stochastic progressive photon mapping" Hachisuka et al., 2009

[2] "Metropolis light transport" Veach et al. 1997

[3] "Primary sample space MLT" Kelemen et al. 2002

[4] "Robust Adaptive Photon Tracing using Photon Path Visibility", Hachisuka et al., 2011

[5] "Practical path guiding for efficient light-transport simulation", Muller et al. 2017

[6] "On-line Learning of Parametric Mixture Models for Light Transport Simulation", Vorba et al. 2014

[7] "Monte Carlo Methods for Physically Based Volume Rendering", Novak et al., 2018


Most, if not all, approaches to GPU programming require the developer to manually identify code that will be executed on the GPU. The code for the GPU is compiled separately from CPU code, so it is clear which binary code is for which processor type. Moving the GPU code to VRAM is often done with DMA, but may be moved to GPU accessible memory using other techniques such as CPU writes across PCIe, or direct CPU memory access by the GPU.

The GPU compiler often resides in userspace, and a system call will hand binary GPU code to the kernel, which can then make the code available to the GPU.


fix over exposure in importance sampling in screen space reflections

After some testing and some research I have found that both the $\sin(\theta)$ and the $\cos(\theta)$ should be removed from the denominator.

the reasoning is that, since we want a probability function of steradians instead of angles, the $\sin(\theta)$ factor (which is essentially used to convert from steradians on a unit sphere to radians on the two polar axis) should be dropped.

Also, $\cos(\theta)$ (as Nathan Reed explains in this blog post) is there to convert from microfacet area to macrosurface area. Since we are generating microsurface samples, this should be removed.

After these corrections we get the formula

$$L_o(p,\omega_o)\approx\frac{1}{N}\sum_{i=1}^{N}\frac{L_i(p,\omega_i)FG}{4(n\cdot\omega_o)}$$

which is free of the artifacts of the formula on the question.


Finding all possible reflection matrices for a given Wythoff construction

I found an answer in the form of the Todd-Coxeter Algorithm from group theory. This still involves iterating over permutations of the group, and removing duplicates. However, the group can be defined using symbols to represent different sequences of operations, and duplicate-matching performed on these symbols rather than elements of the group itself.

For example, if a pair of mirrors A and B have an angle of (PI/3) then they form a cycle of order 3 (since a reflection through a pair of mirrors equals a rotation through twice the angle between them). So the sequence ABABAB maps any element back onto itself. By constructing a path like this for each pair of mirrors, it's possible to define new elements for each step along a path, then eliminate duplicates either through string-matching operations or by constructing a graph and eliminating equivalent nodes. For large sets, this should be much more efficient than multiplying and comparing matrices directly at each step.


Algorithms to "anti-alias" (or somehow improve) binary 1-bit drawings and fonts

PIL's Image.convert function performs dithering by default when you convert the image to 1-bit—not simply thresholding. That's what creates the noise along the edges of the shapes; the antialiasing isn't the problem. Unfortunately, at such a low resolution, the dithering does more harm than good.

I grabbed your original image from your other question and did a straight threshold at 128, and it comes out looking a lot better:

That's probably about as good as you'll be able to do on a low-resolution 1-bit display, absent a bitmap font that's specifically designed for those circumstances.

In PIL this operation can be implemented using Image.point, as follows:

image1Bit = image8Bit.point(lambda x: 0 if x < 128 else 1, mode='1')


The meaning of early ray termination and front-to-back traversal in ray tracing

Front-to-back traversal is the idea that when traversing a ray through the acceleration structure, we want to examine nodes that are in front, i.e. closer to the camera, before other nodes. If you only want the first intersection with a surface (as usual in ray tracing), then if you get an intersection in the front node you don't need to traverse any nodes behind it; as any possible intersection there would be later than the one you already found. This saves time.

For example, in a kd-tree, at each node you have a splitting plane dividing space into two child nodes. You would first descend to the child closer to the camera, and only if no intersections are found there do you descend to the other child.

Early ray termination seems to be a term used in the context of volume rendering. When you're traversing a ray through a medium that attenuates light, such as fog, smoke, cloudy liquids, etc., then beyond a certain distance there is no possibility of the ray having any more effect on the image. If you have thick fog with a visibility of 10 meters, then an object 11 meters away cannot be seen.

You can precalculate the distance at which the attenuation reaches 99% or so (it never quite reaches 100%, but approaches it exponentially). There is no reason to intersect the ray against objects farther than that, so when using the acceleration structure, you can simply discard any nodes beyond that distance. In effect, the ray becomes a line segment terminated at a fixed length. If it doesn't hit anything within that length, there's no point checking further. (Another case where limited-length rays are handy is for shadow rays. You only want to trace them back to the originating light source; any intersections beyond that are irrelevant for the shadow.)


Regarding optimizations:

Could you calculate the min and max corners in modelspace beforehand? Then you could check on the CPU-side the NDCs of these two points and skip the processing of the whole mesh if it is outside of the desired range. This needs only two matrix-vector multiplications on cpu side. And since for rigid bodies these two points do not change over time, you don't have to iterate the meshes in each frame.


Shader sampler2D uniforms instead of attributes

I copied my answer for the same OP at https://github.com/godotengine/godot/issues/19473. But this question differs slightly to the original issue report.
In case you want to manually add support for VERTEX_ID, these are the changes I had to make in 3.2.3-stable for spatial shaders. I wanted to use the vertex ID for creating procedural animations.
In drivers/gles3/shader_compiler_gles3.cpp in the function ShaderCompilerGLES3::ShaderCompilerGLES3() :
actions[VS::SHADER_SPATIAL].renames["VERTEX_ID"] = "gl_VertexID";

and servers/visual/shader_types.cpp in the function ShaderTypes::ShaderTypes():
shader_modes[VS::SHADER_SPATIAL].functions["vertex"].built_ins["VERTEX_ID"] = constt(ShaderLanguage::TYPE_INT);

Then compile the source and it should work. Some other primitives could probably be added in a similar fashion but I have not tested that yet.
One caveat is that the vertex IDs reported by Godot are very different from those reported by other tools such as Blender. This is because vertices are duplicated and in a different order. Godot also only gives you two UV channels to work with. Using UV coordinates also means having a temporal component will require significantly more memory to implement as UVs will typically map sparsely onto a 2D texture, if the data is only used in the vertex shader.
I have a primitive tool working that can map the vertex animation texture generated in blender to one usable by Godot into another grayscale texture (n by 1 texture, where n is number of vertices in Godot). This saves texture memory as no data is duplicated for duplicated vertices. I first need to implement an octree before it will be usable for real-world meshes and not just a set of simple deforming cubes. Once I have that done in my spare time I will release it for free.


Should I prefer an AMD or Nvidia GPU for studying graphics programming, especially with respect to vendor specific debugging tools and ecosystem?

If you want to learn OpenGL X.Y or DirectX A.B then shop for hardware that supports that specification.

If you're using other software (Blender, Maya, 3DSMax, etc) you'll want to make sure it's supported.

If your OS has a history of issues with a particular vendor, adjust accordingly.

Beyond that, it's not going to matter much, especially at the hobbyist level.
Price / performance metrics are still relevant, it's just that the definition of performance is different.


Parametrize set of unordered points in 3d space

You could try natural neighbour interpolation or Sibson coordinates for this. These coordinates are based on a Voronoi tessellation of the points.
The coordinates have been extended to an approach where the points lie on some surface in this publication.


Classify acceleration structure using surface area heuristic

You gave the formula as:

$$Cost(cell) = C_t + \frac{S_L}{S_{Cell}} \times N_L \times t_i + \frac{S_R}{S_{Cell}} \times N_R \times t_i$$

(You put $r_i$ at the right-hand end, but I assume that was a typo.)

There's a more general form of this formula

$$Cost(cell) = C_t + \sum_{k\ \in\ K} \frac{S_k}{S_{Cell}} \times N_k \times t_i$$

where $K$ is the set of all child nodes; that is, to get back to the two-children form you started with, $K = \{L, R\}$.

By generalising the formula like this, you can now calculate the cost for the first possibility the same way as the second two.

You're right that this is a slightly odd question. In practice, your traversal cost is going to be less for a binary tree than a $(> 2)$-ary tree.

  Does this ranking change if we change the spheres with a (scaled) tea-pot that still fits within the bounding volume?

To answer this, you'll need to consider how this affects the surface area, and how this changes each score in turn.

  Can this ranking depend on the position of the camera, when we only    consider visibility rays (so no shadow rays, for example)? You can assume you're working with an orthographic camera.

There's two things to consider here.

1) Does any part of the heuristic depend on the position of the camera? Can moving the camera change any of the numbers that feed into the score?

2) Moving the camera can change how the boxes overlap. Remember that rays that pass through the box need to be tested against the box behind, so needless overlaps make the structure worse.

I think the point is that the heuristic is only a general rule. It's not hard to construct an example with two proposals a and b, where a scores much better but b is actually a better split for a particular camera, and this last part of the question is trying to get you to think about this.


OpenGL task/mesh shader extension NV_mesh_shader not behaving as I expected

I managed to figure it out. A coworker told me that I was misinterpreting the documentation.
layout(max_vertices=3, max_primitives=1) out;
actually means max_vertices and max_primitives per workgroup, not per shader invocation. So I was always only outputting a single triangle, regardless of what I set local_size_x to.

Changing the values to proper numbers resulted in a proper rendering.


How to obtain a parametric equation out of a part designed in Autodesk Inventor?

I guess your question is how to compute the inertia tensor automatically. A parametric equation may not be always possible, but we could write a program to avoid hand calculations as much as possible.
The inertia tensor $\mathbf{I}$ is defined as an integral over the object domain $\Omega$ (see Inertia tensor):
$$
\newcommand{\V}[1]{\mathbf{#1}}
\begin{align}
\V{I} &= \int_\Omega \rho (\V{r}\cdot\V{r} \delta - \V{r}\V{r}^T)\,dV
\end{align}
$$
where $\V{r}$ is the position relative to the center of mass (COM),$\rho$ is the density, and $\delta$ is the identity matrix.
Note that $\text{tr}(\V{r}\V{r}^T) = \V{r}\cdot \V{r}$, the above definition  can be rearranged as
$$
\begin{align}
\V{T} &= \V{r}\V{r}^T \\
\V{I} &= \int_\Omega \rho (\text{tr} \V{T} \delta - \V{T})\,dV
\end{align}
$$
Analytic object
For analytic objects you can substitute their equations into the integral and do symbolic integration. This step can be automated with any symbolic computing software.
Object in tetrahedral mesh
We assume the object is simple and has a closed boundary.
An analytic solution is not possible in this case because we do not have an explicit equation for the object. But we could compute the inertia tensor numerically.
The only integral we have to compute is
$$
\begin{align}
\V{C} &=
\int_\Omega \rho \V{r}\V{r}^T\,dV
\end{align}
$$
To compute this integral, this paper could be helpful.
The general idea is to compute $\V{C}_t$ for each tetrahedron $t$, and sum $\V{C}_t$ to get $\V{C}$.
We first compute the integral for a canonical tetrahedron $t0$ whose vertices are
$(0,0,0)$, $(1,0,0)$, $(0,1,0)$, and $(0,0,1)$.
In maple, this canonical integral $\V{C}_{t0}$ is (assume $\rho$ constant)
r:=<x,y,z>;
                                                    [x]
                                                    [ ]
                                               r := [y]
                                                    [ ]
                                                    [z]

with(LinearAlgebra):
rrt:=r.Transpose(r);
                                             [ 2               ]
                                             [x      x y    x z]
                                             [                 ]
                                      rrt := [        2        ]
                                             [x y    y      y z]
                                             [                 ]
                                             [               2 ]
                                             [x z    y z    z  ]

seq(seq(int(int(int(rrt[i,j],x=0..1-y-z),y=0..1-z),z=0..1),i=1..3),j=1..3);
                      1/60, 1/120, 1/120, 1/120, 1/60, 1/120, 1/120, 1/120, 1/60

Then the integral $\V{C}_t$ can be computed by a change of variable.
Assume the coordinates of a world space tetrahedron can be represented as $\V{r}=\V{F}\V{r}_0 + \V{x}$, where $\V{F}$ and $\V{x}$ are constant.
We have
$$
\begin{align}
\V{C}_t &= \rho \int_{\Omega_t} \V{r}\V{r}^T\,dV \\
&= \rho \int_{\Omega_{t0}} (\V{F}\V{r}_0 + \V{x})(\V{F}\V{r}_0+\V{x})^T \det \V{F} \,dV_0 \\
&= \rho \int_{\Omega_{t0}} (\V{F}\V{r}_0\V{r}_0^T\V{F}^T + \V{x}\V{x}^T + \V{F}\V{r}_0\V{x}^T + \V{x}\V{r}_0^T\V{F}^T) \det \V{F} \,dV_0 \\
&= \det(\V{F}) \V{F}\V{C}_{t0}\V{F}^T + \det(\V{F}) m \V{x}\V{x}^T + \det(\V{F}) m \V{\bar{r}}\V{x}^T + \det(\V{F}) m \V{x}\V{\bar{r}}^T
\end{align}
$$
where $\V{\bar{r}}$ is the COM and $m$ is the mass.
Objects in other form
I did not use Inventor before but I suspect the part produced by it may have other formats, such as Bezier surface. For now I would suggest tetrahedralize the object and apply the above method.


An algorithm to find the area of intersection between a convex polygon and a 3D polyhedron?

There are different numerical approximations you could use:

A simple solution is to use brute-force Monte Carlo integration. Distribute $N$ random points on the polygon and calculate the number of points inside the 3D polyhedron $N_i$ using ray tracing. The area inside the polyhedron is $A_i=A\frac{N_i}{N}$, where $A$ is the area of the polygon.

To improve convergence of the above algorithm, you can distribute the points on grid on the polygon and jitter the sample positions within each grid cell. This is called stratified sampling.

Further improvement would be to perform the MC integration in 1D instead of 2D. Pick a sample line $L$ that's on the polygon plane and a sample direction $D$ that's on the plane perpendicular to $L$. Select sample line segment $L_s$ that bounds the polygon within $D$-directed lines (min & max of the polygon projection on $L$ for tight fit). For each sample on $L_s$ first calculate the line segment $S$ that intersects the polygon in direction $D$. Next calculate the length of $S$ that's inside the polyhedron $S_i$. You can calculate $S_i$ using ray tracing, by starting from one end of $S$ and marching along $S$ to intersection points of the segment with the polyhedron, summing up the segments of $S$ inside the polyhedron. More specifically, you need to first check if the start point on $S$ is within polyhedron and set the "inside" state accordingly. Then for each intersection flip the state and add the length from previous intersection to $S_i$ if the state changes from inside to outside. The approximation of the polygon area inside the 3D polyhedron is:
$$A_i=\frac{|L_s|}{N}\sum_{i=1}^N S_i$$


Rays (~1k) intersecting triangle meshes (~100k) with ca. 100 vertices each

Your procedure could be optimized in two directions:

If your rays are structured you could possibly take advantage of space coherence: if a ray intersects a given triangle, then a nearby rays has high possibility of intersecting the same triangle or an neighboring one. Therefore simple parallelism over the rays may not be optimal, you could at least process them in batches based on adjacency.
If you have adjacency information on the triangles of your meshes, you could possibly exploit it when testing a ray against the mesh: the Möller-Trumbore algorithm can inform you if the intersection point of the ray and the triangle plane lies inside or outside each of the edges, so you may have a hint at which (neighboring) triangle to check next. That said, your meshes are quite small so the extra overhead may not be worth it.