In the last few months I’ve been working on a rendering engine capable of rendering an entire planet. I’ve seen many solutions for procedural planet rendering (generating the height map on the fly), so I was focusing on a different problem: how to render a planet using already existing elevation data, and which is the most efficient way to store that data at large scale.

Sources

• For development/debug/demonstration purposes I’ve been using the elevation data of the Earth available at http://viewfinderpanoramas.org/dem3.html. I use the 3 arc second version which has full coverage.
• For texturing, I use the images available at https://www.h-schmidt.net/map/.
• For loading the textures in jpeg I use Sean Barrett’s stb_image.h from https://github.com/nothings/stb.
• For rendering the atmosphere I implemented the algorithm described in the article https://hal.inria.fr/inria-00288758/document. (I haven’t implemented everything from here, for example there are no light shafts yet.)

Performance

I measured the performance on my (at this point, quite weak) laptop:
CPU: i7-6700HQ 2.6 GHz
GPU: NVIDIA Geforce GTX 950M, 4 GB
HDD: 1TB, 7200 RPM

First, some words about the rendering: the geometry is divided into meshlets (triangle batches containing 64 triangles) to make the pipeline suitable for mesh shaders, and to allow further optimizations like meshlet-culling which is not implemented yet. Currently, there is space allocated for 16384 meshlets, and the “culling” shader doesn’t do much, only decides which meshlet slots are active, and fills the indirect command arguments accordingly.

I measured the GPU performance using PIX. Rendering the above frame (during movement) takes about 11 ms. Most of it (~9.5 ms) is spent in the ExecuteIndirect call that draws the planet with a simple forward rendering method, at the moment. Sky rendering is ~1 ms. The remaining 0.5 ms is for resource barriers and generating new requested meshlets (vertices, normals, uv coords).

The tessellation, lod selection is done by the CPU, on the main thread. I measured the lod update with _rdtsc(). The results were unfortunately quite unstable, but it’s worth noting that the measured cycles were almost always below 1.5 million cycles, even when travelling with high speed close to the ground. When the camera is way above the ground, like in the beginning of the video, the measured values are around 100-200 thousand cycles.

It’s important to note that the loading of the height data from file to gpu memory happens on a different thread, only the requests get issued on the main thread.

Tools

Beside the main program, I made some tools for various purpose, I’ll describe them briefly.

• Shader Compiler: putting together a graphics (even a compute) pipeline in D3D12, allocating descriptor tables, setting the descriptors, binding the descriptor tables, setting root constants could be a pain in the ass. It’s not the worst thing in the world, but it’s better to have something that fprintfs some C++ code based on a short description of the layout, rasterizer settings, render targets, depth-stencil state, struct definitions in hlsl, etc… That’s what this tool does. Also, it compiles the HLSL code using D3DCompile(). Every other program that has shaders uses this tool to compile them.
• Height Generator: this tool converts the raw elevation data organized in HGT files to the format I use directly for rendering. After developing, and debugging it, it took ~24 hours to generate the necessary data with it.
• Sky Generator: this tool is for generating the textures required for sky rendering and aerial perspective.

Design Philosophy, Coding

Working on this project the main focus was on designing algorithms, data structures while taking account the capabilities of the hardware. So the focus was not on coding, in fact I wanted to avoid code-related design issues at all. For this I had two major “rules”:

• Making sure that iteration is as fast as possible, that is, I can change algorithms, data structures, hard-coded constants quickly without thinking about the code and waiting for compilation too much.
• Writing easily understandable code, that is, by looking at the code it should be relatively easy to figure out what the hardware does when runs that code, what data it loads, how it modifies it, and what data it writes.

In order to achieve these, I don’t use templates, STL, the concept of RAII, or any third party library beside the Win32 API, D3D12, XInput, stb_image.h, and some standard C functions (math functions, functions for string manipulations, memcpy). I use lot’s of raw pointers, relative pointers and plane old data structures with global functions. I have only one compilation unit per program which can be built by running a batch file. For text-editor and debugger I use Visual Studio 2017.

