GLSL. Can someone explain why gl_FragCoord.xy / screenSize is performed and for what reason? I'm new to shaders and know that you can color pixels with gl_FragColor but sometimes there is this thing:

vec2 uv = gl_FragCoord.xy / screenSize; // or resolution or smthn depends on the implementation

If gl_FragCoord is like pixel coordinates, what does uv get?

Why is it often done in GLSL?

If someone could even draw an example of which part of the screen will be UV it will be very helpful!

First, gl_FragCoord.xy are screen space coordinates of current pixel based on viewport size. So if viewport size is width=5, height=4 then each fragment contains:

Why are uvs needed? For example I rendered geometry to screen quad and then I need to apply some postprocessing on this quad in another rendering pass. To sample from that quad I need texture coordinates in range [0, 1]. So to calculate them I perform division gl_FragCoord.xy / viewPortSize.


Frustum aligned volume integration - aliasing issue - temporal reprojection Raw volume raymarching especially aligned with frustum seem to produce a lot of aliasing/banding on edges of shadow map. Especially when fragment pos is far because step length is then getting bigger. Making some research how to correct it I found out about temporal reprojection technique.

So basically I have an two 2D textures, current and previous where I store exponential moving average of raymarching result of previous frames like:

$s(n) = w*f(n)+(1-w)*s(n-1)$

where:

$w$ - weight

$f(n)$ - current frame calculation

$s(n-1)$ - result of previous frames calculation

This seem to be working when camera is moving because then bands are moving and blending with each other. But when camera stands still everything returns to normal. I tried to:

add some noise aligned with ray direction to starting position every other frame but banding is still at the same positions
add noise not aligned with ray direction but then bands are shifted ramdomly and raymarching results are wrong

Here (FrostBite engine) it is smooth constantly independent of camera movement.

So is there a way to simulate camera movements? Or any other technique to correct the issue?

Here you can download some content describing subject.

You can look into alias-free volumetric sampling algorithm by Huw Bowles for potential solution to the ray marching aliasing issues. The basic idea is to snap your samples to planes based on the ray direction that's best explained with this Shadertoy demo.


OpenGL with SFML, create an n-pointed star? Completely stuck on how to do this. Yes, it is a homework problem, but its the first we've had in the class so not too familiar with graphics. Also, seeing as its homework, please don't just give the answer, just some pointers in the right direction would be great.

I've tried to do it as a triangle fan, and as a polygon but that just ends up in really weird shapes.

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

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


What is the definition of a Half Edge? I ran into half-edges as part of the result of a Delaunay Triangulation library, but can't find an actual definition of what one is.

I understand vertices, edges and faces, but can't find a concrete definition of a half-edge.

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

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

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

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


Consistent animation speed across systems I recently noticed that in all games, no matter what system you play on or what specs you have, your character will always move at the same speed even if there is  no frame rate limit. I was wondering how this worked because if you were to add the same value to a variable, the animation would be faster on faster processors or graphics cards. I heard somewhere that you could use the time between frames to reduce the speed. I do not see how this would work. So how would you get a consistent speed across all systems?

It works exactly the way you describe. Most 3D games work on keyframe animations, where the human animator sets the important positions (for each joint or whatever is being animated) and the times they should be shown. The times are usually a number of "frames", but these are just intervals of 1/60 s, regardless of what frame rate is actually achieved on the client machine. The computer interpolates between those keyframe positions to decide what to draw at any given time.

For example, let's say an animation has keyframes on frame 0, 6, and 12. These are the positions to show at the beginning of the animation, 6/60 s later (that's 0.1 s), and 6/60 s after that. The first frame is drawn with the pose from the first keyframe.

Let's say the game is running really slowly, and the next frame starts being assembled 0.05 s later (which would give a 20 fps frame rate). The animation system uses interpolation to make a pose that's halfway between the first and second keyframes, because the time is halfway between 0 and 0.1. Next, imagine things are speeding up a bit, and the next frame is drawn 0.04 s later. It's 0.09 s since the start of the animation, so the animation system makes a pose that's 90% from the first keyframe to the second keyframe.

If this frame rate continues, the next frame will be generated at 0.13 s. This time the animation system has to interpolate between the second and third keyframes, taking 70% of the second and 30% of the third. Note that the second keyframe was never drawn in this example, but it contributed to the poses that were drawn, through interpolation.

Now to finish off the example quickly, let's assume the computer really slows down so the next frame is another 0.07 s later, making 0.2 s since the start of the animation. This is exactly the time that the third keyframe was set for, so the animation system just uses that pose. Although the animation was designed to take 13 frames (including the one at the start), we've finished it in just 5 frames of jerkiness, but it still takes the same length of time, 0.2 s.

In summary, it doesn't matter what time each frame actually comes out on the end system, because we can make whatever pose we like via interpolation.


About view matrix and frustum culling I was trying to determine if an object (sphere) is inside a view frustum. My strategy was first to get the view matrix:

glm::lookAt(cameraPosition, lookAtPosition, cameraUp);

Where cameraPosition is the position of the camera, lookAtPosition is calculate by cameraPosition + cameraDirection and cameraUp is pretty much self-explanatory.

After obtaining the view matrix, I calculate the 8 vertices of the frustum in the camera view, then multiply the inverse of the view matrix by each point to get the coordinate of them in the world space.

Using these 8 points in the world space, I can construct the 6 planes by using cross products (so their normals will point inward). Then I take the dot product of the object vector (I get the vector by subtracting the object position with an arbitrary point on that plane) with the normal vector of each plane, if it is positive for all 6, I know it's inside the frustum.

So my questions are:

Is the camera view is always set at (0,0,0) and looking at the positive z direction?
The view matrix transform the world view into the camera view, if I use the inverse of it to transform the camera view to the world view, is it correct?

My frustum is opposite with my cameraDirection which causes nothing to be seen on the screen (but it still drawing stuff in the back of my camera).

It is  common that the camera is pointing towards negative Z, in order to satisfy right hand rule. So you need to know the convention.

As for theinverse matrix of the camera matrix. Sort of yes. You can do it. However the camera matrix is a homogeneous space. While generally our other matrixes are affine spaces. If you take care of the distinction and do the proper camera divide, if tge camera needs one. Then yes you could do that.


The generations of Intel's integrated graphics Which among Intel's integrated GPUs depart from their predecessors in a major way? Which are merely incremental improvements?

BACKGROUND

Intel integrates GPUs with many of the processors they sell. With each generation of processors (Ivy Bridge, Haswell, Broadwell, Skylake, Kaby Lake, etc.), Intel updates the graphics architecture. Much is written in the popular press regarding each of these: does your favorite video game make only 98 frames per second on the latest, or does it achieve a previously unheard-of 99? Bar charts are published, and so on.

None of the articles, or practically none, seems to grasp the underlying technology. I would be interested to learn whether Ivy Bridge, Haswell, Broadwell, Sklyake, Kaby Lake, and so on, are each merely incremental improvements in the same, basic graphics architecture, or whether, say, Broadwell (or Kaby Lake, or whatever) was an entirely new design.

I would likewise be interested to learn if one generation had been upgraded over the previous in a way which made OpenCL-type usage much more practical or capable.

Apparently, Broadwell was the first capable of Vulkan, but I do not know whether this capability represents a fundamental architectural change or if it's just a driver issue.

Most of the popular press seems pretty clueless at this technical level. Can you shed any light?

You want to look up "Intel Graphics Developer Guides", they are quite thorough.  Here's a link to the latest:

https://software.intel.com/sites/default/files/managed/22/6d/6th_gen_graphics_api_dev_guide_dec_2015.pdf

And here's the architectural Cliff Notes for the past five generations of Intel GPUs:

"In several important respects, the 6th gen Intel Core processors' graphics hardware architecture is similar to that of 5th generation Intel® CoreTM processors"

"The graphics hardware architecture of the 5th gen Intel Core processors is similar to that of the 4th gen processors."

"The graphics hardware architecture in the 4th gen Intel Core processor is similar to that in the 3rd gen processor."

"The new generation of microarchitecture, codenamed Ivy Bridge, provides another jump in functionality and performance over Sandy Bridge microarchitecture."

"With the introduction of the Intel microarchitecture codenamed Sandy Bridge, the graphics processor has moved onto the same die as the CPU, and is now referred to as “processor graphics”. In addition, processor graphics has enjoyed numerous architectural improvements that yield significant performance improvements over previous generations of Intel® integrated graphics parts."


PDF when sampling from Rectangular Area Lights So I was reading up on how to sample Rectangular Area lights, and the procedure is pretty straightforward. That is, I can get random points using the equation.

$ p = X_o + \epsilon_1V_1 + \epsilon_2V_2$

Where $X_o$ is my corner point for the rectangle and $V_1$ and $V_2$ are the 2 vectors spanning it And the epsilons are the uniform random numbers.

The thing is almost every where I go, they talk about how the PDF of this sampling is $1/A$ (constant) with the exception of 1 place discussed below, and that when calculating the radiance coming from this area light when sampled by this way we should thus multiply the whole thing by the Area of the rectangle.

However when I check their algorithms, they leave out this Area factor.

Compute Graphics Principles and Practice and Fundamental of Computer Graphics by Peter Shirley describes the calculation of the radiance returned by sampling a point on the area light, (Notice the black underlined factor describing Area of the jth light source)

However both of these books when describing their pseudocode, leave out this $Area$ factor.

However when I read Peter Shirley's research paper regarding Direct illumination techniques, He gives the PDF for rectangular Area lights as

$ 1/ \frac{||V_1 \times V_2||}{2}$

Which is the area of the triangle not a rectangle.

So ultimately I have 2 questions.

1) What is the correct PDF term for rectangular area lights?

2) Do I need to multiply this term when calculating the radiance emitted from this light source.

