Lowering register pressure doesn't necessarily give you any performance boost though. I recently went through this exercise myself on GCN architectures (for a simple ray tracer) and reduced register pressure so that it increased occupancy from 2 to 4, which had no impact on performance. It's generally a good idea to reduce the pressure if you need to hide memory latencies, but it really depends on what the real bottlenecks in your shaders are.

Since you are working on a GPU ray tracer, you might get better improvement in performance by paying attention to the divergence of threads and try to improve this instead. For example by trying to group rays based on location and direction to reduce divergence in KD-tree traversal per wave. I'm not familiar with the paper you are referring to, and you may already have considered this though.


You can look into this: https://cescg.org/wp-content/uploads/2017/03/Dvo%C5%99%C3%A1k-Estimating-Approximate-Plane-of-Symmetry-of-3D-Triangle-Meshes.pdf
There are definitely other papers on the subject though


The lightmap representation used in Valve's stuff, like any finite-storage representation of the surfaces response to light, is only approximate, and the error in the approximation is higher or lower in different scenarios.

Since it doesn't store the light in a direction tangent to the surface (where it would be exactly zero), you don't get an exact zero for this case. However, this isn't that important:

Normal maps generated from height maps will never yield a normal vector exactly perpendicular to the geometry, and tend not to be not that steep anyway.
The 3 generated lightmaps are normalized so that their sum is equal to the light for the "straight up" direction. This reduces directional contrast, but also eliminates light from the tangent direction.


How to handle different data types (uint8, uint16, float) in Volume Datasets?

Ok after searching a whole lot on the internet I found a slide which answered most of my questions and doubts. I'm dumping it here in-case anybody else is curious what these datasets usually contain and what should be (in my opinion) the correct way to handle mappings.

So apparently the CT scanner always records the densities of the substances in Hounsfield Scale. The Hounsfield Scale itself is a little interesting. For body organs and stuff, the scale has the range from -1024 to 3000. This was the scale I talked about earlier as when normalized the values usually lie between 0-4096. However metals usually go past this 3k boundary and some metals like stainless steel and gold having HU number as high as 13k and 30k respectively. This is usually called the "extended" Hounsfield Scale.
Some volumetric datasets contain data for computer simulations mostly fluids etc. For these type of datasets either density or flow velocity is captured. Precision is important for these datasets hence float is used.

Hence the uint8 datasets are probably grey-scale data converted from the original HU units. uint16 datasets contain the original Hounsfield data while float may contain density or velocity.

Now for the mapping part, as others stated one way might be to just divide by the max value of the data type or the max value of the data itself, to get the uint16 or float data in 0-1 range then multiply by 255 to get a grey-scale image. However doing so may lower the contrast and hide minute details. For example, dividing a uint16 dataset that contained only HU units in range 300-2000 with the $2^{16}$ will map the HU units to only a small range of grey values. Dividing by 2000 in this case is surely better. However after reading how CT-scanners work, it's better to let the user set the minimum and maximum range of Hounsfield units they wish to observe. This is usually called the window width. And then map this window width to the whole 0-255 range of grey levels as YardenJ2R mentioned in his second formula. This is more in line with how radiologists visualise the data.

For fancy transfer functions which I haven't implemented yet, we might have to convert the Hounsfield units to opacity or grey scale as mentioned above then map that opacity to color or there might be a way to directly map it to colors using a LUT as defined by the user's control point settings.


OpenCL ray tracer: which is better, geometry code on device side or on host side?

In order to pass data to the device side, in OpenCL you must create the sphere structure on the device as well (that is in the kernel). So you first create your structures as classes or structs on the CPU. The variables in the structs should be of cl_* datatypes. OpenCL provides typedefs for all primitive data types. Then create these structures on the GPU as well.

Care should be taken about data alignment when passing structures from CPU to GPU. If you have passed data in OpenGL, you'll know what this is about.

After you create your structures on the CPU, you encode them in a OpenCL buffer, and pass this buffer to the kernel.


Function to convert HSV angle to RYB angle