With these concepts, I can compile the main program and the tools from scratch in about 2-3 second per program, using -Od (with -O2 it takes only slightly longer).

Definitions

Let be an oriented, differentiable 2-manifold, and its Gauss map (where is the 2-sphere). For such a surface let be the volume form generated by the embedding , that is, for . Let’s assume is invertible. The normal distribution of is defined by the equation

that is, it’s the coordinate of the pull-back of the volume form of by the inverse of the Gauss map, expressed in the basis defined by the volume form of the 2-sphere.

Next, let’s defined the projected area of along a given direction . The projection associated to can be defined by the volume form , . The projected area of along can be defined as

Computing the SGGX normal distribution

Let’s compute the normal distribution for the case where M is an ellipsoid defined by

where is a symmetric, positive definite matrix with eigenvalues , which are the sizes of the main axes of . For a given , let’s choose vectors such that form an orthonormal basis of . Then we have

,

since the Gauss map of the 2-sphere is the identity map. So we got

Now let’s compute the map . This map is the multiplication by the Jacobian of .

This Jacobian can be written as

,

and after some straitforward calculation we get

where we used that is an orthonormal basis. Putting these together, we have

Computing the projected area

Now let’s fix a direction , and compute the projected area . Note that the matrix viewed as a linear transformation, gives a diffeomorphism .

This diffeomorphism can be used to compute the projected area by transforming the integral to the 2-sphere:

Using that the Gauss map is , a simple calculation shows that

Moreover, for , we have

We can express this volume form as for some function . We can compute by, for example, choosing a point and such that form an orthonormal basis. Then , so

By setting , we have , and the projected area is

I’ve been working on an implementation of a mass-spring system using Hamiltonian mechanics which I’d like to apply for hair, fur simulation in the future. I’ll probably write a post about it (it depends on how it goes), but in the meantime, I thought I write something about how I handle resource synchronization between GPU and CPU in my current project (link) (this is the project I use for experimentation, for example the code for my previous posts about normal mapping is in this project).

I use Direct3D 12 which provides the below objects and functionalities that are relevant for us at the moment:

• ID3D12Resource: This object represents any kind of resource on the GPU (meshes, textures, any piece of memory, or even just a “virtual” resource without memory binding to it).
• ID3D12CommandQueue: This object represents a command queue of the GPU, it can be used for dispatching commands (e.g. prerecorded command lists, instruction for presentation) to the GPU.
• ID3D12Fence: This is the main object for synchronization. A fence is just an atomic 64-bit unsigned integer which can be set (both by the CPU and the GPU), can be read (by the CPU), and can be used to block a thread (a thread on the CPU or an execution of a command queue on the GPU) until the fence reaches a certain value.

. Command Queue

I use the struct below to represent command queues.

struct CommandQueue
{
	ID3D12CommandQueue* d12commandQueue;
	u64 lastSignaledFenceValue;
	ID3D12Fence* d12fence;
};

Here’s how it works: after the required commands are dispatched to the d12commandQueue, to be able to tell later whether the GPU has finished the execution of the these commands, a synchronization point can be made in the queue by increasing the lastSignaledFenceValue by one, and dispatching a signal command to the d12commandQueue which sets the d12fence to the lastSignaledFenceValue. This process is done by the function called signalNext.

static void signalNext(CommandQueue* commandQueue)

So in the future, to check whether the GPU has executed the dispatched commands, we have to check whether (the current value of) lastSignaledFenceValue is smaller or equal than (the value in the future of) d12fence.

Of course, in theory, the fence can overflow, but since nanoseconds is about years, overflow is basically impossible in practice.

Tracked Resource

I use the structs below to represent a synchronizable resource:

struct FenceSlot
{
	ID3D12Fence* d12Fence;
	u64 requiredFenceValue;
};

struct TrackedResource
{
	ID3D12Resource* d12Resource;
	FenceSlot useFenceSlots[4];
	FenceSlot modifyFenceSlot;
	D3D12_RESOURCE_STATES stateAfterModification;
};