The correct PDF for an area light is indeed $1/A$. The reason why you see the area term at the top is because you are importance sampling the light, meaning you need to divide by its PDF. This is a Monte Carlo estimator with one sample. Shirley simplifies this by simply writing $1/(1/A) = A$.

I am not sure how he got the $1 / \frac{\|V_1 \times V_2\|}{2}$.


Is there a way for me to record every shading sample for a deep texture Im developing a deep texture processing engine. Right now I have a software renderer to generate the raster samples. Now, for bigger renders it would be nice to do this on hardware.

Is there a standard mechanism on how to capture the data of each shader evaluation in OpenGL? I could do some sort of depth peeling but it seems to me as a bit too brute forcelike approach. I don't mind the time cost just as long as it does not take ages to implement the infrastructure needed.

If you really want to obtain every fragment (shading sample) within the camera frustum, regardless of occlusion, one way to do it would be to disable depth testing, and then use a large SSBO together with an atomic counter to append data into the buffer from the fragment shader.  Each fragment could then store a record of its screen position, depth, color, and any other desired information into the buffer to be read out later.  The buffer has to be large enough to store every fragment rendered, as it will throw away any extra data once it's full.

Using one huge SSBO for the whole screen likely isn't ideal, as there will be huge contention on its counter with every fragment trying to increment it.  A better way might be to slice up the screen in tiles with an SSBO+counter per tile.  There is a limit of 16 bound SSBOs on all current GPUs, and on AMD GPUs, there is a further limit of 8 atomic counters (source: Sascha Willem's OpenGL Hardware Database).  Those limits will determine how many tiles you can create. This will reduce contention and also ensure that even if you run out of memory in one tile, you'll still have valid data in the others.


Turn an enclosed region into a triangle mesh I have a (2d) region defined as a list of vertices that I'm currently drawing as a line loop in opengl. I want to draw this region solidly. I started by using triangle fan, but that doesn't work if the region has vertices positioned at the same angle (relative to the center) but different radii, because then that whole region is draw solid.

So, I think I just need some general method to convert an enclosed spline/list of points into a triangle mesh. Maybe there's a better way?

If the number of vertices is not too large, I would suggest looking into "ear clipping". It's a straightforward method to triangulate possibly-concave polygons, as long as they don't have holes. However, it's $O(n^2)$ in the number of vertices, so may not a great choice if you need to handle a large number of verts. Dave Eberly has a paper on ear clipping and also an open-source implementation in C++. There are many other links available on the web, too.

The basic idea is to find an "ear" of the polygon, which is defined as a triangular subset of the polygon, formed by three consecutive vertices of the polygon. Once you find an ear, you add it to the output triangle list, and remove it from the original polygon (which reduces the number of vertices by one). Repeat this until the polygon is gone.

The reason it's $O(n^2)$ is that you have to repeatedly search the list of vertices to find another ear to remove. There are more efficient triangulation algorithms as well, but they're more difficult to implement.


Encode two R11G11B10F buffers into one using Y'CbCr? We are trying to optimize Lighting in fillrate. Integrated graphics chips are very slow in memory access.

There are two Light Buffers: one for Diffuse and one for Specular. Both use format R11G11B10F. Is there any way to pack two buffers into one? For instance, using color model Y'CbCr:

You might want to look at The Compact YCoCg Framebuffer. It uses a 2-channel buffer to store luminance for every pixel and the two chroma components in half the pixels each, forming a checkerboard. It also uses an edge-aware upsampling filter at the end of the frame to reconstruct the missing chroma components and convert back to RGB.

You could extend this to a single 3-channel buffer to store the diffuse luminance, specular luminance, and checkerboarded chroma (assuming based on the diagram in your question that it suffices to store just one pair of chroma values).

You might also consider using R10G10B10A2 format instead, unless you really need the range of the floating-point format. It's the same amount of bandwidth but will give you a good deal more precision (R11G11B10F only gives you 5–6 bits of mantissa per component). Plus, you can use the 2 bits of alpha to store a couple of flags, such as your diffuse/specular flag.


How can I draw regular polygons along the horizontal axis? I drew simple triangle wave (that you can see It on graph toy).

   f1(x) = frac(x/(PI+5))*(PI+5)
   f2(x) = frac(-x/(PI+5))*(PI+5)
   f3(x) = min(f1(x),f2(x))*1.5

I want to change above shape like Approximating Pi and repeat It periodically along the horizontal axis.

then I want to cut negative parts of the vertical axis:

Update
Note that I know how can I draw It by shader In ShaderToy.I Just need drawing mathematics In desmos or graphtoy,because I want Extract values of wave.I want to create polygon wave form.
I drew Periodic shape that I want to change the number of polygons
Desoms Link.

You can use the equation of a circle to do this. A circle can be drawn using the equation:

$r^2 = x^2 + y^2$

You can solve for the upper half of this by rewriting it in terms of $y$ like this:

$y^2 = r^2 - x^2$

so

$y = ±√(r^2 - x^2)$

But we only want the positive half. So we say just:

$y = √(r^2 - x^2)$

You can then repeat it by using (((x - 2) % 4) - 2) as the input value. In the end the equation is this:

$y = √(1^2 - (((x - 2) mod 4) - 2)^2)$

It looks like this:

In graph toy, it would be:

sqrt(1 - pow((((abs(x) - 2) % 4) - 2),2))


"Anchor" vehicle to 3D path with OpenGL / SDL As a CG assignment with OpenGL and SDL, I want to create a clone of the videogame WipeOut instead of the classical car racing game. Examples of motion in 2D have been given and are rather simple to understand; the complexity in 3D is that the road is not planar - when the vehicle approaches a cliff or a vertical loop, it should be translated on the y-axis and rotated to be mantained parallel to the road. I have no idea on how to get and update such informations (translation and rotation for the road). I thought about creating a bezier curve with informations about its derivative (overlapping the road) and keeping track of the movements of the vehicle with respect to the curve, but I don't know if it's the right (simplest) way of doing this. By the way, being able to track the position of the vehicle relative to the curve would also solve some collision problems by limiting its distance.

Of course, this feature is not mandatory - I could just keep a planar road - but I'd like to understand how to do it. Suggestions?

I think the most common way to have a the "up and down" movement of a vehicle is to actually build a physics simulation.
It does not have to be a very sophisticated one and how much realism you want is up to you.

But just in each frame update the velocity of the car in downward direction (gravity) until it hits the road. Then if the road has an upward slope it will also push the car up (by colliding with the car).

Now you have a new problem and that is the problem that the car remains parallel to the ground-plane. And the physical calculation behind this is a bit more complex.

When detecting if your car hit the ground you would have to detect the exact position of the hit (which part of the car hit the floor?). Then using the gravity, mass and some material properties (all of which you could just set globally and assume as easy values for simplicity) you can calculate how the reaction force to the collision will be. And you know that the force will be exerted on the collision point.
Now when you have the center of mass of the car, and the collision point, and the force created by the ground-collision, you can calculate how to rotate the car.

This is only a very rough sketch of how I learned to handle such things in a "real" way. I hope this was not too confusing, if you have additional questions I can rework the answer later.


Are there LOD Algorithms that work through removing objects? I want to implement a LOD algorithm that reduces the number of objects rather than the number of vertices of the objects.

Imagine a forest, then close up I want to see all the trees but moving away from it there is no need to display every single tree, so I want to remove random trees (up to a threshold).

So far I wasn't able to find current algorithms for this, but that may just be due to a lack of terminology: I have been looking for LOD algorithms. Maybe this technique exists under a different name, so I wanted to ask if there are such algorithms out there already?

Here's how I'd do it:

Assign an sequential id to each tree. This can be an implicit id like an index in an array if you have this structure. The important part is that you can create N stable groups of trees that will be removed based on distance.

Then decide at which distance you want to start removing each group of trees. Each tree group will have a specific distance at which it gets removed from the rendering. You must choose the appropriate granularity based on your use-case.

Then, for each tree, compute its distance from the viewpoint and compare that with its removal distance to decide if it should be rendered or not. To get a nicer effect, use a fade-out distance over which the trees progressively disappear.

Essentially the process is to tag trees with some LOD parameters (removal distance, fade distance, etc.) and remove them based on distance from the viewpoint.


Expressing a reflected ray at a point as a function of the incoming ray and the normal to the surface at that point Said L the incoming ray (at a point $P$).  $R$ is the reflected ray. $N$ is the normal to the surface at the point $P$.  $\alpha$ is the angle between $N$ and $L$ (and $N$ and $R$ also). It is assumed a smooth surface (like a plane). $P_{R1}$ is the projection of the vector $R$ on the plane.  $R_2$ is the vector that joins the vector $R$ and $P_{R1}$.
Now you have,

$R= P_{R1} + R_2$  (sum of vectors)

$P_{R1} = -L+R_2$

Thus:
$R=2R_2-L$

Now, exploiting the property of dot product (namely scalar product)

$R_2 = (N\cdot L)/|N|$          (because $R_2$ is $L\cos(\alpha)$ because is the projection of $L$ on $N$)

where $N\cdot L$ is the dot product and $|N|$ is the magnitude (the module) of vector $N$.
 But the vector $N$ is the normal that typically has $|N|=1$
Thus:

$R = 2(N\cdot L)-L$

But my book (3D Computer Graphics, Alan Watt) tells that:

$R_2 = (N\cdot L)N$

Rather for me is:

$R_2=(N\cdot L)/|N| = N\cdot L$    (because $|N|=1$)

I don't understand this point.

$N \cdot L$ is the scalar length of the projection of $L$ on $N$ (assuming $|N|=1$).

$(N \cdot L)N$ is the actual vector with the length of $N \cdot L$ in the direction of $N$.

There is often a confusion when people say "projection of this on that" whether they mean just its length, or the actual vector resulting from the projection. In this case you want to combine it with another vector, so you need the vector version. The correct reflection formula is thus: $$R = 2(N \cdot L)N - L$$