I just discovered that Adobe color source includes HSV-RYB hue mapping functions (replicated in Ben Knight's Kuler-d3). Apparently Adobe uses uses piecewise linear gradients rather than the polynomial that I was trying to use (and technically it is a CMY wheel not a RBY one, I believe). Here are the relevant stops:

RYBstop HSVstop
60      35
122     60
165     120
218     180
275     240
330     300

Here is the graph of these showing the perceptual mapping:

So there'll be artifacts at the junctions where the mapping function is not smooth (which can be seen in the original Adobe image). And here is an imagemagick script for creating the color wheel.

#!/bin/bash
#assemble gradient pieces (pre-caculated stops), then join them all
convert -size 600x300 gradient:#000000-#181818 -rotate -90 grad1.png
convert -size 600x310 gradient:#181818-#2A2A2A -rotate -90 grad2.png
convert -size 600x215 gradient:#2A2A2A-#555555 -rotate -90 grad3.png
convert -size 600x265 gradient:#555555-#7F7F7F -rotate -90 grad4.png
convert -size 600x285 gradient:#7F7F7F-#AAAAAA -rotate -90 grad5.png
convert -size 600x275 gradient:#AAAAAA-#D4D4D4 -rotate -90 grad6.png
convert -size 600x150 gradient:#D4D4D4-#FFFFFF -rotate -90 grad7.png
convert grad?.png +append grad.png
#create hue
convert grad.png -alpha set -virtual-pixel Transparent -rotate 180 -distort Arc '360 -90 300' +repage -gravity center -crop 600x600+0+0 +repage h.png
#create saturation
convert -size 600x600 -alpha on radial-gradient:white-none s.png
#create value
convert -size 600x600 -alpha on radial-gradient:none-white v.png
#combine h,s,v in hsv colorspace
convert h.png s.png v.png -combine -set colorspace HSV -colorspace sRGB ryb.png


Why would I specify image format in a shader layout qualifier?

From the same wiki

  Image variables can be declared with a format qualifier; this
  specifies the format for any read operations done on the image.
  Therefore, a format qualifier is required if you do not declare the
  variable with the writeonly memory qualifier. Write-only variables
  cannot be used as in any reading operations; this includes calling
  load and atomic (read/modify/write) functions. So if you want to read
  from an image, you must declare the format.

  The format defines how the shader interprets the bits of data that it
  reads from the image. It also defines how it converts the data passed
  for write operations when it writes it into the image. This allows the
  actual Image Format of the image to differ between what the shader
  sees and what is stored in the image, sometimes substantially.


Kaiser windowed sinc filter for mip mapping

I was curious so I tried it. I used the Kaiser Window from Wikipedia: $\frac{I_0\left(\alpha\sqrt{1-x^2}\right)}{I_0\left(\alpha\right)}$

This page has a slightly different formulation using $\pi\alpha$ instead of $\alpha$. Given the comments in the page you linked about the Kaiser window being very similar to Lanczos and how they even overlap in the graph, I suspect they also used the former formula. You can see how close they are if you plot both on Wolfram.

This left me wondering if that's what Mitchell originally meant as it's not "noticably better". I could not see the difference at all with my test texture (it's there, just not noticeable). However, using the second formula, the difference in the plots is much more striking. And there is a noticeable difference in the resulting images. Whether it's better is a matter of interpretation and depends on the filter size.

With $x \in \left[-1,1\right]$ , I used this C++ windowing function:

#include <boost/math/special_functions/bessel.hpp>
float bessel_i0( float x )
{
    return boost::math::cyl_bessel_i( 0, x );
}

float window_kaiser( float x )
{
    const float alpha = 4.0f * M_PI;
    float k =
        bessel_i0( alpha * sqrtf( 1.0f - x * x ) ) /
        bessel_i0( alpha );
    return k;
}

This is the second version. Remove the M_PI to get the first version which is more like Lanczos.

The modified Bessel function of the first kind seems to be a C++17 / boost thing. If you have neither, this quick hack will do:

float bessel_i0( float x )
{
    float r = 1.0f;
    float term = 1.0f;
    for( int k = 1; true; ++k )
    {
        float f = x / float(k);
        term *= 0.25f * f * f;
        float new_r = r + term;
        if( new_r == r )
            break;
        r = new_r;
    }
    return r;
}

Just don't consider that a general replacement for a proper bessel function. However, it's fine for computing the window as far as I could tell.

As the second wikipedia page I linked shows, that window_kaiser can be used on a symmetric filter with $N$ weights by using $x = \frac{2n}{N-1}-1$, where $n$ is $0, 1, 2 \dots N-1$.

Computing the weights of the sinc filter itself is another thing entirely as it depends on filter size, how much you're resizing, etc. I suggest posting a separate question if you don't know how do to that.


A common mistake that causes this behavior is to first apply the translation and then the rotation when rendering the scene. I think that's what's happening. If you apply the rotation first and then the translation it should behave as expected. Give that a try


Why replicating the higher bits of RGB565 when converting to RGBA8888?

There's actually a reasonably good mathematical reason for doing bit replication:

First note that the n-bit string, $N$, actually represents the value $\frac{N}{2^n-1}$ and we want to produce the m-bit string, $M$, where $n<m$ and
$$\frac{N}{2^n-1}\approx\frac{M}{2^m-1}$$

We first scale numerator and denominator with
$$\frac{N.(2^n+1)}{(2^n-1)(2^n+1)}\approx \frac{M}{2^m-1}$$
and this simplifies to
$$\frac{N.(2^n+1)}{2^{2n}-1}\approx \frac{M}{2^m-1}$$

In your case, $n\in \{5,6\}$ and $m=8$ and we can "stop" here, but but the process can be repeated, (ad nauseum), if m >> n.

We next make the approximation...
$$\frac{N.(2^n+1)}{2^{2n}}\approx \frac{M}{2^m}$$ which simplifies to
$$\frac{N.(2^n+1)}{2^{2n-m}}\approx M $$

Note that $N.(2^n+1)$ is equivalent to repeating the n-bit string, to create a 2n-bit string, and the division shifts off the $2n-m$ LSBs to leave an M bit result.

QED

Of course, the  'correct' calculation is $M=\lfloor(\frac{(2^m-1) N}{2^n-1}+\frac{1}{2}\rfloor$ but this approximation, in general, works most of the time. Of course there are times when it's inaccurate, but IIRC only by one bit and relatively infrequently.


Cost of parameter passing in webgl pixel shaders?

My experience working with shader compiler stacks a few years back is that they are extremely aggressive, and I doubt you will see any perf difference, but I would suggest testing as much as you can.

I would generally recommend providing the compiler (and human readers) with more information where the language allows it, marking parameters according to their usage. Treating an in-only parameter as in/out is more error prone for humans than compilers.

Some detail: shaders run almost entirely using registers for variables (I worked with architectures that supported up to 256 32 bit registers) - spilling is hugely expensive. Physical registers are shared between shader invocations - think of this as threads in a HyperThread sense sharing a register pool - and if the shader can be complied to use fewer registers there's greater parallelism. The result is shader compilers work very hard to minimize the number of registers used - without spilling, of course. Thus inlining is common since it's desirable to optimize register allocation across the whole shader anyway.


What is the oldest hardware able to feasibly support Vulkan?

Khronos publishes a list of Vulkan compliant hardware here. The list also shows compliance to different flavors of OpenGL. I don't know all the card models, but for NVIDIA at least the GeForce GT 630 is supposed to work with Vulkan.

To answer your question a bit more, to support Vulkan a video card must support a certain number of programmable shader stages. Shaders is probably the biggest evolution of graphic cards in the last 10 years in my opinion. This already eliminates the graphic card you mentioned in the question.

This early Vulkan presentation from Khronos, states that "Vulkan 1.0 provides access to OpenGL ES 3.1 / OpenGL 4.X-class GPU functionalitybut with increased performance and flexibility." That can give you a basis for finding which individual features a graphic cards must theoretically support to have a Vulkan driver. See the new features introduced in OpenGL 4 at Wikipedia.

Among them:

Tessellation shaders
Indirect draw commands
Transform feedback (more flexibility in transform feedback)
Texture buffer objects (new texture formats only)

See also the sections "New Features" of the OpenGL 4.x specs.


How to get the nearest polygon face and height of intersection point

To get the height intersection once you have found the polygon, you can build a plane of the polygon and calculate the distance to the plane. For example you can displace the point along the planes normal and calculate the length of the normal for when the point is displaced into the plane. Implicit planes are neat for this, they are build with the normal of the plane (which is the normal of the polygon) and one constant $d$. $d$ can be calculated by arbitrary point on the plane (e.g. one of your polygon vertices).

$d = -\vec{n} \cdot q$ and $\vec{n} + d = 0$ (where $\vec{n}$ is the normal and $q$ the point/vertex).

The neat part is: if you calculate the dot product of your point you want the distance of, you have the signed distance toward the plane, i.e. if $p$ is your point and $l$ is the distance between the point and the plane, then

$l = (\vec{n} + d) \cdot p$ (for this you assume the fourth coordinate of a point as $1$)

Now that you have the distance to the plane, you need to check, weather the intersection of normal through point with plane also intersects the polygon. You could for example take the vectors $\vec{v}_i = v_i - p$ where $v_i$ are your polygons, and then take the dot product $\vec{n} \cdot \vec{v}$. If they all have the same sign, your intersection doesn't hit the polygon (assuming the polygon to be convex).

Iff it doesn't hit the vertices you need the distance to the closest two vertices. Since you have already done the dot product test (with the sign) for all polygon vertices, you know that the shortest distance can't be to the edge between the two vertices, and thus it suffices to take the shortest of the two point->vertex distances.

As to which polygon is closest, I'm not sure if there is a better trick, but I'd just transform the point into Model Space of the mesh and calculate the distance to any front facing face.


Update the uniform and draw again.

auto tMVP = glm::translate(MVP, glm::vec3(1,0,0));
glUniformMatrix4fv(loc_MVP, 1, GL_FALSE, &tMVP[0][0]);
gluSphere(...);

tMVP = glm::translate(MVP, glm::vec3(2,0,0));
glUniformMatrix4fv(loc_MVP, 1, GL_FALSE, &MVP[0][0]);
gluSphere(...);

Using glTranslate would have done the exact same thing behind the scenes.

You can instead use instancing instead where you upload all the transforms into a buffer and then tell the gpu to use that.

However the more portable way involves allocating a vertex attribute for the data. This involves using glVertexAttribDivisor.

So a full mat4 would need 4 vertex attributes. You can decrease that to 2 by separating out the translation as a vec3 and the rotation as a quaternion and pass the view-perspective as a uniform and combining it all in the vertex shader.


Understanding distinct vs. not distinct knots in B-splines

So, the number of polynomials that you piece together depends on multiple factors, including the number of control points and the degree of the curve. If the number of points is equal to the order, then you have a single polynomial, if the order is $n$ and you have $n+1$ points, then you are piecing together 2 polynomials, etc...

Distinct points are points that are different from each other.
Take the vector [0,1,2,3,4] where all knots are distinct and contrast it with [0,1,2,2,3,4] where the knot 2 has multiplicity 2 (it's not distinct since it appears twice). These knots lower the degree of the curve locally by 1 for each repetition.

The overall curve can be defined in many ways, the simplest way to look at it is as DeDoor's algorithm, which is a generalisation of DeCastlejeau's.

I am going to do a bit of self promotion and link my own notes of B-Splines. Gitlab's PDF viewer butchers them a bit so I recommend downloading them and reading them with Okular or Adobe reader.

Regardless, that document has more information than what can be put into a stack exchange answer and does indeed try to convey the intuition about how to build B-Splines, both mathematically and algorithmically.


Why does the integral of NDF over a solid angle equals the area where micronormals belong to that angle?

$D(\omega)$ is defined as the area ($m^2$ unit in the numerator) of the microsurface with normals pointing in the direction $\omega$. $\mathcal{M}'$ is defined as the portion of the microsurface with normals point in the direction $\omega \in \Omega'$. So it's natural that the integral of $D(\omega)$ over $\Omega'$ gives the area of $\mathcal{M}'$.

Substituting the definition of $D$ into (4) gives:

$\int_{\Omega'} D(\omega_m) d\omega_m = \int_{\Omega'} \int_{\mathcal{M}} \delta_{\omega}(\omega_m(p_m)) d p_m d\omega_m$

A fundamental property of the Dirac delta function is that it integrates to 1 over its domain. Unfortunately, we are integrating over $\Omega'$, not $\Omega$. Here I will stop pretending that I am actually able to use formal mathematical terms, and I will just note that the integrand has non-zero values only where both $d p_m \in \mathcal{M'}$ and $d\omega_m \in \Omega'$, which is what (3) says. As long as at least one of our integration domains is restricted to the set we care about, we will get the same result. So I'll just move the prime around a tad:

$\int_{\Omega'} D(\omega_m) d\omega_m = \int_{\Omega} \int_{\mathcal{M}'} \delta_{\omega}(\omega_m(p_m)) d p_m d\omega_m$

Now we can apply the Dirac delta identity, since we're integrating over the whole domain:

$\int_{\Omega'} D(\omega_m) d\omega_m = \int_{\mathcal{M}'} d p_m$

Ta-da!


What parts of unbiased path tracing pipeline can benefit from OpenCL?

In short: All of it can benefit from OpenCL.

Path tracing is done by taking many individual ray samples from a scene and combining them into a final image.  As the number of samples increase, so does the quality of the resulting image.

If you are doing regular path tracing (not photon mapping or bidirectional path tracing or anything like that), each sample needs just READ ONLY access to the scene so can work completely independently of the other samples.

That fact makes path tracing astonishing parallelizable.  If you had an output image resolution of 1000x1000, for a total of 1,000,000 pixels, and took 10,000 samples that would mean you needed to take 10,000,000,000 (10 billion) samples.

Each one of those samples is parallelizable, so if you had 10 billion cores, you could make use of all of them.

There are pragmatic concerns when doing this including cache coherency, but the costs of getting that wrong with the naive implementation are much smaller than the benefits of going massively threaded in the first place.

Here is some more information about the details of naive path tracing if you are interested:
http://blog.demofox.org/2016/09/21/path-tracing-getting-started-with-diffuse-and-emissive/


How do I swap the buffers when using a Win32 window with Vulkan?

I have now come to the point where I need to implement SwapBuffers(Window) type of function

Then stop doing that. You need that function in OpenGL because the windowing code is explicitly not part of OpenGL, so to write cross-platform GL code, you need an encapsulation layer.

That is not the case for Vulkan. While you still need a layer to create the window, handling the display, swap-chain images, and so forth is all part of the Vulkan API. The user's code should have complete and total control over this stuff; your abstraction layer need not be involved.

Indeed, since Vulkan has no concept of a "default framebuffer", there's no way for you to abstract away presentable images without user support. Just do what GLFW3 does in its Vulkan support: provide a cross-platform way for the user to create a Vulkan surface.


How to compute A and B in projection matrix

Note: I'll be using column vector notation (eg. $Ax=y$), thus if you are using row notation, transpose the matrix and multiply from the left instead of the right. (eg. $y=x^TA^T$)

A perspective transformation matrix transforms a view frustum into a rectangular view volume. Note that the resulting view volume might not be a cube and might not be centered at 0.

One example of a perspective transformation matrix. There are many more, and yours is included.

Let $l=left,\\ r = right,\\ b = bottom,\\ t = top,\\ n = near,\\ f = far$

$M_{persp}=\begin{bmatrix}
1& & & \\
 &1& & \\
 & &\frac{n+f}{n}&-f\\
 & &\frac{1}{n}&
\end{bmatrix}$

Since the resulting view volume is not the canonical view volume, we need to apply another matrix.

A orthographic transformation matrix transforms a rectangular view volume into a cubic view volume with two corners at (-1, -1, -1) and (1, 1, 1), which happens to be the canonical view volume.

Intuitively it will need one translation and one scaling.

$\text{T}_{ortho} = S\left(\frac{2}{r-l},\frac{2}{t-b},\frac{2}{n-f}\right)T\left(\frac{-(l+r)}{2},\frac{-(b+t)}{2}\frac{-(n+f)}{2}\right)$

$M_{ortho} =
\begin{bmatrix}
    \frac{2}{r-l}& & & \\
    & \frac{2}{t-b} & & \\
    & & \frac{2}{n-f} & \\
    & & & 1
\end{bmatrix}
\begin{bmatrix}
    1 & & & \frac{-(l+r)}{2} \\
    & 1 & & \frac{-(b+t)}{2} \\
    & & 1 & \frac{-(n+f)}{2} \\
    & & & 1
\end{bmatrix}
$

$M_{ortho} =
\begin{bmatrix}
    \frac{2}{r-l}& & & -\frac{l+r}{r-l} \\
    & \frac{2}{t-b} & & -\frac{t+b}{t-b} \\
    & & \frac{2}{n-f} & -\frac{n+f}{n-f}\\
    & & & 1
\end{bmatrix}
$

Now you can combine them together to form a perspective projection matrix.

$M_{proj} = M_{ortho}M_{persp}$

$M_{proj} =
\begin{bmatrix}
    \frac{2}{r-l}& & & -\frac{l+r}{r-l} \\
    & \frac{2}{t-b} & & -\frac{t+b}{t-b} \\
    & & \frac{2}{n-f} & -\frac{n+f}{n-f}\\
    & & & 1
\end{bmatrix}
\begin{bmatrix}
1& & & \\
 &1& & \\
 & &\frac{n+f}{n}&-f\\
 & &\frac{1}{n}&
\end{bmatrix}$

$M_{proj} =
\begin{bmatrix}
    \frac{2n}{r-l}& & \frac{l+r}{l-r} & \\
    & \frac{2n}{t-b} & \frac{b+t}{b-t} & \\
    & & \frac{n+f}{n-f} & \frac{2fn}{f-n}\\
    & & 1 &
\end{bmatrix}$

If you have a free-moving camera you will need a camera matrix too, thus the final projection matrix centered on the camera is.

$M_{proj} = M_{ortho}M_{persp}M_{camera}$


Real GPUs, like CPUs, have a cache hierarchy. So the rasterizer doesn't directly store to memory, but stores into a cache that eventually gets flushed to memory.

Rasterization will usually have some spatial locality (it will generate writes to nearby pixels close together in time), so if those pixels are in the same cache line it will reduce the overall number of memory transactions. This is more effective for filled polygon rasterization than for lines, but even lines have some locality (especially if multiple nearby lines are being rasterized at once).

Another trick is to store the framebuffer in a tiled format instead of linearly left-to-right  and top-to-bottom. You'd pick a tile size such as 4×4 or 8×8 and store the pixels in each tile in a contiguous memory block, then make the whole framebuffer out of an array of tiles. This gives you locality in both dimensions, not just one, which improves cache performance. And tiling is easy to do in hardware, as it comes down to just swizzling the order of some X and Y address bits.


Physically realistic received irradiance values given some BRDF?

I didn't get this part:

  That way each ray would just multiply its solid angle by its radiance
  value, and you'd get an "Irradiance" value of W/m^2.

I'll try to explain how the whole thing works though. You have a light source, it emits light rays. Those light rays bounce around the scene and eventually end up at your camera film. So at the end of the day you have multiple rays hitting the film surface and the radiance of those gets multiplied with the sensitivity function. For efficiency reasons you usually start backwards - from the camera. But if you work with symmetric BRDFs this should not matter. So what arrives at your film is radiance, you multiply that radiance with the sensitivity function and integrate within a "film pixel" to get the result in terms of your screen pixel.

Edit: To clarify what I meant by integrated out:

If you have radiance which is a measure of W/(sr * m^2), you can integrate out sr or m^2 like so:

$$E(x) = \int_{\Omega}{L(x,\omega)\cos\theta\,d\omega}$$

From where you get irradiance at a point $x$, which gives you the energy arriving from all directions at point $x$ (W/m^2). Now you can go further and find out the energy that arrives at some surface with area A, by integrating over all points on the surface:

$$\Phi = \int_{A}{E(x)\,d\mu(x)} = \int_{A}\int_{\Omega}{L(x,\omega)\cos\theta\,d\omega\,d\mu(x)}$$


Could a hexagonal pixel array store an image more efficiently?

Here's my take on it.  A pixel is not a square, and it isn't even a rectangle.  A pixel is a point (infinitely small) that has a color associated with it.

The only way I personally have ever seen pixels viewed (interpreted) by a display is to use "nearest neighbor" sampling where the pixels were on a rectangular grid, which means that the color of any given space on a display is the color of the pixel that it is nearest to.

This is a fancy way of saying "pixels are rectangular and laid out on a grid", but stay with me on this :P

As a result, image formats have their pixels stored in a grid as well, with the assumption that nearest neighbor in a grid will also be used to display it.  For instance, many images will have anti aliasing built into them so that they will look good when displayed on a "nearest neighbor grid".

Interactive applications (games), may use textures that are not meant to be displayed as is on nearest neighbor grids so are kind of an exception to that rule.  They do this because as part of their execution, they do anti aliasing, bilinear texture sampling, etc, so that whatever picture they push out to the display will look good, when the display shows it as a nearest neighbor grid!

Now, getting closer to your question: would a hexagonal grid have any advantages?

I think that yes, it would!

First off, I think nearest neighbor would look better.  I don't have any real proof of that sorry, but the hexagon more closely approximates a circle, and since it isn't a regular grid of data, i think your eye is getting a better distribution of data.  Sorry, that is a little hand wavey.

I think the big part of why it would look better though is that linear filtering would be taking information from 6 neighbors instead of 4, and would be interpolating on 3 axis instead of 2.  More information from less regularly spaced samples than a grid gives you = better resulting image.

Doing cubic interpolation would also be better than cubic interpolation on a grid, so the quality scales up as you scale up the quality of your algorithm too.

As far as whether it stores data more efficiently, the fact that it can do better filtering with less data means to me that yes, it could store data more efficiently.

And i guess lastly... maybe you could use these properties to your advantage.  Maybe you could have an image format stored in a lower res hexagonal format, and then before you needed to display the image at runtime, you could use sampling algorithms to convert it back to a grid.

Maybe there'd even be a good way to do this efficiently in a pixel shader, so it would use less memory at runtime too?

It's an interesting idea (:

PS - how cool would it be to have an analog display (no individual pixels, but a continuous colored surface) that looked at the pixel data you sent not as rectangles, but instead as sample points on a continuous surface.  Maybe a bit out there though....


Rendering moist/wet materials like human anatomies for surgical simulators

The reason why ray-tracing is often used for achieving as-realistic-as-possible wet appearance is their reflective behavior. A thin liquid film on the surface reflects photons like a mirror. If you target real-time rendering, you will have to fake the reflections (or rely on RTX :) There are two ways I can quickly think of: placing reflection probes in the scene using cube-maps or using screen space reflections.
Reflection probes with cube-maps

good for static scenes, you can pre-compute many of them
allow objects outside of the current viewing frustum to be seen in the mirrors
inaccurate due to limited number of probes

There is a great extensive tutorial by Sébastien Lagarde on reflections which continues (what a coincidence) with part B on wet surfaces and there is part C as well. Make sure to check them out!
Screen-space reflections

fine also for dynamic scenes
only objects in the current viewing frustum can be seen in the mirrors
accurate within the frustum

Screen space reflections are explained as part of this awesome series by David Lettier.
Subsurface scattering
Regarding the discussion on subsurface scattering, any organic tissue exposes this property. We just observe it the most on skin or leaves, but the same effect is observed inside the body. Just read section 3.3 of this paper:

One significant casual effect under endoscopic environment is a strong subsurface scattering on the surface of internal organs.

I suggest to use at least some post-processing technique for that.


What are the advantages of having two chroma channels and a luma channel vs. three chroma channels?

From a historical context, the YUV color encoding was originally designed to support the existing black-and-white infrastructure. In addition, as jooja mentioned, humans are generally less perceptive to changes in chroma and more perceptive to changes in luma. Both NTSC and PAL use less bandwidth for the chroma channel than the luma channel. In a digital context, we can use chroma subsampling to achieve the same effect.

YUV/YIQ/YDbDr is used for color PAL/NTSC/SECAM television broadcasts.
YPbPr is used for component video cables.
YCbCr is used for everything else.
YCgCo is like YCbCr but is easier to compute.


What's the meaning of the "dynamic" in "high dynamic range"?

Dynamic, adjective:

  1: pertaining to or characterized by energy or effective action; vigorously active or forceful; energetic:

  ...

  4: of or relating to the range of volume of musical sound.

HDR allows you to deal with high light energies. And such high energies could be characterized as being "vigorously active or forceful; energetic".

But I think the fourth definition is more suitable. Graphics and rendering are forms of signal processing, much like sound. So they share a lot in common. For example, the term "aliasing" originated in digital conversion of sound.

"High light energies" are the visual equivalent of "volume" in sound. So the "dynamic range" that you can record, the "range of volume" is high, as opposed to low.


it is quite easy to measure the local max frequency in an image (at least as a low resolution mask, with some regularization).
Several papers of the MIT graphics group have been around detecting and processing from this kind of clue, with regular or coded aperture cameras.
e.g. Defocus Magnification and Image and Depth from a Conventional Camera with a Coded Aperture.


Let's look at your formula:

pixel_y = (y_world / z_world + 1) / 2 * ScreenHeight();
pixel_y = ScreenHeight() - 1 - pixel_y;

Let's call y_world / z_world for short a, H for ScreenHeight(). So you have:

Y = (a + 1) / 2 *H
Y = H - 1 - Y

Let's replace the first in the second:

   Y = H - 1 - (a + 1) / 2 * H
=> Y = H * (1 - (a + 1) / 2) - 1
=> Y = H * ((2 - a - 1) / 2) - 1
=> Y = H * (1 - a) / 2 - 1

Which is the formula you said every one uses (except with an off-by-one error):

pixel_y = (1 - y_world / z_world) / 2 * ScreenHeight();

Off-by-one error aside, looks like you have the right understanding, just didn't work out that the formula everyone uses is the optimized version of what you came up with.


Legendre Polynomial equation in Spherical Harmonics

The notation Re() and Im() refer to the real and imaginary parts of a complex number.  Mathematicians and physicists are accustomed to using spherical harmonics (and Fourier transforms too) that are complex-valued, due to the factor $e^{im\phi}$.  You would then also have complex coefficients, in general, in the spherical harmonic expansion of a (real or complex) function.

Using Euler's formula, $e^{im\phi} = \cos(m\phi) + i\sin(m\phi)$.  So that factor encodes cosine and sine waves  (that oscillate $m$ times as you move around the equator of the sphere) in its real and imaginary parts, respectively.

When we know we're going to be working strictly with real-valued functions, it may be more convenient to use real-valued variants of the spherical harmonics, where $e^{im\phi}$ is replaced by either $\cos(m\phi)$ or $\sin(m\phi)$, and using real coefficients instead of complex ones.  We trade a single complex coefficient for two real coefficients, so we haven't lost any information or flexibility; it's essentially a change of basis.  This formulation just explicitly ensures that everything always comes out real.


Performance of vector graphics versus bitmap or raster graphics

As TheBuzzSaw said, it does depend on lots of things, including implementations of the rasterized graphics vs the vector graphics.

Here are some high performance vector graphics methods that are rendered using traditionally rasterization methods.

Loop and Blinn show how to render a vector graphics quadratic bezier curve by rendering a single triangle, and using the texture coordinates in a pixel shader to say whether a pixel is above or below the curve:
http://www.msr-waypoint.net/en-us/um/people/cloop/LoopBlinn05.pdf

The basic idea is that you set your triangle corner positions to be the 3 control point positions, and you set the texture coordinates at each corner to be (0,0), (0.5,0) and (1,1) respectively.  In your shader, if the interpolated texture coordinate (x*x-y) is < 0, the pixel is underneath the curve, else it's above the curve.

You can see a faux implementation of it on shadertoy here:
https://www.shadertoy.com/view/4tj3Dy

As for the second method, here's a method from Valve, where distances to a shape are stored in a texture, instead of pixel data, allowing vector graphics to be drawn by using texture sampling.  The decoding is so simple, it could be implemented even on fixed function hardware using only an alpha test!
http://www.valvesoftware.com/publications/2007/SIGGRAPH2007_AlphaTestedMagnification.pdf

To give you an idea of how well the second technique works, this 1024x768 mustache image was generated from a 64x32 source image that had a single color channel! (aka 2KB uncompressed)

I also wrote some stuff about it on my blog:
http://blog.demofox.org/2014/06/30/distance-field-textures/

Here is some sample OpenCL code to show just how simple it is:

float alpha = read_imagef(tex3dIn, g_textureSampler, textureCoords).w;
float3 color = (alpha < 0.5f) ? (float3)(1.0f) : (float3)(0.0f);

Both of these techniques are super fast, and blur the line between vector and rasterized graphics a bit.  They are rendered using rasterization techniques, but have zooming / scaling properties like vector graphics techniques.


If I understand you correctly, the problem could be defined as the (edge) intersection of a cube with a sphere. In that case, you could shade the edges by calculating the distance to the edges of the cube.

//note: from http://iquilezles.org/www/articles/distfunctions/distfunctions.htm
float udBox( vec3 p, vec3 b )
{
  return length(max(abs(p)-b,0.0));
}

Aligning the cube poses another "interesting" issue - you probably want to use the sphere-normal and screen-up to generate a basis.


You basically need to plot the scatter map with varying colors for each point

scatter(x, y, [], c);

The main difficulty is to calculatec so that it’s proportional to the density.

To do so you can divide the map into a grid of a certain resolution, say [res_x, res_y]. Then count for each grid cell how many points are there. This count is an approximation of the density at that grid cell. Finally you can assign to all points within this grid cell the c value equal this count.

Probably you also need to normalize c to [0, 1]


The core calculation seems to be

  sub_v3_v3v3(vec, td->iloc, tc->center_local);

  radius = normalize_v3(vec);

  tratio = ratio * td->factor;

  mul_v3_fl(vec, radius * (1.0f - tratio) + t->val * tratio);

  add_v3_v3v3(td->loc, tc->center_local, vec);

Which in C-ish psuedo code looks like

vec = pos - center; //find vector from center to point
radius = length(vec); // get the length of the vector
vec = vec / radius;  //normalize the vector

tratio = ratio * individual_sensitivity; // adjust sensitivity per vertex

vec = lerp(vec * radius, vec * t->val, tratio); //linearly interpolate between the original and where the point would be on the sphere

pos = center + vec; // store result

So the code projects the point to the sphere and linearly interpolates the point between the original location and the sphere location.


Why does an affine transform work only on three of the corners?

An affine transformation doesn't have enough freedom to do what you want. Affine transforms can be constructed to map any triangle to any other triangle, but they can't map any quadrilateral to any other quadrilateral.

One way to see this is that the matrix for a 2D affine transform has only 6 free coefficients. That's enough to specify what it does to 3 points, but not 4.

For this task you likely want a projective transformation. A 2D projective transformation has enough freedom to map any convex quadrilateral to any other.


Determining shortest vector between a ray and an axis-aligned 3D bounding box (correcting a near-miss)

If you know the rounding error, you could subtract the error from your vectors so that rays would overlap the bounding box. This could be done in your intersection tests.

In 3D, computing shortest distance between a vector and face of a bounding box would be another way to shift the vector to overlap the bounding box. But I doubt this approach would work well as it essentially tries to fix a rounding problem.


There is a trick to shadows like this, but only this exact usecase. While you can certainly draw a true shadow using many techniques. The trick is you can just scale the object flat onto the floor plane.
So all you need to do is draw the object 2 times one time with a matrix that is scaled to zero on the floor axis. You can also use this technique to draw shadows that arent entirely top down by adding skew to counter the now rotated shadow.
This is super easy to implement. No extra cameras. But works only on one plane. This is a special case of shadow volumes without the need to draw the volumes.


Is long term continuous use of GPGPU safe for my GPU?

Running a GPU at full capacity will reduce its lifespan through electromigration; the speed at which the chips are damaged depends on how hot it is.

A desktop computer has enough room that the designer can put in a cooling system to handle the worst-case thermal situation.  Running your GPU at continuous full load won't overheat it too badly, and you can expect it to last several years.

A laptop, on the other hand, is strongly space-constrained.  Typically, the cooling system is designed to handle brief bursts of heat interspersed with long periods of near-idle operation.  You may not be able to run it at full load 24/7: the BIOS or operating system will slow things down to keep from overheating.  In any case, running at full capacity will likely cause things to burn out in a year or less.


Why cubic curves provide the minimum curvature interpolants?

Splines confuse me which is one reason I asked somebody else to write that chapter.

But I like this explanation:
https://www.johndcook.com/blog/2009/02/06/the-smoothest-curve-through-a-set-of-points/


Realistic chainmail reflections for videogames

It seems that you're mixing physics/engine programming with graphics programming. I wouldn't say that dealing with collisions is related to rendering.
Also, I don't really understand your problem with toroids, but in case, I found one paper related to this term - "3D ChainMail: A Fast Algorithm for Deforming Volumetric Objects" by [Gibson1996]. I would say that if you're interested in creating clothing then you'll have to deal with deformations and animations and not only reflections. And, in case if you'd like to just work with reflections/rendering, then you don't need to deal with physics in any way. Static mesh would be enough, but I'm just speaking based on informations that you've provided until this point.

Normal cloth seems too simple, and the reflections look bad (they have to be in the texture, or the cloth would have to be metalic)

I wouldn't agree with that and I don't really understand what you've meant by that. In many cases, reflections are indeed in the texture. Or more specifically - cubemap, local environment map, irradiance volume and many other names, used for explaining GI voxel-based methods. https://www.kostasvardis.com/files/research/PhD_Thesis_Kostas_Vardis_book_online_lowres.pdf - [Vardis2016] PhD thesis is nice place to check for methods that you might be interested in. I wouldn't say that normal cloths are too simple. I was never really interested in that area but there are number of papers and overall, creating realistic clothing in not an easy task.

Normal cloth seems too simple, and the reflections look bad

I can only guess that by "bad", you've meant that there are little to none reflections. That's because most of clothing is rather rough and you've got little materials with big specular term in everyday life. That's why there is little reflection to it.
I don't know what is your current level of knowledge. Also, I don't know if you're going to creating real-time interactive application or if you're going to create offline renderer. In any case - if you have Turing graphics card or newer, I would suggest using Ray Tracing for most realisting reflections. You can start with cwyman tutorials - http://cwyman.org/code/dxrTutors/dxr_tutors.md.html also you should read microsoft's and nvidia's introduction to DXR. Also, here is a relatively easy path-tracer, created by Matt Petineo - https://github.com/TheRealMJP/DXRPathTracer
If you get through above tutorials, then you can start reading about denoising and state of the arts methods introduced in recent years - https://alain.xyz/blog/raytracing-denoising Raytracing is a subject, which exists since the begining of computer graphics, so there are a lot of papers related to that. But I think that MJP's example is good enough for start. Also, read related papers that he's stating in github and his code. They're mostly related to avoiding artifacts and fireflies.


Is this smaller grid for Poisson disc sampling still correct?

You approach will work, and in general any square of size <  will work fine because the invariant "at most 1 point per square" is valid. Extrapolating from your idea, it means that one should get the best reduction in area checks for infinitesimally small squares. But we don't want to do that, right?

There are two reasons why  is usually better than :

How much state do you want to manage? Each square can either be "occupied" or "unoccupied" and contains a pointer to its sample point (if any).

For very large domains with small r, storage costs can be hard to manage. The problem is much worse if you're working with many dimensions. See Gamito et al.'s paper, Figure 12. Even in 4D you can be severely storage limited.
In your proposed grid size, each candidate points checks 24 squares, as opposed to 20 in the usual method. There is a per-square cost (e.g. checking whether square is occupied and fetching sample point if it is) which increases the overall cost in your method.

Basically, your method works, but it is usually not more efficient.

However, the state of the art in unbiased Poisson-disk sampling techniques actually successively reduces the grid size as the uncovered area reduces. In addition to Gamito et al.'s paper, also see Ebeida et al.'s paper.

P.S. Sorry, I don't know how to do math formatting in Stack Exchange.


The values for each pixel in the Z-buffer are interpolated from the values at the corners of the triangle during rasterization. To make this work, the projection matrix produces depth values that are a function of the reciprocal of the original depth value in camera-space.

That is, the $Z'$ value in each transformed vertex is of the form $aZ + b$, where $a$ and  $b$ together determine the near and far planes, and the transformed $W'$ value is just the original $Z$. After doing the homogeneous divide by $W'$, the transformed $Z'$ becomes $a + b/Z$. By using the reciprocal of $Z$ in the depth value like this, you can then just interpolate linearly in screen space from the three corners to each pixel to get correct results. Since the last step produces floating point results, and most Z-buffers use unsigned integer values, you'll need to bias, scale and quantize before writing to the  actual buffer.

I seriously recommend Eric Lengyel's book Mathematics for 3D Game Programming for the nitty-gritty details of all this kind of stuff, he also explains perspective-correct interpolation of vertex attributes and all the other math you'll need to know to write a rasterizer.


I believe a common solution is to split the camera transform used to project the grid from the camera transform that is used to render the grid. At perspectives close to top-down, the two cameras coincide, but as the viewing camera gets close to a horizontal perspective, the projection camera deviates and tries to keep a minimum inclination, i.e. it hovers somewhere above the view camera and looks down slightly.

The tricky bit is making sure that the field of view of the projection camera always covers the region of the scene seen from the render camera. I don't have a resource at hand that details how to compute the appropriate transforms, and it might be tedious to derive by hand.

A different solution is to grab the signal processing toolbox: The artifacts seen in your image are essentially aliasing, caused by insufficient sampling of the wave heightfield by the projected grid. Therefore, one solution is to filter the heightfield appropriately, depending on the projected area of a grid cell. I believe this is used in offline rendering of oceans, and it essentially makes sure that waves at the horizon go flat. However, I'm not sure how feasible this is in real-time rendering, since you would need high-quality anisotropic filtering to make this approach look reasonable.


The original use for the stencil buffer was for non-3D content. Say you're making a 3D game. Most of the screen is the 3D content, but you also have some GUI around the outside which is not drawn through OpenGL but through a 2D canvas API. You can mask off the shape of the 2D GUI using the stencil buffer, and then make sure that your 3D/OpenGL content won't affect it and won't cause overdraw (which wastes fragment shading time).

As you've seen, it's also used for techniques like shadow volumes and some reflection techniques. In general, it can be used for techniques where you draw an object and then want to run another fragment shader (or some more draw calls) to alter that object's appearance, by adding reflections, lighting, or screen-space effects. In that case, when you draw the object the first time, it's very cheap to fill the stencil buffer at the same time. Then, your later pass (to add the reflections or whatever) can use the stencil test to only run on those fragments that eventually did get drawn on the screen.

Another common use is when rendering a scene in multiple parts. For example, in a flight simulator, it's common for the cockpit (and aircraft exterior) and the terrain (and other aircraft) to be done in completely separate passes. This lets you use a lower-precision depth buffer without causing z fighting, because you can use different depth ranges for the near content and the far content. But if you do this, you can't reuse the same depth values for both passes. In this case, you can fill the stencil buffer when you're drawing the cockpit, and use the stencil test to avoid overdrawing the cockpit when you draw the terrain. (You need to do more to handle glass effects on the cockpit windows: this description is an oversimplification.) You can use similar techniques for drawing the player's weapon and/or hands in first-person shooters.

The stencil buffer is used in many techniques: too many to list them all individually. But in general, if you want to have some render passes only affect some fragments, and the fragments you want to affect are those corresponding to particular objects or shaders, then the stencil buffer is a cheap way to enable that.


In any modern desktop hardware, there shouldn't be: uniform buffers will just be memory blobs read by the shader. (See also: Uniform Buffer Objects)


Dealing with wasted space when using SIMD on CPU for vector maths

There's a third option you've neglected. Because of how graphics APIs work, typically you have to copy the data between the CPU and GPU anyway - even if you're on a mobile SoC where they share the same physical memory. This copy is the ideal time to take the hit of changing format, since you have to touch all of the memory anyway. So the copy in user memory, owned by the CPU, can be stored padded for SIMD alignment. Then, when you need to send that data to GPU memory, instead of doing a memcpy, just read in the 4-vectors and write out 3-vectors. It's less memory bandwidth than a normal copy (you're reading just as much, but writing less) and negligible CPU work.

If you've got that fourth float doing nothing, you might as well use it and make your co-ordinates homogeneous instead of just xyz. Depending on what compute you're doing on the CPU, that might make some things easier (transforms are the obvious example, but also reducing divisions if you were renormalising your normals a lot before).

This is similar but not identical to a trick I used when implementing OpenGL texture upload. On ARM, most Neon instructions operate at memory bandwidth rate, so if you're already doing a memcpy, it's free to do a pipeline-worth of instructions at the same time. You can use those free instructions to do pixel format conversion and reordering to Z-order or U-order.


How do you implement perspective correct Gouraud shading across a triangle?

(Is gouraud shading even supposed to be perspective correct?)

Originally, it wouldn't have been perspective correct, but on (hardware) systems these days it will be. FWIW Dreamcast had perspective correct Gouraud shading because, once you are doing perspective correct texturing, it is relatively little additional cost to do Gouraud "correctly".

  I have flat shading working and I have perspective correct texture mapping working. Is the way your raster Gouraud shading much different from the perspective correct texture mapping?

I'm assuming that you've started with homogeneous 'world space' vertices, i.e. something like,

$[x,y,z, 1, u,v, R,G,B]$,

applied a projection matrix to get

$[X_p, Y_p, Z_p, 1, u, v, R, G, B, w]$,

and then, after clipping, divided through by w to obtain screen coordinates

$[X_p/w, Y_p/w, Z_p/w, 1/w, u/w, v/w, R/w,... ]$.

For convenience, let's rename $u/w$ and $v/w$ as $u'$ & $v'$, $1/w$ as $w'$, $R/w$ as $R'$ etc.

To do perspective correct texturing I assume you are linearly interpolating, per pixel, the $u'$, $v'$, and $w'$ values, computing the reciprocal of the interpolated  $w'$ and multiplying that by $u'$ and $v'$ to obtain the per-pixel U & V coordinates.

To do perspective correct Gouraud, simply use the same approach with the interpolated $R'$ , $G'$ and $B'$ values. The most expensive part is typically the division (i.e. reciprocal) but that can be shared with the texturing.


OpenGL Instanced Rendering with Level of Detail

IMO, a variation of your second approach is going to be the simplest, and could be quite fast.

Have you profiled anything yet? A CPU should chew through 10000 matrix multiplications extremely quickly. If you wanted to get even more performance, you could SIMD / thread the update. That said, the slowest part is probably going to be the bandwidth of sending 10000 matrices over the PCIE bus.

But, again, this is all supposition. Try it out, and profile.

I was thinking something like this:

(I apologize for the DX functions. I'm not super familiar with the corresponding OpenGL functions. That said, OpenGL should be quite similar)

// Variable definitions

vector<matrix> lod0;
vector<matrix> lod1;
vector<matrix> lod2;

ID3D11Buffer lod0VertexBuffer;
uint lod0VertexCount;

ID3D11Buffer lod1VertexBuffer;
uint lod1VertexCount;

ID3D11Buffer lod2VertexBuffer;
uint lod2VertexCount;

ID3D11Buffer instanceBuffer;

/// ..................

/// Per frame...
lod0.clear();
lod1.clear();
lod2.clear();

for (uint i = 0; i < bushes.size(); ++i) {
    float distance = distance(bushes[i].pos, camera.origin);
    if (distance < 50.0f) {
        lod0.push_back(bushes[i].matrix);
    } else if (distance < 100.0f) {
        lod1.push_back(bushes[i].matrix);
    } else {
        lod2.push_back(bushes[i].matrix);
    }
}

D3D11_MAPPED_SUBRESOURCE mappedResource;
ID3D11Buffer vertexBuffers[] = {
    lod0VertexBuffer,
    instanceBuffer
};

// Update lod0 matrices
context->Map(instanceBuffer, 0, D3D11_MAP_WRITE_DISCARD, 0, &mappedResource);
memcpy(mappedResource.pData, &lod0[0], sizeof(matrix) * lod0.size());
context->Unmap(instanceBuffer, 0);

// Draw lod0
context->IASetVertexBuffers(0, 2, vertexBuffers, nullptr, nullptr);
context->DrawInstanced(lod0VertexCount, lod0.size(), 0, 0);

// Update lod1 matrices
context->Map(instanceBuffer, 0, D3D11_MAP_WRITE_DISCARD, 0, &mappedResource);
memcpy(mappedResource.pData, &lod1[0], sizeof(matrix) * lod1.size());
context->Unmap(instanceBuffer, 0);

// Draw lod1
vertexBuffers[0] = lod1VertexBuffer;
context->IASetVertexBuffers(0, 2, vertexBuffers, nullptr, nullptr);
context->DrawInstanced(lod1VertexCount, lod1.size(), 0, 0);

// Update lod2 matrices
context->Map(instanceBuffer, 0, D3D11_MAP_WRITE_DISCARD, 0, &mappedResource);
memcpy(mappedResource.pData, &lod2[0], sizeof(matrix) * lod1.size());
context->Unmap(instanceBuffer, 0);

// Draw lod2
vertexBuffers[0] = lod2VertexBuffer;
context->IASetVertexBuffers(0, 2, vertexBuffers, nullptr, nullptr);
context->DrawInstanced(lod1VertexCount, lod2.size(), 0, 0);

If this ends up being too slow, you could try to do this calculation in a compute shader, then draw the result using an indirect draw call. See AZDO.


Is there an efficient way to achieve color blending in compute shaders?

If there is no such function, why have noone made an extension for this yet?

Just because something might be "useful" doesn't mean it would be possible or cheap.
GPU hardware achieves blending by using an entirely separate piece of hardware, one that preserves the order of the fragment values generated by unordered fragment shader executions. These execution units, commonly called "ROPs", are distinct hardware from fragment shaders.
CS's can't just do that. They have a fundamentally unordered execution model, and what you're talking about requires at least some ordering between them. That's not simple or cheap to implement, even if it's possible at all.
The closest you can get to that is fragment_shader_interlock, which can force parts of FS execution to be ordered by primitive order and to execute critical sections that prevent overlapping executions. But as the name suggests, that's a feature of fragment shaders (and not one universally supported), not a feature of compute shaders.
Also, it doesn't really make sense to use a compute shader for post-processing effects, since FS's are perfectly capable of dealing with that sort of thing, and actually have access to ROPs.


Estimating the position of vertexes in a 3D model

If you can construct a mesh that is representative of the tongue, like a spring mass system, then you can attach 'hooks' at your two control points and get a close approximation of the deformation of the tongue.


Forward ray-tracing means that rays starts from the light source and are traced towards the eye (light-surface-eye), this is a highly inefficient but accurate method. Reverse raytracing is simply starting in the eye and then trace them towards the light source. These are two classes of ray tracers.

  Ray casting is the use of ray-surface intersection tests to solve a variety of problems in computer graphics and computational geometry. The distinction is that ray casting is a rendering algorithm that never recursively traces secondary rays, whereas other ray tracing-based rendering algorithms may do so.

Ray-marching is an implementation of a ray tracer, often used to visualize volumes. Normally ray casting uses a formulae to determine the point of intersection with the surface. In ray marching such a function is not required, instead the intersection is evaluated at several points along the ray. Either these points have a fixed distance between them or every iteration the distance to the next point is based on the distance to the closest surface.

The difference between those are that ray-tracing is a concept to solve the rendering equation, in the end resulting in a certain color for each pixel. Forward and reverse ray-tracing are classes of ray-tracers, ray-casting is a simple algorithm to do ray-tracing and ray-marching is another algorithm to do ray-tracing. As by the comment many other algorithms have been developed to solve the problem with different advantages and disadvantages.


UV coordinates are defined such that (0, 0) is your lower left (or upper left, depending on what you're working with) corner of the image and (1, 1) the opposite corner (technically (1, 1) is outside the UV range, but we ignore this for easiness sake). Thus, if you want the exact middle, you obviously get (0.5, 0.5)

Now with your image, the first pixel lower left is the pixel at (0, 0) as well. The pixel in the opposite corner is (4095, 2047).
To get to the (1, 1) corner, you need to take that corner pixel and divide it's coordinates by the image's width and height respectively. here you will see, that the (1, 1) does not fit exactly:

$u = 4095 / 4096 = 0.99976$
$v = 2047 / 2048 = 0.99951$

The same principle applies to any other pixel of your image, you just have to divide by the image's width in pixels.
Your center pixel is at (2047, 1023) (or the image borders on this pixel and three of its neighbours). The resulting UVs are:

$u = 2047 / 4096 = 0.49975$
$v = 1023 / 2048 = 0.49951$


How to calculate the RGB intensity values at center of face using flat shading?

You shouldn't solve for $\cos^{-1}(\theta)$.

The ambient component is as follows (we are assuming the ambient light is white):
$$ambient\ color = (ambient\ light\ color)(ambient\ intensity)(diffuse\ color) = (1,1,1)(0.1)(0.6,0.7,0.8) = (0.06,0.07,0.08)$$

The diffuse component is as follows:
$$diffuse\ color = (light\ color)(light\ intensity)(diffuse\ color)(\cos(\theta)) = (1,1,1)(1)(0.6,0.7,0.8)(0.8) = (0.48,0.56,0.64)$$

Now we just need to add them up:
$$final\ intensity = (diffuse\ color) + (ambient\ color) = (0.06,0.07,0.08) + (0.48,0.56,0.64) = (0.54,0.63,0.72)$$

For more information, here's a pretty decent summary of various lighting calculations.


Combination of Triangles and Quads in a mesh

Topologically correct is a very vague.

I believe you think of Polygon-Meshes when you say meshes, which is to represent the surface of the object by a patch of polygons, in most cases triangles or quads. There are other options to model meshes, one idea would be to use signed distance functions as for example in constructive solid geometry.

So back to the question. Is it correct to mix the type of polygon? Yes. There is, mathematically and theoretically nothing stopping you from doing it.

The only real 'correctness' criteria that comes to my mind is the one of the mesh being a manifold. This is in (overly) simplficated the property that the mesh has a defined outside and inside and no strange places where you can't decide. (This happens for example if the mesh intersects with itself, but it get's freaky). And having different types of primitive does not violate the manifoldness of the mesh.

But most algorithms can actually deal with non-manifolds. There is no problem in rasterizing or raycasting a non-manifold as these algorithms process each primitive separately.

The more interesting question here is: Why DONT we mix the type of primitive?

And here the answer is very simple: Efficiency. You want to design an algorithm that handles your mesh as efficient as possible, and so you assume it's only one type of primitive. And at some point in time, because triangle-meshes are really handy, hardware and hardware-api began to adpot mainly to processing triangles. And now they are good at that and even other primitive types get split into triangles since splitting and rendering triangles is still faster than rendering the polygons directly.

So in short: You can mix polygon types, but if you want to be fast, use only triangles since they are simple and widely used.


How does vsync affect fps exactly when not at full vsync fps?

It depends on how a missed frame is handled by the driver.

One option is to just wait until the next vsync, causing a hitch of 32 ms and if the application is just at the limit of 16 ms can cause fluctuations.

The next option is to queue the frame for display next frame but don't wait on it. This will still cause a visual hitch but the application can then immediately start on the next frame instead of being forced to wait 16 ms. If the next frame is done faster then the late frame may even never be shown.

The final option is to not block and push the late frame to the display immediately which can cause tearing.


How to translate the center of an equirectangular projection?

An equirectangular projection treats the x coordinate as the angle theta around the vertical axis going from 0 to 360 degrees. These angles match the longitudinal angles of a globe. The y coordinate is treated as the angle phi that represents the latitudes of a globe. They typically go from +90° at the zenith (or north pole) to -90° at the nadir (or south pole).

So to translate the center of projection of an image by some amount (x0, y0), you  can simply add x0 to the x value of every coordinate (wrapping around to the other side at the edges). The y coordinate is similar - add y0 to every coordinate.  However, instead of wrapping around, you need to reflect the value. So if your new coordinate is less than 0, simply take the absolute value. If it's greater than the height, take the difference between the new value and the height of the image and subtract that from the height of the image. If it was less than 0 or greater than height then you need to also add half the width to the x coordinate.

So putting it in pseudocode:

newY = oldY + translateY;
if (newY < 0)
{
    newY = fabs(newY);
    translateX = fmod((translateX + (width / 2.0)), width);
}
else if (newY >= height)
{
    newY = height - (newY - height);
    translateX = fmod((translateX + (width / 2.0)), width);
}
newX = oldX + translateX;
if (newX < 0)
{
    newX = width + newX;
}
else if (newX >= width)
{
    newX = newX - width;
}


How to properly direct light sample GGX in path tracing

This may simply be a typo or oversight. Note that in the previous section they've been discussing the microsurface BxDFs $f^m_r(\mathbf{i},\mathbf{o},\mathbf{m})$ and $f^m_t(\mathbf{i},\mathbf{o},\mathbf{m})$, which are naturally evaluated at a microsurface normal, but in this section they've moved on to talking about macrosurface BxDFs, which have the microsurface integrated out.

AFAICT, the intent of equation 19 is just to say that the total BSDF is a sum of reflection and transmission parts, so it doesn't really matter ($\mathbf{m}$ is a free variable there and could be named anything), but to be less confusing it should probably have been $\mathbf{n}$. Maybe they just forgot to update it during a round of edits, or something.

Yes, you can evaluate a BRDF/BSDF with any combination of input and output directions and it should give the correct answer.


Why don't you just use Perlin noise twice on the same grid, or volume? Each with slightly different parameters (a phase shift, or different pseudo-random vectors). In this case both component of your float2 are smoothly defined by a Perlin Noise field. A float2 $v$ could be defined as $v = \{P(u,v), P(u+0.5, v+0.5)\},$ where $P(u,v)$ is the Perlin noise function. However, this will make the $x$ and $y$ component of $v$ very close depending on the scaling of the grid values.


How can I reproduce the MatchMoving/MotionTracking technique from Just Dance?

The Kinect SDK has a skeletal tracking demo in the c++ api (BodyBasics). The skeleton you get is the same "topology" unless it fails to acquire the skeleton. So you can easily compare the left shoulder of two captures as you always know which point it is.

You'll have to experience with metrics. A simple point to point metric would probably be sufficient though. Maybe with some weights on different body parts though.


Downscaling by skipping rows and columns – any example images?

A basic nearest neighbour scale function would work like so. Here I also use 24:8 fixed point arithmetic to make it more micro friendly.

const int originalWidth = 320, originalHeight = 240;
int targetWidth = 50, targetHeight = 50;

char *pSource = sourceBuffer;
char *pDest = destBuffer;    // you can use RGB instead of char

int yStep = (originalWidth << 8) / targetWidth;
int xStep = (originalHeight << 8) / targetHeight;

for(int y = 0, sampleY = 0; y < targetHeight; y++, sampleY += yStep) {
    // Point to start of scaled Y source line
    char *pSourceLine = pSource + ((sampleY >> 8) * originalWidth);

    // Output the line, scaling only in X (we already know scaled-Y)
    for(int x = 0, sampleX = 0; x < targetWidth; x++, sampleX += xStep)
        *pDest++ = pSourceLine[sampleX >> 8];
}

This will also work to scale image up as well as down, it should work on any input vs output size.


How can I offset/shrink a triangular polygon in GLSL?

Given triangle ▲ABC, we bisect angle ∠BAC with line AD, derived with Angle Bisector Theorem:

BA / BD = CA / CD

Point E represents our objective refined position on the desired resulting inset triangle. As it lies upon angle bisector AD, it is equidistant from sides BA & CA, forming identical right triangles ▲AFE & ▲AGE. We can now use Sine for Right Triangles to find the length of AE:

AE = EG / Sin(∠EAG)

That's all the math we need, so let's cook up some GLSL!

We start off with all the typical attributes: position, normal, and transformation matrices, but since the vertex shader only works on a single vertex, we need to add the neighboring vertices as additional attributes. This way, each vertex will find its own "point E", creating the resulting inset triangle. (Note: I don't call them "B" & "C" here, because they are not yet in screen space.)

    attribute vec3 left; //vertex to the left of this vertex
    attribute vec3 right; //vertex to the right of this vertex

Speaking of screen space, I'm also including the aspect ratio of the display, (and making it a uniform, in case the window is resized.)

After preparing varying normals for the fragment shader, and transforming the face into clipping-space, we can get down to the business of applying the above math:

        attribute vec3 left; //vertex to the left of this vertex
        attribute vec3 right; //vertex to the right of this vertex
        uniform float aspect;
        varying vec3 vNormal;
        varying vec2 vUv;

        void main() {
            vNormal = normal;
            vUv = uv;

            mat4 xform= projectionMatrix * modelViewMatrix;
            vec4 A = xform * vec4( position, 1.0 );
            vec4 B = xform * vec4( left, 1.0 );
            vec4 C = xform * vec4( right, 1.0 );

            vec3 CB = C.xyz - B.xyz;
            vec2 BA = B.xy - A.xy;
            vec2 CA = C.xy - A.xy;
            float lengthBA = length(BA);
            float lengthCA = length(CA);
            float ratio = lengthBA / ( lengthBA + lengthCA );
            vec3 D = B.xyz + ratio * CB.xyz;
            vec3 AD = D - A.xyz;
            vec3 bisect = normalize(AD);

            float theta = acos( dot(BA, CA) / (lengthBA * lengthCA) ) / 2.0;
            float AE = 1.0/(sin(theta)*aspect);
            newPos.z += AE/length(AD) * (D.z - A.z);
            newPos.x += bisect.x*AE;
            newPos.y += bisect.y*AE;

            gl_Position = newPos;
        }

This code gives us the results below.

Note, there are a few edge-cases having to do with almost backface-culled triangles being flipped by this process, and I began to address this in code, yet decided to simply avoid these cases for now. Perhaps I'll revisit it when I get this project done.


Are the terms screen and the monitor in computer graphics used interchangeably?

The terms screen, monitor, and display are used more or less interchangeably to refer to the real-world physical device that shows the image you're looking at. (They're not 100% interchangeable: monitor means a discrete display device such as you'd use with a desktop PC. However, a built-in display such as the one in a laptop, tablet, or phone is never called a monitor.)

The image plane means the abstract geometric plane that identifies the surface of the display, but extends infinitely in all directions, i.e. isn't bounded by the edges of the display. The viewport or viewing window is the rectangle within this plane that identifies the edges of the display (or more generally, the edges of the image you're rendering, in case your app is running in a window and not fullscreen).

The term screen might also be used more colloquially to mean various other things, such as: various coordinate systems on the image plane ("screen space"), the viewport or view frustum ("the object is off-screen"), the back buffer ("render it to the screen"), and probably more.


Phong shading in OpenGL: line of light from the center of the world

I forgot to set the wEye uniform, so it was always (0, 0, 0).


GL_DEPTH_CLAMP has no effect (reverse z)

It is compatible, but instead of clipping it now clamps (hence the name). If the depth of your primitives is outside the [0,1] range it ends up behaving as if you have no depth-test at all. The only way to mitigate that is to either scale down the depth range so it fits in [0,1] or to use GL_NV_depth_buffer_float to disable clipping and clamping completely.


How can I perform a triangle inside test in polygon meshes?

There is an alternative method that relies on flood filling. First arrange your edge data into a loop where the edges are forming a counterclockwise loop. Then start at an arbitrary point on the loop and pick edges joining that point. Use the outbound boundary edge and cross it with the other outbound edge, if it points in the direction of the face normal then it's an edge to be included, if not discard it. From this edge continue until you hit a boundary edge, at which point you terminate the fill. Continue at a yet to be visited boundary edge vertex.


What is the "dipole approximation" for subsurface scattering?

The assumption underlying such model is the same as lots of other models for skin rendering; the subsurface scattering can be approximated as a diffusion phenomenon. This is good because in highly scattering media, the distribution of light loses dependency from the angle and tends to isotropy.

The dipole approximation is a formulation for the resolution of such diffusion problem in an analytical fashion.

Basically they start by approximating the BSSRDF as a multiple scattering and single scattering component. The multiple scattering is then defined as:

Where $F_t$ are Fresnel terms and $R$ is the diffusion profile expressed as function of the distance between the entry and exit point.

This $R$ is referred to as diffusion profile and they formulate this profile via a dipole approximation. The contribution of the incoming light ray is considered to be the one of two virtual sources: one negative beneath the surface and one positive above it (that's why dipole)

Here in the picture r is the $\|x_i - x_o\|$ above.  The contribution of those light sources is dependent on various factors such as the distance of the light from the surface, scattering coefficient etc. (See below for a more detailed description of the formula itself).

This model only account for multiple scattering events, but that's fine enough for skin. It must be noticed though that for some translucent materials (e.g. smoke and marble) the single scattering is fundamental. That paper propose a single scattering formulation, but is expensive.

The diffusion profile is usually approximated for real-time application as a series of gaussian blurs (like in the seminal works of D'Eon et al. in GPU Gems 3 then used for the Jimenez's SSSSS) so to make it practical for real time scenarios. In this wonderful paper there are details on such approximation.
A picture from that paper show actually how good is this formulation:

As side note the dipole approximation assumes that the material is semi-infinite, however this assumption doesn’t hold with thin slabs and multi-layered material such as the skin. Building on the dipole work, Donner and Jensen [2005] proposed the multi-pole approximation that accounts for the dipole problems.
With this model instead of a single dipole, the authors use a set of them to describe the scattering phenomenon. In such formulation the reflectance and transmittance profiles can be obtained by summing up the contribution of the different dipoles involved

EDIT:
I am putting here the answers to a couple of @NathanReed 's questions in the comment section:

  Even with the diffusion profile approximation, the BSSRDF model still requires integrating over a radius of nearby points on the surface to gather incoming light, correct? How is that accomplished in, say, a path tracer? Do you have to build some data structure so you can sample points on the surface nearby a given point?

The BSSRDF approximation still need to be integrated over a certain area, yes.

In the paper linked they used a Montecarlo ray-tracer randomly sampling around a point with a density defined as:

$\sigma_{tr}e^{-\sigma_{tr}d}$

Where that sigma value is the effective extinction coefficient defined below ( it is dependent on scattering and absorbing coefficient, which are properties of the material) and d is the distance to the point we are evaluating. This density is this way defined because the diffusion term has exponential fall-off.

In [Jensen and Buhler 2002] they proposed an acceleration technique. One of the main concepts was to decouple the sampling from the evaluation of the diffusion term. This way they perform a hierarchical evaluation of the information computed during the sampling phase to cluster together distant samples when it comes to evaluating the diffusion. The implementation described in the paper uses an octree as structure.
This technique, according to the paper, is order of magnitude faster than the full Monte Carlo integration.
Unfortunately I never got myself into an off-line implementation, so I can't help more than this.

In the real-time sum-of-Gaussians approximations the correct radius is implicitly set when defining the variance of the Gaussian blurs that need to be applied.

  Why one positive and one negative light? Is the goal for them to cancel each other in some way?

Yes, the dipole source method (which dates way before Jensen's paper) is such defined to satisfy some boundary condition. Specifically the fluence must be zero at a certain extrapolated boundary that has a distance from the surface of $2AD$ where

Being $F_{dr}$ the fresnel reflectivity of the slab considered and that sigma value is the reduced extinction coefficient described below.

EDIT2: I have expanded (a tiny bit) some of the concepts in this answer in a blog post: http://bit.ly/1Q82rqT

For those who are not scared by lots of greek letters in a formula, here's an extract from my thesis where the reflectance profile is briefly described in each term:


Modulo vs Scaling vs Capping when performing filters

When doing a filter operation such as a blur, in many cases the filter itself will be normalized so its values sum to 1.0, precisely to avoid this problem.

In cases where out-of-range outputs do need to be handled, the most common way is what you called capping. In graphics it's more commonly known as clamping or saturating. This is the default behavior if, for instance, you run a pixel shader on a GPU and output out-of-range values.

Renormalizing the output image might be necessary in some special cases depending on what you're doing. I can't think of anytime that I've seen the modulo behavior used—it would be a pretty weird operation to do on a visual image.


What does it *mean* to scale in an arbitrary direction?

When you scale along the X-axis, the X-coordinate (parallel to the axis) gets stretched, while the Y-coordinate (perpendicular to the axis) remains the same. You can think of scaling along an arbitrary axis as stretching along some diagonal.

Here's a pic of a square being scaled along the main diagonal (the axis pointing to <1, 1> ) by factors of 2 and 0.5.

You can see that the distance along the axis (green line) gets scaled while the perpendicular distance (red line) stays the same.

Mathematically, this is the same as rotating until your chosen axis lies along the X-axis, scaling by X, then rotating back again.


Relationship between spatial frequency of an image and Aliasing

Like in any other kind of signal processing, the relationship is Nyquist's theorem. An image is a discrete sequence of samples of a continuous signal. If the original signal has frequency components higher than half the sampling rate, then there will be aliasing. To put that another way, if you look at the real-world size of a pixel, any details smaller than two pixels wide will be aliased. This applies to synthetic images as well as to photographs. If you define a procedural texture via a mathematical function, the function is continuous but you only evaluate it at certain points.

The problem is, of course, if you only have the sampled image and not the original signal, the high-frequency components have already been aliased, so you can't measure them directly. You have to use statistical techniques to guess which details in your aliased image were originally higher-frequency signals.


There is a small issue in the original code. Indeed, real array values ranging between [0, 255]. However, the mask values ranging between [0, 3*3] (as dn = 3 in the example above).

After fixing the issue:

The correct code (with image reading and all data normalized between [0, 1.0]):

import numpy as np
from numpy.random import default_rng
import matplotlib.pyplot as plt
from PIL import Image

# Read an image from the disk
real = Image.open("example.jpeg")
real = np.array(real) / 255.0
# If it is a RGB image, transform it to a gray image
# TODO: Could use any function (like luminance)
real = (real[:,:,0] + real[:,:,1] + real[:,:,2]) / 3.0

# Create the mask by shuffling
dn = 6
rng = default_rng()
# In the original code, there was a problem here.
# Indeed, the mask values ranging between [0, dn*dn]
# So the generated mask needs to be properly normalized
dither = rng.choice(dn*dn, size=(dn,dn), replace=False) / float(dn*dn)
print(dither)

# Generate the results
res = np.zeros_like(real)
width, height = res.shape[1], res.shape[0]
for y in range (0,width):
    for x in range (0,height):
        i = x % dn
        j = y % dn
        res[x][y] = 1.0 if (real[x][y] > dither[i][j]) else 0.0

# Display the result
fig = plt.figure()
ax1 = plt.subplot(1, 3, 1)
ax1.title.set_text('Real Image(300x300)')
plt.imshow(real, cmap="gray")
ax2 = plt.subplot(1, 3, 2)
ax2.title.set_text('Dithering Matrix(3x3)')
plt.imshow(dither, cmap="gray")
ax3 = plt.subplot(1, 3, 3)
ax3.title.set_text('Resultant Image(300x300)')
plt.imshow(res, cmap="gray")
plt.show()


Sphere intersection occlusion (for hybrid raytracing)

Given that I didn't miss anything, you can probably cut this down to a problem in the 2D space. Viewing onto the plane defined by the center points of the spheres and your camera origin, the scene looks like this:

The spheres become circles with the center points $C_1$ and $C_2$, and the intersection circle is now only 2 points with only the closer one $P$ being interesting. The camera/eye is arbitrarily set to the point $E$.

Calculating if one point on the spheres is visible or not is easy: Simply check whether or not the angles at point $P$ between $E$ and $C_1$ respectively $E$ and $C_2$ are both greater (or equal to) 90 degree1.

If $P$ is visible, some part (e.g. at least that point) of the intersection circle is visible. Otherwise the whole intersection circle must be occluded by one of your spheres, namely the one which creates an angle of less than 90 degree.

Here is how it looks if $P$ is not visible from $E$:

You can clearly see how that point is occluded by the circle around $C_2$ and that the angle between $E$ and $C_2$ in $P$ is less than 90 degree.

1 Having an angle of exactly 90 degree means that the line between $E$ and $P$ just touches the respective circle/sphere in point $P$ as a tangent.


It almost certainly depends on implementation. I have used CPU only keyers that had no problems keying full hd in real time.

But back to the question would a shader be faster? Most likely, especially if your chroma key is a classic color transform. Then the GPU would shine as it would get to use it massive multicore architecture to full potential.


There are multiple different ways to do tone-mapping, including both local and global methods. I won't discuss local methods (i.e. methods based on using surrounding information as part of the tone mapping process) here. Global methods can be broken down into two components:

A formula for mapping from $[0,\infty)$ to $[0,1]$.
A method of applying this to the color.

The simplest formula out there is just $out = in / (in + 1)$, but the results are also fairly boring. A fancier curve is the filmic tone mapping function developed by John Hable, which is popular in the video gaming industry:

float hable(float x)
{
    float A = 0.15;
    float B = 0.50;
    float C = 0.10;
    float D = 0.20;
    float E = 0.02;
    float F = 0.30;

    return ((x*(A*x+C*B)+D*E)/(x*(A*x+B)+D*F))-E/F;
}

You can also get even fancier if you can constrain the input domain to a known, fixed peak value.

As for how to apply this function to an actual color value, there are again multiple techniques. A very simple way is to do it per-channel, i.e. color.r = map(color.r); color.g = map(color.g); color.b = map(color.b);. But this has the downside of warping the appearance of bright colors towards white (e.g. orange becomes yellowish under tone-mapping). Some in the movie industry tend to like this look, for whatever reason.

To do color-preserving tone mapping, you can perform the calculation on a single value only and then apply the same scale to all three channels. A good choice for the value to use is to choose the “brightest” of the three components:

float sig = max(color.r, max(color.g, color.b));
color.rgb *= vec3(map(sig) / sig);

However, this approach also has downsides, namely that it takes very bright objects to unnaturally saturated versions of
themselves. For example, a bright sun might suddenly appear various shades of yellow and orange - which also looks unnatural on photographic content.

So ultimately, what I do is first desaturate the color based on how bright it is, and then apply tone mapping after the fact, leading to an algorithm that looks something like this:

// Calculate the desaturation coefficient based on the brightness
float sig = max(color.r, max(color.g, color.b));
float luma = dot(color.rgb, vec3(0.2126, 0.7152, 0.0722));
float coeff = max(sig - 0.18, 1e-6) / max(sig, 1e-6);
coeff = pow(coeff, 20.0);

// Update the original color and signal
color.rgb = mix(color.rgb, vec3(luma), coeff);
sig = mix(sig, luma, coeff);

// Perform tone-mapping
color.rgb *= vec3(map(sig) / sig);

But YMMV.

References:

https://github.com/haasn/libplacebo/blob/master/src/shaders/colorspace.c#L664
You can also poke at ITU-R Report BT.2390-2 for some more ideas and guidelines.


how to draw mesh equal distance from irregular shape geometry bounding box

One way to visualize an object's bounding box would be like this:

$\hskip2.5in$

As Rahul mentioned in a comment, you can achieve this by adding a total of 12 cylinders at the edges of the bounding box, connecting the spheres at the corners.

If you want to render the bounding box transparently, I don't see an easy way to hide the intersections of the spheres and cylinders without knowing what program you are using. You could try to use partial spheres for the corners and see if that helps, i.e.

$\hskip2.5in$

Maybe this is already sufficient for your purpose. (This is exactly 1/8th of a full sphere.)


Should ideal specular multiply light colour with material colour?

What "ideal" means in this context is that there is no divergence in the direction of light reflection vectors (i.e. no roughness) but that they are all considered to be perfect reflections from an optical flat surface. However, even for optical flat surfaces Fresnel equation is still applied in BRDF evaluation that changes specular reflection and is potentially wavelength dependent, which is a common characteristic with metals.

However, in smallPt there's no Fresnel approximation (e.g. Schlick's approximation) but the surface reflectance is simply multiplied with constant term which is wrong. Also 0.999 specular albedo is very high for any real material. Below is a table of specular reflectange of various metals & dielectrics at normal incidence for reference:


Can pixels on the image plane be mapped to world space coordinates?

What you are asking about is 'mouse picking'. If a user clicks a pixel on the screen how can you get corresponding world space coordinates of an object being clicked on.
Imagine a ray from that pixel, through the view volume all the way to the far plane. That ray is a set of all possible points the user is clicking on. So we want to perform ray/poly ntersections with every object in the scene, and return the first intersect.

To transform the pixel into world space, you first have to turn your x,y point in screen space into NDC coordinates (-1 to 1) use -1 for your Z component and 1 for W to create a 4 float component vector and then multiply that vector by the inverse projection and then inverse view matrices.
There are many examples of this process if you want to see code, but this is the basic process.
Now we have a ray from the camera/eye position into the scene, we want to do intersection tests with our objects. There are many ways of achieving this, take your pick.

Hope this helps, good luck!


Why is spherical harmonics used in low frequency graphics data instead of a spherical Fourier transform?

Spherical harmonics really are the "spherical Fourier transform" you're looking for. The kind of hack you mention in comments, of doing a 2D Fourier transform on a lat-long projection, suffers from all the problems you usually have when you try to project a sphere onto a plane: not all spatial relations in the sphere are well-represented in the plane. If you take stretching at the poles as an example, that clearly results in a frequency shift. High-frequency components in the sphere get shifted down to low frequencies in the plane, and when you do the inverse transform, any inaccuracy shows up as high-frequency components at the pole.

Such a transform is useless for analysing the spectrum of the original signal, because it distorts the frequencies (and distorts them differently in different parts of the sphere). It's also poor for representation, because of the introduction of these high-frequency artefacts.

Spherical harmonics use different basis functions so they can work directly on the lat-long without treating it like a plane. If you plot the response of the spherical harmonics on a plane projection of the sphere, they look kooky (like a map of the earth does), but they look much more natural on a sphere. They don't suffer from artefacts at the poles, and they don't frequency-shift signals. They're also rotation-invariant, so it doesn't even matter where you put the poles. They don't even have poles, in the sense of discontinuities, which is a property you can't achieve from a plane projection however many ad hoc co-ordinate transforms you do. Rotation-invariance is useful for applications like environment maps, where you'd like to rotate the spherical-harmonic transform of the environment (light probe, or whatever) into the local co-ordinate system before feeding it to the shader. There isn't a projection onto the plane that has this property.


Your ambient lighting contribution is missing from the second one (:


Math for adaptive supersampling for raytracing

I found a paper from 1997 that addresses this subject: http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.29.6799


For light sources with larger solid angle and where the shadow caster is relatively closer to the light than the receiver, you get notable soft shadowing effect. So if you render larger light sources closer to the shadow receiver it's important to handle soft shadows properly for realistic lighting. Even with the Sun which has quite small solid angle, you can still see fairly large penumbras from tall buildings.

PCSS is fairly simple algorithm for implementing contact hardening shadows, but it has some notable issues. PCSS fails when you have soft and hard shadows intersecting (blocker search fails), so it's not a good algorithm to use when large penumbras are needed. See in the below image where shadows from the box & triangle intersect. This artifact is particularly disturbing when animated.

There's a family of "back-projection" shadow techniques though, which is able to handle larger penumbras. However, there's a little research done how to implement this efficiently for real-time applications that I'm aware of.

In back-projection the principle idea is to project shadow map texels to the light source for every pixel to calculate light occlusion, which when implemented naively is obviously extremely slow and furthermore introduces shadow gaps & overlaps. However, there are some techniques to try to avoid the artifacts, such as adjusting the microquad extents as shown below.

Another improvement to handle the gaps and overlaps is to use bitmask soft shadowing algorithm, where occlusion bitmask determines the occluded areas of the light source.

The most recent paper I have seen on the topic is "Real Time Area Lighting. Now and Next" by Sam Martin from 2012, where the cost was several milliseconds for a single light, but it's definitely worth the read if you are interested of the topic, to get general idea about the state-of-the-art real-time shadow techniques.


How to ensure a fragment shader only runs once

I don't think your approach is viable for a variety of reasons.

a GPU is heavily cache-dependent for performance. Pulling the same pixels over and over will just measure how fast a cache hit is. Sampling random pixels will only measure a cache miss. Neither are useful statistics.
each texture2d call will pull in multiple pixels to leverage that, in a fragment shader, neighbouring fragments will want neighbouring pixels.
The driver could optimize your shader to pull out the loop if it sees that each iteration does the exact same computation or if the computed result is not used.


Normal vector in Phong Illumination Model should be normalized?

Short answer: yes.

Longer answer: yes, because the vectors you’re using are meant to represent directions, not directions-and-distances. Think of it in terms of light: it doesn’t matter how far a photon’s traveled, whether it’s from the sun or from a lamp on your desk—once it arrives at a surface, it’s going to get reflected in exactly the same way.

From a mathematical standpoint, you’re making lighting calculations based on the angles between the incoming light, the view point, and the surface normal. The dot product is defined as the length of the two input vectors multiplied by each other and by the cosine of the angle between them. To get just the angle—which you need for the lighting model—you need to factor out the length of the two vectors, i.e. normalize them.


Procedural terrain: What's the best approach to calculate noise in the GPU?

r a similar problem (a tree of combined noise functions, evaluated on the GPU), I found a good method is to generate a shader from the expression tree.  Each predefined node corresponds to a single shader function, e.g.

float simplexNode(vec3 pos) {
  // ... implementation of simplex noise
}

or

float sumNode(float val1, float val2) {
  return val1 + val2;  // ... possibly clipped or whatever
}

You'd then create a shader like

#version 300 core
// ... include node definitions here
in vec3 vertexPos;
out float vertexHeight;
void main() {
  float f1 = simplexNode(vertexPos);
  float f2 = worleyNode(vertexPos);
  float f3 = sumNode(f1,f2);
  // .. a bunch more
  vertexHeight = f76;
}

Constant parameters could either be uniforms, or just pasted straight into the code.

One limitation of this method is that the noise functions and methods have to be local, i.e. the evaluation can be done pointwise.  In my case, some nodes were blends and convolutions, so I had to evaluate pieces of the expression tree into textures, blend or convolve the textures into new textures, then use these textures as leaves for the rest of the expression tree.  All of the combinators you mentioned seem to be local, so this may not be a problem for you.