A fence slot is just a fence associated with a fixed fence value. A fence slot is considered valid if the fence hasn’t reached the given value, and it is invalid otherwise. A fence slot can be waited on by waiting until the fence reaches the given value. A valid fence slot of a resource represents that the resource is being used by a queue.

A resource (in general, not only here) should be modified (write) by only one thing at once, and not used by any other thing during modification, and can be used (read) by multiple thing at once, while it’s not being modified. Keeping this in mind, a TrackedResource object has one fence slot which indicates that the resource is being modified, and four fence slots which indicate that the resource is being used. (Of course the number four is completely arbitrary here, and in fact, the fence slots don’t have to be stored inside the resource as an array, but as a forward list, allocating them dynamically from a memory pool for example, if necessary.)

Additionally, the TrackedResource struct has a D3D12_RESOURCE_STATES member, which indicates the current state of the resource. More precisely, it indicates the state of the resource after a modification finished, that is, after the fence slot used for modification becomes invalid.

Functions for resource synchronization

The two key functions for assuring proper synchronization of the resources are the markUse and markModify functions.

static void markUse(
	TrackedResource* resource,
	CommandQueue* commandQueue
);

static void markModify(
	TrackedResource* resource,
	CommandQueue* commandQueue,
	ID3D12GraphicsCommandList* commandList,
	D3D12_RESOURCE_STATES stateToModify
);

To mark that a resource will be used (but not modified) by a command queue, we have to call the markUse function. First, this function checks whether the modifyFenceSlot of the resource is valid. If it is, it will make the queue wait for the modifyFenceSlot. After that, the function records the commandQueue->d12fence and (commandQueue->lastSignaledFenceValue + 1) values to a (currently invalid) useFenceSlot of the resource. (So the function assumes that the queue will use the resource until its fence reaches the next fence value to be signaled. This is also assumed by the markModify function.)

To mark that a resource will be modified by a command queue, we have to call the markModify function. First, this function makes the queue wait for all the valid fence slots (modify and use) of the resource. Then it records the commandQueue->d12fence and (commandQueue->lastSignaledFenceValue + 1) values to the modifyFenceSlot of the resource. Also, if the stateToModify parameter differs from the value of resource->stateAfterModification, the function records a resource barrier to the commandList which makes the proper state transition of the resource and sets resource->stateAfterModification to stateToModify.

Issues

In my project this synchronization model is pretty useful. It’s also simple, requires only little code and during debugging it’s relatively easy to see what’s going on. Of course, it has to be used cautiously, otherwise deadlocks can easily occur. For example, after calling markUse and/or markModify on some resources, one has to call the signalNext function on the command queue, otherwise the fence slots set by markUse and/or markModify will stay valid forever.

However, I encountered a pretty nasty deadlock that I considered a bug, not a misuse of the model. Perhaps you’ve already found it… Let’s consider for example what happens if we call markModify on a resource, and after that, in the same “pass” (before calling signalNext on the command queue, so the lastSignaledFenceValue of the queue is unchanged) we call markUse or markModify again on the same resource. In this case the queue is forced to wait for a fence value that should be set by it later… But the solution to this problem is pretty easy: the markModify function just has to check that all the fence slots that we make the queue wait on is not the same as the new fence slot that is being recorded.

This post is the continuation of the previous post. (link)

In this post, I will continue from where I finished the previous post. First, I will talk about how to compute the normal vectors for a surface defined by a height map, and how they get transformed when we transform the underlying height map. After that, I will derive to proper formula for the normal vector of a surface whose vertices are offset with a height map. I used this formula for computing the normal in the project I talked about in the previous post. (For more info, please check the previous post.)

Connection between the normal and the height map

Let’s take a height map, i.e. a function. Let’s assume that this function is smooth. The graph of this function, , defines a smooth surface in .

What is the normal vector of this surface? Or in other words, what is the normal map, corresponding to the height map ? To calculate that, we can take the cross product of the tangent and bitangent vectors, and normalize the result. So, for the unnormalized normal vector (denoted by ) we have