What are Spherical Harmonics & Light Probes? What are Spherical Harmonics & Light Probes? How useful are they in computer graphics? What exactly do they do? I've heard the word spherical harmonics & Light Probes everywhere, from siggraph presentations to blog posts.

Recently Matt Pettineo posted a 6 part blog series on them but I still don't understand what they are.

Is it another way to improve ambient lighting?

Spherical harmonics

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

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

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

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

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

Light probes

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

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

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

Using them together

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

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


Does a point have a length? I am going though a functional source that has the following class:

 class Point3D
 {
       int x ,y,z;

       int vSize2() const
       {
           return x*x +
                  y*y +
                  z*z;
       }

       int vSize() const
       {
         return sqrt(vSize2());
       }
 };

It is hard to understand the functionalities of the above two functions (vSize2() and vSize()). It seems to me that it is calculating the length between the point and the origin.

Any thoughts?

No, a point does not have a length. A point is only a location - it has no extent in any direction.

You are correct in guessing that the function vSize() returns the distance from the origin to the point. The other function, vSize2() returns the square of that distance. This is used in calculating the distance, and in some cases it may be useful to work directly with the square of the distance instead of the distance. This avoids calculating square roots, so can lead to improvements in speed and accuracy.

Note that as Nero points out, the functions are using only ints (integer variables and return values). The distance between a point and the origin is not in general an integer value, even if the coordinates of the point are all integers. The distance squared returned by vSize2() will be accurate, but the distance itself returned by vSize() will be accurate only to the nearest integer.


What method is used for baking grayscale curvature maps I'm looking for a way to generate curvature maps like this:

The map stores the Convexity and Concavity of each pixel as a value from 0 to 1, where .5 is no curvature and 1 and 0 are the extremes.
I searched for how to calculate the curvature of meshes and found these questions on StackExchange:
What is the simplest way to compute principal curvature for a mesh triangle?
https://gamedev.stackexchange.com/questions/92621/cavity-map-generation
https://blender.stackexchange.com/questions/146819/is-there-a-way-to-calculate-mean-curvature-of-a-triangular-mesh
But I'm not sure if these are about the method actually used in the example.
What algorithm is used to bake these maps, and where can I find more information about it?

This is how I ended up doing it:
Instead of using vertex colors like most of the solutions I found online, I calculated the angle of each edge, and then used a shader to draw the lines in UV space. This Also gives me control over how thick the lines should be.
In some cases the angle calculated isn't right, but that should be easy to fix.

Here is the project: https://github.com/Jummit/curvature-baker


Derivation of Wikipedia's path tracing diffuse BRDF Wikipedia's path tracing algorithm has the following line for the diffuse BRDF:

Color BRDF = 2 * m.reflectance * cos_theta;

I'm a little confused as I've been reading that Lambertian surface BRDF is $\frac{\rho}{\pi}$ where $\rho$ is the reflectivity (albedo?) of the surface. Where does the $2$ come from?

When you perform regular Monte Carlo integration over a hemisphere using $N$ samples, each sample represents $\frac{2\pi}{N}$ steradians. So the Monte Carlo integration for Lambertian BRDF is:

$$\frac{2\pi}{N}\sum_{i=1}^N\frac{\rho}{\pi}L_i*Cos\theta_i$$

For path tracing, you only take one sample per path segment, so because $N$=1, the above sum becomes:

$$2\pi\frac{\rho}{\pi}L*Cos\theta = 2\rho L*Cos\theta$$


Parallax/relief mapping with gl_FragDepth and conservative depth New to writing shaders, please forgive any naivety in the question!

Parallax mapping in any of its ray-marchy variants (e.g. steep parallax, cone step mapping, POM) generates fake geometry on a material, generally nearly indistinguishable from real geometry and at a much higher resolution than is practical using GPU tessellation. However, the geometry created does not extend outside the initial bounds of the triangle being rendered (unlike with GPU tessellation) - so e.g. a parallax mapped sphere still has circular edges.

I'm want to try adapting a parallax shader (preferably OpenGL, any high quality variation of parallax mapping will do though I gather cone stepping is pretty much the state of the art) to write outside of the poly bounds using gl_FragDepth (and enabling conservative depth output).

This seems to be theoretically possible (see this Unity example Modifying depth when parallax mapping | Unity Community). However that post did not include any shader code and got no responses.

Is anyone aware of an existing implementation of this (with source available), or else could give me a good idea of how I should proceed with implementing it?

Many thanks.

Unfortunately, gl_FragDepth won't enable you to write outside of the rasterized triangle's screen-space bounds. All it enables is to modify the values written to the depth buffer within the triangle, replacing the triangle's true depth with a shader-calculated value.

What you're seeing in the screenshot in the Unity forum thread you linked to is two intersecting parallax-mapped surfaces, where the depth output enables occlusion between the two surfaces to be correctly resolved in the screen-space region they both cover. But if you look at the right side of the screenshot, you can see the straight-line edge of the bottom brick polygon.

However! You might still be able to accomplish the goal of getting correct silhouette edges by using discard in the fragment shader when the parallax ray gets outside the bounds of the original polygon. This would involve modifying the ray-march code to check if the march has crossed a triangle edge, and exit with discard if so, rather than continuing to march.

You'd have to find a way to get the information about triangle edges into the fragment shader, though, as this isn't provided by OpenGL. Maybe make your vertex/index buffers available via SSBOs; then you can use gl_PrimitiveID in the fragment shader to index into them and lookup the vertex data.

Unfortunately I don't know offhand of any implementation of something like this, although I'm sure someone out there must have tried it before.


Where does the cosine factor comes from in the GGX PDF? The GGX NDF, as it appears on the paper where it is presented is:

$$D(m)=\frac{\alpha_g^2\space\chi^+(m\cdot n)}{\pi\cos^4(\theta_m)(\alpha_g^2+\tan^2(\theta_m))^2}$$

It is equivalent, in the range $[0,\frac{\pi}{2}]$, to the following formulation (notation simplified by me):

$$D(m)=\frac{\alpha^2}{\pi(\cos^2(\theta)(\alpha^2-1)+1)^2}$$

by finding the derivative of the inverse of the expresions used for importance sampling ($\theta_m=\arctan(\frac{\alpha\sqrt{\xi_1}}{\sqrt{1-\xi_1}})$ and $\phi_m=2\pi\xi_2$) presented in the same paper we can derive the following PDF:

$$p(m)=\frac{\alpha^2\cos(\theta)\sin(\theta)}{\pi(\cos^2(\theta)(\alpha^2-1)+1)^2}$$

When finding the corresponding PDF for, for instance, the Blinn NDF the only difference is a $\sin(\theta)$ in the numerator that comes from the fact that the NDF is a distribution over a hemisphere. So where does the $\cos(\theta)$ factor in the GGX PDF come from?

I think it is explained after formula (37) in the paper but i don't understand the reasoning behind it.

PD:here are the Blinn NDF and PDF as i understand them for reference

NDF: $\frac{n+1}{2\pi}\cos^n(\theta)$  --- PDF: $\frac{n+1}{2\pi}\cos^n(\theta)\sin(\theta)$

EDIT: here is a link to the paper in case it makes the question easier to answer

Normal distribution functions are defined a bit differently than you might expect. They're not strictly a probability distribution over solid angle; they have to do with the density of microfacets with respect to macro-surface area. The upshot is that they're normalized with an extra cosine factor:
$$ \int_\Omega D(m) \cos(\theta_m)\, d\omega_m = 1 $$
This cosine factor accounts for the projected area of microfacets onto the macrosurface. When importance-sampling the NDF this cosine factor must be accounted for as well.

This is actually hinted at in the text of the paper, as it says

  The equations for sampling $D(m) | m \cdot n |$ are:

(just above equations 35–36). The $|m \cdot n|$ factor there is equal to $|cos(\theta_m)|$. So the quantity being sampled is not just the NDF itself, but includes this extra factor. Then, as you note, the sine shows up as part of the area element for spherical coordinates.

I wrote a blog post a few years ago that goes into more detail about this.


What algorithm(s) are behind Google Street View's 3d mapping mode? I recently noticed that in some locations you could switch to a 3d mapping mode in Google Street View.
Which algorithm(s) did they use to generate the 3d models of the streets? I assume it was more than just capture, given the huge amount of data that would be necessary to produce an accurate representation of the streets.

Google said how in 2012, avoiding explaining too many details

The never-ending quest for the perfect map

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

and here it's a video of the feature.

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

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

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

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

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


What's the justification of multiplying colors together? In Phong shading model (and other shading models as well), the light color is multiplied together with the surface color component-wise. I want to know the justification for this. Is this some what related to the physics of light?

When you multiply colors together the surface color acts as a filter. Remember that colors are scaled from 0 to 1.

That way with white incoming light a green object will remain looking green.

This does bring a few issues because physical light isn't as orthogonal (where the 3 primary colors are actually overlapping frequency response functions) but for the vast majority of cases you can tweak the values enough for it to stay looking realistic.


Ray vs AABB algorithm that also gives which side was hit? There are plenty of well-known algorithms for determining if a ray hit an Axis-Aligned Bounding Box (AABB), like Andrew Woo's covered in a response here.

Using exit = origin + (dir * tMax) will give you the exact position of exit using distance through the ray. (and same for tMin and the entry)

However, I need a way to not only do that, but also get if we hit the left, right, top, bottom, front, or rear face of the AABB.

The naive approach would be to do some kind of floating-point comparison given the minpoint/maxpoint and try to see what side (if any) the point exactly equals in each dimension, or exactly equals within a tolerance. That seems error-prone, I was wondering if there was a better way.

Is there an algorithm for this?

Given an intersection algorithm like the one you linked which first generates a $t$ value for every plane of the AABB, then computes $t_{min}$ and $t_{max}$ from those, you simply need to figure out which of the six initial values correspond to your final $t_{min}$. There will be two matches if your ray hits an edge, three if it hits a corner. In that case, this may include a back face (see the fix below).

