Merge multiple Point Clouds into a single 3D object

CloudCompare has statistical outlier removal (called SOR) function that can eliminate all lonely points for you. I would filter each *.ply first and then try to carry out a cloud to cloud registration. For registration, you can use MeshLab or use 30 days trails of commercially available software such as Faro Scene. Scene won't accept *.ply file though, so you need to convert your files to *.e57 first. Good luck.


Why is the transposed inverse of the model view matrix used to transform the normal vectors?

Here's a simple proof that the inverse transpose is required. Suppose we have a plane, defined by a plane equation $n \cdot x + d = 0$, where $n$ is the normal. Now I want to transform this plane by some matrix $M$. In other words, I want to find a new plane equation $n' \cdot Mx + d' = 0$ that is satisfied for exactly the same $x$ values that satisfy the previous plane equation.

To do this, it suffices to set the two plane equations equal. (This gives up the ability to rescale the plane equations arbitrarily, but that's not important to the argument.) Then we can set $d' = d$ and subtract it out. What we have left is:

$$n' \cdot Mx = n \cdot x$$

I'll rewrite this with the dot products expressed in matrix notation (thinking of the vectors as 1-column matrices):

$${n'}^T Mx = n^T x$$

Now to satisfy this for all $x$, we must have:

$${n'}^T M = n^T$$

Now solving for $n'$ in terms of $n$,

$$\begin{aligned}{n'}^T &= n^T M^{-1} \\
n' &= (n^T M^{-1})^T\\
n' &= (M^{-1})^T n\end{aligned}$$

Presto! If points $x$ are transformed by a matrix $M$, then plane normals must transform by the inverse transpose of $M$ in order to preserve the plane equation.

This is basically a property of the dot product. In order for the dot product to remain invariant when a transformation is applied, the two vectors being dotted have to transform in corresponding but different ways.

Mathematically, this can be described by saying that the normal vector isn't an ordinary vector, but a thing called a covector (aka covariant vector, dual vector, or linear form). A covector is basically defined as "a thing that can be dotted with a vector to produce an invariant scalar". In order to achieve that, it has to transform using the inverse transpose of whatever matrix is operating on ordinary vectors. This holds in any number of dimensions.

Note that in 3D specifically, a bivector is similar to a covector. They're not quite the same since they have different units: a covector has units of inverse length while a bivector has units of length squared (area), so they behave differently under scaling. However, they do transform the same way with respect to their orientation, which is what matters for normals. We usually don't care about the magnitude of a normal (we always normalize them to unit length anyway), so we usually don't need to worry about the difference between a bivector and a covector.


Can i compile my shaders to SPIR-V without using glslangValidator?

Typically compiling GLSL to SPIRV would be something you do as part of the build process, so that you end up distributing SPIRV instead of GLSL.  Then you don't need to compile it at runtime.  Assuming you're using the Vulkan SDK, the binary for glslangValidator (or a shell script to build it) should be part of the installed SDK.

That said, if you want to compile GLSL to SPIRV at runtime, you can link against the glslang library and do that.  Take a look at the Programatic Interface section in the readme.


How can I implement next event estimation?

You can provide an intersect_any(ray, scene, min_t, max_t) function, that returns true at the first intersection. So it has the option to exit early compared to the usual intersect function. Additionally you do not care about normals, texture coordinates, positions, so you can avoid computing those in many cases. You need max_t though to limit the ray to the segment between the two points you are testing.

In the worst case scenario (there is nothing between the two points) you will still need to iterate over all objects unless you are using acceleration structures.


The interpolation looks fine. The main problem here is that the hash function you're using isn't very good. If I look at just one octave, and visualize the hash result by outputting hash(PT).x, I get something like this:

This is supposed to be completely random per grid square, but you can see that it has a lot of diagonal line patterns in it (it almost looks like a checkerboard), so it's not a very random hash, and those patterns will show up in the noise produced by it.

The other problem is that your hash only returns gradient vectors in [0, 1], while they should be in [−1, 1] to get gradients in all directions. That part's easy to fix by remapping.

To fix those problems, I switched the code to use this hash function (which I learned from Mikkel Gjoel, and is probably due to a paper by W.J.J. Rey):

vec2 hash(vec2 co) {
    float m = dot(co, vec2(12.9898, 78.233));
    return fract(vec2(sin(m),cos(m))* 43758.5453) * 2. - 1.;
}

Note that due to the trig functions it's going to be a bit more expensive than your version. However, it considerably improves the appearance of the resulting noise:


Short answer: Yes, but it won't be pretty.

Explaination: While there are algorithms to automagically UV map polygon soup, the mapping probably wont be ideal. UV mapping is an art really. Choosing where to hide the seams and where to put them to limit stretching and optimize texture space.

That said, there are lots of tools that can make manual UV mapping a heck of a lot easier. For example: http://knowledge.autodesk.com/support/maya-lt/learn-explore/caas/CloudHelp/cloudhelp/2015/ENU/MayaLT/files/Create-UVs--Automatic-Mapping-htm.html


Your attribute locations could be wrong, because you're not querying them:

glVertexAttribPointer(0, 3, GL_FLOAT, false, 7*4, 0);
glVertexAttribPointer(3, 4, GL_FLOAT, false, 7*4, 3*4);

The first parameter should be queried using glGetAttribLocation:

glVertexAttribPointer(glGetAttribLocation( program, "position" ), 3, GL_FLOAT, false, 7*4, 0);
glVertexAttribPointer(glGetAttribLocation( program, "i_color" ), 4, GL_FLOAT, false, 7*4, 3*4);

If using more modern OpenGL, the locations can be set manually in the shader, eliminating the need for glGetAttribLocation:

layout(location=0) in vec3 position;
layout(location=3) in vec4 i_color;


Do you have to group and optimize your draw and shader calls in Vulkan?

Vulkan will not have the CPU-time performance costs associated with state changes (validation overhead, etc). But Vulkan command buffers do have state. They have two very major pieces of state: descriptor sets and the pipeline object.

Both of these objects are monolithic, dealing with potentially lots of pieces of information. The (potentially) huge quantity of stuff these objects control should strongly suggest to you that binding new objects is not cheap. And that's in terms of GPU execution speed, not CPU validation or whatever else.

Vulkan in general tries to make slow things in the API look slow and fast things look fast. Render passes are big, huge and complicated, so clearly you shouldn't just be starting and stopping them whenever you feel like it. Pipeline objects control a huge amount of information, so binding new ones feels like a heavyweight operation.

By contrast, changing what buffers get used for vertex input feels like a short, quick, focused operation. So it's probably cheaper than changing the pipeline that uses those buffers.

So yes, you should avoid state changes in Vulkan, just as much if not moreso than in OpenGL (since many state changes in Vulkan require entirely new pipelines).


Only show opaque objects behind transparent objects

The simplest solution is to do three passes:

Render opaque meshes to a buffer (front-to-back, depth read/write on)
Render translucent meshes to another buffer (front-to-back, depth read/write on)
This makes sure that only the closest translucent mesh is rendered.
Alpha blend translucent mesh on top of opaque buffer using opaque depth buffer. (depth read on)
This will take care of opaque meshes blocking translucent meshes and preserving translucent opacity levels. Depth test for values less or equal for small early z test (potential) optimisation.

EDIT: As mentioned below you could use the same depth buffer for the translucent pass and do the alpha blend then, the key is rendering the translucent pass front-to-back and depth write/testing.


The most performant way to organize vertex data on modern GPUs

Apologies, I was going to add a comment to your question thread but found I ended up elaborating too much.  My development experience is from a DX 11 perspective so some of this may not hold in OpenGL

The memory locality of the data does certainly play a major part.  But there is a few other items that do play into this, the width of the data as you know.  I've had a few GPUs and have had a stark performance hit on some based upon certain tipping points.  That is the for example, on my old AMD r290 if I recall, you could basically pass up to 4 float4s through to the vertex buffer with no additional cost over 1, but once I posted > 4 then there was a measurable drop (in frame rate).   This is purely off my recollection and so is anecdotal as best.  But it holds that GPU architecture is changing all the time, and what niche techniques that give you a performance edge today maybe the bane of your performance tomorrow.  I've tried myself to stick to the advised known best on the assumption that hardware vendors will build their GPUs to the most common techniques of passing payloads via VB (obviously with some looks from constant buffers).

Having said that, you have asked the question about ordering the vertices, and certainly that will help marginally.  The performance gain is when you combine this with index buffers which then allows the hardware to optimise and cache already calculated vertices.    You certainly may get even more gains from Triangle strips etc which specifically ordered as per your thoughts.  Most of model rendering I do is indice/vertex based optimised models with instancing, I use a small amount of lookups for cyclic motion effects (tree branches for example), in these cases the entire branch of the tree is looking up the same value.  So caching can be taken advantage of here also.

All I can say in summary from myself is:

Try not to overthink optimisation at the start, while your still developing your game/application try not to prematurely optimise.  If you need to come back and add new capabilities, you might lose all the effort.  You have my habit maybe ><, love trying to get the best performance and techniques
GPU architectures are varied, and strengths of one can reflect weaknesses in another.  AMD and Nvidia are famous for getting devs to optimise games to their architecture, for a reason, each has advantages/disadvantages.  Taking the middle road on development might be the best place, don't use any hardware vendor based features (opinion here).  (add to that, some packed formats only exist in AMD).

These are just a few thoughts and experiences I have had.  There are a lots books out there that you should get your hands on around these very topics.  I've not seen to many prescribe what you're proposing, but that doesn't mean its wrong.  Good luck.


An effect like that is most likely done as a second pass, after storing a copy of the color buffer.

When looking carefully, you can see some of the vertices from the shape of the effect above the character. So this scene in particular probably does something like this:

First render the scene into a color buffer, without the distortion effect.
Because it's not possible to both read from and write to a same frame buffer at the same time, copy it as is (blit) into another frame buffer. From now on, new draws will go into that new frame buffer.
In the new frame buffer, draw a mesh with the shape of the effect, using a shader that will fetch the texture of the first color buffer, with an offset depending on the normal. Also note that:

The mesh is opaque (no blending): it's the texture that makes it look transparent.
Although we don't see this in the screenshot, it's also perfectly possible to use the depth buffer as usual and have the effect intersect with the scene.

Finally a simpler effect, like rain drops on the screen, could be done this way:

First render the scene into a frame buffer, without the distortion effect.
Copy the entire frame buffer into another frame buffer using a custom shader. The shader will apply some offset to the texture coordinates (either directly or maybe using a distortion texture), thus distorting the image.


Given a picture of a unit square at an arbitrary angle and position, how can I reconstruct the distance and orientation?

This problem is known as camera calibration. The goal is estimating the various parameters of a "camera":

The intrinsic parameters: Things like focal length, distortion of the lens, sensor parameters etc.
Extrinsic parameters: Orientation of position of the camera.

To robustly determine those from an image that contains labelled points you usually need more than one image. A typical way to do this is recording many images of a chessboard pattern. But as soon as you have determined all the intrinsic parameters, you could save those for future applications as they do not have to be estimated over and over again.

I recommend reading the introduction to camera calibration of openCV:

https://docs.opencv.org/3.1.0/dc/dbb/tutorial_py_calibration.html

https://docs.opencv.org/3.1.0/d9/db7/tutorial_py_table_of_contents_calib3d.html


Texture prefiltering for Linearly Transformed Cosines

I have ported their WebGL demo in DX shader model 4, but I have not implemented such texture preprocessing. Instead, I use some abstract pattern (like a comuflage) that is twice wider and higher. So, in a reflected area only internal part is rendered if roughness is low, and with a greater roughness, mip-map with lower resolution is taken. That spreads pattern in zone of reflection realistically.
It could be some starting point for your implementation.

It is possible to make some "texture expansion" like in article, but some special code to do this is needed. This is "special", it modifies Gaussian.

But in many cases, reflected pattern could be very simple, and such expansion can be eailly created by artist and mip chain generated. Again -good for start point.


Am I passing the correct matrices for my shader? openGL/C++

I wasn't placing the lightposition into view space before passing it in.

vec4 LightPosition_view=view*LightPosition; //put lp into view space with the rest

When it comes to these things, you have to think in terms of spaces, and I was multiplying the view matrix and the model matrix, into view space, but not the lightposition itself.

I hope that helps anyone else facing the same/similar problem.


Light coordinate in Phong shading and Normal interpolation

Typically, you calculate lighting in either world-space or view-space (or camera-space). The benefit to calculating it in camera-space is the increased precision, although this may involve introduce additional math operations in shaders depending on your setup. So it wouldn't be in screenspace, or you would need to reconstruct the positions (which could happen in deferred shading depending on your G-buffer layout).

I assume you are talking about a punctual light type (such as a point light). If so, if your light position equals your camera position, then the view and light angles are the same, which can be factored out in the BRDF that you use, so you save a dot product.


What is the purpose of this cube in color quantization?

The cube looks like just an RGB cube, i.e. X,Y,Z map to R,G,B. Since there are 3 color channels, the 3D cube is a convenient way to visualize the space of all possible RGB colors. However, I don't know what the pie-wedge-looking shape outlined inside the cube is.


Sharing code between multiple GLSL shaders

There's a bunch of a approaches, but none is perfect.

It's possible to share code by using glAttachShader to combine shaders, but this doesn't make it possible to share things like struct declarations or #define-d constants. It does work for sharing functions.

Some people like to use the array of strings passed to glShaderSource as a way to prepend common definitions before your code, but this has some disadvantages:

It's harder to control what needs to be included from within the shader (you need a separate system for this.)
It means the shader author cannot specify the GLSL #version, due to the following statement in the GLSL spec:

  The #version directive must occur in a shader before anything else, except for comments and white space.

Due to this statement, glShaderSource cannot be used to prepend text before the #version declarations. This means that the #version line needs to be included in your glShaderSource arguments, which means that your GLSL compiler interface needs to somehow be told what version of GLSL is expected to be used. Additionally, not specifying a #version will make the GLSL compiler default to using GLSL version 1.10. If you want to let shader authors specify the #version within the script in a standard way, then you need to somehow insert #include-s after the #version statement. This could be done by explicitly parsing the GLSL shader to find the #version string (if present) and make your inclusions after it, but having access to an #include directive might be preferable to control more easily when those inclusions need to be made. On the other hand, since GLSL ignores comments before the #version line, you could add metadata for includes within comments at the top of your file (yuck.)

The question now is: Is there a standard solution for #include, or do you need to roll your own preprocessor extension?

There is the GL_ARB_shading_language_include extension, but it has some drawbacks:

It is only supported by NVIDIA (http://delphigl.de/glcapsviewer/listreports2.php?listreportsbyextension=GL_ARB_shading_language_include)
It works by specifying the include strings ahead of time. Therefore, before compiling, you need to specify that the string "/buffers.glsl" (as used in #include "/buffers.glsl") corresponds to the contents of the file buffer.glsl (which you have loaded previously).
As you may have noticed in point (2), your paths need to start with "/", like Linux-style absolute paths. This notation is generally unfamiliar to C programmers, and means you can't specify relative paths.

A common design is to implement your own #include mechanism, but this can be tricky since you also need to parse (and evaluate) other preprocessor instructions like #if in order to properly handle conditional compilation (like header guards.)

If you implement your own #include, you also have some liberties in how you want to implement it:

You could pass strings ahead of time (like GL_ARB_shading_language_include).
You could specify an include callback (this is done by DirectX's D3DCompiler library.)
You could implement a system that always reads directly from the filesystem, as done in typical C applications.

As a simplification, you can automatically insert header guards for each include in your preprocessing layer, so your processor layer looks like:

if (#include and not_included_yet) include_file();

(Credit to Trent Reed for showing me the above technique.)

In conclusion, there exists no automatic, standard, and simple solution. In a future solution, you could use some SPIR-V OpenGL interface, in which case the GLSL to SPIR-V compiler could be outside of the GL API. Having the compiler outside the OpenGL runtime greatly simplifies implementing things like #include since it's a more appropriate place to interface with the filesystem. I believe the current widespread method is to just implement a custom preprocessor that works in a way any C programmer should be familiar with.


Why shouldn't Bump, Normal and Displacement maps be gamma corrected?

The fact that a photo is stored in gamma space doesn't have to do with PNG or JPEG, it has to do with the fact that it's a photo.

The camera detects a 50% intensity but is going to save it as a 187 instead of 127, because that is more efficient to keep the color information that matters to the human eye when you have only 256 different possible values. And so does Photoshop if you tell it so.

For normals or displacement maps, you don't have this need so there is no reason to use gamma space for them (using gamma for them would actually be a mistake, since it would affect how precision is distributed). Thus when you decode a normal texture or a displacement texture, you don't need to apply gamma.

If the normal map is encoded in a different space though, as can be the case for a G-Buffer for example, then you have to apply the corresponding correction.


Problems rendering a plane via `drawElements`

This looks like a case of 16bit truncation in your indice buffer. For a 257x257 grid you will be addressing more than 65536 vertices so you will need to use a 32bit indice buffer.

gl.drawElements(gl.TRIANGLE_STRIP, buffers.numIndices, gl.UNSIGNED_SHORT, 0);

change to :

gl.drawElements(gl.TRIANGLE_STRIP, buffers.numIndices, gl.UNSIGNED_INT, 0);

Also change your buffer allocation code to match this type (Uint16Array -> Uint32Array for Javascript).


Yes and no. The RGB values are most likely gamma compressed.

Consider gamma calibration images like the one below.

For a well calibrated monitor the gray value should look the same as the pattern of black and white lines.

If you load in the image you will likely find RGB values of $186$.

$(186/255)^{2.2} \approx 0.5$, which means that a halftone pattern of black and white lines of equal width correspond to RGB values of 186.

However, the gamma compression applied was more likely the sRGB gamma compression than an overall gamma of $2.2$.

$$
\frac{w_\text{white}}{w_\text{white}+w_\text{black}} =
\begin{cases}
12.92V  &\text{ if } V \leq 0.0031308 \\
1.055 V^{1/2.4} + 0.055 &\text{ if } V > 0.0031308
\end{cases}
$$

where $V = 0.299R+0.587G+0.114B$


Standard ways to create xml scenes for research purposes

(In case someone encounters similar issue)

I end-up in Way 1, and it seems the standard way in academia.

Once I am familar with Blender, the scene creation becomes less challenging; still, it takes some time to create a simple scene. Fortunately for academic papers the number of scenes are a few (< 10).

Unfortunately the old Mitsuba parser does not work for the latest Blender version. So I wrote a custom 2.8x Blender add-on for Mitsuba, as well as creating new, custom material nodes for my own renderer.


What are the performance implications of the optional flags used when mapping a buffer?

It seems like the right thing to do would be to always specify that flag yet it is optional.

You can map buffers for reading; invalidation is negatively useful for that. Also, you can map a buffer and only overwrite part of it. Invalidation is negatively useful for that too.

glMapBufferRange has to be able to work outside of the narrow use case of mapping for overwriting an entire range.

  Are there performance reasons (or other reasons) why you would not want to invalidate the data that is in a buffer that you are writing over?

Invalidation is a mechanism to allow the implementation to not delay the act of mapping the buffer (due to prior GPU commands that read from it) if you're going to overwrite all of the area in the buffer. The implementation can allocate a chunk of memory, and you can write to that. Then, it'll schedule the DMA from that chunk into the actual buffer storage when you unmap it, handling the synchronization operations as well. For whole buffer invalidation, it can allocate GPU-accessible memory and completely replace the buffer's storage.

Of course, this only matters if you're mapping frequently. If you're using persistent mapped buffers for streaming (and you should be, if you're using DSA), then you don't care about invalidation. You're going to handle synchronization manually.


Understanding the irradiance map generation

Yes, the convolution is effectively evaluating the following integral:

$$
\text{irradianceMap}(\vec n) = \frac{1}{\pi} \int \text{environmentMap}(\vec \omega_i) \, \max(\vec n \cdot \vec \omega_i, 0) \, d\omega_i
$$

This convolution isn't about whether the environment map is equirectangular, cubemap, or some other projection. We still have to do the convolution regardless, because we want to get the total diffuse light (irradiance) coming to a surface with a given normal, from all directions in the hemisphere around that normal. The environment map to start with just has light from one direction (radiance).


@PaulHK explained it very well in the comment section:

  Ambient is a (rough) approximation of indirect lighting. For example, an outdoor object would be lit by the sun but the surfaces facing away from sun would be black if using only lambert shading, we can use dark gray for the ambient colour to mimic the effect of indirect lighting from the ground and sky

If you were to visualize this with a shading model, say Phong, the Ka, Kd and Ks values would correspond to the following lighting calculations.


Is it possible to decode and play e.g. 60fps video as 30fps, if computer is too slow for 60fps?

It is likely somewhat possible, depending on how the video was encoded. The player can easily skip displaying half the frames but that won't gain us much. What needs to happen is to skip decoding some of the frames.

Modern video codecs nearly all use a mix of I, P and B frames. Most are quite flexible in which frame types are used and only specify what a decoder should be able to handle, leaving the encoder free to pick the most appropriate way to encode a video.

Because of this, even if we want to show only half the frames, we may still end up having to decode more than that. It is possible to build a video stream where every other frame can be skipped because it is not referenced by other frames. But I suspect that most encoders will prioritize encoding efficiency over this rigid requirement. So in practice, we will only be able to skip decoding a smaller number of frames, unless the video was specifically encoded to support this scheme.

The other issue is player support. The player would need to recognize that computing power is insufficient and specifically target half the frame rate. In practice, most that I've seen will skip erratically instead. They expect the lack of computing power to be temporary as that is the most common scenario.

The expected fix for a permanent lack of computing power is to use a lower quality video stream. For devices with well known computing power (eg. phones, tablets, etc), this can be handled automatically by not presenting the high quality alternative. I don't know specifically if YouTube's "auto" quality setting handles this the same way it handles insufficient bandwidth but it could certainly be done and would be a much better solution.

TL;DR: It's possible but usually not worth doing.


Why is there a difference between the CIE XYZ colour gamut vs CIE RGB?

Correct, the missing curvy section in the green-cyan-blue area represents where the red component would have to go negative to express those colors in CIE RGB coordinates.
RGB and XYZ are, at one level, just different coordinate systems covering the same color space—the space of all colors visible to typical human vision. In a mathematical sense, when used as coordinates, there's nothing wrong with negative RGB values (as long as the overall luminance of the color remains positive). But it does present a problem for storing or transmitting such values, as conventional image formats and display signal protocols like HDMI etc only allow for positive values.
On another level, various RGB color spaces are used because they more or less directly represent the actual red, green, and blue subpixels on the display. Those can't emit negative light of course, and so the RGB triangle in color space represents the gamut of colors that can be produced by the display.
Unfortunately because the spectral locus is curved, there's no way that 3 primaries can cover the whole thing. All RGB spaces inevitably cut off a big chunk of the highly saturated green/blue colors.
The XYZ space sort of has the opposite problem. All visible colors can be represented using only positive values in XYZ, but the XYZ primaries themselves are not physically possible colors—they're well outside the visible gamut. So, there's a big chunk of XYZ space that is not valid as a color. And it's not trivial to determine exactly what values are or aren't valid, as you have to test whether they fall inside or outside of the curved spectral locus. This also means that you need more bits per component to get good precision, if you store/transmit images expressed in XYZ—8 bits won't do, probably not 10 bits either, maybe 12 bits would do the job.
If we want to make displays that cover more of the visible gamut, we will eventually need to move to 4, 5, or more primaries. However, that doesn't mean we need 4- or 5-dimensional color spaces. A futuristic display device could be fed by images in XYZ space, for instance, and the device would decide how best to generate each color using the primaries it has available.


OpenGL shader introspection/interface queries for images

Those are not accessible via introspection. It should also be pointed out that some of them don't have to "match".
The shader-defined format doesn't have to be the same as the image-bound format. The image-bound format defines how writing from the image interprets the data provided by the shader. The shader-defined format determines how reading from the image interprets the data in the image. So the two values are fundamentally doing different things and need not match either each other or the texture's internal format. There are rules for how conversion works between these and the image's actual format.
For the access type, you could just pass GL_READ_WRITE for everything and let the shader decide how it's going to access the storage. Given that Vulkan has no equivalent to this value, it probably doesn't matter performance-wise. Also, note that the multibind function glBindImageTextures, which binds multiple images in one call, always uses GL_READ_WRITE (and doesn't have a specified format). Since the whole point of that function is performance, it almost certainly doesn't matter if you always use read-write from the binding side.


Software to compose separate models into an articulated model?

Is there any software which would let me import the existing models,
  define the bone hierarchy and the model which corresponds to each
  bone, and export a single articulated model in T-pose?

I found that it is relatively easy to create articulations using Blender.
This tutorial might be helpful: Blender 3D: Noob to Pro/Bones | Wikibooks


What is a 0 "null" normal map? Or Unlit shader

In the most common normal map encoding, an RGB value of (0.5,0.5,1) would be uniformly flat. If the shader is meant to be unlit, you should be able to remove anything in the code that’s using the surface normal, light / view direction, etc.


DirectX 12: Uploading multiple buffers to 1 default heap using the UpdateSubresources API without explicit copies

UpdateSubresources is a helper function. Inside it is another memcpy. I would recommend simply rewriting it - source is in d3dx12.h. Then you will be able to copy your data no more than necessary.


How to convert a 3D model into an OpenGL fixed pipeline representation

You can compile it into a

Vertex vertices[] = {{posx, posy, posz, texu, texv, ..},
                     {posx, posy, posz, texu, texv, ..},
                     //...
                    }

and then loop over them:

glBegin(GL_TRIANGLES);
for (int i = 0; i < arraycount(vertices); i++){

    glTexCoord2f(vertices[i].texu, vertices[i].texv);
    //...
    glVertex3f(vertices[i].posx, vertices[i].posy, vertices[i].posz);
}
glEnd();

Another solution would be to mimic what the modern vertex buffer based solutions do by creating a VertexAttributeLayout struct and setting it up. Then in the loop you interpret the layout struct and interpret the data as needed.


How was computer graphics libraries or gtk,qt programmed?

You're asking about two different types of library. OpenGL and Vulkan aren't just libraries: they're hardware interfaces. They are implemented by GPU vendors so that application developers can access the GPU, like a special kind of hardware driver. They're implemented in C or C++ and the core functionality is implemented by writing some structures in memory that the GPU itself can understand, then writing to registers on the GPU to make it read those structures and do the work. There isn't any C code in OpenGL to transform triangles, compute their projection onto the screen, and shade them, because that work is all done by the GPU.

OTOH, Qt and Gtk are libraries for making applications. A lot of their work isn't about graphics: they provide a framework for applications to manage themes and interactivity, letting the application developer wire together different parts of the application with callbacks. This is just normal library code. They also provide GUI classes such as buttons and checkboxes. These will each have something like a "draw" function which draws the graphics on the screen. A few years ago this would commonly have been implemented using the OS's 2D drawing API (e.g. GDI on Windows), which offers operations like "draw a filled rectangle", "blit an image", but nowadays they're more likely to use OpenGL or Vulkan operations. The UI components are defined in images and they're used as textures to generate the image of the whole window via a sequence of OpenGL draw calls.

You wouldn't implement OpenGL or Vulkan: this is something that GPU vendors create. You can't really do it yourself without detailed knowledge of the GPU hardware. (However, Mesa is an example of a software implementation of OpenGL, which doesn't depend on any GPU.)

If you wanted to implement a GUI library like Qt, you might start by defining an abstract Widget class with a draw method, and then implementing concrete subclasses of it to draw particular GUI components. Add extra methods and functionality as you need it.

To be honest, I don't think that would be a valuable project to work on. There are many smaller GUI frameworks already, and another toy one is not going to improve the world. You might treat it as just a learning project, but I think you'd learn more by studying one of the existing frameworks and writing an interesting GUI application on top of it.


There is the whole derivation of it but I'll be discussing a brief overview.
This is for the perspective projection  where the line joining the eye and the center of the projection/image plane is perpendicular to it. Like here

As we can easily see, by similar triangles $\triangle ABC$ and $\triangle AEF$ we have

$Y_p / Y = D/-Z$

where $Y_p$ is the projected $Y$ coordinate on to the image plane. $AB = D$ is the total distance from eye to the image plane and is usually set to 1.

Hence we have,

$Y_p = (Y*D)/-Z$

Most APIs also set that the camera is facing towards the negative Z axis following the convention for right handed coordinate system, hence the negative sign. The projection plane is also assumed to be a rectangle instead of a square. It's 2 units long in the $Y$ direction and centered at the origin i.e in the range [-1,1] so it's 1 unit long in the positive Y axis and 1 unit long in the negative. However horizontally, its in the range $[-A, A]$ so when calculating $X_p$ we have to correct for this factor.

$A = Width\ of\ rectangle/Height\ of\ rectangle$

Similarly for $X_p$ and $Z_p$ we have,

$A*X_p = (X*D)/(-Z)$

$X_p = (X*D)/(-Z*A)$

$Z_p = D$

We can easily see see that

$tan(Y_{fov}/2) = BC/D = 1/D$

$D = 1/tan(Y_{fov}/2)$

Here, $Y_{fov}$ is the angle of field of view in the Y axis. Since we are only talking about the upper half its $Y_{fov}/2$.

Hence one way to write $X_p$ and $Y_p$ is

$X_p = \frac{X}{-Z} * \frac{1}{A*tan(Y_{fov}/2)}$

$Y_p = \frac{Y}{-Z} * \frac{1}{tan(Y_{fov}/2)}$

This is what you are seeing in the matrix. The $X$ and $Y$ come from the vector at which the transformation matrix is applied and the $-Z$ term is accounted for by the $-1$ in the 4th row. This $-Z$ is stored as the $w$ coordinate which is then used in perspective division in which we divide by $w$.

For more information you can check this site here

Note:- Haven't discussed properly about NDC space, perspective divide and about near and far planes in $Z_p$ as that would take too much time. Just wanted to show you that it can be derived.


Jagged texture's edge on perfectly straight UV of sphere mesh

It looks like distortion because the trapezoids are rendered as two triangles by the GPU. This leads to incorrect UV interpolation.

Try this: draw a horizontally symmetrical trapezoid with the top being much shorter than the bottom. Then draw one of its diagonals to split it into two triangles. Place a mark at the midpoint of that diagonal. It should be obvious the mark is not at the horizontal center of the trapezoid. Yet that mark is where the middle of the interpolated texture coordinates will be. Your border will go through the mark but also through the center of the top and bottom edges.

You should not have that problem when producing a final image with any decent renderer.


for (i = 0, i < 3, ++i)
    glTexCoord3d(texturePoint[i].x, texturePoint[i].y, 0);
    glVertex3d(v[i]->x(), v[i]->y(), v[i]->z());

This code calls glTexCoord3d three times and glVertex3d only once. You need more braces around your for loop.


In texture mapping, should we set texture coordinate to every vertex or can I skip vertices in between?

Typically, texture coordinates are interpolated from the vertices of a triangle during rendering. This can be seen in two ways. evaluating your texture at vertices and the interpolating the result or interpolating the texture coordinates and evaluating your texture at the fragment level.

If you do not supply vertices with corresponding texture coordinates then there is no guarantee that the interpolated texture coordinates or colours for the fragments of the rendered triangle will behave nicely. By skipping a few vertices, the vertices would need some default texture coordinate such as $(0,0)$, which would show up as a weird artefact in your rendering. The vertices B and E in your example are required to have texture coordinates in order to have a continuous set of interpolated texture coordinates or colour values.


1
For each texel of the generated diffuse cube map (which has $6$ faces times $32 \times 32$ texels), the diffuse equation combines all texels of the environment map (which has $6$ faces times $64 \times 64$ texels). Thus $32 \times 32 \times 6 \times 64 \times 64 \times 6$ as stated.
Since the article considers SH with $9$ coefficients, and each coefficient needs to combine all texels of the environment map, that's $9 \times 64 \times 64 \times 6$.
2
"Bruteforce" doesn't necessarily refer to graph theory. In this article it refers to the naive approach of going through every single texel of the diffuse map and each time combine every single texel of the environment map.

Appendix
Let's describe the problem at hand.
As an input, we have an environment map: it's a set of six textures forming a cube, that represents the light coming from all directions. The light is assumed to be coming from sources far enough that, at our scale, it's as if they were infinitely far. An environment map typically looks like the left figure of this illustration from the GPU Gems article:

As an output, for a any point of a given object that we're rendering, we want to know the diffuse irradiance affecting it. What does that mean?
A parenthesis
Models describing how the perceived color of a surface varies depending the light reaching it, typically consist in integrating all incoming light over an hemisphere, with a weight that depends both on the direction of the light and of the observer.
For opaque materials, it is common to consider two components contributing to that color: the diffuse and the specular. The diffuse is due to light slightly penetrating the surface, getting partly absorbed, then getting out in random directions. The specular is due to light bouncing off the surface without entering it.

Back to the problem
So we want to know the color of the points of an object, we know that to compute it we need to integrate light in all directions, and we have a cube map that represents light coming from all directions. So far so good.
A problem though is integrating basically means we need to go over every single texel of the hemisphere (so half the environment map: $(64 \times 64 \times 6) / 2$ texels), do some math, and add that to the result. That's a lot of computation for each point, which we'd like to avoid if possible.
We know that the contribution of one light depends on the light direction. If we consider the environment to be static (the lighting doesn't change), then we can isolate the part of the computation that only depends on the light and surface normal (and not on the material or the observer), pre-compute it and store it to use later. For the diffuse term, that's the diffuse irradiance mentioned earlier, and it typically looks like the center figure of the illustration. Each pixel represents the irradiance term for a given surface normal. The right figure is the specular environment map, computed in a similar way but with a different integral.

With that texture, instead of painfully integrating, we just fetch the precomputed result, combine it with the rest of the equation (the material diffuse color, etc.) and voilà.
Finding the 32x32x6x64x64x6
Since the result looks very blurry, it probably doesn't need to be $64 \times 64$ for storing, and the article assumes $32 \times 32$. But still, each of those $32 \times 32 \times 6$ texels' color is computed separately, by integrating over the entire corresponding hemisphere. That part never disappeared: we just moved it from the rendering stage to a pre-computation stage.
So that's still $(32 \times 32 \times 6) \times (64 \times 64 \times 6) / 2$. That's half the figure stated in the article. I suppose the author was either assuming a model that integrates of the entire sphere, or just dropped the $0.5$ factor for the sake of clarity.
Finding the 9x64x64x6
We can't help but notice this irradiance diffuse environment map is really blurry: surely there must be a way to store its information using a lot less than $32 \times 32 \times 6$ terms, and hopefully reduce the amount of work as well. That's what spherical harmonics do, by compressing the irradiance in just a few coefficients.
I am not familiar enough with the math to engage in an explanation, but suffice to say it is generally considered that for diffuse, order three SH are sufficient. Order three SH involve nine coefficient, and again, for each of these we need to integrate over the entire hemisphere.
So that's $(9 \times 64 \times 64 \times 6) / 2$ texels to fetch and process. Here too a $2$ factor is missing from the article.


Filling pixels under or above some function

Compute the tangent to the line at that point, Use it to make a tangent vector, then use the point on the line, the point being tested and the new point to create two vectors. Take the dot product of the vectors. Positive values are on one side of the line, negative are on the other side. (above and below).


Using Software Rendering to make an External GPU

Theoraticaly this can work, but I am not sure if it would be very efficient. Modern GPU's are tailored for doing low level graphics operations and the drivers that manage them are also tailored for efficiently overseeing these operations to render your scene.
Now you can ask all the hardware questions related to using Raspberry Pi, how to access its cpu etc to Raspberry Pi SO.
As per software side, I think P. Shirley's books on
ray tracing should get you started on rendering on cpu's rather than gpus. Using multiple threads is extensively discussed in the following issue of the source code of the book. So once you had written your single threaded cpu ray tracer, you can adapt it easily to execute in multiple threads.

Update 21-01-2021
I noticed that there is already a driver for Raspberry Pi 4 provided by mesa, which gives you access to OpenGLES 3 functionality. I think it is a decent trade off, I was even surprised to see that there is actually a graphics driver for Pi.
Unless you insist on ray tracing your scenes, I would suggest to use that for your experiments.


How to use trackball for rotating an camera matrix in eye frame

I'm assuming the rotation matrix generated from the arcball is 2D, so I think there are a few clear top level steps. But I'm not sure about the details of some of them.
First you'll need to extend the 2D rotation to 3D which should be easy.
Next you'll need to rotate the rotation matrix itself so its axis aligns with your Forward vector. I'm not certain how to do this step.
Then, use the resulting matrix to transform the Up vector and reconstruct the camera matrix using the new vector and the other values and vectors.


How to use 32 bit Integers for Element Indices in WebGL 1.0?

Enabling the extension is as simple as:

var uints_for_indices = gl.getExtension("OES_element_index_uint");

if uints_for_indices is null then the extension is not available and you'll have to split your data anyway.

If it's not null then you can now pass gl.UNSIGNED_INT for the index type into gl.drawElements.


How do I create my own bitmap using a text editor?

In addition to the pixel values, the width and height of the image are also required, and the colour depth (how many colours are available per pixel). Some formats also specify an alpha value, for transparency, which may be defined per pixel in addition to the red, green and blue components.

There are many different image file formats but a simple one to experiment with to gain understanding is ppm, one of the Netpbm formats.

Netpbm defines the following file types:

PBM: portable bitmap (pixels are either on or off)
PGM: portable graymap (grayscale pixels)
PPM: portable pixmap (full colour pixels)

These have the advantages of being simple to understand, and easy to generate in a text editor for working with simple examples.

For example, a PPM file is simply the characters "P3" to identify it as PPM, followed by the width, height, maximum colour value and the pixel values in English reading order (left to right then top to bottom). Each pixel colour is a red value, a green value, and a blue value. All these values are space separated and human readable.

P3                         # Indicates an ASCII PPM file
2 2                        # Width and height of the image
255                        # Maximum colour value
255 0 0     0   255 0      # A red pixel (255 0 0) followed by a green pixel
0   0 255   255 255 255    # A blue pixel followed by a white pixel

While this makes the format very easy to work with for small examples, it will also make the file size much larger than other formats which use bytes to represent colour components instead of writing out the human readable numbers. Many formats also include compression (either lossless or lossy).

Once you are accustomed to writing PPM files in ASCII characters, indicated by the file starting with "P3", you can move on to PPM files using bytes, indicated by the file starting with "P6". This is a similar format to what many programs use internally, but is no longer human readable in a text editor as not every byte will correspond to a printable ASCII character.

You also mention that a pixmap can have pixel values that refer to colours stored elsewhere in the file. I have avoided this format as it is not necessary for understanding how a pixmap works. This is known as using a palette (a list of predefined colours so each pixel is represented by an index of a colour in the palette). This is necessary for formats that only allow a small number of distinct colours, but for full colour images there is no need for a palette - each pixel is defined as red, green and blue components without having to refer to anything else.


Differential or intersection of islands of 2D polygons in different layers

In CGAL, there is the Arrangement package that allows to build a topologically valid planar partition given a set of segments, and the Regularized Boolean Set-Operations that provides boolean operations between polygons.


It's a variant to projective texture mapping.

The simple way to do this in the fragment shader by using the position of the fragment to decide whether it is close enough to the plane of the laser to light up.


Could we dispense the near clipping plane?

The near clipping plane in rasterizing setups with perspective projection is there to avoid a divide by 0 and bound the possible depths for orthogonal projection.

With a bit of extra math you can make the depth stored in the depth buffer linear again. However you still need the guard against the divide by 0.

In raytracing the near plane can be at 0 no problem. There is nothing technical stopping you from putting it at -10 however (besides that it would be very odd).


How to randomly draw Quaternions within a specific range of Euler angles for rotation?

The simplest way to do this would be to generate the Euler angles randomly, then use them to create the quaternion. Most existing quaternion implementations in game engines and 3D math libraries will have a FromEuler ( Vector3 ) function.

Since you seem to be rolling your own quaternions rather than using a library, here's the conversion formula.

As a small optimization, since quaternions use half-angles, you could generate random values in the range ±15 degrees (pi/12) and feed these to the sin and cos functions directly.


Extending an orthogonal set of vectors (graphics application?)

The most common use is almost certainly computing normal vectors for a surface, using a cross product. Although strictly speaking, this does not quite qualify as the two initial vectors are usually not orthogonal. This is used for shading calculations and various other algorithms. I know of at least one where the initial surface is explicitly a rectangle so it would qualify.

The sampling of BRDFs is another application. In many cases, the math is done for a canonical space with the surface normal pointing up and must be transformed to the actual surface. This is done by building a full 3D basis from the surface normal vector.

Anisotropic BRDFs typically have two vectors as input: the surface normal and one vector to define the orientation of anisotropy. Depending on the implementation, it may be assumed to be orthogonal to the normal. Otherwise, it is made orthogonal. Then, again, a full basis is built from those two.

There is a class of ray intersection algorithms which work by projecting the surface in two dimensions. This is done by starting with the ray and generating two more orthogonal vectors which define the projection plane.

This is certainly not an exhaustive list but they are the main ones that came to mind. I am not aware of any uses for higher dimensions but I don't know about every corner of computer graphics either.


Interpreting Image Processing Math Equation

It sounds exactly like that!

$\Omega$ is the image domain, it can be either discrete or continuous. No idea how it is defined in that particular book.
$x=(x_1, x_2)$ are the 2D image coordinates
$u_0, u, n$ are functions that return intensity values at the given location $x$.


VAOs are essentially 'plumbing' objects that help get data from your buffers into your vertex shader ready for drawing. On the GPU this is handled by the vertex fetch stage, so a VAO is basically just a bundle of vertex fetch state corresponding to a particular model.

You can think of it as a bundle of pointers, but with somewhat more complexity than CPU-side pointers. Since OpenGL buffers are just a blob of untyped bytes by default (that's why you get a void* when you map them), the VAO handles type-casting (to vec3 or whatever) when fetching the data. It also deals with the size of the data not necessarily matching the stride between consecutive elements, and makes sure that each data item (position, normal, uv, etc.) is transmitted to the correct 'channel' in the vertex shader input.

You can render without a VAO in legacy OpenGL, but this is generally more of a PITA, since you need to bind each vertex attribute separately each time you draw the object. Encapsulating all these bindings in a VAO lets you set them up once, then you have a 'plug and play' object that simply needs to be bound before a draw call, and all your vertex data will be ready to go.


Is there a performance implication if the frame buffer resolution and the physical screen resolution are different?

For PCs with discrete displays, the resampling to the screen resolution may happen either on the computer prior to the signal being sent to the display, as part of window compositing, or on the display side. I don't think it's anything to worry about in terms of FPS, in either case, but it possibly could affect latency.

It depends whether you're in "true" fullscreen (exclusive) mode or not. If the app renders to a window, including a borderless fullscreen window, then the window compositor (part of the OS) will be running each frame, and any upsampling needed can happen "for free" as part of that operation. In other words, I wouldn't expect an upsample during compositing to cost noticeably more than it does to do compositing in general (which is very cheap).

On the other hand, if you're in true fullscreen and your app has exclusive control of the display, it's possible to set the video mode to something non-native and send your rendered image directly to the display, in which case the display's scaler chip will handle it, again without performance cost to you, since it's happening in parallel to anything going on with the GPU.

You mentioned Android, so I assume you're working on a mobile device. In that case, I'm not sure if there is a scaler chip, or if the mobile GPU is just expected to always output at native resolution, in which case scaling would always be done on the GPU by the compositor software. (By the way, the reason the scaling doesn't show up in profiling is likely because it's being done by the OS rather than by your app.)

That being said, though upscaling is very cheap, it does imply an extra pass that could increase input latency. With a borderless fullscreen window, if no other windows are visible on top of your app, and no scaling is required, the compositor could be short-circuited so that your app writes directly to the display output rather than passing through the compositor. (This is an optimization that the OS would have to perform when it detects that all the conditions are right.) Similarly, on the display, if no extra image processing is required including scaling, color correction, etc. then some stages of the display hardware might be bypassed. In either case, you could see lower latency when you render directly to native resolution than when you require scaling. This is a matter of concern for VR in particular.


What is the relation between Aliasing and Flickering?

Flickering can be a form of temporal aliasing. It's a similar phenomenon to spatial aliasing such as jaggies, but it occurs in time instead of space.

For instance, a common cause of image flickering in graphics is when the camera or geometry is in motion, and geometric features fluctuate in pixel size as they move. For example, imagine a railing with thin vertical bars. Depending where a bar appears relative to the pixel grid, it might get rendered as 2 pixels wide, only 1 pixel wide, or it might not appear at all. And in motion, it may rapidly fluctuate between these states, creating a visually objectionable flicker.

Another common cause of image flickering is specular surfaces with a bumpy normal map and a high specular power. The specular highlights can flicker in motion, due to their alignment with pixels changing from frame to frame.

Antialiasing strategies that address only spatial aliasing will often produce an image that looks good in a static screenshot, but turns into a flickery mess as soon as things start moving. This is one reason why temporal antialiasing has become popular in games in recent years.


You're trying to map all three UAVs to the same descriptor slot (slot 0).

ranges_madhu[0].Init(D3D12_DESCRIPTOR_RANGE_TYPE_UAV, 1, 0, 0, D3D12_DESCRIPTOR_RANGE_FLAG_DATA_VOLATILE);
ranges_madhu[1].Init(D3D12_DESCRIPTOR_RANGE_TYPE_UAV, 1, 0, 0, D3D12_DESCRIPTOR_RANGE_FLAG_DATA_VOLATILE);
ranges_madhu[2].Init(D3D12_DESCRIPTOR_RANGE_TYPE_UAV, 1, 0, 0, D3D12_DESCRIPTOR_RANGE_FLAG_DATA_VOLATILE);

The baseShaderRegister value (the third parameter to these calls) indicates which descriptor slot to use. You have all of them set to 0, which won't work. They need to be distinct. For example, they could be 0, 1, and 2.


How to use GetViewProjMatrix().TransformVector(LineDirection) in UE4?

I don't think you can simply use TransformVector when the matrix you're transforming by involves a projection matrix. To fully project to screen space, you have to divide by W, which TransformVector doesn't do (it simply multiplies by the matrix without translation).

Also, because of the divide by W, transforming vectors (as opposed to points) by a projection matrix is probably not what you actually want. If what you want is the screen-space vector between two points, you should transform and divide-by-W each point, then subtract them to get a screen-space vector. Otherwise the vector will likely be inaccurate.


What software was used to make KDA - Villain?

The video was created by a studio called Trizz; you can see their page about the project here. The site doesn’t have any information about what software they used, but you can achieve this look (with enough work) in pretty much any modern renderer. If you use Blender, Cycles is a good one; with other 3D software, Redshift, Octane, and Arnold are popular options.
If you’re also looking to make the kind of physics-based animation—cloth, fire, snakes, etc.—that appears in the video, you’ll probably also want to check out Houdini.


Guidance with Sphere-Map calculation in GLSL

You probably need to calculate the reflection vector per pixel, not per vertex. The calculation is non-linear, so doing it per vertex and interpolating will not give good results unless the model is finely tesselated.

Send norm and toCameraVector down to the pixel shader, normalize them both and move all your calculations under // Reflection Map to the pixel shader. That should give you some better results.


total internal reflection implementation in pbrt-v3

The code does seem a bit ugly since it doesn't assign pdf and instead returns 0.

By returning 0 it's basically saying the same thing as pdf = 0 because f.IsBlack() returns true.

This is acceptable, it's basically saying the BSDF doesn't allow light go through at the given $w_o$


Apply Projective Transformation 2D Rectangle

Since homography is not my topic, I can't really tell if I am missing some important edge cases, but after some research, I think I got the basics. I will use the same variable names, the guy in the video uses on slide 15 (slide numbers are in the bottom right corner on the red background).

You start by defining 4 points in the plane $X$ where you say you know their projected points in plane $x$. You assume, that the transformation between $X$ and $x$ is linear and can be described by a matrix $H$. Means
$$x = HX$$

Now you want to calculate the Matrix $H$ that perfoms the projection:

$$H = \begin{bmatrix}h_1&h_2&h_3\\h_4&h_5&h_6\\h_7&h_8&h_9\end{bmatrix}$$

$h_1$ to $h_9$ are unknowns. So how do you get them? Well, you have 4 times:

$p_i' = Hp_i$

where $p_i'$ and $p_i$ are the projected and the original points. This full equation looks like this:

$$\begin{bmatrix}x_i'\\y_i'\\1\end{bmatrix}= \begin{bmatrix}h_1&h_2&h_3\\h_4&h_5&h_6\\h_7&h_8&h_9\end{bmatrix}\begin{bmatrix}x_i\\y_i\\1\end{bmatrix}$$

Perform the multiplication and you will get this:

$$\begin{bmatrix}x_i'\\y_i'\\1\end{bmatrix}= \begin{bmatrix}x_ih_1+y_ih_2+h_3\\x_ih_4+y_ih_5+h_6\\x_ih_7+y_ih_8+h_9\end{bmatrix}$$

You can also write this in some kind of nested vector fashion:

$$\begin{bmatrix}x_i'\\y_i'\\1\end{bmatrix}= \begin{bmatrix}\bar{h_1}^T p\\\bar{h_2}^T p\\\bar{h_3}^T p\end{bmatrix} = \begin{bmatrix}p^T\bar{h_1}\\p^T\bar{h_2}\\p^T\bar{h_3}\end{bmatrix}$$

where
$$\begin{matrix}\bar{h_1} = \begin{bmatrix}h_1\\h_2\\h_3\end{bmatrix}&\bar{h_2} = \begin{bmatrix}h_4\\h_5\\h_6\end{bmatrix}&\bar{h_1} = \begin{bmatrix}h_7\\h_8\\h_9\end{bmatrix}\end{matrix}$$

Now here comes the trick (got it from this link). If  $p'$ is equal to the vector $Hp$ their cross product must be the zero vector $\bar{0}$:

$$p' \times Hp = \bar{0}$$

Applying the cross product will give you three equations:

$$\begin{matrix}y_i'p^T\bar{h_3} -  p^T\bar{h_2} &= 0 \\ p^T\bar{h_1} - x_i' p^T\bar{h_3} &= 0\\x_i' p^T\bar{h_2}-y_i' p^T\bar{h_1}&=0\end{matrix}$$

Write this as a matrix equation and you get:

$$\begin{bmatrix}\bar{0}^T &-p^T&y_i'p^T\\p^T&\bar{0}^T&-x_i'p^T\\-y_i' p^T&x_i' p^T&\bar{0}^T\end{bmatrix} \begin{bmatrix}\bar{h_1}\\\bar{h_2}\\\bar{h_3}\end{bmatrix}= \begin{bmatrix}0\\0\\0\end{bmatrix}$$

Now, if you look at the matrix, you might notice, that the 3 columns are not linear independent. $-y_i'$ times the first column minus $x_i'$ times the second one, yields the third column. So you can drop one.

If you drop the last column and replace all the vectors with their corresponding entries, you get:

$$\begin{bmatrix}0&0&0 &-x_i&-y_i&-1&y_i'x_i&y_i'y_i&y_i'\\x_i&y_i&1&0&0&0&-x_ix_i'&-y_ix_i'&-x_i'\end{bmatrix} \begin{bmatrix}h_1\\h_2\\h_3\\h_4\\h_5\\h_6\\h_7\\h_8\\h_9\end{bmatrix}= \begin{bmatrix}0\\0\end{bmatrix}$$

You can multiply the second row by $-1$ and swap it with the first one (that does not change the solution of the system) to get the exact same formula as presented in the video on the slide 17.

If you apply the formula to all 4 points you get the system:

$$\begin{bmatrix}
0&0&0 &-x_1&-y_1&-1&y_1'x_1&y_1'y_1&y_1'\\
x_1&y_1&1&0&0&0&-x_1x_1'&-y_1x_1'&-x_1' \\
0&0&0 &-x_2&-y_2&-1&y_2'x_2&y_2'y_2&y_2'\\
x_2&y_2&1&0&0&0&-x_2x_2'&-y_2x_2'&-x_2' \\
0&0&0 &-x_3&-y_3&-1&y_3'x_3&y_3'y_3&y_3'\\
x_3&y_3&1&0&0&0&-x_3x_3'&-y_3x_3'&-x_3' \\
0&0&0 &-x_4&-y_4&-1&y_4'x_4&y_4'y_4&y_4'\\
x_4&y_4&1&0&0&0&-x_4x_4'&-y_4x_4'&-x_4' \\
\end{bmatrix} \begin{bmatrix}h_1\\h_2\\h_3\\h_4\\h_5\\h_6\\h_7\\h_8\\h_9\end{bmatrix}= \begin{bmatrix}0\\0\end{bmatrix}$$

The problem is now, that you have 9 unknowns but just 8 formulas. According to the answers in this question, one possible solution is to constraint $h_9$ to be 1 and solve the resulting linear system. However, another way is to solve it as an underdetermined least-squares problem. This might be done using the SVD. However, as I said in the comments, this is a rather complicated algorithm. An easier one is to use QR decomposition, which can be a substep of the SVD, but this is also not too easy to code yourself. Since you don't seem to be that experienced in advanced linear algebra, you might try the constraint technique first, before trying the SVD.

You can also use more than 4 points, but then you get an overdetermined system, which can be solved with the SVD, QR decomposition or any other least-squares solver.


Help me find out what this texture mapping technique is called

What you see in the image called a UV map. That is, it is simply texture coordinates to be looked up encoded in a image. Same thing happens in all texture lookup in 3D there is a underlying sampler that picks where to pick texture color from.

Image 1: Image showing UV map of two overlapped triangles and sampled texture with same UV coords

Here are the sources for those images, please do not overwrite the sources.

Demo showing the UV map
Demo showing the texture


Each cube face is made up of a quad, or two triangles.  Let's take 3 points from either the quad, or one of the triangles, and call them $A,B,C$.

You can get the normal to this surface by normalizing the cross product of $\overline{AB}$ and $\overline{BC}$.  You can then use that surface normal to calculate lighting.

You can do similar calculations to get the tangent and bitangent to then calculate $u,v$ values for texture mapping.

I'm not sure why the edges not having normals makes you think you can't calculate lighting for the cube faces, but you totally can.

You could even use the vertices to calculate a normal for the edges if you wanted, but I'm not sure why you would want to do that :P


Why is the approximation valid, in the formula provided by Brian Karis?

"the information that $f(l_k,v)$ carries is low-frequency enough": As IneQuation explains, low-frequency was used to refer to the detail of the brdf function. I did actually mean that $f_r$ was low frequency though (which is the case with diffuse lighting), not $L_i$.

"with respect to $L_i$": what this means it that, since there aren't any large peaks in $f_r$, no part of $L_i$ will be wheighted much more heavily than any other in the integral.

If each function was, say, reciprocal of the other (meaning in this case that the peaks in each one would cancel out when multiplied by the other) the correct integral would be very different to the separated integral where this cancelation does not occurr. but since the two functions are not correlated (or they shouldn't be) this difference is likely to be fairly small.

Also, if you use a really high frequency BRDF that takes into acount only incident light around the normal and you separate it, you will get something that's way off. in this example the estimate using separation is 4x larger that the correct result: https://www.desmos.com/calculator/c4zgypthrc

Experimental results can be generated easily. here is an example: https://www.desmos.com/calculator/poyswojjsa it shows a low frequency brdf $f_r$ and a high frequency incident light $L_i$ function being multiplied before and after integration. the ratio between the correct result and the one obtained using separation is of about 1.6

Addendum regarding your latest comment on my other answer:

$$\frac{1}{N}\sum^{N}_{k=1}f(L_k) = \sum^{N}_{k=1}(f(L_k)\cdot \frac{1}{N})$$

Here we have "moved" or distributed $\frac{1}{N}$ into the sum. This formulation is equivalent but looks a bit more similar to an integral. if we set $\Delta L = \frac{1}{N}$ then, we can get

$$\sum^{N}_{k=1}f(L_k)\Delta L\approx \int_{L}f(L)dL$$

keeping in mind that the separation that occurs is, conceptually, separation of integrals, and that every integral "needs" its own $dL$. it follows that every montecarlo estimator needs its own $1/N$

I hope this helps you understand where the extra$\frac{1}{N}$ comes from


How does the Lane Riesenfeld algorithm work?

What you're looking for is called de Boor's algorithm. It lets you compute a point on a b-spline curve by doing a series of linear interpolation (LERP) calculations. So, it works very much like the de Casteljau algorithm for Bezier curves. In fact, the de Casteljau algorithm is a special case of de Boor's algorithm.

A link
Another one
And another


Ray Transformation to Object Space for Motion Blur

Lerping the ray positions/directions between keyframes should be equivalent to lerping the inverse matrices between keyframes and transforming by the lerped matrix. Trouble is, if the keyframes have different rotations, that lerped matrix will in general be something "weird", with shearing, nonuniform scale, etc.

I wouldn't be surprised if some of your intersection and shading routines don't work properly in such a distorted coordinate system, unless you've specifically tested and hardened them against such cases. (For example, taking the dot product of two unit vectors doesn't give the same answer in a sheared coordinate system as in an orthonormal one.)

This is just a guess, but it might work better if you choose an interpolation method where the translation, rotation, and scale (if applicable) are lerped separately (using quaternions for the rotation part) and re-combined.


pseudoinverse Jacobian and adding more control in computer animation

I think I found your misunderstanding, but it's IMHO based on a little inconsistency (or at least lack of clarity) in the book.

  But, the θ is always zero vector because J+J=I, so θ=0z.

This is not true, $J^+J$ is not necessarily $I$. $JJ^+=I$ but the multiplication by the pseudoinverse is not commutative. Let's look at this in more detail:

$JJ^+ = JJ^T(JJ^T)^{-1}= (JJ^T)(JJ^T)^{-1} = I$

$J^+J = J^T(JJ^T)^{-1}J = ?$

The underlying misunderstanding comes from the assumption that the matrices $J^T(JJ^T)^{-1}$ and $(J^TJ)^{-1}J^T$ are equal (since the latter would produce $I$ if right-multiplied by $J$), which is not true. They are both pseudoinverses of $J$ but they are not necessarily the same matrices, in fact they don't necessarily both exist.

As the book says, the $(JJ^T)^{-1}$ only exists if the rows of $J$ are linearly independent, and in turn Wikipedia says that $(J^TJ)^{-1}$ exists if the columns of $J$ are linearly independent. However, they can't both exist if the matrix is not square, since one of the dimensions would necessarily have more vectors than the dimension, which can't all be linearly independent. Specifically, for a normal redundant joint system you probably have more columns than rows, making only the pseudoinverse $J(JJ^T)^{-1}$ exist.

We can, however, see that for a square matrix they are indeed the same and actually the real inverse:

$J^T(JJ^T)^{-1} = J^TJ^{-T}J^{-1} = (J^{-1}J)^TJ^{-1} = I^TJ^{-1} = J^{-1}$

$(J^TJ)^{-1}J^T = J^{-1}J^{-T}J^T = J^{-1}(JJ^{-1})^T = J^{-1}I^T = J^{-1}$

(But in that case the system would have a unique solution and there's no point to angle control anyway.) But you cannot resolve the brackets with a non-square matrix, since the individual matrices aren't square and don't have an actual inverse. The only situation when both matrices exist are square $J$s and then they are indeed equal, but then there's also no point to angle control.

But I understand why you are confused and I am a little, too, since Parent actually seems to use the alleged identity $J^T(JJ^T)^{-1}=(J^TJ)^{-1}J^T$ in equation 5.20, where he derives the pseudoinverse's relevance for solving the inverse kinematics problem. He computes everything with the left inverse but than concludes with using the right inverse. I guess that works conceptually because they're both pseudoinverses, even if it doesn't make complete mathematical sense (to me and you at least). He seemed to have used the assumption of a square matrix as a trick for transforming the problem and then used the conclusion for a non-square matrix, which might be conceptually correct.

I guess it simply doesn't matter which of the two pseudoinverses we use, since only one of them ever exists (for a non-square Jacobian) anyway and if the left inverse would exist (making your above assumption true), that means you have more constraints than DoF and thus your system is overconstrained and there's no point in angle control.


What does "region to be shaded" mean in the PCF paper?

Texture mapping is the process of providing a mapping from the surface to a location on the texture. That is, given a particular point on a surface, you can find a specific point in a texture that corresponds to it. Giving positions a UV texture coordinate is a common means of defining a texture mapping, but you can use a function of arbitrary complexity.

Texture mapping is at least theoretically reversible. That is, given a particular texel, you can use the texture mapping to determine where on an object that texel appears. Of course, not every texture mapping is reversible. Some mapping cause a texel to be used in multiple locations on a surface, and not every mapping uses every texel in the texture.

But the important part of the reverse mapping operation is this: a texel maps to some area of the surface, not a single point.

Shadow mapping is a form of texture mapping. It defines a mapping from the scene to a texture by performing perspective projection of the points in the scene onto the surface of a texture, based on a light position.

If you reverse this mapping, a texel in the shadow map maps to an area in the scene. This is the "region in the scene".


Creating a vanishing point perspective shader

It should be possible to do this properly by using a modified projection matrix for your entire scene that has its far plane at infinity, as detailed here. That should allow you to properly render points infinitely far away from the camera. The presentation says that such a point is represented by a direction (a vector with .w = 0), so it should be as simple as replacing the position of your vertex with the direction vector pointing towards the desired vanishing point.

One caveat is that any interpolation between the vanishing point and any other vertex would be incorrect, because the perspective correct hardware interpolator can't properly handle infinite values. This shouldn't affect lighting assuming normals are the same, but it will make texturing considerably harder. Probably the easiest way here would be to do an inverse projection of the point back to object space in the fragment shader and define the texture coordinates to be based on the object space position of the fragment.


Most mobile GPUs use tiled rendering, which can make MSAA surprisingly cheap relative to desktop GPUs. If the app is doing simple forward rendering without filtering/blurring passes, the GPU may be able to do MSAA rasterization and resolve directly on chip (per tile) without ever storing the MSAA samples out to memory. This may make MSAA both faster and less power-hungry than a technique like FXAA that relies on filtering the framebuffer (which necessitates storing it out to memory and reading it back in when sampling).

As for image quality, it's difficult to compare (and is subjective anyway). FXAA might be better at reducing jaggies depending on the MSAA level, but at the cost of over-blurring edges. FXAA can also address shader aliasing, for example aliasing due to specular on a rough normal map and such things, which MSAA won't help with. Ultimately you'll have to make a choice whether the aesthetic benefits (if any) outweigh the increased cost.


since that is the case for window-space coordinates in the vertex stage.

There are no "window-space coordinates in the vertex stage". Or at least, not in the vertex shader.
Window-space is not relative to the viewport; it's relative to the window (or framebuffer, in this case). The viewport transform is used to convert from NDC-coordinates (which are relative to the viewport) to window coordinates (which are not).
If you need your fragment shader to have viewport-relative coordinates, then you need to either compute them in the fragment shader (by providing the viewport as a uniform) or pass NDC coordinates from the vertex shader as an in/out variable pair (just compute gl_Position.xy/gl_Position.w).


ShadowMap peter-panning and wrong direction

I found out that simply removing minus sign in shader in converting coords to [0, 1] solved direction problem. But scale is still off - here is an example:

void LightClass::GenerateViewMatrix()
{
    XMVECTOR eyePos = { m_position.x, m_position.y, m_position.z };
    XMVECTOR focusPos = { m_lookAt.x, m_lookAt.y, m_lookAt.z };
    XMVECTOR upVec = { 0, 1, 0 };

    m_viewMatrix = XMMatrixLookAtLH(eyePos, focusPos, upVec);
}

void LightClass::GenerateProjectionMatrix(float screenDepth, float screenNear)
{
    float fieldOfView = (float)3.14f / 2.0f;
    float screenAspect = 16.0f / 9.0f;

    m_projectionMatrix = XMMatrixPerspectiveFovLH(fieldOfView, screenAspect, screenNear, screenDepth);
}

Edit [Actual solution - ignore previous posts]:

I followed PaulHK solution (thanks, it's second time you've helped me) about checking once again View/Projection matrix. Basically my calculations were wrong. I debugged my convertions from model to lightspace and changed Vertex Shader from:

output.lightViewPosition = mul(input.position, worldMatrix);
output.lightViewPosition = mul(output.lightViewPosition, g_lightViewMatrix);
output.lightViewPosition = mul(output.lightViewPosition, g_lightProjectionMatrix);

To correct:

matrix lightViewProj = mul(g_lightProjectionMatrix, g_lightViewMatrix);
output.lightViewPosition = mul(input.position, lightViewProj);

That simple changed fixed my results. With simple Pixel shader code below I am getting this results:

input.lightViewPosition.xyz /= input.lightViewPosition.w;

tex.x = input.lightViewPosition.x / 2.0f + 0.5f;
tex.y = -input.lightViewPosition.y / 2.0f + 0.5f;

depthValue = shadowMap.Sample(SampleType, tex).r;
return color = float4(depthValue, depthValue, depthValue, 1.0f);

However, even if it seems correct, my suspicion is that I need to compare my shadowmap to light Z-buffer so I modified my Pixel Shader and added little tint:

input.lightViewPosition.xyz /= input.lightViewPosition.w;

tex.x = input.lightViewPosition.x / 2.0f + 0.5f;
tex.y = -input.lightViewPosition.y / 2.0f + 0.5f;

depthValue = shadowMap.Sample(SampleType, tex).r;
color = float4(depthValue, depthValue, depthValue, 1.0f);

lightDepthValue = input.lightViewPosition.z - bias;
if(lightDepthValue < depthValue)
{
    lightIntensity = saturate(dot(input.normal, input.lightPos));

    if(lightIntensity > 0.0f)
    {
        // Determine the final diffuse color based on the diffuse color and the amount of light intensity.
        color += (ambientColor * lightIntensity);

        // Saturate the final light color.
        color = saturate(color);
    }
}
return color;

I may not be right, but without this checking, you would get errors with more complicated shadows. Because my previous PixelShader code basically just project shadowmap on plane without any actual depth checking. Therefore it would lay shadowmap on everything it can without even checking if that's correct.


What kind of transformation when the aspect is changed for a rotated shape?

Let's talk about linear geometric transformations in homogeneous coordinates.  In your question there are mainly two kind of geometric transformations involved: planar rotation and planar scaling, the computation-friendly analytical re-definitions of which are found hereUnified frameworks of elementary geometric transformations.

The 2D homogenous rotation matrix that rotates 2D objects clockwise by $\theta$ radians is as below:

$$R=\left[
\begin{array}{ccc}
 \cos \theta & \sin \theta& 0 \\
 -\sin\theta & \cos\theta & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]$$

A 2D scaling due to unequal aspect ratio $\rho$ in this case is a stereohomology with $x$ axis $(0,1,0)^T$ as the stereohomology interface (a line of invariant points), and the normal direction $(0,1,0)^T$ as its stereohomology center per definition 3.9 in page 12 and line No. 5 in Table 1, page 9 of the article Unified frameworks of elementary geometric transformations. Then its scaling matrix can be obtained as below:

$$S=\color{gray}{\left[
\begin{array}{ccc}
 1 & 0 & 0 \\
 0 & 1 & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]+\left(\rho-1\right)\dfrac{\left[
\begin{array}{c}
 0 \\
 1 \\
 0 \\
\end{array}
\right]\cdot \left[
\begin{array}{ccc}
 0 & 1 & 0 \\
\end{array}
\right]}{\left[
\begin{array}{ccc}
 0 & 1 & 0 \\
\end{array}
\right]\cdot \left[
\begin{array}{c}
 0 \\
 1 \\
 0 \\
\end{array}
\right]}=}\left[
\begin{array}{ccc}
 1 & 0 & 0 \\
 0 & \rho  & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]$$

When aspect ratio normally is $\rho=1$, suppose under this circumstance the original object has the homogeneous coordinates of $n$ points in $3\times n$ matrix $X$, and the transformed object in $3\times n$ matrix $Y$, then the geometric transformation between them is the rotation you have mentioned:

$$Y=R\cdot X$$

When the aspect ratio abnormally is $\rho\ne 1$, both $X$ and the $Y$ after rotation are as a matter of fact scaled along the $y$ axis:

$$X'=S\cdot X$$
$$Y'=S\cdot Y$$

Now you want to find the geometric transformation $T$  which satisfies $Y'=T\cdot X'$. Since:

$$Y'=S\cdot Y=S\cdot R\cdot X =S\cdot R\cdot S^{-1}\cdot X'$$

then

$$T=S\cdot R\cdot S^{-1}$$
$$=\color{gray}{\left[
\begin{array}{ccc}
 1 & 0 & 0 \\
 0 & \rho  & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]
\cdot
\left[
\begin{array}{ccc}
 \cos \theta & \sin \theta& 0 \\
 -\sin\theta & \cos\theta & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]\cdot
\left[
\begin{array}{ccc}
 1 & 0 & 0 \\
 0 & \dfrac{1}\rho  & 0 \\
 0 & 0 & 1 \\
\end{array}
\right]=}\left[
\begin{array}{ccc}
 \cos \theta & \dfrac{\sin \theta}{\rho } & 0 \\
 -\rho  \sin\theta & \cos \theta& 0 \\
 0 & 0 & 1 \\
\end{array}
\right]$$

It is a rotation under the generalized meaning.


GPU Architecture: What are the theoretical limits on triangle throughput in modern GPUs?

Your question is not that easy to answer, as there is no standard definition of what "triangles/second" actually means.

What shading is used? What calculations are done? How are the triangles represented? This and much more has to be taken into account when looking for benchmarks on GPU-performance.

To actually compare the performance of GPUs you would have to look at a range of different measurements. GFlop/s and memory-speed as some basic ones, but also look at full benchmark scores (of which there are I am sure many).

There are websites that compile these measurements and display them in a comparison mode. A quick search brought me to this page which seems to have everything needed (I entered an example search).

This helps to compare the performance of two GPUs but if you want to know if a specific GPU is able to do the specific thing you want it to do, that may be very hard to answer and impossible without much more knowledge about what you want to do.

Calculating the maximum performance is also very hard because I believe the GPU vendors do not share everything they do on their GPUs so it can be hard to estimate the actual performance.


What would be the correct way to calculate saturation in this case?

There is unfortunately no good answer to this question. Simply it wont work. There is no good way to define colorful, it this context. Cie is trying to capture the physical measurement. It however does not succeed very well in relating the colors to each other.

Colors on the very outer arc represent spectral distributions of close to Dirac delta function. So one could construct a model that says that a color is very colorful when it is a Dirac delta.

There is a unforeseen consequence of this definition though. Namely the magenta colors do not exist as Dirac Deltas. As these colors do not exist in the spectrum. So they consist of mixture of 2 wavelengths only. This would mean they are less colorful than most other colors.

Other problems

Unfortunately, xyY is not perceptually uniform. So a straight line on the xyY does not represent interpolations between 2 color mixtures. Therefore making a polar transformation means you will have different color bases on same coordinates. Also precieved color do not really move over to your model. To do this properly you would need to do a extremely sophisticated transformation.

There are many problems with converting color to a polar coordinates in that that is exactly contrary how vision works. White is also a bit problematic in this context. The distance to full saturated signal is different for each of the 3 different cones in the eye. Hell, even what is while depends on the surrounding colors and ambient color conditions. So aim afraid your trying to force a worldview that does not exist.

Finally

What would this be useful for?


Interpolate vertex attributes with $z$ AFTER homogeneous divide

The $Z$ of NDC space is related to $1/Z_\text{view}$ but not the same. With a typical projection matrix,  they're related by an affine remapping,
$$
Z_\text{NDC} = a \, \frac{1}{Z_\text{view}} + b
$$
where $a, b$ are some constants related to the near and far planes. So, in general you wouldn't be able to substitute $Z_\text{NDC}$ for $1/Z_\text{view}$ in the interpolation expressions. The $a$ factor isn't a problem, as it would factor out and cancel, but if $b \neq 0$ then it will produce some extra terms and throw off the interpolation.
If I'm hearing you right, it sounds like you're concerned about the idea of needing to prepare every single vertex attribute value (multiplying by $1/Z_\text{view}$) before interpolation and then reverse that afterward. Well, there is an easier way: we can create perspective-corrected barycentric coordinates $w'_0, w'_1, w'_2$, defined by
$$
w'_i \equiv \frac{w_i / Z_i}{w_0/Z_0 + w_1/Z_1 + w_2/Z_2}
$$
Then we can interpolate any attributes by $w'_0 A_0 + w'_1 A_1 + w'_2 A_2$, and it will be perspective-correct interpolation. You can verify this is equivalent to your original interpolation formula in your question. If we have more than three scalar attributes to interpolate (which we usually do), this would be the cheaper way to do it.


Naive Path tracer produces much darker images than Explicit light sampling

It looks like you're clamping all your samples to 0-1 in line 215. Apply clamping only when displaying the image, not when accumulating samples.


Way to estimate the amount of overlap between a rectangular prism and a frustum?

The typical prism/obb/aabb is made up of 6 planes, and the frustum is also typically made up of 6 planes.
To generate a new volume that is the intersection of the two original volumes, use plane clipping to clip one against the other. This technique is generic enough that it can be used with any object like a mesh.
That is take all the planes from the frustum then clip them one by one against the prism. After each clip operation, the new volume becomes the input to the next clipping operation.
What is left will be the clipped volume representing the intersection of the two. Compute the volume of the new shape and divide it by the original volume and you will get a value that is somewhere between [0,1].
This technique is often used to generate decals in games (though the volume computation is not done). Searching the web for decal generation alogorithms/tutorials should help you with the implementation details.
And if you have specific questions I am happy to add more details to this answer.


What are the best ways to start graphics programming?

There are many possible answers depending on your OS, goals and wanted level of abstraction. For low-level programming you can start by learning OpenGL, Metal or Direct3D. Most used shader languages are HLSL based. For higher level there are things like ShaderToy and engines like Unity (C#) or Unreal (C++) where you can experiment with graphics coding.


You are right to be confused. What I think they should have written:
$$
L( x \leftarrow \Psi ) = L_{in} \delta(\Psi - \alpha)
$$
using $ \alpha $ instead of $ \Theta $, which is already used as a dummy variable in the integral.

You should look up the Dirac delta function to learn its meaning and properties. In this context, you can imagine the $ L$ above as representing a very concentrated beam (a laser) coming from the angle $ \alpha $. Practically, to do the integral over $ \Psi $ with $ \delta (\Psi - \alpha) $ present in the integrand, remove the integration and replace all occurrences of $ \Psi $ with $ \alpha $. Then it should be clear how they arrive at the next line, which should read, for all $ \alpha $: $$ \int f_r(x, \alpha \rightarrow \Theta) \; \cos(N_x, \Theta)\; d \omega_\Theta \leq 1 $$.

The fact that this condition is sufficient follows from two facts: 1. that any function (e.g. $ L $) can be approximated by a sum of many $ \delta $ functions, and 2. that everything is linear. In other words, if I write $ N(L) $ and $ D(L) $ for the numerator and denominator of the left hand side of (2.21), then you can see that if $ L = L_1 + L_2 $, then $ N(L_1+L_2) = N(L_1) + N(L_2) $ and $ D(L_1+L_2) = D(L_1) + D(L_2) $. So if I know $ N(L_1) \leq D(L_1) $ and $ N(L_2) \leq D(L_2) $, then I know $ N(L) \leq D(L) $ since $$ N(L) = N(L_1) + N(L_2) \leq D(L_1) + D(L_2) = D(L) $$.


The damping force you mentioned $f=-k \frac{\dot{l} \cdot l}{|l|} \frac{l}{|l|}$ is a special case of $f=-k \dot{C} \frac{\partial C}{\partial \mathbf{x}}$.

Let
$$
\begin{align}
C(\mathbf{x}) &= \lVert\mathbf{l}\rVert-l_0 \\
\mathbf{l} &= \mathbf{x}_1 - \mathbf{x}_2 \\
\mathbf{\dot{l}} &= \mathbf{v}_1 - \mathbf{v}_2 \\
\mathbf{v} &= \dot{\mathbf{x}} = \begin{pmatrix}\mathbf{v}_1 \\ \mathbf{v}_2\end{pmatrix}
\end{align}
$$
where $l_0$ is the rest length of the spring, and $\mathbf{x}_1$, $\mathbf{x}_2$ are the positions of two end points.
The damping force is
$$
\begin{align}
\mathbf{f} &= -k \dot{C} \frac{\partial C}{\partial \mathbf{x}} \\
&= -k (\frac{\partial C}{\partial \mathbf{x}})^T \dot{\mathbf{x}} \frac{\partial C}{\partial \mathbf{x}}\\
&= -k (\mathbf{l}^T\mathbf{l})^{-1/2} \mathbf{l}^T \begin{pmatrix}\mathbf{I}& -\mathbf{I} \end{pmatrix} \dot{\mathbf{x}} (\mathbf{l}^T\mathbf{l})^{-1/2} \begin{pmatrix}\mathbf{I}& -\mathbf{I} \end{pmatrix}^T \mathbf{l} \\
&= -k (\mathbf{l}^T\mathbf{l})^{-1} \begin{pmatrix}\mathbf{l}^T& -\mathbf{l}^T \end{pmatrix} \mathbf{v} \begin{pmatrix}\mathbf{l} \\ -\mathbf{l} \end{pmatrix}\\
&= -k \frac{(\mathbf{l}\cdot \mathbf{\dot{l}}) }{\lVert \mathbf{l} \rVert^2}
\begin{pmatrix}\mathbf{l} \\ -\mathbf{l} \end{pmatrix}
\end{align}
$$

This matches your damping force equation. The $-\mathbf{l}$ part is just the force for the other particle, which has the opposite direction as expected. Bold letters are used to denote vector quantities.

As mentioned in [baraff98], this damping force is along the direction $\frac{\partial C}{\partial \mathbf{x}}$,
and its magnitude is the projected component of velocity $\frac{\partial C}{\partial \mathbf{x}}^T \mathbf{v}=\frac{\partial C}{\partial \mathbf{x}}^T \mathbf{\dot{x}}=\dot{C}$.


Artefacts on top and bottom of utah teapot

I'm fairly certain your problem lies in not handling so-called degenerate Bezier patches correctly, in particular (as Joojaa noted), the computation of the surface normal.  I say "so-called degenerate" because, geometrically, the surfaces are often perfectly well behaved. It's just that some assumptions people frequently make regarding the parametric equations may not hold.

Books such as Gerald Farin's, snappily titled, "Curves and Surfaces for CAGD. A Practical Guide" will give more details, but I'll try to summarise two simple cases. Now assuming your Bezier is defined as $\bar{B}(u,v)$ the usual two causes of problems are:

Zero derivatives:

To compute the normal at $(a,b)$ one normally (pardon the pun) computes the two derivatives,

$\frac{\partial }{\partial v}\bar{B}(a,b)$ and $\frac{\partial }{\partial u}\bar{B}(a,b)$

or scaled versions thereof, to obtain two tangents and then take the cross product.

(Implementation note: Since we can use a scaled version of the tangent in the calculation, we really don't have to calculate the actual derivative. For example, especially at the corners, differences of control points can yield a scaled tangent. However for brevity in this discussion we will assume the actual derivatives).

A common occurrence, at least at the corners, is that the 1st partial derivatives at the location can be zero, which leads to an incorrect normal, i.e. a zero vector.

In the case of the tops and bottoms of the teapot, one whole edge of a number of the (bicubic) Bezier patches has been collapsed to a point, i.e. all 4 control points are the same, and thus, say

$\bar{B}(0,v)==\bar{B}(1,v)$   and   $\frac{\partial }{\partial u}\bar{B}(0,0)=\bar{0}$.

The surface, however, is completely well behaved so you can simply choose another derivative that starts at that collapsed point. In this case, say, choosing

$\frac{\partial }{\partial v}\bar{B}(1,0)$

for the second tangent.

Having said this, you still have to check that your first derivatives are not zero for another reason (e.g 2 or 3 coincident control points), in which case, you can fall back (thanks to L'Hopital's rule) on the second (or if that's zero, even the third!) derivative(s) to obtain valid tangents.

Parallel tangents:

Another similar problem can arise if your two tangents, are parallel. - Farin has a good example in his book.

In this case, I think you may need to look at using something like $\frac{\partial^2 }{\partial u \partial v}\bar{B}$ or, possibly, just fall back to using a small UV offset to approximate a vector.


How to interpolate vertex normals for Phong lighting when barycentric coordinates are negative?

I have figured out the solution to this, and want to share it. What I am doing is I calculate different barycentric coordinate for each of the step I need to do: z-buffer, interpolation of normal, etc.

The correct way is to use the barycentric coordinate that is calculated when checking z-buffer to test whether a pixel inside a triangle. Because that barycentric coordinate is guaranteed to be inside the [0, 1] range (otherwise the pixel wouldnt have been drawn), it is the correct barycentric coordinate to use.


How is Raster and Random Scan Display different?

A raster display draws every pixel on the screen in every frame whether there is something to show or not.

A random scan display activates only those pixel which are occupied by an geometric primitive.

So yes, they both use pixels, but the difference is in how they draw the pixels onto the screen.


BC7 texture with alpha has the same size as BC7 without alpha

BC7 always effectively encodes 4 channels. That's part of the nature of the format. Turning on alpha means that the compressor will read and use the alpha value in your source image, but that's all.
All BC7 compressed images of a given dimension have the same size, since BC7 compresses each 4x4 group of texels into 128-bits. But that doesn't mean they have the same quality.
BC7 is an adaptive encoding format. While each block has a fixed size (both spatially and bit-wise), each block stores how many channels that block actually compresses and the way that they correlate with one another. The specifics are complex, but the general upside is that the compressor figures out on a per-block basis what manner of compression works best for a particular block. If the compressor can look at a block and see that the alpha for all texels is 1.0, then it can spend more bits of the 128 bits on the colors and fewer on the alpha.
As such, a texture "without alpha" is exactly like a texture with alpha where every texel has an alpha of 1. These blocks will spend the majority of their bits on the non-alpha parts of the texels, only wasting a few bits on a "mode" field that lets the decompressor know that the alpha is a constant 1.0.
This also means that, if your texture has a varying alpha, but parts of it have a constant alpha, BC7 can encode these constant alpha areas with fewer bits, so that more bits can be spent on improved color encoding.
So you shouldn't feel that there is a problem because the two compression switches came out with the same number of bytes.


PBR - How to handle multiple BSDFs and material type

Although the specular term of the Cook-Torrance BRDF is widely used in Path Tracing, the original BRDF was meant initially for rasterization (it was proposed before the advent of Path Tracing). The idea of summing weighted diffuse and specular components (i.e. $R=dR_d+sR_s, d + s \leq 1$) is not capable of representing the complex results that may emerge from the interaction between diffuse and specular material elements (e.g. refraction, Fresnel for plastics, etc).

Actually, a better physically based model for material representation would be that of layers. For instance, smooth and rough conductors could be modeled through a unique layer represented by the specular term of the Cook-Torrance BRDF. Smooth and rough dielectrics could be modeled through a unique BSDF. Plastics, on the other hand, could be modeled through the use two layers: the topmost layer representing a dielectric (a BSDF) and the second layer, underneath the dieletric, representing a diffuse material (e.g. a lambertian or Oren-Nayar BRDF). That two-layer setup for plastics is capable of simulating the Fresnel reflectance of plastics (through the dielectric layer), the plastic color (through the light reflected by the colored diffuse layer underneath the dielectric), the light absorbed by the dielectric (through the implementation of Beer's Law), etc.

The layer-based material model approximates quite well what actually happens in real life. I have already written a related answer, with more details, and that may help you (I've included some references in that answer): https://computergraphics.stackexchange.com/a/5761/5681

The original Cook-Torrance does not properly conserve energy (I've already posted a question about that here). In 2007, Walter et al. have published the paper Microfacet Models for Refraction through Rough Surfaces, where they propose a Cook-Torrance based model amenable for Path Tracing.

Regarding layered material, I've implemented the system presented in the paper by Weidlich and Wilkie,
Arbitrarily layered micro-facet surfaces, where a quite general layered material model is presented. After that, more sophisticated layered material models have been proposed. One of them is
A comprehensive framework for rendering layered materials
, by Wenzel Jakob. There is a new paper (published in 2017, I think), that supports subsurface scattering for top layers (I will add it here when I find it).

I think that PBRT uses layer-based materials too, but I do not know how they have actually implemented it.