and for the normal vector, we have . Now, what happens if we scale the height map by some number , that is, we consider the height map ? The partial derivatives will get multiplied by as well, so to compute the new normal from the old one, we have to do the following steps:

1. Divide the old normal by its z coordinate to get the old unnormalized normal back.
2. Multiply the first two coordinates of the old unnormalized normal by to get the new unnormalized normal.
3. Normalize the new unnormalized normal.

As you can see, the normal vector does not behave so nicely during the transformation of the height map, the partial derivatives (gradient) behave in a much simpler fashion. In fact, if we apply a more complicated transformation to the height map (for example multiply, add, divide two height maps together, take the exponential, logarithm, sine, square root of a height map, or whatever you like), the gradient of the new height map can be computed relatively easily and intuitively from the old gradient map, using the chain rule for derivation. Then the transformed normal map can be computed from the transformed gradient map.

In my opinion, it’s always more straightforward to work with the gradient map directly instead of the normal map. It’s not only that gradient maps can be transformed in a more naturally way, but in case of linear sampling and mip level generation they produce more correct results. This easily follows from the fact that differentiation is a linear operation, but the operation of taking the normal vector is not. For example, in this project, to compute the mip levels for a normal map I generate all mip levels of the underlying height map, and generate each level of the normal map from the corresponding level of the height map, instead of generating the biggest level of the normal map from the height map, and generating the mip levels of it by downsampling.

Also, as we will see below, if we want to apply normal mapping properly, what we need directly is not the normal map of the height map, it’s rather the gradient map of it.

Normal mapping

So let’s see, what we have and what we want. We have a surface . As I wrote at the end of the previous post, we can assume that the texture coordinates of a point is . We also have a height map . Let’s define the function in a way that is the normal vector of the surface at the point . (This map is called the Gauss map of the surface .) From and we can define a new surface by offsetting the surface along its normal by the amount , that is,

.

(This is the equation I use in the vertex shader to offset the vertices.) What we want is to find the normal vector of the surface . For that, we have to compute the tangent and bitangent vectors of this surface, and take the cross product of them. For the tangent vector, we have . The terms look familiar, except for one: we haven’t seen the term yet. This term represents how the normal vector changes if we move along the tangent direction on the surface. It can be viewed as a second order term for describing the surface (the first order terms are the tangent and bitangent vectors), and can be computed from the first and second partial derivatives of the surface . Without going into too much details, I say that in differential geometry there is a thing called the shape operator (or Weingarten map) for a given surface. This operator tells us how the normal vector changes along a vector from the tangent plane on the surface, and this is exactly what we need. It can be represented by a 2×2 matrix which for the surface can be computed as

where denotes the usual scalar product on .

I mention that the shape operator is related to the curvature of the surface: its eigenvalues are the principal curvatures, its determinant (the product of the eigenvalues) are the Gaussian curvature, the half of its trace (the sum of the element in the main diagonal, or the sum of the eigenvalues) is the mean curvature of the surface.

For us, the important thing is that we have . (Don’t care about the negative sign, let’s just accept it as a convention.) So for the tangent vector of the offset surface we have

,

where is the usual (but not orthogonalized!) TBN matrix. Similarly, for the bitangent vector of the offset surface we have

.

And actually, that’s it. As usually, the cross product of the tangent and bitangent vectors gives the unnormalized normal , and normalizing that we get the required normal vector .

Computation in practice

In order to do the above computation in the pixel shader, we need the following data

• TBN matrix (at least the tangent and bitangent vectors),
• the height map,
• the gradient map (or indirectly, the normal map) of the height map,
• the shape operator of the mesh.

The shape operator can be placed in the vertex data, next to the TBN matrix. Note, that even if the TBN matrix is not orthogonal, the normal vector is always perpendicular to the tangent vectors, and has unit length, so it can be computed on the fly, doesn’t have to be stored.