This will not require a tolerance as $t_{min}$ is always one of the six values. In fact, it is always one of three intermediate values and you can use that to make your check simpler, faster and better. Given code, from the link above, such as:

float tmin = max(max(min(t1, t2), min(t3, t4)), min(t5, t6));

you could rewrite as:

float tminx = min(t1, t2);
float tminy = min(t3, t4);
float tmin = max(max(tminx, tminy), min(t5, t6));
int axis = 2;
if (tmin == tminx) axis = 0;
if (tmin == tminy) axis = 1;

As an added bonus, unlike testing all six values, this will not match the back face in the event that a far edge is hit.


Sphere Texture Mapping Algorithm I am trying to generate UV coordinates automatically to map a texture on a sphere.

The problem is with seams when going from 0.8 to 0.2 there is a discontinuity and interpolation in fragment shader is incorrect. It should go from 0.8 to 1 and then from 1 to 0.2. But instead it will interpolate directly between 0.8 to 0.2.

https://en.wikipedia.org/wiki/UV_mapping

How to handle this discontinuity from 1 -> 0.

Here is what is currently happening at the seams

A few options:

Make the 0-1 discontinuity explicit in the mesh. That is create a 2 sets of vertices that lie exactly on the line where the value would be 1 or 0 (one set gets 1 and one gets 0) then connect the vertices up like they would makes sense.
Switch to a different texture mapping projection, I prefer a cubemap because it minimizes distortion at the poles and you don't need to have seams in the mesh.


Camera rotation effect on pitch, yaw, roll I'm trying to calculate what the effective pitch, yaw, and roll are for an image.
I have the image's original pitch, yaw, and roll.  Then the image (camera) is rolled at some arbitrary angle.  I want to calculate the new relative pitch, yaw, and roll.

I know if pitch and yaw are 0, the new roll is just original roll + camera rotation angle. Likewise, if yaw and roll are 0, the new pitch is just the original pitch + camera rotation angle. For an arbitrary angle, I know both pitch and roll are affected. I imagine there is a simple equation for this, but I'm having trouble finding a reference.

Thanks for the suggestion to look up Euler angle rotations.  https://www.learnopencv.com/rotation-matrix-to-euler-angles/ is a good reference to answer this question.


Strange artifacts while rendering a dielectric material I am trying to render a dielectric sphere inside the Cornell Box. Some circles are generated on the sphere. Here is the photo generated by path tracer:

And below is my main function:

function trace(surfaces::Array{Surfaces, 1}, ray::Ray, depth::Int64, maxDepth::Int64)
if depth >= maxDepth
    return Vec3(0,0,0)
end
t, material = nearest(surfaces, ray, Inf)
if typeof(material) == Empty
    return Vec3(0,0,0)
end
if material.isLight == true
    return material.emittance
end
ρ = material.reflectance
n = material.normal
if material.isDielectric == true
    cosθ = -dot(n, ray.d)
    ref = material.ni/material.nt
    if cosθ < 0
        cosθ = -cosθ
        n = -1.0*n
        ref = 1/ref
    end
    disc = 1 - ref*ref*(1-cosθ*cosθ)
    if disc < 0
        R = ray.d - 2*dot(ray.d, n)*n
        In = trace(surfaces, Ray(ray.s + t*ray.d, R), depth+1, maxDepth)
        return In
    elseif disc > 0
        R0 = (ref-1)/(ref+1)
        R0 = R0*R0
        Rθ = R0 + (1-R0)*(1-cosθ)*(1-cosθ)*(1-cosθ)*(1-cosθ)*(1-cosθ)
        u = rand()
        if u < Rθ
            R = ray.d - 2*dot(ray.d, n)*n
            In = trace(surfaces, Ray(ray.s + t*ray.d, R), depth+1, maxDepth)
            return In
        else
            R = ref*ray.d + (ref*cosθ - sqrt(disc))*n
            In = trace(surfaces, Ray(ray.s + t*ray.d, R), depth+1, maxDepth)
            return In
        end
    end
end
BRDF = ρ/Π
R = hemiRand(n)
In = trace(surfaces, Ray(ray.s + t*ray.d, R), depth+1, maxDepth)
return Π*BRDF*In
end

And I think the fresenl effect of the area light on the sphere looks strange and wrong! Am I correct?

EDIT: I searched a bit and I guess it's the moire effect. But I don't know how to get rid of it.

I solved my problem by adding a small value to t value when the ray is getting refracted.
But I still doubt the result is true. Below is the new generated photo:

Why is the area light reflected twice on top of the sphere? Is there still a bug in my path tracer? The glass sphere is closer to the camera than the area light.


Transformation between aspect ratio I am trying to measure width of an object in image. The object size is 2D object 15cm X 15cm.

When I capture a frame with resolution 1280 X 720 (16:9) I am getting that the object width is 50px.

But when I capture a frame resolution 640 X 480 (4:3) I am getting that the object width is 30px.

The distance between the object and the camera did not changed during the test.
I didn't find a transformation that can convert the width in pixels from 1280 X 720 to 640 X 480. I think that the aspect ratio is the issue here but I didn't succeed to solve it out.

How can I know that the width in pixels is both resolutions is actually the same object size?

It looks to me like the vertical dimension is the one to use for scaling:

$$\frac{50\text{ px}}{720\text{ px}} \approx \frac{30\text{ px}}{480\text{ px}}$$

There's a slight error there, but it's only a few pixels. I'm not sure if the object sizes you mentioned were intended to be precise to the nearest pixel (30 and 50 are roundish numbers).

Typically when switching between aspect ratios, the vertical FOV would remain constant while the horizontal FOV widens or narrows to match the aspect ratio.


Optimal 2D graphics resource I am a some-what experienced programer, and I have been working on creating a 3D graphics engine utilizing the OpenGL system, primarily focused on use in games. However, I also want to do some work with functional programming, and wanted to create a clean UI for more functional programs, like a task manager. An example of the kind of effect that I want to create could be somewhat like the Visual Studio Code UI. My question is: is there a specific library that I should use for creating desktop programs, or should I continue to use OpenGL's 3D aspects, but keep everything in a single Z axis plane, effectively creating a 2D program? Or is there a better way to create normal program visuals, that wont be cheating out of a 3D graphical space?

Thanks for your help.

EDIT:
I am currently using C++/C, but I am also happy to expand my knowledge base to include new languages.

As you said you want to create something like this
As you know these are not graphic intense rendering.
Thus, (i guess but not 100%) Windows do these rendering on CPU not on GPU.
GPU is for intense rendering with lighting, textures, 3D models, etc. Where you can take the benefit of algorithm built in GPU circuit, Which is what OpenGL does.
So i would recommend to use built-in libraries that don't use GPU and use CPU for rendering, Although you can use OpenGL or Direct3D as well.
Some libraries that are built in different programming languages are :-

Java - Jframe, Graphics, Buffered image,  Canvas. (Note : all these libraries have different functionality and you have to use them together for software rendering See this for more detail)
C# - System.drawing,  (i don't have much experience of C#, So i can't name any more libraries), also C# form application.
JavaScript - Canvas check this out
C++ - you need third party libraries like QT, GTX etc. Qt learning resource, related resource
Python - Pygame (No personnal eperience) See this for detail


Rotate a model or a view with mouse input to produce intuitive results I am writing 3D graphics software in Java using the LibGDX framework (which is a front-end for OpenGL), and am baffled by a transformation problem.  How do I rotate the orientation of the camera or models being displayed by the dx and dy inputs from a mouse?

Click here to view the full source code of a concise demo in github.com

But, here is the method where I'm certain the problem and solution will be found:

@Override public boolean touchDragged( int screenX, int screenY, int pointer )
{
    lastX = screenX-lastX;
    lastY = screenY-lastY;

    // distance of mouse movement
    screenAng = (float) Math.sqrt( lastX*lastX + lastY*lastY );
    // direction vector of the AOR
    screenAOR.set( lastY/screenAng, lastX/screenAng, 0f );
    if ( touchedButton == 0 )
    {
        cubes[ selectedCube ].modelInstance.transform.rotate( screenAOR, screenAng );
    }
    else
    {
        camera.rotateAround( Vector3.Zero, screenAOR, -screenAng/5.5f );
        camera.update();
    }

    lastX = screenX;
    lastY = screenY;
    Gdx.graphics.requestRendering();
    return true;
}

EDIT: Thanks to Andrew Wilson for his comment below.  I edited the camera section and it works properly now.  Still working on the model section...

    {

        // transform the screen AOR to a world AOR

        // get the camera transformation matrix
        tempMat = new Matrix4( camera.view );
        tempMat.translate( camera.position );
        tempMat.inv();

        // transform the screen AOR to a world AOR
        worldAOR = transformThings( tempMat, screenAOR, new Vector3() ).nor();

        // apply the rotation of the angle about the world AOR to the camera
        camera.rotateAround( Vector3.Zero, worldAOR, screenAng/5.5f );
        camera.update();
    }

EDIT: I edited the model section and it still is giving me some grief - it looks like it's orbiting the world origin instead of rotating around the model origin.  Here's what it looks like so far:

    {

        // transform the screen AOR to a model AOR

        // get inverse of the camera transformation matrix
        tempMat.set( camera.view );
        tempMat.translate( camera.position );
        tempMat.inv();

        // get inverse of the model transformation matrix
        Vector3 tempPos = cubes[ selectedCube ].position;
        cubes[ selectedCube ].modelInstance.transform.translate( tempPos.cpy().scl( -1 ) );
        tempMat2.set( cubes[ selectedCube ].modelInstance.transform );
        tempMat2.inv();

        // multiply the camera and model matrices together
        tempMat.mul( tempMat2 );

        // transform the screen AOR to a model AOR
        modelAOR = transformThings( tempMat, screenAOR, modelAOR ).nor();

        // apply the rotation of the angle about the model AOR to the model
        cubes[ selectedCube ].modelInstance.transform.rotate( modelAOR, -screenAng );
        cubes[ selectedCube ].modelInstance.transform.translate( tempPos );
    }