The shape operator is not large, it’s four floats. But anyway, is this the shape operator really necessary? As I said above, it can be seen as the second order terms of the surface, so in some way, they are less significant as the first order terms (the tangent and bitangent vectors). For the shapes in the video in the previous post I didn’t see too big difference between using the shape operator in the computation, or assuming it is zero. I saw only a little change in the rendered image when I switched on and off the shape operator in the shader (and hot-reloaded it). I want to do more experiments with it in the future, but it’s obvious from the above calculation, that the bigger the vertex offset (the values in the height map), the more significantly the shape operator influences the normal. Also, of course, the matrix of the shape operator depends on the actual surface, the more it curves, the bigger the values in the matrix are. For example, the shape operator of the plane (not surprisingly) is the zero matrix, and the shape operator of a sphere with radius is

.

The case when TBN is orthogonal

Let’s finish with assuming that the TBN matrix is orthogonal, the shape operator is zero, and compute the normal in this case. Hopefully, we will get back the formula that is usually used for normal mapping.

If the shape operator is zero, we have

, and

.

And if TBN is orthogonal, we have

.

Does this look familiar? I hope it does. The last term in the above equation is the unnormalized normal of the height map that we’ve seen before. That is what is usually stored in the normal map, and that’s what we multiply by the TBN matrix to get the normal vector of the normal mapped surface.

I started a project recently in which I’m experimenting with some random stuff. I think it can be called a “handmade” project, currently I use only the C runtime library, Windows API, Direct3D12 and the AVX instruction set. (I don’t use COM objects or the d3dx12.h helper header.)

Anyway, I made a simple Perlin noise generator, and used it for creating height maps, and some images with “infinite detail” (that is, they can be zoomed into arbitrary many times). Let me call the latter ones simply just fractals in the following. I also mixed the two types, that is, I made a height map which is a fractal as well. Here is the link to the github repo:

https://github.com/SandorDeak/Experiments

I guess the code is in a very messy state, and it will definitely stay in this state for a while, since I don’t have any reason to make a proper interface and code clean-up, and anyway, it’s just the code. Below, there is a demonstrating video. Unfortunately the recording caused some frame rate issues and popping.

Currently the images are generated by the CPU, and for the fractals, they are constantly recomputed asynchronously and uploaded to the GPU when they are ready. The non-fractal, static images (height maps) and the meshes are generated at startup.

About the height mapping: after generating the height maps with Perlin noise, I computed the normal maps from them. I used the height map for offsetting the vertices in the vertex shader (yep, no nice enough tessellation yet, perhaps in the future), and used the normal map in the pixel shader to compute the proper normal. As you can see in the video, the height offsetting is relatively large, so it’s important for the normal mapping to be precise enough, that is, the tangent, bitangent vectors, and the scaling should work well. In the remaining of this post I describe briefly what normal mapping is, how it is usually done (more precisely, I describe what I’ve seen in other tutorials about how it is usually done). In the next post, I will derive the precise formula for how it should be done, at least in theory.

Preliminaries

To make the geometry of a 3D object look more detailed, than the geometry defined by the vertices, we usually use normal maps. A normal map is just a texture whose value at a point is the more detailed version of the normal vector, expressed in the basis defined by the tangent, bitangent and the less detailed normal vector. The tangent vector encodes the direction and rate in which the first texture coordinate changes along the surface. The bitangent vector is the same but for the second texture coordinate. More precisely, if we have an function, which is the parametrization of our surface, where is some 2-dimensional subset of , and also, we have a function, which gives the texture coordinates (uv coordinates) of the surface; then the tangent and bitangent vectors are the first and second colums of the matrix , respectively. (For a function , let’s denote its Jacobi matrix by .) Note, that .

Before going further, let’s take a well-known example. Let’s say that our surface is just a triangle, with vertices , and with correspondig texture coordinates . Then both and are just parametrizations of triangles, namely is the triangle whose vertices are , and ,

, and

.

Hence and , so the tangent and bitangent vectors are the columns of the matrix . This expression is usually used to calculate the tangent and bitangent vectors for a triangulated mesh at the vertices . Of course, a vertex can be a vertex of more than one triangles (and usually is). In that case, for a vertex we can simply take the average of the tangent and bitangent vectors computed for the triangles containing that vertex. It’s important to note that the normal vector is always perpendicular to both tangent vector (after all, this is how the normal vector is defined), so we have to make the averaged tangent and bitangent vectors perpendicular to the normal vector.

It’s important to notice that there is no reason for the tangent and bitangent vectors to be orthonormal (have unit length and perpendicular to each other). But we should usually aim for that, since that would mean that the texture is mapped to the surface without any distortion. Since the mesh is usually fixed, this means that we should define the texture coordinate function in a way that the columns of the matrix are orthonormal. Intuitively, this just means that we have to unfold the surface, flatten it without any stretching. If that can be done, it’s great, but it’s usually not the case. Without going into details, I mention that for a smooth surface, such exists only if the surface has zero Gaussian curvature, and that’s just simply not the case for most surfaces.

Despite that, (at least what I’ve seen in every tutorial about normal mapping) it is assumed that the TBN matrix (the matrix consisting of the tangent, bitangent, and normal vectors) is orthogonal. I guess the reason for this is that it’s usually a good enough approximation, but as we seen above, it can be good enough only if the texture coordinate function is chosen well enough. In this project, this was not the case, and I got horrible results when I orthogonalized the TBN matrix.

Another problem occurs when we want to scale the height map. For example, how should we change the normal vector, if we offset the vertices by the value of the height map multiplied by two? What I’ve seen so far, people usually use an empirical intensity value for the normal map, which is used to scale the coordinates of the normal map corresponding to the tangent and bitangent vectors (the first two coordinates). This, of course, makes sense, but only if there is no actual height map being used for the normal map, otherwise the scaling of the height map and the normal map has to be synchronized properly. Actually, synchronizing them is pretty easy, but I keep the explanation for the next post.

One thing before wrapping up: notice that we can assume that the parametrization of a given surface can be selected in a way that its domain is the texture coordinates, that is, we have , and the texture coordinates of the point is just . This can be achieved by replacing with . This simplifies things a little bit, since then, the tangent vector is , and the bitangent vector is . In the next post, I will assume that this holds for .

Reading this very useful series about the graphics pipeline I realized I had never thought about what is the proper way to interpolate between vertex attributes (outputted by the vertex, domain or geometry shader to the pixel shader) using the vertices of the triangle after perspective divison. I spent some time thinking on that, and although it’s not hard to figure out and the result is described in the link above, I thought I describe one possible way to get the solution.

Let’s assume we have a triangle before perspective divison with vertices , , with some corresponding vertex attributes , , (e.g. texture coordinates, color, normal). After perspective divison we have , , . Let’s have . The question is, how should we compute the proper vertex attribute corresponding to ?

Let’s think backward: assume we have . Then, of course, the vertex attribute corresponding to is simply . Now, what is the point we get from after perspective divison? Using the notation we have

Well, it’s easy to check that the sum of the coefficients in the above interpolation between , and is . So by setting and , we can write .

And basically, we are done: for a point , the corresponding vertex attribute is , where and can be computed from the two equations and , where .

Of course, we can solve for and explicitly (it’s just a system a linear equations with two variable), and we get

, and ,

where . We can go further by substituting and into the equation for

Therefore, to do the interpolation in practice, we have two main choices:

• For given and values (that is, basically, for a given pixel inside the triangle) we calculate and , and use them to calculate all vertex attributes by a simple linear interpolation.
• For all vertex attributes , we compute the corresponding values, and use them to calculate the vertex attributes directly from given and .

Of course, in both cases we have to compute as well.