Here is an image that the demonstrates the problem:

I want the object to rotate in the direction it is dragged by the mouse, no matter which direction it happens to be orientated at the time. As it is now, when I first drag the mouse to the right, the object rotates to the right about the screen Y axis as expected; but then when I drag the mouse upward I want the object to rotate upward about the screen X axis, but instead it spins to the left about the screen Z axis. Think of it like a floating ball in a bowl of water - whichever way you swipe at it, it ought to rotate in that direction.

It seems to me that the mouse movement is transforming the objects directly in their local coordinate system; but instead I think I need to transform the axis of rotation itself from the Screen Coordinate System into the Object Coordinate System before applying it to the object. I just can't figure out how to do it.

I would really appreciate any insight or help to resolve this; I'm running out of hair to pull out... Thanks in advance.

Here's a more complete answer from my comment. The end product we want is:

Projection * ViewTranslate * Model * ViewRotate

All are 4x4 matrices. Where model is your models transformations. To implement you can just separate out the view's translate and rotate during your draw.

If you're using an abstraction that requires you to use Proj * View * Model. Or you just really don't want to ruin your own abstraction that follows that pattern. You'll need to modify model to cancel out views rotational part and put it on the right side:

Model = ModelTranslate * View^-1 * ViewTranslate * ModelRotate * ViewRotate
      = ModelTranslate * ViewRotate^-1 * ModelRotate * ViewRotate

Here's some of my code that does the above. It does not matter if model or view is translated, objects will still rotate around their local origins:

glm::vec2 diff = (prevMousePos - mousePos) * rotSpeed;
glm::mat4 rotX = MathHelp::matrixRotateY(diff.x);
glm::mat4 rotY = MathHelp::matrixRotateX(diff.y);

// Decompose the translate and rotate part of view
glm::mat4 viewT, viewR;
Decompose(cam->view, viewT, viewR);

// Decompose the translate and rotate part of model
glm::mat4 modelT, modelR;
Decompose(exampleObj->model, modelT, modelR);

exampleObj->model = modelT * glm::inverse(viewR) * rotX * rotY * viewR * modelR;

This preserves the existing rotation and translation in exampleObj->model as the function gets called continually applying small rotations. And assumes you'll be using Proj * View * Model.

Lastly decompose just separates out the 3x3 rotation/scaling from the translational part:

void Decompose(glm::mat4 m, glm::mat4& t, glm::mat4& r)
{
    glm::mat4 I4 = glm::mat4(1.0f); // Identity
    for (unsigned int i = 0; i < 4; i++)
    {
        for (unsigned int j = 0; j < 4; j++)
        {
            if (i == 3 || j == 3)
            {
                r[i][j] = I4[i][j];
                t[i][j] = m[i][j];
            }
            else
            {
                r[i][j] = m[i][j];
                t[i][j] = I4[i][j];
            }
        }
    }
}


Barycentric Coordinates of a Point on a Line Given a point and a line segment, how do you calculate the barycentric coordinates of the point to the line?

bool BarycentricCoordinates(const Vector3& point, const Vector3& a, const Vector3& b, float& u, float& v, float epsilon)

If you have a point $C$ on a line defined by a point $A$ and a point $B$, here is how you calculate the barycentric coordinates $(t_0,t_1)$:

$t_0 = \frac{sqrt((C_x-A_x)^2+(C_y-A_y)^2)}{sqrt((B_x-A_x)^2+(B_y-A_y)^2)}$

$t_1 = 1.0 - t_0$

In other words, you calculate the length of the vector from $A$ to $C$, and divide that by the length of the vector from $A$ to $B$.  This gives you $t_0$, and you can subtract that from 1.0 to get $t_1$ since barycentric coordinates always add up to 1.0.

A more mathish style version of the above is this:

$t_0 = \frac{\left\| C-A\right\|}{\left\| B-A\right\|}$

$t_1 = 1.0 - t_0$

But it's important to note that $C$ must be on the line defined by $A$ and $B$.  If that isn't the case, let us know and let us know what it is you are trying to do specifically.  For instance, maybe you want to project some other point $D$ onto the line $\overline{AB}$ and then find the barycentric coordinates of that projected point?


How do graphic engine developers debug their 3D graphics code? I'm currently trying to create an amateur game engine, and I often run into very basic bugs (eg I change the node hierarchy code and suddenly nothing is displayed). Some of the methods I will use include to track down these bugs include:

Use my IDE's debugger and lots of printf calls to look at what's in memory.
Make lots of tests with different parameters, placing the camera in different positions, to get a better feel of what the code is doing.
Just brute force my way through by re-reading each line of code until I find the one that's misbehaving.

These methods all seem sluggish and inefficient to me. The problem is, that 3D rendering involves a lot of unintuitive operations with not much to visualize except the final "product"; when I'm running into a bug (especially "nothing is displayed" types of bugs), I often feel like I'm fumbling around blindly.

So.

Let's say you're a developer working on Unity, Godot, or some other 3D engine, and you just changed some steps in your rendering pipeline, and now, for whatever reason, the wireframe renderer displays incoherent results when entity batching is activated. The "wires" are crisscrossing all over the screen instead of representing the edges of the batched entity's triangles.

How do you track down and fix that bug? Do you have a specific workflow? Do you use some kind of GPU analysis tool? Or do you just fumble around blindly?

In other words, what methods, tools and best practices can engine developers (not engine users) use to analyze the data and behavior of their 3D graphics code (not just its performance) to efficiently address defects and bugs?

I'm programming mostly in OpenGL. In this case, the very first thing I do is setting up KHR_debug extension. It provides a way to subscribe to messages generated by GL drivers - from errors to warnings to detailed information about what the driver is doing. Especially open source Linux drivers (Mesa) generate lots of useful messages.

When it comes to tools, my first choice is usually ApiTrace. It allows you record a trace of all rendering API calls and then replay it and inspect state at any point in time (this includes contents of buffers, render targets, etc.). It also allows you to share a trace with other people, which is useful when working in a team. Then, there are other tools such as RenderDoc, NVIDIA Tegra Graphics Debugger, Xcode's graphics debugger and some others. Whichever you chose, it will improve your work a lot.

Another thing is testing. One thing is building a collection of simple scenes excercising single feature of the engine (such as different material settings, different animation techniques, etc.).
When a bug is introduced into the codebase, this helps debugging it (since you have simple to analyse test cases readily available) and verifying that a fix doesn't accidentally break something else. I guess it's quite common thing to do; Unity3D developers do it, for example.

You could extend this and automate it by creating test harness, which would build and run all test on regular basis, and verify that no regressions happened. If you have lots of test scenarios, this might be feasible, although obviously requires additional investments into infrastructure.

Then there are more "standard" unit tests. One area I've applied them is node-based material system, which can target different rendering techniques (forward, deferred, whatever) and APIs. I approached it by writing unit tests which generate shader code from material and other settings, compiler it and issues a simple draw calls with given inputs (pixel shader inputs). Outputs are written to a render target and fetched for comparison with expected values. IMHO this is a nice  thing to do when you have complex material system, esp. with shader code generation.

Finally, more "traditional" approaches can be applied to areas which can be separated from rendering. Things like spatial data structures, asset loading and similar can be unit/integration tested, as applicable in specific case.

Edit: Here's a post from Unity3D developer describing how they test graphics code.


GLSL shapes signed distance field implementation explanation? I'm trying to understand all primitives In this article.After hard work I just learn Sphere and Cube , I wrote description here.

I could draw 2D shapes but I have trouble understanding 3D below codes:

Torus - signed - exact

float sdTorus( vec3 p, vec2 t )
{
  vec2 q = vec2(length(p.xz)-t.x,p.y);
  return length(q)-t.y;
}

I can't imagine how this shape drew and I don't know why two length used here.

Cone - signed - exact

float sdCone( vec3 p, vec2 c )
{
    // c must be normalized
    float q = length(p.xy);
    return dot(c,vec2(q,p.z));
}

I don't know why dot product used here?

I'm confused please help me to understand this shapes.I hope someone explain each shape.

Torus

A torus is defined by two parameters: the major radius, and the minor radius. The major radius (t.x) is the radius of the big ring (in red in the diagram), and the minor radius (t.y) is the radius of the circular cross-section. The x and y here are just indices into the vector: they're unrelated to the x and y axes.

The left side of the diagram shows the xz plane, and the right side shows a cross-section with y going horizontally.

p.xz is the projection of the point onto the plane of the torus, so length(p.xz)-t.x gives the distance from that projected point to the inside of the torus. (I've labelled it $r$ on the diagram.) p.y is the distance from the point to the plane of the torus, so the length of the vector ($r$, p.y) is the distance from the point to the inside of the torus. (Computing the length of a two-element vector is a handy way to use Pythagoras' theorem to find the hypotenuse of the triangle.)

Cone

Like the torus, the cone code uses the fact that the shape is symmetric about an axis: the z axis in this case. The cone is defined by a unit vector c that's normal to the surface of the cone, and the apex of the cone is at the origin. Similar to the above, length(p.xy) gets the distance from the point to the z axis, so that we can work in the plane where the cone is a triangle. vec2(q,p.z) is the projection of p into this plane. The dot product projects this vector onto the direction c: how far the point is along this direction tells us how far it is from the surface of the triangle.


Color space of MERL BRDF database? I know that by now it might be considered as a kind of sad recurrent joke, but by chance, does anyone here has the least information about the colorspace of RGB values in the MERL BRDF measurement database ?