Here’s a story I think worth sharing. At the beginning I learned C++ programming mostly on my own by watching tutorials, CppCon videos reading posts and practicing. At some point (in my graphics engine project) I needed dynamic polymorphism, but I didn’t want to use OOP style inheritance, I just wanted to have a piece of memory that knows its type that can be changed dynamically, in a relatively convenient way. So I had this genius idea to make a struct containing an enum value for storing the current type, and a buffer for storing the data. What should the alignment and the size of this buffer be? To calculate that in “the most general way” in compile time I did the following. The possible types were specializations of a struct template, where the template parameter was the enum for the type. Then I used constexpr functions which iterated through the possible types, and calculated the required alignment and size for the buffer that could store them. Now it seems completely ridiculous.

Later I learned that there is this thing called union in C++… It could be a coincidence or it could be completely my fault, but it seems very unfortunate that I hadn’t met the union keyword in the sources I had learned C++ from, and learned it only pretty late . Namely because I think unions are very useful and I’ve been using them a lot.

I have found some sources about how to transform normal vectors and why they transform in a way the do. These sources seemed to rely on the fact that the transformations we use are affine transformations (multiplying with a matrix and adding a vector). Using linear algebra and vector calculus, I would like to give a simple derivation of how the normal and tangent vectors are transformed by a general transformation.

Let be a differentiable function, and is a smooth surface. Let’s denote the image of under by . If is “nice enough” is another “well-behaved” surface. (For example, you can keep in mind the usual case where , where is a vector, and is an invertible matrix. Then, is “nice enough”.)

Let’s fix a point . Let’s denote the tangent plane of at by . Similarly, let and be the tangent plane of at . Let’s define first how the tangent vectors in are transformed by to , so let’s make a function from in some meaningful way… It turns out that the most meaningful way for that is to do the following. Let’s fix a vector. We can represent this vector with a small piece of curve in M, for which and . Now let’s apply to . Obviously, , so by taking the derivative of at , we get a vector . And that’s what we wanted! That is . Let’s calculate it, using the chain rule for derivation

,

where is the Jacobi matrix of at .

So we got that is a linear transformation, it is the multiplication by . For the case where we have , therefore, the Jacobi matrix is constant, it is for every .

Now let’s see how the normal vector is transformed. For simplicity, let’s use in the following. First, what is the normal vector of at ? It is a (nonzero) vector perpendicular to . Similarly, the normal vector of at is a (nonzero) vector perpendicular to .

First, let’s assume that is invertible. In this case, for every vector we can write

.

Since is invertible, every vector can be written in the form of for some . (Namely, for ). So we got that the vector is perpendicular to (and nonzero), so it’s a normal vector of of . Therefore, we can say that the way transforms the normal vector at a point is by multiplying it with the inverse of the transpose of the Jacobi matrix of at , if this matrix exists. If is on orthogonal matrix, then of course, we have .

What if is not invertible? This roughly means that as a function from to collapses some (one or more) directions into zero at . If these are not the direction in , we can still define the normal vector, but I won’t go into details now. Otherwise, we can’t, because actually, in this case is ill-defined, and there is probably a crack in at .

There seem to be many sources about premultiplied alpha out there, so here’s another one. 🙂 This one makes the most sense to me, and I couldn’t find this explanation in this form elsewhere. (Although I haven’t done too much research.)

When blending images together it’s useful to think about them as filters, functions, applied to the colors already in the framebuffer. An image
can be viewed as a transparent light source where is the emittance and is the transparency. Then applied to a color results in the color , that is, it absorbs some light, and adds to the remaining light its own emittance value . This is just the blend equation using premultiplied alpha. Because of that, premultiplied alpha seems a more physically plausible blending, than non-premultiplied alpha.

Let’s see how should we blend together two images. it’s actually quite simple and natural. Let’s have two images and an apply them to the color . We get

So for the image we have . This shows how to blend images together. Expressing the transparencies with the alpha channels, we got the blend equation for the alpha channel

Associativity is usually a required property of blending equations, which is the property that for each images we have . The blending above is by definition associative, since it was defined simply by function composition, which is always associative.

Welcome! There will be some hopefully interesting things here, but until then feel free to explore my blog, read about my earlier projects and ask questions.