Was interested in this myself and did some digging, I was wondering if they where spectral, but it appears not (although some interesting work on spectral BRDF is here http://web.cs.wpi.edu/~emmanuel/courses/cs563/write_ups/cliffl/cliffl_spectral_brdf.html )

This chat seems to indicate they are camera native linear, probably with sRGB primaries.

https://groups.google.com/forum/#!msg/pbrt/p9zMvZlaodA/4RD7-We8h4gJ


What is performance better - add four materials to mesh or duplicate every vertex and use packed texture? I am making ski jumping game (in Unity) with procedurally generated hill mesh. I want to add lines on hill like this (mesh is created by script, but I've manually assigned materials in Blender)

Because this is procedurally generated (on runtime) (hills can have different sizes) I cannot bake texture in Blender. Mesh for large hill has around 800 vertices. There is my question - is it better to duplicate vertices and then use one material or use four different materials without duplicating vertices?

Hello and welcome Jonek2208!

Unity specific questions probably are better suited for the Unity Graphics Forum. However, I'll try to give you an answer as I do have some experience with the engine.

In my opinion, the only way of knowing for sure which method is the best would be to implement both and profile them on your target device. Do not fool yourself in thinking that some basic algorithmic analysis is all you need to do. Always profile and compare your solutions, no matter the circumstances.

One thing I can do is to tell you the trade-offs between the two methods, and also provide a hint at another possible solution. Let's go!

When using the built-in rendering pipeline, Unity batches objects by materials. A batch is considered to be a change in GPU state followed by one or more draw calls. In other words, if you have four materials, you increase the number of batches when compared to having just one.

If you are targeting mobile devices, reducing the batch count is one of the keys for achieving good performance, and most mobile developers would choose to use a single material. On the other hand, duplicating vertices incurs a cost in memory, and hence you get better batching at the expense of using more memory. If you are targeting consoles or PCs, you can have a much higher batch count and memory budget, so the solution that you think is simpler to implement, expand, and maintain should be your best bet.

Another solution would be to use textures instead of vertex colors. For this, you would need to compute the uv-map procedurally when generating the vertices. You could place all textures into an atlas and try to use a single material. This is the usual solution for mobile. If you are targeting PCs/consoles, you could even afford procedural textures in your fragment shader.

If you are using the latest version of Unity, you can use the SRP Batcher, which batches per shader instead of per material. Hence, you can use 4 materials, but have a single batch if all of them use the same shader. Note some devices may not have this capability, especially on mobile.

All in all, there are a lot of things to consider in any kind of application. That is why I started talking about profiling on your target device. Honestly, it is the only way of knowing for sure which method is the best. You should build a test scene with thousands of generated meshes and then profile each method. Select the one that fits best in your budget (FPS/memory/etc) and the long-term goals of your project.

Good luck! :)


How can I control memory for command allocators and command queues in D3D12? The D3D12 programming guide suggests specific allocation strategies for command queues and allocators, but the creation methods for those objects operate directly on queue/allocator pointers without taking a heap and an offset as an argument. How do I explicitly mark the memory I want to use for command processing?

You don't.

These objects typically work with CPU memory, not device memory. And to the extent that they involve device memory, such allocations tend to be rare and/or fixed in size (a queue may have a small spot of device memory that the hardware queue reads commands from or something, but even that is implementation-dependent).

Vulkan allows you to give objects which may allocate non-device memory allocation parameters to allow you to allocate said memory for them. But even then, there is no guarantee that they will always allocate memory through your functions. Such systems can bypass your allocation routines, and only tell you allocator that it is allocating X bytes of memory (this is typically done for memory pages that have to be allocated in a special way). It doesn't let you actually provide it.

That documentation does not seem to be talking about controlling how those objects allocate memory, but when and how you create those objects. Command allocators and queues are critical, so you create X allocators and Y queues at the start and keep them around forever.


Fast approximation of spherical texture coordinates I'm working on a Software Rasterizer on a sub-30 MHz RISC CPU My current focus is on zooming-in on a slowly rotating planet. For maximum quality, no 3D polygons are involved with the planet - the planet texturing is purely image-based. This is what I have:

the planet rotates only around Y axis and is always viewed straight
on with camera at the equator (so you can't see poles, hence the
default pole distortion artifact is not an issue here at all)

midpoint Bresenham circle algorithm to compute 8 endpoints at each iteration, effectively giving me 4 horizontal scanlines
simple linear remapping between texture row and planet's scanline - e.g. if the planet is 64 pixels tall on screen and texture is 128
pixels tall, I simply skip every other row
same linear remapping within each scanline - e.g. if the current scanline is 64 pixels wide, and texture is 128 pixels wide, I simply
skip every other texel in current texture row

The above works and is great as a first working and reasonably fast  prototype. However, when the planet rotates, despite the low screen resolution, due to the simple linear remapping (within each scanline), it does not feel and look like planet rotation, as the perspective/distortion/stretching at the sphere edges (left and right edge of the circle) is obviously missing.

What is the equation that can give me the exact texture coordinate, for a current texture row, of each on-screen pixel of that circle's scanline ?

Please note that I can't afford to compute sin/cos (or use precomputed tables). Division is expensive, but I could perhaps hide its latency via interleaving with other instructions (while division is being computed). So, ideally, I could compute it just via add/sub/mul/bitshift.
Obviously, due to symmetry, that computed coordinate will be reused for 3 other points, so we're in reality computing that for only 25% points of the circle.

Incidentally I implemented something similar for similar hardware some 20 years ago (:

IIRC, I calculated 1D look-up-table using acos for parallel projection of a cylinder and just scaled it for each scanline to reduce computation cost. The LUT just gives you x-coordinate offset to the texture and you add constant offset for each scanline for the rotation.

For y-axis the input texture can be preprocessed to have proper spherical distortion applied and you can just linearly scale the texture vertically. Though it's not that much extra computation to do it run-time either.

This doesn't give you proper perspective projection but might be enough for your purposes.

Edit: Didn't notice you said you can't afford LUT. You could do the table look-up every N pixels and linearly interpolate inbetween. This was a common method in the past to do affordable perspective correct texture mapping for software rasterizers to avoid div/pixel


How could I check the correctness of my result of fluid simulation? I wrote a particle based fluid simulating program. It's hard to tell if I get the right result. The visualized result seems reasonable, but some part of it looks weird. I don't know wether it's a feature of fluid. Is there some accurate method to verify if my program is right?

Amending some details：

My program is a 2D particle-based simulating program. The fluid is compressible. The implementation is almost based on a classical paper:

  Müller, Matthias, David Charypar, and Markus Gross. "Particle-based fluid simulation for interactive applications." Proceedings of the 2003 ACM SIGGRAPH

I solved the Navier-Stokes equation with iteration method. It only considered pressure, gravity, viscosity and surface tension.

Compare it with someone else's software. Run some standardized test and find out if you get roughly the same answer as others. If you get the same answer, than the probability of having your code right is quite high.

Some tests:

Flow past cylinder. In 2d take rectangular domain, cylinder in the middle, inflow on the left, outflow on the fight and calculate the force on the cylinder. Here is benchmark comparing handful of codes.
Buoyancy flow. Closed box, hot plate on bottom, cold plate on top, hot fluid starts to rise due to the buoyancy force. Here is benchmark.
Rising bubble, benchmark.

But unfortunately it might be quite difficult to compare your code to scientific codes in those benchmarks. I guess you implemented something as SPH or stable fluids which are not made for accuracy but for stability.

Take for example the flow past a cylinder. I would start the test with very small Reynolds number and then measure the force on cylinder as you increase accuracy of your simulation(lower time step, increase subdivision or increase number of particles). Does the force converges to some number? If no, than you have a problem, if yes, than have a look at the benchmark paper and compare your result with others.

This method is very similar technique to one I use for testing my raytracer. I just render test scene with someone else's renderer and compare it with my result. Do they converge to the same result? If yes than I have it right, if no, than I have it wrong.


Does atomic functions on same memory location cause an performance issue? is it worse to do for each invocation an atomicAdd on the same memory location than doing for each invocation an atomicAdd to a unique memory location? (I know, atomicAdd does not make sense when each invocation writes to another location... but this is theoretically)
I am asking, because I wrote the marching cubes algorithm and while creating the mesh I additionally want to calculate the center of mass as well.
I am doing this by having four integers in the SSBO (centerX, centerY, centerZ and countTriangles)
each of these 4 integers will be increased by the atomicAdd function within the GeometryShader.
When everything is done, I only need to divide the centerX/Y/Z variable by countTriangles and I am done.
if these atomics on one memory location would slow down my shader, I would outsource it to another maybe CPU function.
Thanks in advance

Will having all invocations incrementing the same atomic counter be slower than not doing that? Sure. How much slower? That depends on many factors, but it's hardly going to be a crippling performance issue.
The whole point of atomic memory operations is that there is some hardware mechanism that allows them to work reasonably efficiently. Not as efficiently as bumping a local variable of course, but they shouldn't be assumed to cause a massive performance penalty.
On some hardware, atomic counter operations aren't even memory operations. AMD implements them as specialized hardware within their shader system, such that atomic counter operations don't touch actual memory storage until after the rendering process that uses them is complete.
Obviously, the best thing you can do is profile on your hardware-of-choice, but failing that, I would see no reason not to use this mechanism as you see fit.


How to organize Lights and Models? I'm interested in how to organize in code (I prefer C++) lighting.
Problem is: I have Scene which contains Model where model contains Shader, someTextures and geometry Mesh. Now I want to add Light with the following fields:

struct Light
{
    Vector3f Color;
    float AmbientIntensity;
    Vector3f Direction;
    float DiffuseIntensity;
    float SpecularIntensity;
    Vector3f Position;
};

To apply light I need to know the values of this field and the position of the fragment pixel of geometry. This means I need to make a function like Light::lightModel(Model mdl) to get data from Model which I want to light. But also I want to make Light debuggable - show it in the scene like a Model (cube). How can I light all Models and do not light the light itself?  Make two Scenes? Scene<Model> models and Scene<Light> lights and for loop for each like this?

for (const auto& light : lights->data()) {
    for (const auto& model : models->data())
    {
        light->lightModel(model);
        model->draw(camera);
    }
}

If I understand you correctly then you are probably looking for a way in which lights in the scene interact with models in the scene. Also, you want each light to be visually debuggable in scene, maybe by drawing a cube. Here's one of the simplest way to do it. There could be more better & optimized ways of doing this.

You can have a class for Each light type eg> PointLight, DirectionalLight etc each with their own properties. You need to have a class which then keeps track of all the lights in the scene. This class will have list of light types in the scene. This class will also have functions for adding active lights from the scene to it's list of lights. Usual get-set functions for accessing individual lights will be required too.

class LightsManager
{
// maybe create singleton object here..
public:
   void GatherPointLights(PointLight* light);
   PointLight GetPointLight(int id);

private:
   int numPointLights;
   std::vector<PointLight*> vecPointLights;
}

void LightsManager::GatherPointLights(PointLight* light)
{
    vecPointLights.push_back(light);
    numPointLights = vecPointLights.size();
}

In your class where you create objects in the scene, you can then add light objects from the scene to this LightManager using the function you created earlier.

Assuming scene class where all the objects including Lights are created...

void Scene::Init()
{
    // load mesh
    ...
    ...
    // add lights
    PointLight light1;
    light1.color = Color::Red;
    light1.position = vec3(0,0,0);
    LightManager::getInstance()->GatherPointLights(&light1);
}

If LightManager creates singleton object then you can then access all the lights in the scene by their type during actual rendering to extract their properties in a class where you render actual scene geometry by passing light color, position etc to a shader which draws scene geometry.

Now drawing the actual debug geometry for any light object, say cube for point light, I would make sure my PointLight type object has a way of rendering cube at a Point light position, using it's color property to draw cube of that same color with a flat color shader.

As per making sure lights are not lighting actual light debug geometry, if your cube drawing shader is a flat color shader representing color of a light, then it will not be affected by any light in the scene.

Hope this is what you asked for.


Multiple VAO share a EBO (OpenGL 3.3) As mentioned here, multiple VAOs can share a single VBO. However, can multiple VAOs share a single EBO?

For example, I want to create 2 rectangles. So, I do the following:

const GLuint indices[6] = {
    0, 1, 2,
    0, 3, 2
};

// First
glGenVertexArrays (1, &VAO1);
glBindVertexArray (VAO1);
glGenBuffers (1, &EBO1);
glBindBuffer (GL_ELEMENT_ARRAY_BUFFER, EBO1);
glBufferData (GL_ELEMENT_ARRAY_BUFFER, sizeof(indices), indices, GL_STATIC_DRAW);
// creating VBO, setting attributes etc
// Second
glGenVertexArrays (1, &VAO2);
glBindVertexArray (VAO2);
glGenBuffers (1, &EBO2);
glBindBuffer (GL_ELEMENT_ARRAY_BUFFER, EBO2);
glBufferData (GL_ELEMENT_ARRAY_BUFFER, sizeof(indices), indices, GL_STATIC_DRAW);

glBufferData(GL_ELEMENT_ARRAY_BUFFER) is called twice for different EBOs with same data. As I understand, this leads to double GPU memory consumption and takes twice loading time.

So, would the following be correct?

const GLuint indices[6] = {
    0, 1, 2,
    0, 3, 2
};
glGenBuffers (1, &EBO);
glBindBuffer (GL_ELEMENT_ARRAY_BUFFER, EBO);
glBufferData (GL_ELEMENT_ARRAY_BUFFER, sizeof(indices), indices, GL_STATIC_DRAW);

// First
glGenVertexArrays (1, &VAO1);
glBindVertexArray (VAO1);
glBindBuffer (GL_ELEMENT_ARRAY_BUFFER, EBO);
// creating VBO, setting attributes etc
// Second
glGenVertexArrays (1, &VAO2);
glBindVertexArray (VAO2);
glBindBuffer (GL_ELEMENT_ARRAY_BUFFER, EBO);

Your code is non-functional, but only because (in the core profile) you cannot bind a buffer to GL_ELEMENT_ARRAY_BUFFER unless a VAO is also bound. So the part that creates the buffer would need to use a different binding point, which is totally fine:

const GLuint indices[6] = {
    0, 1, 2,
    0, 3, 2
};
glGenBuffers (1, &EBO);
glBindBuffer (GL_ARRAY_BUFFER, EBO);
glBufferData (GL_ARRAY_BUFFER, sizeof(indices), indices, GL_STATIC_DRAW);

Everything else is perfectly valid.


Appropriate filter steps for resolving a droplet pixel I have a series of images obtained from a microscope as follows

I need to create a series of filters that would isolate the centre section    on the yellow droplet.(This is highlighted by the bright yellow circle area), and remove all the noise and other pictures.

I tried using just RGB  channel threshold but its not working. Also ambient light may also change so I am not getting a consistent output. Can anyone provide a MWE for the workflow of the filter required in MATLAB to achieve this?

First simple thresholding is usually not a good idea. However, if you want to do that, I suggest to also consider other color spaces. In your case, HSV might be a better choice; You can convert your image into HSV-image then you can define a range on Hue channel to select yellow.

Also, take a look at this article; it is about implementing Adaptive Thresholding.


How to convert from XYZ or sRGB to ACEScg/AP1 I'm confused about how to properly convert from sRGB or XYZ color space to ACEScg. Specifically whether I have to account for differences in white point explicitly.

The aces github page gives the following matrix to convert from XYZ to AP1:

 1.6410233797, -0.3248032942, -0.2364246952,
-0.6636628587,  1.6153315917,  0.0167563477,
 0.0117218943, -0.0082844420,  0.9883948585

FGED2 gives the following matrix to convert from XYZ to sRGB with white point D65:

 3.240970, -1.537383, -0.498611,
-0.969244,  1.875968,  0.041555,
 0.055630, -0.203977,  1.056972

And this shadertoy gives the following matrix to convert from sRGB to AP1:

0.613097, 0.339523, 0.047379,
0.070194, 0.916354, 0.013452,
0.020616, 0.109570, 0.869815

So I figured if I combine the XYZ->sRGB and sRGB->AP1 matrices it should be the same result as XYZ->AP1. However there is a small, but noticeable difference:

 1.66058,    -0.315295,  -0.24151,
-0.659926,    1.60839,    0.017298,
 0.00900358, -0.00356713, 0.913644

I assume the difference comes from the way the white point is handled? My question is, which is "correct". I have use cases for both XYZ->AP1 and sRGB->AP1, but now I'm irritated that my matrices don't allow me to round-trip.

The discrepancy comes from chromatic adaptation, which is included in the matrix from shadertoy. If you combine XYZ->AP1 with D65->D60 from the aces page you end up with the same matrix that I got by chaining XYZ->sRGB and sRGB->AP1 (with chromatic adaptation) together.


Nice way to make noise seamless (without tampering with the noise generator) Here's the problem :
I have multiple noise generators in my program I Can't change.
Perlin Noise, Diamond Square, Voronoi, etc.
Unfortunately for me, I want the noise to be tileable, no matter the generator :

The idea is to create a torus-like strategy game with the user choosing a noise and then having it tiled to play on it.

I'm very well aware it's a rather poor design to make it tileable afterwards, but it's a time saver for the coding part.

Thank you in advance for your help and ideas.

If you can't change the textures, I see 2 possibilities:

Use mirrored repeat on the noise textures so they tile seamlessly.  (taken from here)
Do like Gimp does (or many other alternatives) and create tilable versions of the noise textures. The idea is to blend the left and right, and bottom and top edges together in pairs. Here's how to do it in an image editor: translate the image by 1/2 in both axis so the corners are in the middle. Blur, smudge, etc. to remove the seams. Revert the translation.

Keep in mind that the randomness and quality of the noise will of course be affected by either technique.


Determine paint color from a photo I ran into this problem while trying to find the paint color of the doors in my house: is it in principe possible to automatically determine the paint color from a photo?

Different types of colour standards (such as RAL, https://en.wikipedia.org/wiki/List_of_RAL_colors) mention RGB values, but of course lighting conditions, material, camera and probably other factors highly affect the RGB value that is measured.

Does a solution to this problem already exist?

An idea I had is to hold an official color reference chart in the photo next to the target object (e.g. my door), so the transformation due to lighting conditions and camera between the known RGB values of the colors on the reference chart and those in the photo could somehow be established. Then the inverse transformation could be applied to the color of the target object and the closest matching paint color could be found.

Could this idea work and if so, what kind of transformation between reference and measurement could be expected?

In theory yes, but practically speaking, I doubt it would be very easy. First, the color may not be representable by the color space that your camera of choice produces images in. We've recently started seeing wide color gamut (P3) cameras and displays. There are colors these cameras record and displays produce that fall outside of the sRGB color space that most cameras and displays support. Likewise, there are visible colors that are outside of the P3 color space. So the first hurdle is to have a device which records and stores the colors you're dealing with.

Next you have to do the lighting correction once you get the image into the computer. This is actually not as hard as it sounds. There are products that allow you to do this using a system similar to what you describe. Here's a video for one that works with video rather than stills.

Once you've done all that, you need some way of translating that color from some RGB space into something that your local paint store can understand. I'm not sure how to do that. Most paint store employees are not trained in color science and have to use a limited toolset provided by the paint companies.

Is there any way you could get a small, inconspicuous piece of the actual object? Most paint stores have a machine that will take a physical object and create a paint that's a fairly close match to that object. Whether their matches are close enough for your needs is something you'd have to determine. I've used them before and found that while they get very close, I'm still able to see the difference between the original and the new paint, unless I can paint up to corners.


Rendering equation in terms of paths rather than directions Path tracing seems to perform Monte Carlo integration where the space being integrated over is the space of paths (e.g. chains of rays) but the rendering equation is stated in terms of outgoing rays. Naively you would have to perform the integration at each bounce to model the rendering equation but this would take exponential time. It seems that there should be a proof that either a) integration over the space of paths is equivalent to recursive integration over the space of rays or b) sampling paths is a valid numeric technique for evaluating the rendering equation. While somewhat intuitive neither of those facts is obvious to me. Is there a proof of some such claim?

This is a very good question. There is a common misconception that Monte Carlo, or integration is applied "recursively" on the rendering equation. That is not what's happening. Numerical integration methods are tailored to problems of the form:

$$I = \int_{\Omega}{f(x)d\mu(x)} \approx \sum_{k=0}^{N-1}w(x_k)f(x_k)$$

Note that this is not the case for the rendering equation. Namely the fact that it is an equation, and the unknown radiance function $L$ shows up on both sides. In order to apply numerical integration techniques (such as Monte Carlo), one has to bring the rendering equation to the above form. Fortunately this is a well studied problem. The rendering equation is a specific case of a Fredholm integral equation of the second kind. If one can show that the scattering operator in that equation is a contraction $||T||<1$ (which holds if brdfs are energy conserving, there are cases where this is not required but that goes outside the scope of my answer), then we may perform a Neumann expansion (https://en.wikipedia.org/wiki/Liouville%E2%80%93Neumann_series).

The commonly used form of the rendering equation is the solid angle formulation (here $\sigma$ is the solid angle measure, $\sigma(\omega) = \sin\theta\,d\theta d\phi$):

$$L(x,\omega_o) = L_e(x,\omega_o) + \int_{\Omega}{f(\omega_o,x,\omega_i)L_i(x,\omega_i)\cos\theta_i\,d\sigma(\omega_i)}$$

In order to perform the expansion, I will rewrite this in the area formulation (here $\mu(x)$ is the Lebesgue area measure):

$$L(x_1 \rightarrow x_0) = L_e(x1 \rightarrow x_0) + \int_{M}{f(x_2 \rightarrow x_1 \rightarrow x_0)L(x_2 \rightarrow x_1)\cos\theta_{x_1}\frac{\cos\theta_{x_2}}{||x_2-x_1||^2}V(x_2,x_1)\,d\mu(x_2)}$$

Where $M$ is the set of all surface points in the scene, $V(x,y)$ is the visibility function which is $1$ if there is nothing between $x$ and $y$ and $0$ otherwise. If the normals of the surfaces at point $x_1$ and $x_2$ are respectively $N_{x_1}$ and $N_{x_2}$, then $\cos\theta_{x_1} = N_{x_1} \cdot \frac{x_2-x1}{||x_2-x_1||}$ and $\cos\theta_{x_2} = N_{x_2} \cdot \frac{x_1-x2}{||x_2-x_1||}$. As for the radiance, $L(x_2 \rightarrow x_1)$ gives the radiance arriving at $x_1$ from the direction of $x_2$. And the brdf notation relationship is: $f(x_2 \rightarrow x_1 \rightarrow x_0) = f(x_1 \rightarrow x_0, x_1, x_1 \rightarrow x_2)$ (that means $\omega_o = x_1 \rightarrow x_0$ and $\omega_i = x_1 \rightarrow x_2$).

I will rewrite the above for simplicity and conciseness as:

$$L(x_1 \rightarrow x_0) = L_e(x1 \rightarrow x_0) + L_r(x_1 \rightarrow x_0)$$

Let us now split the incoming radiance into direct illumination arriving at $x_1$ and indirect (at least one bounce) illumination. The direct illumination is obviously due to direct light source rays arriving at $x_1$:

$$L(x_1 \rightarrow x_0) = L_e(x1 \rightarrow x_0) + \int_{M}{f(x_2 \rightarrow x_1 \rightarrow x_0)(L_e(x_2 \rightarrow x_1) + L_r(x_2 \rightarrow x_1))\cos\theta_{x_1}\frac{\cos\theta_{x_2}}{||x_2-x_1||^2}V(x_2,x_1)\,d\mu(x_2)} =
L_e(x1 \rightarrow x_0) + \int_{M}{f(x_2 \rightarrow x_1 \rightarrow x_0)L_e(x_2 \rightarrow x_1)\cos\theta_{x_1}\frac{\cos\theta_{x_2}}{||x_2-x_1||^2}V(x_2,x_1)\,d\mu(x_2)} + $$

$$\int_{M}f(x_2 \rightarrow x_1 \rightarrow x_0)\Bigg[\int_{M}f(x_3 \rightarrow x_2 \rightarrow x_1)L(x_3 \rightarrow x_2)\cos\theta_{x_2}$$

$$\frac{\cos\theta_{x_3}}{||x_3-x_2||^2}V(x_3,x_2)\,d\mu(x_3)\Bigg]\cos\theta_{x_1}\frac{\cos\theta_{x_2}}{||x_2-x_1||^2}V(x_2,x_1)\,d\mu(x_2)$$

All I did here was separate the sum $L_e(x_2 \rightarrow x_1) + L_r(x_2 \rightarrow x_1)$ into two integrals, and then I additionally expanded $L_r$ into the corresponding integral using the recursive definition. You can expand this until infinity. For conciseness, I will rewrite the integration as an operator: $L_r = TL$. Now we can see the relationship to the Neumann expansion. Since the rendering equation can be written as:

$$L = L_e + TL$$

The solution is formally given as:

$$(I-T)L = L_e$$
$$L = (I-T)^{-1}L_e$$

Applying the Neumann expansion yields:

$$L = \sum_{i=0}^{\infty}T^iL_e$$

Note that each term of the sum is an increasingly dimensional integral. The first term is obviously just $L_e(x_1 \rightarrow x_0)$, the second term $TL_e$ is the integral for $L_e(x_2 \rightarrow x_2)$ that I wrote when I split $L_e$ and $L_r$, and so on. What each term gives you is the energy coming from a light source after $i$ bounces. The first term gives you the radiance emitted from $x_1$ towards $x_0$ (paths of length 0), the second term gives you the radiance that is due to direct illumination of $x_1$ scattered towards $x_0$ (paths of length 1), then you have the radiance due to lights two bounces away (paths of length 2) etc. With each bounce each new integral in the sum gains two dimensions (if we are just integrating over incoming directions). I would also like to emphasize that the integration variable $x_2$ in $TL_e$ is not the same $x_2$ as the one in $T^iL_e$ (this means that in the example with integration, $x_2$ in the first integral is not the same as in the second one - they are just integration variables). Another important fact is that within each sum each path starts at the film $x_0$ and ends at a light source ($L_e(x_2\rightarrow x_1)$ is non-zero only for points $x_2$ lying on light sources).

Now that we have a sum of integrals, we can apply our numerical techniques to each integral to estimate the sum. An obvious optimization, is reusing the samples used to compute $T^iL_e$ to compute $T^{i+1}L_e = T(T^iL_e)$. This yields a formulation where we can formally write out the sum of integrals as one integral with integration over all paths that start at the camera film and end at a light source. This also illustrates why Monte Carlo is a preferred technique: since we want to estimate an infinite dimensional integral, and Monte Carlo's convergence doesn't depend on the dimensionality unlike standard quadrature rules (there's also the point that it doesn't care about the smoothness of the integrand, but that's a double edged sword).

As you may have noticed, with the path integration formulation an exponential path count growth is not required. Splitting a path into multiple new ones is a technique known as splitting. While it is beneficial for high energy paths, on average with each bounce the energy gets lower due to attenuation, so this is in most cases counterproductive. On the other hand early path termination (for low energy paths) is more often than not beneficial in terms of efficiency (that's simply cutting off the sum at some point and not computing the remaining infinitely many terms, because you deem they have a contribution that is too low), and that's where Russian roulette comes in.

All of this is actually contained in Kajiya's paper on the rendering equation, even though he doesn't go into the details I went. He just refers to Rubinstein's book: https://dl.acm.org/citation.cfm?id=539488 (note - the first edition, later editions do not have the part he is referring to). In this book how to solve Fredholm integrals of the second kind is described (what I explained above, but more formally).

I hope my explanation was useful, and will continue being useful for future readers.


How to calculate curvature of isoline of image I've got a discrete grayscale image (a.k.a. 2D rectangular array of floats). Its continuous representation (whether received via reconstruction with sinc, with a cubic kernel, with a triangle kernel...) is a continous 2D function, so it has isolines. I'm interested in calculating the magnitude of the derivative of the isoline at a given point.

Pic:

Red line is part of an isoline, and green arrows are example derivatives of that isoline at 2 points.

Here's my approximate solution:

float getCurvature(vec2 p) {
    vec2 grad = getBilinear(gradients, p);
    if(grad = vec2(0, 0)) return 0;
    else grad = normalize(grad);

    vec2 gradP = perpendicularLeft(grad);

    float val = getBilinear(img, p);
    float valLeft = getBilinear(img, p+gradP);
    float valRight = getBilinear(img, p-gradP);
    return (val - (valLeft + valRight) * .5f);
}

It basically computes the second directional derivative in a direction perpendicular to the gradient. It works fine for my purposes, except that it's somewhat slow in the CPU implementation (and I have reasons not to port it to GPU).

Any faster (and/or more accurate) method?

Edit: I found a twice-faster method, still calculating the second derivative, but in a simpler way. It's equation (5) on this page. It works, but the result is more crude. I think I should switch to a higher-quality differentiation kernel than [+1, 0, -1].

Gradient magnitude and angle can be found using a Sobel operator. It's more sophisticated than the [1, 0, -1]. This is commonly used either directly for edge detection, or as part of a more complex algorithm.

Or, if I've misunderstood, and you want only the curvature of the 2D isoline, this answer on Math Overflow might be of use.