A cute little trick to running classic IIR filters on the GPU

In my previous behemoth post I mentioned the struggle of implementing a comb + notch filter that looked passable for analog composite video. The blurring effect of adding the notch with FIR filtering was far too great, and there’s no way TVs at the time did it that way.

Nagging issues like these haunt me until I solve it, and I started experimenting with IIR filtering since it’s one of those topics I never really explored on GPUs, due to it being very GPU unfriendly. A notch IIR biquad filter is extremely effective on paper, and there has to be a way to do it well in parallel.

A biquad IIR

This is the building block of a ton of audio filters out there, since they can easily model classic analog filters and are easy to tweak parametrically. In the Z plane they have two zeroes and two poles:

// Z is a complex number representing where frequency response hits 0.
// P is a complex number representing a pole
// where the frequency response hits infinity (or pure resonance).
// Frequency response at a given frequency is measured by setting
// z = exp(2 * j * w)

// I probably screwed this formula up somehow,
// don't read too much into it.
H(z) = (z - Z) * (z - conj(Z)) / ((z - P) * (z - conj(P));

Having conjugated poles and zeroes ensures the filter coefficients are real, so the actual filter looks more like this in code:

out[x] =
  in[x] * b0 + in[x - 1] * b1 + in[x - 2] * b2 + // FIR portion
  out[x - 1] * a1 + out[x - 2] * a2 // IIR portion

The Z plane is often a little abstract, but in this case, we will design the notch directly in the Z domain. For example, we place zeroes where we want to remove frequency response, e.g. at the chroma subcarrier, then place poles close to it at the same frequency. The intuition is that close to the carrier, the response will fall to 0, and further away the pole and zero mostly negate each other, leaving a unity frequency response. The response can be made as sharp as desired, but sharper notches lead to more ringing.

E.g. for a pole set at 0.9 times the zero. Quite sharp. Getting this kind of response out of a normal convolution requires a comical number of taps. The response at DC here isn’t exactly 0 dB, but the FIR portion can be rescaled to ensure unity gain.

struct
{
  float b0, b1, b2, a1, a2; // a0 is implicitly 1.0.
} push = {};

// Compute a biquad IIR notch filter.
double w = 2.0 * muglm::pi<double>() / 27.0;
constexpr double fir_q = 0.99;
constexpr double iir_q = 0.85;

// Don't do a complete notch since that leads to really bad ringing.
// We just need a fairly narrow band-stop though ...

if (options.system == System::NTSC)
  w *= 315.0 / 88.0;
else
  w *= 4.43361875;

// Preconvolves two zeroes at exp(i * w) and exp(-i * w).
double b0 = 1.0;
double b1 = -2.0 * fir_q * cos(w);
double b2 = fir_q * fir_q;

// Preconvolves two poles at exp(i * w * Q) and exp(-i * w * Q).
// Flip sign here since we're designing the filter as B(z) / A(z) and
// when evaluating the filter we flip the signs.
double a1 = 2.0 * iir_q * cos(w);
double a2 = -iir_q * iir_q;

// Normalize the FIR gain to obtain a 0 dB gain at DC.
double fir_gain = b0 + b1 + b2;
double target_fir_gain = 1.0 - a1 - a2;
double fir_scale = target_fir_gain / fir_gain;

b0 *= fir_scale;
b1 *= fir_scale;
b2 *= fir_scale;

push.b0 = float(b0);
push.b1 = float(b1);
push.b2 = float(b2);
push.a1 = float(a1);
push.a2 = float(a2);

cmd.push_constants(&push, 0, sizeof(push));

The trick – a spicy parallel scan

While IIRs have infinite impulse response, we’re only really interested in the output of the active signal area, which is ~1440 pixels. A 1440-tap convolution is fine if we want to use FFT convolution, but that’s quite overkill for what’s literally a 5-tap filter. Surely there are ways to do IIRs effectively without falling back to serial algorithms.

Many parallel GPU algorithms boil down to scan operations in some way, or have scans as a primary ingredient of their madness. I found an older paper from 1990 (Prefix Sums and Their Applications, Guy E. Blelloch, section 1.4) which goes through this algorithm in a way that was easy enough to digest. It’s not the state of the art I think, but it gets the job more than adequately for my needs. For the IIR filter, the FIR portion of it is handled trivially with a filter kernel. We only need to concern ourselves with the feedback portion of the filter.

Usually we think of scans as just simple binary operators.

for i in range(len(data)):
  data[i] += data[i - 1]

This is common enough that we have special subgroup operations for them, like subgroupInclusiveAdd. This is not enough to implement IIRs since it cannot capture the “decay” of an impulse response. The two intuitions from the paper is that the binary operator can be as complex as we want it to be as long as it follows some properties that are expected of binary operators, and the data doesn’t have to be a scalar either. For an IIR with two histories, we need a 2×2 matrix and 2 element vector, then define a binary operator which is effectively:

struct Biquad
{
  mat2 A;
  vec2 B;
};

Biquad init(float v)
{
  // Note that GLSL is column major here.
  return Biquad(mat2(vec2(a1, a2), vec2(1, 0)), vec2(v, 0));
}

Biquad scan(Biquad prev, Biquad current)
{
  return Biquad(prev.A * current.A, prev.B * current.A + current.B);
}

The intuition is that the A matrix keeps track of impulse response decay over time, while the B vector serves as the history buffer. Once the scan is done, the result is extracted from B.x and written out to memory. The performance is excellent and processes my frame at like sub 10 microseconds. Can’t complain about that.

For the full shader example, see https://github.com/Arntzen-Software/parallel-gs/blob/main/analog/shaders/iir.comp, but the general gist of it is:

256 thread workgroup – work on N pixels per thread

The horizontal resolution is expected to be reasonable, so N is about 6 in my use case. Not having to split the scan over multiple kernels is good.

This means we can load 8 inputs and compute the FIR portion of it for 6 inputs all in registers. The initial scan can also be done in registers:

float xs[PER_THREAD + 2];
Biquad quad[PER_THREAD];

[[unroll]]
for (int i = -2; i < PER_THREAD; i++)
  xs[i + 2] = LoadInputInRow(PER_THREAD * int(index) + i);

// Do the FIR portion and produce 6 values.
[[unroll]]
for (int i = 0; i < PER_THREAD; i++)
  quad[i] = init(xs[i + 2] * b0 + xs[i + 1] * b1 + xs[i + 0] * b2);

[[unroll]]
for (int i = 1; i < PER_THREAD; i++)
  quad[i] = scan(quad[i - 1], quad[i]);

Custom scan with magic biquad operator

Then we do the custom scan over the subgroup using ShuffleUp.

for (uint step = 1; step < gl_SubgroupSize; step *= 2)
{
  Biquad prev;
  prev = subgroupShuffleUp(quad[PER_THREAD - 1], step);
  if (gl_SubgroupInvocationID >= step)
  {
    [[unroll]]
    for (int i = 0; i < PER_THREAD; i++)
      quad[i] = scan(prev, quad[i]);
  }
}

… and then complete the scan via groupshared memory with same strategy. Groupshared memory is kept at a minimum since we only need to store the last value per subgroup. That’s about it, and with that I have a decoded composite signal I’m happy with.

For very large IIR filters, this approach gets out of hand quickly since we need a big N x N matrix to hold the recurring state per element, but chaining together biquads is done trivially by just scanning multiple times. Since we can do the horizontal filter in one compute workgroup, the intermediate steps can easily live in shared memory if desired to avoid excessive memory bandwidth.

Emulating old junk from yesteryear – or my obsession making native resolution PS2 emulation look good

Lately I’ve been on a kick to tackle the latter part of PS2 graphics emulation that never seems to come up, analog TV emulation. 99% of the PS2 emulation space is all about upscaling to the max, cranking out 4K with polygons so sharp they can cleave mountains in half, but that’s not particularly to my taste, which is why paraLLEl-GS focuses more on the super-sampling aspects rather than raw upscaling. Earlier blog post on the paraLLEl-GS here if you have no idea what I’m referring to.

For example, a raw native render with progressive scan (640×448) with paraLLEl-GS backend:

With Vulkan HW renderer in PCSX2, a basic 2x upscale would look like:

which is obviously higher resolution, but does not attempt to resolve aliasing either, and as upscaling factors go up, the mismatch between texture resolution, polygon counts and output resolution create a jarring effect for me. The approach I tend to favor is super-sampling and keeping the resolution native, even if it means a more blurred resolve. Here’s how it’d look with 16x SSAA. This is quite overkill for most cases, but why not.

There is a special mode in paraLLEl-GS that can scanout the 16x super-samples at double the resolution, which is effectively 4x SSAA at double the resolution. It can look great when it works well, but the key word is when. This approach is not playing to the strengths of paraLLEl-GS at all, but it’s a thing.

The main problem is that 2D images like HUD elements remain at native resolution with integer scale, which can create a jarring look, and post processing passes can mess up things in some cases. My main focus is on the native resolution, super sampled output, since it has very few gotchas compared to other upscaling styles.

The stopgap solution – FSR1

For a while, I’ve relied on FSR1 to do the upscaling, which is just a temporary hack. While properly anti-aliased 3D can look surprisingly good with FSR1 when blown up to large resolutions (what it is designed for), 2D game elements still create questionable artifacts. Text upscaling starts looking like those old HQ2x, SuperEagle and xBR filters that used to be popular with SNES emulators. That’s likely because FSR1 is “just” edge-aware Lanczos filtering at its core. Still, at SD resolutions there’s only so much you can do to blow it up to a modern display. The result here is of course quite blurry, but looking at it from a reasonable distance, it can look sort-of okay on a good day.

The only proper way to display SD content is in my opinion to use a CRT, or in lieu of getting a nerd tan, CRT shaders. Simulating CRTs has been done to death to the point Hel is blushing, so I feel a bit uncomfortable even trying to write about it, but this post wouldn’t be complete without it. RetroArch has like 23428323 CRT shader presets already, and there’s nothing novel about any of this. However, there are some considerations for PS2 that most CRT shaders don’t target:

  • Analog TV input path
  • 480i instead of 240p
  • HDR output
  • High refresh rate simulation

A lot of CRT shaders assume a high quality VGA-style monitor. Nothing wrong with that of course, but I find most of them a bit “too” good for PS2. Where we’re going, we’ll need that analog fuzz too.

Analog signal paths

I’ve been deep down the rabbit hole looking up old specifications to get a deeper appreciation for the brilliant engineering involved in making color TV work all those years ago, and the PS2 era was the last hurrah for SD CRTs, 480i warts and all. I’ve sort of gone through all of this stuff before back in my early days of programming, but I’d like to think I’m a little smarter than I was back then, and I learned far more details than I knew before.

From a graphics programming PoV this is hardly “difficult” stuff like debugging random GPU hangs at 10pm on a Friday coming from the latest and greatest AAA game, but hey, sometimes you just gotta relax with some good old DSP coding to stay sane. The signal processing for NTSC and PAL is fairly straight forward and it’s actually a good entry point into signal processing for graphics programmers since it combines very visual things with actual real world problems.

The cables

Component is the highest quality analog cable (actually, 3 cables!) available to consumers. It supports progressive scan which very few PS2 games natively supported. Simulating this cable is quite trivial.

Composite is the yellow cable “everyone” used. This one is the hardest since separating luma from chroma in the same single-channel signal is actually kinda complicated to do well. Most of the time I spent was trying various strategies for this problem and try to come up with a look that seems authentic and “tastefully shitty” 😀

S-Video is very similar to composite, except that luma and chroma are separate signals. Chroma is packed into one cable and chroma decoding is basically same as for composite signals, where phase and amplitude dictate the hue and saturation. Main difference to component is that bandwidth for chroma should be quite low, so more color smearing should be present.

Us Euro-bros technically had RGB SCART too, but I’ve never seen those cables myself for a Sony console. I believe our old N64 and GameCube had those, but I’d have to double check next time I visit the family …

Reference literature

While it is very easy to find a ton of content online about old TV standards, your random YouTube video is not going to have the minute detail needed to implement much. I found “Video Demystified – A Handbook for the Digital Engineer (2005)” by Keith Jack which has ton of great detail that is often left out to make sure my implementation stays grounded.

The bridge between digital and analog – BT.601

BT.601 defined how to take NTSC and PAL signals and turn them into the digital domain. It is the foundation for all digital video today. The primary interest for us is that the standard defines a 13.5 MHz sampling rate. This was supposedly chosen since it was convenient for both NTSC and PAL and filtering requirements. This scheme works out to 720 horizontal pixels for NTSC and PAL when you account for 525 lines at 29.97 Hz and 625 lines at 30 Hz. The H-sync part of the analog signal pads out to a bit over 800 pixels, but that’s irrelevant here.

I also learned just recently why it’s called YUV444, 422, 420, etc:

The “4:2:2” notation now commonly used
originally applied to NTSC and PAL video, 480i and 480p Systems
implying that Y, U and V were sampled at 4×,
2× and 2× the color subcarrier frequency,
respectively. The “4:2:2” notation was then
adapted to BT.601 digital component video,
implying that the sampling frequencies of Y,
Cb and Cr were 4×, 2× and 2× 3.375 MHz,
respectively.

Now you know! The Nyquist frequency of SD video is thus 6.75 MHz. Talking about video like this is a little weird, but it will make more sense later since analog signals like NTSC and PAL have a certain bandwidth, which limits how much horizontal resolution we can cram into it. The number of video lines is hardcoded. The video signal is really just a 1D signal if you look at it in an oscilloscope (haven’t seen those since my Bachelors …).

The next question was to determine the sampling rate of the PS2 CRTC. Given it generates analog signals, it’s not obvious that PS2 would use a standard frequency here. However, several old forum threads I found did indeed suggest that PS2 is based around the 13.5 MHz rate. This is likely because it could use off-the-shelf video DACs at the time. The nominal maximum horizontal resolution of PS2 is 640 pixels, but it seems like overscan is supposed to stretch the 640 pixels out to fill the screen anyway. Even if 480 lines are “visible” in NTSC, games typically just render 448 lines since the top and bottom is eaten by overscan by most TVs.

The actual CRTC clock seems to be 54 MHz, because when programming the CRTC, you’re supposed to set some dividers which ends up determining the resolution. E.g. 640 pixel width is a divider of 4, 512 pixels a divider of 5 and so on. It conveniently supports most relevant horizontal resolutions like this.

My armchair theory for how this works is that the CRTC runs at 54 MHz and uses the divider to do a zero order hold (aka nearest filtering) which is then fed into the DAC. Either that or the console is really doing horizontal linear filtering up to 640 pixels and does composite video generation at 13.5 MHz. I’m not sure what really happens, so I just have to guess. I’m not quite obsessed enough to have an oscilloscope to suss out micro-details like these from a real PS2. The PS2 is technically able to emit 1080i and 720p signals, so there’s no reason why it couldn’t do analog video processing at the full 54 MHz.

Step 0 – Prepare a test image

To debug any of this we need test images. I don’t exactly have test signal generators that TV engineers had back in the day, so I had to synthesize something. The use for these is to validate various edge cases in the NTSC and PAL decoding process.

The basic idea is to have some 75% color bars at the top to validate the color pipeline. The middle section is a sweeping increase of horizontal frequency up to Nyquist of 6.75 MHz. Each 1 MHz section is delimited by a single line color bar.

The lower portion is a sweep with increasing frequencies. This is used to test the comb filter.

Step 0.5 – Install GNU Octave

Doing FIR filter design by hand is uh … not something I have time for. Anything beyond the basic windowed sinc design process gets annoying very quickly. GNU Octave is a free alternative to Matlab that does what we need for these tasks.

Step 1 – Generate a bandlimited signal

While we can do composite video generation at the 13.5 MHz rate, it is not easy to avoid aliasing, especially since we will be modulating with a chroma subcarrier that will shift the spectrum all over the place, potentially creating aliases out of thin air. The handbook calls for a 2x oversampling to avoid this.

First, the image is padded out to fill 720 horizontal pixels for BT.601 reasons. Then, integer scale the image up to 2880xN (to match with how I understand the CRTC to work), and downsample that with a low-pass to 1440xN which completes a clean 2x oversampling. A polyphase upsampling filter is of course possible too to avoid the intermediate upscale, but that’s needless complexity for something that literally takes 10 microseconds on the GPU 😀 The immediate upscale is not written to memory at least, the filter just assumes a 2880 pixel input and it’s just sampling the input texture redundantly instead.

Since these are fundamentally analog signals, we only filter horizontally. The natural image type for these is 1D Array. All the processing shaders are 1D as well.

(The 1504 width is just extra padding for convolutions.)

Step 2 – Convert to YUV / YIQ / whatever

It’s all the same thing really. NTSC is the oddball one with YIQ but in practice this difference is completely irrelevant. The original idea of IQ was to take the blue and red difference signals and rotate them by 33 degrees to make it so that I would align with skin tones better, and give I more bandwidth than Q. However, the handbook doesn’t seem to give it too much consideration. Especially not for composite signals which are not meant to be sent over the air.

The implementation of this is really just doing a 3×3 matrix multiply with RGB, nothing special.

Step 3 – Low-pass Y and low-pass the hell out of chroma

NTSC and PAL have fairly narrow bandwidth defined for their broadcast signals, but composite signals don’t really have those strict limits. However, since the digital input has Nyquist at 6.75 MHz the handbook calls for bandlimiting the signal to about 6 MHz anyway. Super sharp falloff filters just lead to ringing.

Component video is specified by BT.1358 and doubles the sampling rate to 27 MHz. Y should fall off at about 11 MHz and chroma half that. Interpreting that in interlaced terms, the falloff should start at 5.5 MHz, getting close-ish to the Nyquist limit for SD video. After filtering luma and chroma, just decode back to RGB and we have a nice signal, done.

Chroma is nuked down to ~1.3 MHz or so for composite / S-Video. Sometimes, 0.6 MHz is called for it seems, but it’s quite unclear … S-Video can ignore the modulation + demodulation part if we just want to simulate the smear, but it was easier to let it go through the full chain for completeness. Despite the chroma bandwidth being so horrible, it sort of looks good. Amazing how terrible our eyes are at seeing color detail. I wonder if this is where the esoteric 4:1:1 YCbCr subsampling mode comes from now that I think about it …

E.g. NTSC luma filter, which starts falling off at about 4.2 MHz and stops at 6.75 MHz-ish. Broadcast NTSC is capped well below this bandwidth.

Chroma encode, with ~1.3 MHz passband:

Main difference for PAL is that luma passband is a bit higher. PAL-B (which Norway used) seems to specify 5 MHz rather than 4.2 MHz. Generating these filters can be done easily with firls and friends in Octave.

The handbook calls for a gaussian kernel for chroma to avoid any ringing, but I missed that memo 😀 Either way, these are implemented with trivial convolutions which GPUs eat up like butter.

Step 4 – Modulate chroma and combine into composite signal

This is the first point where NTSC and PAL differ significantly. While NTSC and PAL have different color subcarrier frequencies, PAL is also a bit more sophisticated.

// NTSC
C = Y + I * cos(wt + offset) + Q * sin(wt + offset)

// PAL
C = Y + U * sin(wt + offset) + sign * V * cos(wt + offset)

The “sign” of V is where Phase Alternate Line comes in. It flips every scanline. In broadcasting, this was meant to fix bad colors being introduced during broadcasting through complicated terrain (and boy do we have that over here in Norway). Bad phase shifts being introduced by NTSC will manifest as hue shifts. The basic idea behind PAL is that if phase shifts are introduced by signal reflections the sign flipping of V every scanline ensures that the decoded hue is complementary from scanline to scanline. This manifests as the Hanover Bar artifact. By averaging chroma from scanlines, the errors cancel out and we recover the correct color with slightly less saturation. The cost is of course reduced vertical chroma precision, but given how comically smeared chroma is horizontally, I’m not sure this matters, and digital video uses 4:2:0 subsampling anyway.

Now, broadcasting considerations are kind of irrelevant for something like composite video (I would think), but I’m not sure if PAL TVs skip the filter for anything not coming from an antenna. I kept the vertical chroma filter in my implementation because it’s neat to have.

The NTSC chroma subcarrier is constructed such that every scanline completes 227.5 cycles. Every line flips chroma phase which is very convenient and makes luma and chroma separation less complicated. The NTSC chroma pattern is a checkerboard as a result.

PAL is more annoying here. A half cycle method would not work since V is already flipping every line, so PAL chose 287.75 cycles. On top of that a tiny 1 / 625 cycle offset per line is added for … reasons. The V flipping leads to a chroma pattern where U takes a diagonal pattern while V has a pattern along the other diagonal.

Frame progression is also a concern. As fields are scanned, each time the same field is drawn, the chroma subcarrier should have opposite phases. This follows naturally from how fields are drawn. The period of NTSC is 4 fields:

  • Field 0: 262 lines, chroma carrier starts in + phase, ends in + phase due to even number of lines
  • Field 1: 263 lines, chroma carrier starts in + phase, ends in – phase due to odd number of lines
  • Field 2: 262 lines, chroma carrier starts in – phase, …

As we can see, things repeat after 4 fields and this point is easy to miss. The artifacts introduced by composite video should ideally cancel themselves out over time which manifests itself as flickery noise rather than a horribly glitched image.

PAL is annoying and has a longer cycle of 8 fields due to the three-quarter of a cycle setup. After completing 2500 lines (625 * 4), the chroma subcarrier has completed an integer number of cycles, and the sign of V is back to where it started.

After modulation, the NTSC color bar signal looks more like:

and the next line flips phase as expected:

A note on really old consoles

Some old Nintendo consoles (and likely others) emit NTSC and PAL in non-standard ways. E.g. NES is infamous for shifting the chroma carrier by 120 degrees instead of 180 which leads to very particular artifacts. See NESdev Wiki for more detail.

Step 5 – Validate composite output against color bar reference charts

It’s easy to mess up RGB to YUV conversions and the compositing process. The handbook had reference outputs for RGB inputs in NTSC and PAL where I could confirm that the math was indeed correct. The things to check are that minimum and maximum of the signal are what they should be. NTSC, at least the US variant of it maps black to 7.5 IRE (just think of it as some abstract voltage) and white to 100.0 IRE, and it tripped me up a bit at first since the NTSC color bars were defined in terms of this shifted and scaled IRE value.

Looking at the peaks and valleys of the generated composite signal in RenderDoc is enough since we just need to eyeball it. Close enough is good enough.

Step 6 – Separate chroma and luma – NTSC

This is non-trivial and a source of endless head scratching. Chroma information lives in the frequency spectrum around the carrier, but so does higher frequency luma detail.

The basic theory for comb filters is to take advantage of the opposing phase of the chroma carrier. For NTSC, I used this basic structure from the handbook:

In code, this translates to:

void chroma_separation_ntsc()
{
  float line0 = sample_previous_line();
  float line1 = sample_current_line();
  float line2 = sample_next_line();
  float estimate = 0.5 * (line1 - 0.5 * (line2 + line0));
  store(estimate); // Bandpass the result.
}

Ideally, by adding two neighbor lines, chroma should cancel out and only luma remains. Subtracting, we cancel luma and only chroma remains.

We know that chroma won’t (or at least shouldn’t) exist outside its bandwidth, so the result is run through a bandpass filter that centers around the carrier frequency, and we have an estimate for the modulated chroma signal. Since the composite signal is Y + C, we subtract the chroma estimate from composite signal to form a Y estimate. Chroma can now be demodulated and low-passed to remove the harmonics introduced by demodulation.

This filter works “perfectly” for regions where the chroma is constant, but not so much where there are discontinuities. This results in “chroma dots” where the color subcarrier bleeds into decoded luminance.

Notice the dot pattern on the bottom of the image.

Thus, different colors modulate the
luminance intensity differently, creating a
“dot” pattern on the scan line between two col-
ors. To eliminate these “hanging dots,” a
chroma trap filter is sometimes used after the
comb filter.

In the real world of analog circuitry, having a perfectly locked signal like this is probably also not too realistic to assume.

The literature also calls out for notch filtering as an approach.

I attempted combining a comb filter with a notch filter on top to reduce the artifacts, but it is quite tricky to create a notch filter that works well. A simple FIR notch filter with zeroes is easy enough to make:

% 2x oversampled, so 27 MHz.
ntsc_subcarrier_w = ntsc_subcarrier / (2 * bt601_rate);

% Hope you know your complex numbers.
subcarrier_zero = exp(2 * pi * j * ntsc_subcarrier_w);

% To get a real-valued signal
% we need zeroes at positive and negative freqs.
% use real() here just a precaution.
% It should not be needed with perfect FP math accuracy.
subcarrier_notch =
  real(conv([1, -subcarrier_zero], [1, -conj(subcarrier_zero)]));

% Normalize DC gain at 0 dB
subcarrier_notch = subcarrier_notch' / sum(subcarrier_notch);

freqz(subcarrier_notch);

This filter is convolved with a simple low-pass to complete the luma decoding filter.

This approach just leads to severe blurring for NTSC and a band-stop filter approach just lead to less severe blurring and ringing instead, so I’m not sure what should be done. Unfortunately, the handbook isn’t clear on what kind of filtering is called for here. IIR notch filter designs can be super sharp to surgically carve out the carrier, but IIR filtering is also a massive pain in the ass on GPUs. It’s also likely to ring heavily too, which I found rather annoying in my testing.

E.g. 3-line comb NTSC without the notch (integer nearest upscale from 640×240):

and with notch:

Yikes. There’s no way this notch approach is correct. It’s like we’re getting double vision here. It does clean up the chroma dots though, so … yay?

Going beyond these base techniques there’s adaptive filtering where the filtering strategy changes based on which kind of case we’re dealing with. And even more sophisticated is taking advantage of temporal information (look ma’, TAA has been a thing since forever :D) since N fields in the past we have complementary chroma phase perfectly aligned to our pixel grid.

Very cool stuff, but I doubt consumer TVs at the time would have those. The added latency for doing this kind of analysis doesn’t sound like something you’d want for games at least …

Either way, I’m not designing high-end TV circuitry in the late 90s/early 2000s here. We can just flip on S-Video to simulate the perfect Y/C separator, so at some point I have to decide that I’ve done enough filter masturbation and move on.

Step 7 – Separating chroma and luma in PAL

This was way harder than expected and I had to bang my head against the wall for a while to come up with a good solution. The ~90 degree shift every scanline means the basic comb filter for NTSC won’t work at all.

The handbook has two main strategies here. Either a delay line which is slightly longer than a scanline to align the phases:

Or use a highly magical “PAL modifier”:

The function of this modifier is esoteric as all hell, but I think the purpose of it is to phase shift the signal by 90 degrees to “realign” the carrier somehow (it will still be off by 0.6 degrees). This filter path with two bandpass filters just got so messy, and I couldn’t figure out how to debug the thing effectively (were the inevitable visual artifacts my bugs or just the filter being bad?) that I eventually gave up and designed my own filter. That’s more fun anyway.

The magic cross filter

I started from first principles and designed a 3×3 kernel that should be able to perfectly pass a chrominance signal and 100% reject any luminance signal that is DC in either horizontal or vertical direction.

To make things simple I started with the assumption of 4 samples per subcarrier cycle to make the examples easier. Given a constant value of 1.0 for U, a signal would look something like:

 0  1  0 -1  0  1  0 -1  0  1  0 -1 ...
-1  0  1  0 -1  0  1  0 -1  0  1  0 ...
 0 -1  0  1  0 -1  0  1  0 -1  0  1 ...
 1  0 -1  0  1  0 -1  0  1  0 -1  0 ...
 0  1  0 -1  0  1  0 -1  0  1  0 -1 ...
-1  0  1  0 -1  0  1  0 -1  0  1  0 ...
 0 -1  0  1  0 -1  0  1  0 -1  0  1 ...
 1  0 -1  0  1  0 -1  0  1  0 -1  0 ...

U = sin(wt) with N + 0.75 cycles per line here. A kernel that satisfies the criteria is:

|  0.25  -0.25   0.0  |
| -0.25   0.5   -0.25 |
|  0.0   -0.25   0.25 |

The sum of all rows and column is 0, meaning that if the signal is DC either horizontally or vertically, the result is completely filtered out. The filter also rejects V signals.

V = +/- cos(wt) and looks like

 1  0 -1  0  1  0 -1  0  1 ...
 0 -1  0  1  0 -1  0  1  0 ...
-1  0  1  0 -1  0  1  0 -1 ...
 0  1  0 -1  0  1  0 -1  0 ...
 1  0 -1  0  1  0 -1  0  1 ...
 0 -1  0  1  0 -1  0  1  0 ...
-1  0  1  0 -1  0  1  0 -1 ...
 0  1  0 -1  0  1  0 -1  0 ...

This is just the same signal flipped horizontally, so:

|  0.0  -0.25  0.25 |
| -0.25  0.5  -0.25 |
|  0.25 -0.25  0.0  |

Then a combined filter is made that accepts U and V signals together. U and V can be perfectly split later during demodulation so that is okay.

|  0.25 -0.50  0.25 |
| -0.50  1.00 -0.50 |
|  0.25 -0.50  0.25 |

This just boils down to a simple diagonal edge detection filter (high pass vertically and horizontally), but actually works quite well.

To deal with the actual 2x oversampled rate of 27 MHz and the PAL subcarrier at ~4.433 MHz, the ~90 degree shift per line is about 1.51 samples, so to make this sort of work, I stretched out the horizontal kernel to a 5-tap filter:

[ -0.2602, -0.2398, 1.0000, -0.2398, -0.2602 ]

The vertical kernel remains the same of

[ -0.5, 1.0, -0.5 ]

Some error is introduced since we’re not sampling the signal 100% correctly anymore (theoretically we need a sinc to reconstruct the signal perfectly), but I measured the reconstructed error to be -40 dB, which is good enough I think. The measured error for U and V were also similar, which indicates no weird artifacts from the V flips.

With this 15-tap kernel, we get a pretty good chroma estimate even in PAL. From here the same ideas as NTSC apply, bandpass the estimate and subtract it from composite signal to get luma.

Notch filtering to cleanup the chroma dots worked way, way better for PAL than NTSC, likely because the carrier has a much higher frequency on PAL, so the low pass behavior of the notch isn’t as devastating to image quality as it is on NTSC.

In the end I think I prefer 3-line comb + notch for PAL and just plain 3-line comb for NTSC. These screenshots are just one still frame (or rather, field). The color fringing will cancel out the next field and it’s hard to show the effect without seeing it at full 60 fields per second.

Bonus hack – PAL60

While PlayStation 2 didn’t support this hack mode, GameCube did back in the day. It’s a non-standard video mode that has same refresh rate as NTSC and vertical resolution, but retains the bandwidth and chroma encoding system of PAL, the best of both worlds! My implementation can trivially implement this by just enabling PAL on 60 Hz games. Only thing I’m not quite clear on is how the 1 / 625 subcarrier offset per line is supposed to work, but it’s a non-standard mode anyway, so eeeeeh.

Step 8 – Validate the color cards in different configurations

With comb filter and notch on NTSC:

As expected, luma detail is murdered around ~3.58 MHz carrier. Also serious color fringing due to the extreme high frequency diagonal patterns. In the pattern at the bottom, no fringing is observed since the comb filter did its job as expected.

PAL is similar, but the carrier and notch moves to ~4.43 MHz instead:

Step 9 – Validate the Hanover Bar filter for PAL

The main feature of PAL is being robust against phase shifts during analog broadcasts. It’s a little unclear if composite inputs cared about this case, but for completeness sake, I implemented it. This path is naturally skipped for S-Video and Component outputs since I can’t imagine a TV caring about that for the more luxurious inputs. It doesn’t take many degrees of error in the phase to get quite different colors for NTSC.

For PAL, the phase error manifests as complementary errors every line.

However, by averaging out chroma vertically, we can recover the original chroma almost perfectly. At worst, a little less saturation.

 

Writing a semi-competent CRT shader

This topic is kind of unfortunate in that it’s done to death already, and it’s a 100% subjective topic meaning that everyone has some kind of opinion, none which agree with each other. Holy wars have been fought over less.

Trying to write anything fresh about this topic is futile in 2026 – the heyday of CRT hobbyist shader development was in the early 2010s – but I felt the need to explain what I did at the very least. If anything, it’s a useful intro to writing your own shader.

https://nyanpasu64.gitlab.io/blog/crt-appearance-tv-monitor/ is a good read too for more background information.

Scanlines – and what not to do at 480i

The most obvious part of a CRT filter is scanlines, however, the idea that CRT images should have clearly visible scanlines is actually an artifact of 240p. For PS2 games, we’re operating at either 480i or 480p for NTSC.

For interlaced video, we expect each individual field to have clearly visible scanlines, but the complete image (fused together by our brains) should not. The beam profile should be tuned as such. What most shaders do is for each scanline to take a gaussian distribution in the Y direction, sampled for the neighbor lines to cover the useful portion of the kernel.

vec3 inv_stddev = ...;
vec3 inv_variance = inv_stddev * inv_stddev;
vec3 gaussian = exp(-0.5 * inv_variance * phase * phase);

// Scaling factors to ensure integral of the distribution remains 1.
const float one_over_sqrt_two_pi = 0.3989422;
gaussian *= one_over_sqrt_two_pi * inv_stddev;

Breathing effect

Another common effect is that very bright scanlines are smeared out, supposedly due to the electron guns not being as stable when they’re driven at high voltages. This can be simulated by varying the standard deviation. It can be subtle, but creates a neat effect I think. Exactly how to come up with the beam profile for a given input voltage is purely up to taste I suppose, I doubt there is a linear relationship between R’G’B’ value and standard deviation of the beam 🙂

// Breathing effect.
vec3 inv_stddev = mix(
  vec3(registers.scan_factor_narrow),
  vec3(registers.scan_factor_wide), sampled);

The sampled RGB value is in gamma-space still, since the CRT gamma curve is due to the phosphor response, not the CRT itself adjusting the gamma curve.

Applying gamma in the appropriate place

The BT.1886 standard calls for a 2.4 gamma for SD content, which is the default and looks good. I also added options to use 2.2 (NTSC legacy) and 2.8 (?!, PAL legacy) for fun. Most CRT shaders I’ve seen apply gamma in this way:

// A little unclear if we should do gamma before or after.
// Before makes a little more sense I think.
color += pow(sampled, vec3(registers.gamma)) * gaussian;

I think it would depend on whether or not the phosphor’s response is a function of how many electrons hit it, or individual “particles” respond to the energy of the electrons hitting them, where the gaussian beam profile is just a distribution of how many particles light up. In the latter interpretation, the code as-is makes sense, while the first interpretation would call for apply a gamma function on the gaussian profile. The visual output as-is looks good to me at any rate.

After this point, all color math happens in linear light, so floating point render targets is a must.

Applying a colored dots mask

Color CRTs get their colors by having colored phosphors that are arranged in some kind of grid. The venerable Trinitrons use vertical stripes of RGB, and I like that look. While CRTs don’t really have a horizontal resolution, there is a “dot pitch”, which sort of dictates resolution. This part is the key to create the “texture” of a CRT.

From what I read online a typical dot pitch for consumer TVs was 0.5 – 1.0 mm, and for a typical 20″ CRT, I estimated a reasonable number of RGB triads to be ~640 or so. Close enough to BT.601 standard horizontal resolution, neat. From what I understand, this value is also referred to as “TVL”, and these values seem ballpark reasonable. When looked at a distance from the screen, the dots blend together nicely as we’d expect. It’s basically just LCD subpixels, just larger.

The dot layout I used was mostly lifted from Lotte’s CRT shader, but the approach can probably be found in a million shaders already. Just alternating stripes of R, G and B.

vec3 grille(vec3 color, vec2 pos)
{
  vec3 mask = vec3(0.25);
  pos.x = fract(pos.x / 3.0);

  if (pos.x < 0.333)
    mask.r = 1.0;
  else if (pos.x < 0.666)
    mask.g = 1.0;
  else
    mask.b = 1.0;

  return color * mask;
}

I suppose a perfect mask of 0.0 should be used here, but it doesn’t end up looking as good as I’d like, even after adding glow effects, so I think the intent behind passing through a portion of other colors is more of a pragmatic decision. At least I cannot think of a physical interpretation of why we’d want to do it.

Avoiding aliasing

The signal we are creating has a ton of high frequency information and we need to be very careful sampling it such that obvious aliasing is avoided. The common mistake is to just render this shader at output resolution and hoping for the best. This will almost surely lead to terrible aliasing patterns in the image. Bad aliasing of the scanlines in Y direction leads to a low frequency pumping pattern in the image which is extremely distracting, and bad aliasing in the X direction leads to a horribly noisy pattern due to uneven sampling of the dot mask.

The easy fix is to render the effect at an integer multiple. E.g. if input image is 240 lines, render to a height that is an integer multiple of that. For the color dot mask, make sure the horizontal resolution is e.g. 3 times the dot resolution (one pixel for R, G, B dots). I ended up with something like width = 640 * 3, and height = 240 * 6 (3x sampling for progressive).

Glow / Bloom effect

A nebulous effect of the CRT is the glow aspect to it. Anyone can tell it’s there, but it’s not entirely clear to me why this happens. Google searches don’t turn up anything useful either. Without knowing the physical reason for it, it’s hard to emulate accurately. Could it be scattering effects inside the thick CRT glass perhaps?

Either way, the common way to emulate this effect is to compute a gaussian blur (lots of those around here) and composite it over the original image. Very similar to the usual HDR + bloom effect that games in the 2010s loved to overuse. The main effect here is that the phosphor dots end up blending together nicely, yet retain the added “texture” that the aperture grille pattern gives. Humans like to see some high frequency detail, even if that detail is completely bogus. That’s the common trick behind video compression after all.

The glow component, boosted up a ton to make it very visible:

With the typical HDR effect in HD games, only very bright pixels participate in the effect, effectively spreading excess light energy over a larger area of pixels. It makes little sense to do anything like this for a CRT shader unless there is a physical threshold where phosphors just randomly start to glow more than they should, but all of this is purely up to taste anyway.

Here’s from Soul Calibur II with progressive scan, component cable emulation and 16x SSAA, without any glow added. The look is very harsh to me. (The full-screen image is needed to see it without the extreme aliasing caused by thumbnailing.)

Some glow on top and it looks like:

Purely up to taste how much to add of course. I like a decent amount of glow.

Feedback / Phosphor persistence

Phosphors don’t turn off right away when they’re lit. It’s very quick though, but adding a few percent of feedback between frames seems to help a bit with making 480i games look better with less flicker. This is not really how things work in the real world I think, but a reasonable approximation.

Scaling the image to screen

We need a high quality rescaler to get the integer sampled CRT to the screen without introducing significant aliasing from the aperture grille or scanlines. The way to go here is simply to use a proper windowed sinc or something like that.

Curvature?

I don’t like it, so I don’t implement it. If you do, make sure to consider resampling it properly to not introduce more aliasing.

Using the correct color space

A point many shaders miss is that the RGB of a modern monitor is not the same as RGB on an old CRT TV. What we think of RGB today is usually BT.709 sRGB which defines a set of color primaries and white point. Old SDTV era video uses BT.601 which is a bit narrower than BT.709. In linear RGB space, this transform is a trivial 3×3 matrix multiply.

I actually learned that the very old NTSC 1953 standard defines a set of primaries that were extremely saturated compared to the standards of today. While the primaries of 1953 were aspirational, it was clearly way too early.

SMPTE refined NTSC to use more reasonable primaries as part of BT.601. PAL primaries are almost exactly the same (TVs tended to use the same phosphor formulation across the world I suppose), but there’s a theoretical difference so I added both for good measure.

Supposedly, Japan kept the use of legacy NTSC 1953 primaries, so that opens an interesting question if the same games actually looked vastly different in Japan compared to the rest of the world? I’ve never heard anything mention this before, so who knows. Either way, I support enabling those primaries for fun. The look of it is quite … something? It would need a HDR monitor with solid gamut to do justice.

Here’s with standard BT.601 primaries: (From Legaia 2, which is a purely field rendered game, hence the scanlines)

and with NTSC 1953 primaries. Almost like the “Interpret sRGB as Display P3” bullshit that phones do these days to “pop”.

Bonus round – HDR

Given all the masking we’re doing which lower average brightness, it’s beneficial to support HDR10 rendering. Now that HDR is widely available on Linux too, enjoying some HDR CRT shaders is a good time 😀

I added a few modes where I can target a specified maximum nit level. KDE at least respects MaxCLL HDR metadata and disables tonemapping if MaxCLL falls within bounds of the display. I also added a no-tonemap option where MaxCLL = 0 (unknown), which makes KDE tonemap how it wants to.

Experimental round – High refresh rate

Black Frame Insertion and their friends have been a thing in emulation for ages, and it works quite well with CRT simulation, especially to sell the de-interlacing effect properly. In my implementation, I query EXT_present_timing and decide how many frames I should insert in-between. There’s a gentler falloff between the frames. The overall screen brightness decreases a lot as expected, but with HDR, we can crank the brightness of the proper frame up to compensate.

It’s still very experimental and any missed frame leads to horrrrrrrible flicker at the moment (big epilepsy warning), so it’s not something I actually recommend, but it’s fun to experiment with.

Getting de-interlacing for “free”

While simulating each field independently with scanlines, we get de-interlacing the same way a CRT would in theory. This is not free from flicker of course, but most games had mitigation strategies for this.

The two CRTCs of the PS2 – FRAME mode

The PS2 supported blending two frames being sent to the video output. Most games render internally at e.g. 448p @ 30 FPS, but since they cannot output that resolution without component cables, the frame is output interlaced over two fields where the CRTC scans every other line rather than every line. That tends to look quite flickery if done as-is given how aliased PS2 graphics are, so what pretty much every game did was using the two CRTCs to blend the two frames vertically before sending it by programming a 1 pixel offset with duplicate inputs. By shifting the offset every field, the 30 FPS progressive image could be scanned out nicely into a 60 FPS interlaced image.

This is the “flicker filter” that some games allowed toggling. Here’s a RenderDoc capture showing that CRTC 1 and 2 are configured with one pixel offset in Y:

After merging and blending the frames together, a smoothed image is sent.

PCSX2 GSdx and paraLLEl-GS have modes to detect this pattern and just scan out the full 480p of course, without the added blurring. Most games fall into this pattern, which is fortunate if we want to avoid interlacing shenanigans, but not all games are so nice to deal with. The “Anti-Blur” option is designed precisely for disabling this filter.

I also added an option to force-disable the automatic progressive scan, mostly to test the video output that the games would actually have output back in the day, which is interlaced video.

FIELD mode

Some games decided that they wanted to render at 60 FPS and sacrifice half the vertical resolution to do so, jittering the rendering to stay in sync with the interlaced output. These are painful to deal with to this day since they absolutely require some kind of de-interlacing solution to look good.

I never got satisfactory results with a typical de-interlacer, but the CRT simulation does a quite good job at it I think. It’s not perfect (interlacing wasn’t exactly perfect on CRTs either), but it’s usable for me to the point now that I can play interlaced games as intended. It’s not really possible to demonstrate this with screenshots.

The really cursed shit

Some games break if you try to promote them to progressive scan, because the games might decide for some stupid reason to use the SCANMSK feature to discard pixels every other line, and rely on the FRAME scanout mode to exactly scan out the pixels that were not masked. Kings’s Field IV is an example of this absolute insanity.

PCSX2 integration

I maintain a patch for PCSX2-git which supports parallel-gs and now this CRT/analog thing. This is super niche stuff that I don’t really expect many people to actually use, but it’s there for those who are interested.

It does what I need at least, and that’s what’s important to me.

Putting it all together

I put together some test video clips in HEVC/PQ/4:4:4/1440p at ridiculously high bitrates. I tried AV1 but my CPU could not decode it in real time, so it is what it is ._.

Clip 1 is: Raw RGB passed into the CRT. It uses the game’s native 480p and widescreen support. Super sampling is 4x SSAA.

Clip 2 is: “PAL60” with 3-line comb + notch. It also uses the more default 4:3 and interlaced video. It has some frame drops which completely botch the interlace and even at comically large bitrates, the aperture grill effect doesn’t translate well, so it is what it is, but it’s a rough approximation of what it should look like.

Anyway, I’m happy with the results. Time to actually play something instead of debugging stuff 😀

A case for learning GPU programming with a compute-first mindset

Beginners coming into our little corner of the programming world have it rough. Normal CPU-centric programming tends to start out with a “Hello World” sample, which can be completed in mere minutes. It takes far longer to simply download the toolchains and set them up. If you’re on a developer friendly OS, this can be completed in seconds as well.

However, in the graphics world, young whippersnappers cut their teeth at rendering the elusive “Hello Triangle” to demonstrate that yes, we can indeed do what our forebears accomplished 40 years ago, except with 20x the effort and 100000x the performance.

There’s no shortage of examples of beginners rendering a simple triangle (or a cube), and with the new APIs having completely displaced the oxygen of older APIs, there is a certain expectation of ridiculous complexity and raw grit required to tickle some pixels on the display. 1000 lines of code, two weeks of grinding, debugging black screens etc, etc. Something is obviously wrong here, and it’s not going to get easier.

I would argue that trying to hammer through the brick wall of graphics is the wrong approach in 2025. Graphics itself is less and less relevant for any hopeful new GPU programmer. Notice I wrote “GPU programmer”, not graphics programmer, because most interesting work these days happens with compute shaders, not traditional “graphics” rendering.

Instead, I would argue we should start teaching compute with a debugger/profiler first mindset, building up the understanding of how GPUs execute code, and eventually introduce the fixed-function rasterization pipeline as a specialization once all the fundamentals are already in place. The raster pipeline was simple enough to teach 20 years ago, but those days are long gone, and unless you plan to hack on pre-historic games as a hobby project, it’s an extremely large API surface to learn.

When compute is the focus, there’s a lot of APIs we could ponder, like CUDA and OpenCL, but I think Vulkan compute is the best compute focused API to start out with. I’m totally not biased obviously 😀 The end goal is of course to also understand the graphics pipeline, and pure compute APIs will not help you there.

Goal for this blog

I don’t intend to write a big book here that has all the answers on how to become a competent GPU programmer. Instead, I want to try outlining some kind of “meta-tutorial” that could be fleshed out further. I’ve been writing compute shaders since the release of OpenGL 4.3 ages ago and I still learn new things.

To abstract or not to abstract

For this exercise, I will rely on a mid-level API abstraction like my own Granite. I don’t think throwing developers into the raw API is the best idea to begin with, but there must be some connection with the underlying API, i.e., no multi-API abstraction that barely resembles the actual API which you’ll typically find in large production engines. Granite is a pure Vulkan abstraction I’ve been chiseling away at for years for my own needs (I ship it in tons of side projects and stuff), but it’s not really an API I’ve intended others to actively use and ship software on. Migrating away from the training wheels quickly is important though and compute makes that fairly painless. Granite is just one of many approaches to tackling Vulkan, and the intent is not to present this as the one true way.

The debugging first approach

Getting something to show up on screen is important to keep the dopamine juices flowing. Fortunately, we can actually do this without messing around with graphics directly. With RenderDoc captures we get debugging + something visual in the same package, and learning tooling early is critical to be effective. Debugging sufficiently interesting GPU code is impossible without this.

Shading language

The debug flow I propose with RenderDoc will rely on a lot of shader replacements and roundtrips via SPIRV-Cross’ GLSL backend, so Vulkan GLSL is the appropriate language to start with. It’s more or less a dead language at this point, but it’s also the most documented language that has support for buffer device addresses, which I will introduce right away to avoid having the brick wall of descriptors and binding models. This is a very compute-centric move, but makes other parts of the API easier to grasp later.

HLSL from the Direct3D ecosystem is a popular option, but as a compute language, HLSL is weaker than Vulkan GLSL in my experience, due to lack of a lot of features that come up in more interesting compute workloads, but being bilingual in this area is unavoidable these days. No matter which language you use, someone will call you a filthy degenerate anyway. :v

Deferring synchronization

Being debugger-centric we can avoid poking at explicit synchronization for a very long time and once we get there, we can simplify a ton. You can do a lot of interesting things with a single dispatch after all.

Writing the first program

Here’s a very basic program that copies some data around. It should trivially build on Linux or Windows on the usual compilers. Make sure to clone or symlink Granite so that it can be picked up by the CMake build.

If you try to run this, the output might look something like this:

...
[INFO]: Enabling VK_LAYER_KHRONOS_validation.
[INFO]: Enabling instance extension: VK_EXT_debug_utils.
[INFO]: Found Vulkan GPU: AMD Radeon RX 6800 (RADV NAVI21)
[INFO]: API: 1.4.328
[INFO]: Driver: 25.2.99
[INFO]: Found Vulkan GPU: NVIDIA GeForce RTX 4070
[INFO]: API: 1.4.312
[INFO]: Driver: 580.328.576
[INFO]: Using Vulkan GPU: AMD Radeon RX 6800 (RADV NAVI21)
[INFO]: Enabling device extension: VK_KHR_external_semaphore_fd.
[INFO]: Enabling device extension: VK_KHR_external_memory_fd.
[INFO]: Enabling device extension: VK_EXT_external_memory_dma_buf.
[INFO]: Enabling device extension: VK_EXT_image_drm_format_modifier.
[INFO]: Enabling device extension: VK_KHR_calibrated_timestamps.
[INFO]: Enabling device extension: VK_EXT_conservative_rasterization.
[INFO]: Enabling device extension: VK_KHR_compute_shader_derivatives.
[INFO]: Enabling device extension: VK_KHR_performance_query.
[INFO]: Enabling device extension: VK_EXT_memory_priority.
[INFO]: Enabling device extension: VK_EXT_memory_budget.
[INFO]: Enabling device extension: VK_EXT_device_generated_commands.
[INFO]: Enabling device extension: VK_EXT_mesh_shader.
[INFO]: Enabling device extension: VK_EXT_external_memory_host.
[INFO]: Enabling device extension: VK_KHR_fragment_shader_barycentric.
[INFO]: Enabling device extension: VK_EXT_image_compression_control.
[INFO]: Graphics queue: family 0, index 0.
[INFO]: Compute queue: family 1, index 0.
[INFO]: Transfer queue: family 1, index 1.
[INFO]: Detected attached tool:
[INFO]: Name: Khronos Validation Layer
[INFO]: Description: Khronos Validation Layer
[INFO]: Version: 327
[INFO]: Detected tool which cares about debug markers.
[INFO]: Allocating 64.0 MiB on heap #1 (mode #3), before allocating budget: (0.0 MiB / 15096.3 MiB) [0.0 / 16368.0].
[INFO]: Allocating 64.0 MiB on heap #0 (mode #0), before allocating budget: (0.0 MiB / 64216.8 MiB) [0.0 / 64372.4].
[ERROR]: Failed to load RenderDoc, make sure RenderDoc started the application in capture mode.
[INFO]: Allocating 64.0 MiB on heap #1 (mode #1), before allocating budget: (64.6 MiB / 15096.3 MiB) [64.0 / 16368.0].

Capturing Vulkan code in RenderDoc

The code we just wrote executes on the GPU, but we have no easy way to observe the code actually running on the device. This is where RenderDoc comes in. Point to the executable we built.

After launching, the capture happens automatically, and when the process terminates, the capture should appear.

Clicking on the copy command and double-clicking the destination buffer, we can see the raw contents:

The zero-initialization flag we passed into buffer creation was technically not needed, but it helped make the capture a little easier to understand. That clear happened automagically inside Granite. Normally, memory is not assumed to be zero-cleared on allocation.

Running some actual shaders

Instead of using copies, we can create our own little memcpy. Here’s an updated sample gist. To keep things simple, we can use shaderc‘s method of compiling GLSL into a C header file.

glslc -o copy.h -mfmt=c --target-env=vulkan1.2 copy.comp

Vulkan 1.2 is used here since that introduced buffer device addresses in core. Building and capturing again gives us:

To inspect the push constants, it’s under uniform buffers:

RenderDoc understands how to resolve pointers in buffers into links that open the relevant buffer.

If you click on an event before the dispatch you’ll see the writes disappear.

Shader replacement workflow with SPIRV-Cross and debug prints

This workflow is extremely powerful for difficult debugging scenarios and I cannot do my job without this. It’s imperative to learn this early. Debug gives you a more traditional step debugger. Depending on the bug, this may be the correct tool (e.g. inspecting a particular broken pixel), but in my experience, when working with a ton of GPU threads in parallel, it’s often required to study it in aggregate to see what is going on, since you might not even know which thread is at fault to begin with.

For now, select Edit -> Decompile with SPIRV-Cross

The Vulkan API uses the SPIR-V intermediate representation and SPIRV-Cross converts this back to equivalent GLSL. Fortunately, this result looks very similar to our original shader. This is one of the main reasons I prefer working with GLSL since the translation back and forth to SPIR-V is the least lossy compared to alternatives.

E.g. we can hack in some debug prints, hit Apply and the dispatch will have messages attached to them. RenderDoc implements debugPrintfEXT by rewriting the SPIR-V to write the printf values back to host. The Vulkan driver itself does not understand how to printf. It will just ignore the SPIR-V command to printf.

Shader replacement like this is not just for debug prints, but you can modify the code and see the results without having to recompile the entire program and recapturing.

For fun, debug print “Hello World” instead, and you have your checkbox ticked off.

Inspecting ISA

If you’re on a driver that exposes it, you can study the machine code. For this purpose, I highly recommend AMD GPUs with RADV driver on Linux. The ISA is arguably the most straight forward to read. All Mesa drivers should give you ISA no matter the graphics card if you don’t have an AMD card lying around.

; s6 = WorkgroupID.x
; v0 = LocalInvocationID.x
; GlobalInvocationID.x =
;   (WorkgroupID.x << 4) + LocalInvocationID.x
v_lshl_add_u32 v0, s6, 4, v0

; v1 = int(v0) >> 31, i.e. a simple sign extend
; s4, s5: Holds dst GPU pointer
; This is a slower 64-bit address computation path which
; only shows up for raw device addresses
v_ashrrev_i32_e32 v1, 31, v0
v_lshlrev_b64 v[0:1], 2, v[0:1]
v_add_co_u32 v2, vcc_lo, s4, v0
v_add_co_ci_u32_e32 v3, vcc_lo, s5, v1, vcc_lo

; Load 32-bits from pointer
global_load_dword v2, v[2:3], off

; In parallel, compute dest address with 64-bit math.
v_add_co_u32 v4, vcc_lo, s2, v0
v_add_co_ci_u32_e32 v5, vcc_lo, s3, v1, vcc_lo

; Wait until the read request completes
s_waitcnt vmcnt(0)

; Store
global_store_dword v[4:5], v2, off
s_endpgm

Understanding how a compute dispatch is organized

Before trying to make sense of this, we need a mental model for how GPU compute executes on the device. This model was more or less introduced by CUDA in 2007 and has remained effectively unchanged since, neat!

At the highest level, in the CPU side, we dispatch a 3D grid of workgroups. In this sample, we just have a 1x1x1 cube or dispatches.

 cmd->dispatch(1, 1, 1);

For every workgroup, there is another 3D grid of invocations. Multiple invocations work together and are able to efficiently communicate with each other. Communicating across workgroups is possible in some situations, but requires some scissor juggling.

Why is there a two level hierarchy?

GPUs are extremely parallel machines. To get optimal performance, we have to map the very scalar-looking code to SIMT. The model employed for essentially all modern GPUs is that one lane in the vector units maps to one thread.

Inside the workgroup, the invocations are split up into subgroups. The mental model to understand the distinction is:

  • Workgroup -> runs concurrently on the same shader core
  • Subgroup -> runs in lock-step in a SIMT fashion

For a workgroup to be running well, the number of invocations in it should be an integer multiple of the subgroup size, otherwise there will be lanes doing nothing, and that’s no fun.

The subgroup sizes in the wild vary quite a lot, but there is an upper legal limit of 128. In practice, these are the values you can expect to find in the wild:

  • 4: Really old Mali Bifrost, old iPhones
  • 8: Intel Arc Alchemist
  • 16: Intel Arc Battlemage, Mali Valhall, Intel (runs slower if not Battlemage)
  • 32: AMD RDNA, NVIDIA, Intel upper limit (runs even slower)
  • 64: Adreno, AMD GCN + RDNA
  • 128: Adreno

Some vendors support multiple subgroup sizes. Usually we don’t have to care too much about this until you graduate to the more hardcore level of compute shader programming, but Vulkan gives you control to force subgroup sizes when need be.

For desktop use cases, catering to the range of 16 to 64 is reasonable. In the example we’ve been looking at, the workgroup size is just 16, so this is not optimal. On mobile GPUs, you might need to consider a wider range of hardware.

The rule of thumb (for desktop) is to use one of three constellations for local_size:

  • (64, 1, 1) for 1D
  • (8, 8, 1) for 2D
  • (4, 4, 4) for 3D

Integer multiples of these are fine too. This should make almost any GPU happy. The maximum limit is 1024 invocations, but I never recommend going that high unless there are very good reasons to.

AMD specifics

For the AMD case, v_ instructions are vector instructions, meaning while it looks like a simple register, there are multiple instances of it, one for every invocation in the subgroup.

s_ instructions are scalar. It runs once per subgroup and runs in parallel with the vector units. Taking advantage of subgroup uniform code can be very powerful.

One way to think of this is that v_ instructions look more like SIMD, except that the SIMD width is much larger than CPU instruction sets:

addps xmm0, xmm1

where scalar instructions look more like regular CPU instructions:

add eax, ebx

Vector load-store looks more like gather/scatter, and scalar load-store is more like normal CPU load-store.

Shader replacement by modifying SPIR-V assembly

In rare cases, it’s useful to do minor in-place modifications to SPIR-V itself. As a very ad-hoc sample, we can attempt to clean up the awkward 64-bit pointer math by using OpInBoundsAccessChain instead of OpAccessChain. Edit -> Decompile with spirv-dis

Now replace OpAccessChain with OpInBoundsAccessChain and apply. While keeping the shader tab open, go back to Pipeline State viewer and look at GPU ISA:

OpInBoundsAccessChain tells the compiler that we cannot index outside the array, which means negative indices and massively large indices are not allowed, and this allows the compiler to emit u64 + u32 addressing format. Sure looks much nicer now. This is way beyond what a beginner should care about, but the point is to demonstrate that we can replace raw SPIR-V too, and also demonstrates how you can inspect the SPIR-V of shaders easily.

Introducing descriptors

We can do a lot with simple buffer device addresses to process data, but there is a limit to how far we can get with that approach if the end goal is game rendering. There are things GPUs can do that raw pointers cannot:

  • Efficiently sample textures
  • Automatic format conversions
  • “Free” bounds checking

CPU ISAs do none of these. In this updated memcpy sample, I introduce two descriptor types, STORAGE_BUFFER and UNIFORM_TEXEL_BUFFER. Using descriptors like this is the “normal” way to use Vulkan, and should be preferred when feasible.

For pragmatic reasons, it’s easier to debug and validate descriptors compared to raw pointers. Raw pointers are also prone to GPU crashes, which are very painful and annoying to debug. Unlike CPU debugging, we cannot just capture a SIGSEGV in a debugger and call it a day.

Already, the resources show up in a more convenient way. Even though the shader didn’t specify 8-bit inputs anywhere, it just works. Typed formats have up to 4 components. This is enough to cover RGB and Alpha channels, which of course has its origins in graphics concepts. texelFetch cannot know if we have R8_UINT or R16G16B16A16_UINT for example, so we have to select the .x component in the shader.

s_mov_b32 s0, s3

; On RADV, descriptors live in a fixed 32-bit VA region:
; 0xffff8000'xxxxxxxx
; Since this is known by compiler, we only need to pass down 32-bits
; and synthesize the upper half with s_movk_i32 (0x8000 is sign-extended).
s_movk_i32 s3, 0x8000

; STORAGE_BUFFER and TEXEL_BUFFER are both 16 bytes.
; Both are loaded here in one go into scalar registers.
; Notice the s_load. These loads go into the constant cache,
; and is almost "free". The same load is shared for all threads
; in the subgroup.
s_load_dwordx8 s[8:15], s[2:3], null
v_lshl_add_u32 v0, s0, 4, v0

; Wait until the scalar load completes so we can use the descriptor.
s_waitcnt lgkmcnt(0)

; Typed load. The descriptor holds information about e.g. R8_UINT.
; Bounds checking is also free.
buffer_load_format_x v1, v0, s[12:15], 0 idxen
v_lshlrev_b32_e32 v0, 2, v0
s_waitcnt vmcnt(0)

; Store, but uses a descriptor instead.
; This is automatically bounds checked and is free on AMD.
buffer_store_dword v1, v0, s[8:11], 0 offen
s_endpgm

The basic Vulkan binding model

Granite implements a rather old school binding model, but I think this model is overall very easy to understand and use. More modern bindless design is introduced later when it becomes relevant.

In the shader, I declare things like layout(set = 0, binding = 1). In Granite, this simply means that we have to bind a resource of the appropriate type before dispatching, e.g.:

 cmd->set_buffer_view(/* set */ 0, /* binding */ 1, *buffer_view);

This papers over a ton of concepts, and makes it very easy and convenient to program. In reality, there are a lot of API objects in play here. When the compiler sees layout(set = 0, binding = 1) for example, it needs to check the provided VkPipelineLayout given to pipeline creation. A set denotes a group of resources that are bound together as one contiguous entity. In the ISA, set = 0 is determined to be initialized in a certain scalar register:

s_load_dwordx8 s[8:15], s[2:3], null

The VkPipelineLayout also contains information about what e.g. binding = 1 means. In this case, the driver happened to decide that binding = 0 is at offset 0, and binding = 1 is at offset 16. Since these descriptors are adjacent in memory we got a lucky optimization where we load 32 bytes at once.

On the API side, we need a compatible VkPipelineLayout object when recording the command buffer to ensure that everything lines up. Granite does this automatically, through shader reflection, which synthesizes a working layout for us.

Based on the contained VkDescriptorSetLayout inside the pipeline layout, it knows how to allocate a VkDescriptorSet from a VkDescriptorPool and write descriptors to it. Then it can bind the descriptor set to the command buffer before dispatching. We can see all of this in effect in the capture. Turn off the Filter and we get:

The descriptor set is updated, then later bound. In reality vkCmdBindDescriptorSets is just a 32-bit push constant, which the shader ends up reading in s3 register.

Deeper understanding with VK_EXT_descriptor_buffer

Managing descriptors is always a point of contention in Vulkan programming if you’re writing the raw API. There’s a ton of concepts to juggle and it’s mostly pretty dull stuff.

As an extension to the original old school model I outlined about, it’s possible to treat a descriptor set as raw memory which gets rid of a ton of jank. Granite supports using this model by opting in to it. Change the sample to use and recapture:

 if (!ctx.init_instance_and_device(
   nullptr, 0, nullptr, 0,
   Vulkan::CONTEXT_CREATION_ENABLE_ROBUSTNESS_2_BIT |
   Vulkan::CONTEXT_CREATION_ENABLE_DESCRIPTOR_BUFFER_BIT))
   return EXIT_FAILURE;

Make sure to make a release build of the test and not a debug build, otherwise descriptor buffers are disabled. Note that this requires a recent build of RenderDoc. The latest stable v1.40 release supports descriptor buffers.

Now we explicitly tell the driver that the descriptor set lives at offset 0 from the bound buffer. If we then inspect the bound descriptor buffer …

Now we can see the raw guts of the storage buffer and texel buffer being encoded. You can even see the 0x40 and 0x10 being encoded there which corresponds to the sizes of the descriptors.

Porting ShaderToy shaders to compute

To get something interesting on screen to end this bringup exercise, we could port some shadertoy shaders. These are super convenient since many of them don’t require anything fancy to run like external textures or anything. I picked some shadertoy arbitrarily. Store this to mandelbulb.glsl, and then we replace our shader with a mandelbulb.comp that calls the shadertoy code:

#version 460

// We're writing 2D images now, so this makes more sense.
layout(local_size_x = 8, local_size_y = 8) in;

layout(set = 0, binding = 0) writeonly uniform image2D output_image;

// Constants used by shadertoys.
layout(push_constant) uniform Registers
{
  vec2 iResolution;
  float iTime;
};

#include "mandelbulb.glsl"

void main()
{
  vec4 color;
  mainImage(color, vec2(gl_GlobalInvocationID.xy) + 0.5);

  // Stores the result to texture.
  imageStore(output_image, ivec2(gl_GlobalInvocationID.xy), color);
}

On the API side, we simply need to create a storage texture and bind it to the shader.

Just with this simple setup, you can go completely nuts and play around with the more math-heavy side of graphics if you want.

Where to go from here?

From here, I think the natural evolution is to learn about:

  • Atomics
    • Living in a world without mutexes: lockless programming with millions of threads
  • Shared memory
  • Subgroup operations
    • Case study: Scalarization
    • Case study: Bindless and non-uniform indexing of descriptors
  • Texture sampling and mip-mapping
    • The bread and butter of graphics
    • Case study: do some simple image processing with simple filters
  • Memory coherency and how to communicate with other workgroups
    • Case study: Single pass down-sampling
  • If relevant, start porting over the code to more shading languages
  • API synchronization and how to keep CPU and GPU pipelined
  • … and maybe only then start looking at getting some images on screen (with compute)
  • Bring up your own Vulkan code from scratch to get rid of the training wheels and make sure you understand how everything comes together

After that, it’s a matter of learning the common algorithms that show up all the time, like parallel scans, classification, binning, etc, etc. This naturally leads to indirect dispatches, and once these concepts are in place, we can design a very simple compute shader rasterizer that renders some simple glTF models. Only when those concepts land do we consider the graphics pipeline.

I designed my own ridiculously fast game streaming video codec – PyroWave

Streaming gameplay from one machine to another over a network is a reasonably popular use case these days. These use cases demand very, very low latency. Every millisecond counts here. We need to:

  • Send controller input from machine A to B over network
  • B renders a frame on the GPU
  • B encodes the frame into a bitstream
  • B sends the result over a network to A
  • A decodes the bitstream
  • A displays image on screen
  • Dopamine is released in target brain

Every step in this chain adds latency and we want to minimize this as much as possible. The go-to solution here is GPU accelerated video compression using whatever codec you have available, usually H.264, HEVC or if you’re really fancy, AV1. Ideally, we want all of this to complete in roughly ~20 ms.

To make this use case work well, we have to strangle the codecs quite a bit. Modern video codecs love latency since it allows tricks like flexible rate control and B-frames. That is what allows these codecs to operate at ridiculous compression ratios, but we have to throw away most of these tricks now. Since the codec cannot add latency, and we’re working on a fixed bit-rate budget (the ethernet cable or WiFi antenna), we’re left with:

  • Hard-capped constant bit rate
    • There is no buffer that can soak up variable bit rate
  • Infinite GOP P-frames or intra-refresh
    • Either choice deals with packet loss differently

When game streaming, the expectation is that we have a lot of bandwidth available. Streaming locally on a LAN in particular, bandwidth is basically free. Gigabit ethernet is ancient technology and hundreds of megabits over WiFi is no problem either. This shifts priorities a little bit for me at least.

Back in my student days, I designed a very simple low-complexity video codec for my master thesis, and after fiddling with Vulkan video and PyroFling for a while, that old itch was scratched again. I wanted to see what would happen if I designed a codec with laser focus on local streaming with the absolute lowest possible latency, what could go wrong?

Throwing out motion prediction – intra-only

This is the grug-brained approach to video, but it’s not as silly as it sounds. Bit-rates explode of course, but we gain:

  • Excellent error resilience
    • Even on local WiFi streaming to a handheld device, this does matter quite a lot, at least in my experience
  • Simplicity
    • Duh
  • Consistent quality
    • With CBR, video quality is heavily dependent on how good of a job motion estimation can do

Intra-only has use cases in digital cinema (motion JPEG2000) and more professionally oriented applications where these concerns are likely more important than squeezing bandwidth. We’re now working at 100+ Mbits/s instead of ~10-20 Mbit/s, so streaming over the internet is no longer feasible outside of peer-to-peer with fiber links. For reference, raw 1080p60 with 420 chroma subsampling is in the range of 1.5 Gbit/s, and it only gets worse from there.

Throwing out entropy coding

Entropy coding is an absolute nightmare for parallelization, which means encoding solely on the GPU with compute shaders becomes an extremely painful affair. Let’s just throw that out and see how far we get. Gotta go fast!

There are codecs in this domain too, but it’s getting very specialized at this point. In the professional broadcasting space, there are codecs designed to squeeze more video through existing infrastructure with “zero” lag and minimal hardware cost. My master thesis was about this, for example. A more consumer oriented example is VESA display stream compression (I’m not sure if it does entropy coding, but the compression ratios are small enough I doubt it). There isn’t much readily available software in this domain, it’s generally all implemented in tiny ASICs. If FFmpeg doesn’t support it, it doesn’t exist for mere mortals.

Discrete Wavelet Transforms

While modern codecs are all block-based Discrete Cosine Transform (DCT) + a million hacks on top to make it good, there is an alternative that tried its best in the 90s to replace DCTs, but kinda failed in the end. https://www.youtube.com/watch?v=UGXeRx0Tic4 is a nice video explaining some of the lore. DWT-based compression has a niche today, but only in intra video compression. It’s a shame, because it’s quite elegant. https://en.wikipedia.org/wiki/Discrete_wavelet_transform

A graphics programmer will be familiar with this structure immediately, because this is just good old mip-maps with spice. Effectively, we downsample images, and also compute the “error” between the high-res picture and low-res picture. With signal processing lenses on, we can say it’s a critically sampled filter-bank. After processing N pixels, we obtain N / 2 low-pass and N / 2 high-pass pixels. The filters designed to do this are very particular (I really don’t know or care how they were made), but it’s basically just a basic convolution kernel, nothing too wild. The number of levels can vary but I chose 5 levels of decomposition.

Once the image is filtered into different bands, the values are quantized. Quantizing wavelets is a little tricky since we need to consider that during reconstruction, the filters have different gains. For the CDF 9/7 filter, high-pass is attenuated by 6 dB, and there are other effects when upsampling the lower resolution bands (zero-insertion). Rather than sweating out new graphs, I’ll just copy paste from my thesis. CDF 9/7 has very similar looking spectrum to the 5/3 I used here.

After normalizing the noise power, higher frequency bands can be quantized much harder than low-frequency bands. This exploits human psychovisual effects. This effect is used during rate control, which is another interesting problem. In the end, the higher frequency bands quantize to zero for most of the frame, with bits being allocated to critical regions of the image.

Classic artifacts of wavelets

The JPEG blocking artifact is infamous. Wavelet’s typical failure mode is that all high-pass information is quantized to 0, even where it shouldn’t be. This leads to a blurring – and if severe – ringing artifact. Given how blurry games these days can be with TAA, maybe this simply isn’t all that noticeable? 😀 Modern problems require modern solutions.

Packing coefficients into blocks really fast

Fiddling with this part of the codec was the thing that took the longest, but I think I landed on something alright eventually.

The basic block is 32×32 coefficients. This forms a standalone unit of the bitstream that can be decoded independently. If there is packet loss, we can error correct by simply assuming all coefficients are zero. This leads to a tiny blur somewhere random in the frame which is likely not even going to be perceptible.

The 32×32 block is further broken down into 8×8 blocks, which are then broken down into 4×2 blocks. This design is optimized for GPUs hierarchy of threads:

  • 1 thread: 4×2 coefficients
  • Cluster of threads (subgroup): 8×8 block
  • Workgroup: 32×32 block (128 threads)

8 coefficients per thread is deliberately chosen so that we can be byte oriented. Vulkan widely supports 8-bit storage of SSBOs, so I rely on that. We absolutely cannot be in a situation where we do bit fiddling on memory. That makes GPUs sad.

Like most wavelet codecs, I went with bit-plane encoding, but rather than employing a highly complicated (and terribly slow) entropy coding scheme, the bit-planes are just emitted raw as-is. I did this in my master thesis project and I found it surprisingly effective. The number of bits per coefficient are signaled at a 4×2 block level. I did some experiments on these block sizes and 4×2 was the right tradeoff. Using subgroup operations and some light prefix sums across the workgroup, it’s very efficient to decode and encode data this way. For non-zero coefficients, sign bits are tightly packed to the end of the 32×32 block. This was mildly tricky to implement, but not too bad.

The details are in my draft of the bitstream.

Accurate and fast rate control

In this style of compression, rate control is extremely important. We have a fixed (but huge) budget we have to meet. Most video codecs struggle with this requirement since the number of bits we get out of entropy coding is not easily knowable before we have actually done the compression. There is usually a lot of slack available to codecs when operating under normal VBR constraints. If a frame overshoots by 30%, we can amortize that over a few frames no problem, but that slack does not exist here since we’re assuming zero buffering latency.

Without entropy coding, we can trivialize the problem. For every 32×32 block, I test what happens if I throw away 1, 2, 3, … bits. I measure psychovisually weighted error (MSE) and bit-cost from that and store it in a buffer for later.

During RD analysis I can loosely sort the decisions to throw away bits by order of least distortion per bit saved. After the required number of bits have been saved through some prefix summing, we have achieved a roughly optimal rate distortion across the entire image in fixed time.

In the final pass, every 32×32 block checks how many bits to throw away and packs out the final bit-stream to a buffer. The result is guaranteed to be within the rate limit we set, usually ~10-20 bytes under the target.

Being able to rate limit like this is a common strength of wavelet codecs. Most of them end up iterating from most significant to least significant bit-plane and can stop encoding when rate limit is met, which is pretty cool, but also horribly slow …

Gotta go ridiculously fast 😀

So … is it fast? I think so. Here’s a 1080p 4:2:0 encode and decode of Expedition 33, which I found to be on the “brutal” end for image compression. Lots of green foliage and a lot of TAA noise is quite hard to encode.

0.13 ms on a RX 9070 XT on RADV. Decoding is also very fast. Under 100 microseconds. I don’t think anything else even comes close. The DWT pass was quite heavily optimized. It’s one of the few times where I found that packed FP16 math actually helped a lot, even on a beast like RDNA4. The quantizer pass does the most work of all the passes, and it took some effort to optimize too. Doing DWT in FP16 does have a knock-on effect on the maximum quality metrics we can achieve though.

Encoding more “normal” games, the quant + analysis pass has an easier time. 80 microsecond encode is pretty good.

Here’s a 4K 4:2:0 encode of the infamous ParkJoy scene.

0.25 ms, showing that 1080p struggles a bit to keep the GPU fully occupied. An interesting data point is that transferring a 4K RGBA8 image over the PCI-e bus is far slower than compressing it on the GPU like this, followed by copying over the compressed payload. Maybe there is a really cursed use case here …

I think this is an order of magnitude faster than even dedicated hardware codecs on GPUs. This performance improvement translates directly to lower latency (less time to encode and less time to decode), so maybe I’m onto something here.

Power consumption on the Steam Deck when decoding is also barely measurable. Not sure it’s less than the hardware video decoder, but I actually wouldn’t be surprised if it were.

Quality comparisons

Given how niche and esoteric this codec is, it’s hard to find any actual competing codecs to compare against. Given the domain is game streaming the only real alternative is to test against the GPU vendors’ encoders with H.264/HEVC/AV1 codecs in FFmpeg. NVENC is the obvious one to test here. VAAPI is also an option but at least FFmpeg’s implementation of VAAPI fails to meet CBR targets which is cheating and broken for this particular use case. It’s possible I held it wrong, but I’m not going to try debugging that.

Over 200mbit/s at 60 fps, I find it hard to tell any compression artifacts without side-by-side comparisons + zooming in, which is about 1.5 bpp. For something as simple as this codec, that’s quite neat. For objective metrics, I made use of https://github.com/psy-ex/metrics.

To even begin to compare a trivial codec like this against these codecs is a little silly, but we can level the playing field a bit by putting these codecs under the same harsh restriction that we have in PyroWave:

  • Intra-only
  • CBR with hard cap
    • Encoder is not allowed any slack, which should make the rate control really sweat
  • Fastest modes (not sure it matters that much to intra-only)

Example command line here:

ffmpeg -y -i input.y4m -b:v 100000k -c:v hevc_nvenc -preset p1 -tune ull -g 1 -rc cbr -bufsize 1667k out.mkv

No one in their right mind would stream like this, but let’s try anyway. 🙂

The video clips from the games are 5s clips I captured myself to raw NV12 video. I don’t think it’s super useful to upload those. My hacked up scripts to generate the graphs can be found here for reference.

I ran the NVENC tests on an RTX 4070 on 575 drivers.

ParkJoy

I included this as a baseline since this sequence has been seared into the mind of every video engineer out there (I think?) … This clip is 50 fps, but since my test script is hard-coded for 1080p60, I hex-edited the .y4m :’)

O_O. One thing to note here is that AV1/HEVC rate control kinda fails in this scenario. It ends up using less than the allotted budget, probably because it has to be conservative to meet the ridiculous hard-capped CBR. The graphs are done using the final encoded size however.

… VMAF, are you drunk?

Back to reality.

The more typical metrics look more like what I would expect. XPSNR is supposed to be a weighted PSNR that takes psychovisual effects into account, but I have no idea if it’s a good objective metric.

Expedition 33

Quite hard to encode. It was this game that actually gave me the last push to make this codec since even at 50 mbit/s with motion estimation, I recall some sections giving the encoder real trouble. Bumping bit rates just never cleaned things up properly for whatever reason …

I don’t know why, but VMAF really likes PyroWave.

Stellar Blade

Surprisingly easy to code. Must be the blurred backgrounds in play.

The use of FP16 kinda limits how high the PSNR can go. This is way beyond transparency, so, whatever.

Street Fighter 6

This scene is arguably a good argument for 4:4:4 …

FF VII Rebirth

More foliage, which I expected to be kinda hard. The game’s presentation is very soft and it shows in the compression rates 🙂

VMAF really seems like a joke metric for these use cases …

Another example of PSNR flattening off due to lowered internal precision.

Conclusion

I’m quite happy with this as-is, and having a 100% DIY streaming solution for my own use is fun.

PlayStation 2 GS emulation – the final frontier of Vulkan compute emulation

As you may, or may not know, I wrote paraLLEl-RDP back in 2020. It aimed at implementing the N64 RDP in Vulkan compute. Lightning fast, and extremely accurate, plus the added support of up-scaling on top. I’m quite happy how it turned out. Of course, the extreme accuracy was due to Angrylion being used as reference and I could aim for bit-exactness against that implementation.

Since then, there’s been the lingering idea of doing the same thing, but for PlayStation 2. Until now, there’s really only been one implementation in town, GSdx, which has remained the state-of-the-art for 20 years.

paraLLEl-GS is actually not the first compute implementation of the PS2 GS. An attempt was made back in 2014 for OpenCL as far as I recall, but it was never completed. At the very least, I cannot find it in the current upstream repo anymore.

The argument for doing compute shader raster on PS2 is certainly weaker than on N64. Angrylion was – and is – extremely slow, and N64 is extremely sensitive to accuracy where hardware acceleration with graphics APIs is impossible without serious compromises. PCSX2 on the other hand has a well-optimized software renderer, and a pretty solid graphics-based renderer, but that doesn’t mean there aren’t issues. The software renderer does not support up-scaling for example, and there are a myriad bugs and glitches with the graphics-based renderer, especially with up-scaling. As we’ll see, the PS2 GS is quite the nightmare to emulate in its own way.

My main motivation here is basically “because I can”. I already had a project lying around that did “generic” compute shader rasterization. I figured that maybe we could retro-fit this to support PS2 rendering.

I didn’t work on this project alone. My colleague, Runar Heyer, helped out a great deal in the beginning to get this started, doing all the leg-work to study the PS2 from various resources, doing the initial prototype implementation and fleshing out the Vulkan GLSL to emulate PS2 shading. Eventually, we hit some serious roadblocks in debugging various games, and the project was put on ice for a while since I was too drained dealing with horrible D3D12 game debugging day in and day out. The last months haven’t been a constant fire fight, so I’ve finally had the mental energy to finish it.

My understanding of the GS is mostly based on what Runar figured out, and what I’ve seen by debugging games. The GSdx software renderer does not seem like it’s hardware bit-accurate, so we were constantly second-guessing things when trying to compare output. This caused a major problem when we had the idea of writing detailed tests that directly compared against GSdx software renderer, and the test-driven approach fell flat very quickly. As a result, paraLLEl-GS isn’t really aiming for bit-accuracy against hardware, but it tries hard to avoid obvious accuracy issues at the very least.

Basic GS overview

Again, this is based on my understanding, and it might not be correct. 😀

GS is a pixel-pushing monster

The GS is infamous for its insane fill-rate and bandwidth. It could push over a billion pixels per second (in theory at least) back in 2000 which was nuts. While the VRAM is quite small (4 MiB), it was designed to be continuously streamed into using the various DMA engines.

Given the extreme fill-rate requirements, we have to design our renderer accordingly.

GS pixel pipeline is very basic, but quirky

In many ways, the GS is actually simpler than N64 RDP. Single texture, and a single cycle combiner, where N64 had a two stage combiner + two stage blender. Whatever AA support is there is extremely basic as well, where N64 is delightfully esoteric. The parts of the pixel pipeline that is painful to implement with traditional graphics APIs is:

Blending goes beyond 1.0

Inherited from PS1, 0x80 is treated as 1.0, and it can go all the way up to 0xff (almost 2). Shifting by 7 is easier than dividing by 255 I suppose. I’ve seen some extremely ugly workarounds in PCSX2 before to try working around this since UNORM formats cannot support this as is. Textures are similar, where alpha > 1.0 is representable.

There is also wrapping logic that can be used for when colors or alpha goes above 0xFF.

Destination alpha testing

The destination alpha can be used as a pseudo-stencil of sorts, and this is extremely painful without programmable blending. I suspect this was added as PS1 compatibility, since PS1 also had this strange feature.

Conditional blending

Based on the alpha, it’s possible to conditionally disable blending. Quite awkward without programmable blending … This is another PS1 compat feature. With PS1, it can be emulated by rendering every primitive twice with state changes in-between, but this quickly gets impractical with PS2.

Alpha correction

Before alpha is written out, it’s possible to OR in the MSB. Essentially forcing alpha to 1. It is not equivalent to alphaToOne however, since it’s a bit-wise OR of the MSB.

Alpha test can partially discard

A fun thing alpha tests can do is to partially discard. E.g. you can discard just color, but keep the depth write. Quite nutty.

AA1 – coverage-to-alpha – can control depth write per pixel

This is also kinda awkward. The only anti-alias PS2 has is AA1 which is a coverage-to-alpha feature. Supposedly, less than 100% coverage should disable depth writes (and blending is enabled), but the GSdx software renderer behavior here is extremely puzzling. I don’t really understand it yet.

32-bit fixed-point Z

I’ve still yet to see any games actually using this, but technically, it has D32_UINT support. Fun! From what I could grasp, GSdx software renderer implements this with FP64 (one of the many reasons I refuse to believe GSdx is bit-accurate), but FP64 is completely impractical on GPUs. When I have to, I’ll implement this with fixed-point math. 24-bit Z and 16-bit should be fine with FP32 interpolation I think.

Pray you have programmable blending

If you’re on a pure TBDR GPU most of this is quite doable, but immediate mode desktop GPUs quickly degenerates into ROV or per-pixel barriers after every primitive to emulate programmable blending, both which are horrifying for performance. Of course, with compute we can make our own TBDR to bypass all this. 🙂

D3D9-style raster rules

Primitives are fortunately provided in a plain form in clip-space. No awkward N64 edge equations here. The VU1 unit is supposed to do transforms and clipping, and emit various per-vertex attributes:

X/Y: 12.4 unsigned fixed-point
Z: 24-bit or 32-bit uint
FOG: 8-bit uint
RGBA: 8-bit, for per-vertex lighting
STQ: For perspective correct texturing with normalized coordinates. Q = 1 / w, S = s * Q, T = t * Q. Apparently the lower 8-bits of the mantissa are clipped away, so bfloat24? Q can be negative, which is always fun. No idea how this interacts with Inf and NaN …
UV: For non-perspective correct texturing. 12.4 fixed-point un-normalized.

  • Triangles are top-left raster, just like modern GPUs.
  • Pixel center is on integer coordinate, just like D3D9. (This is a common design mistake that D3D10+ and GL/Vulkan avoids).
  • Lines use Bresenham’s algorithm, which is not really feasible to upscale, so we have to fudge it with rect or parallelogram.
  • Points snap to nearest pixel. Unsure which rounding is used though … There is no interpolation ala gl_PointCoord.
  • Sprites are simple quads with two coordinates. STQ or UV can be interpolated and it seems to assume non-rotated coordinates. To support rotation, you’d need 3 coordinates to disambiguate.

All of this can be implemented fairly easily in normal graphics APIs, as long as we don’t consider upscaling. We have to rely on implementation details in GL and Vulkan, since these APIs don’t technically guarantee top-left raster rules.

Since X/Y is unsigned, there is an XY offset that can be applied to center the viewport where you want. This means the effective range of X/Y is +/- 4k pixels, a healthy guard band for 640×448 resolutions.

Vertex queue

The GS feels very much like old school OpenGL 1.0 with glVertex3f and friends. It even supports TRIANGLE_FAN! Amazing … RGBA, STQ and various registers are set, and every XYZ register write forms a vertex “kick” which latches vertex state and advances the queue. An XYZ register write may also be a drawing kick, which draws a primitive if the vertex queue is sufficiently filled. The vertex queue is managed differently depending on the topology. The semantics here seem to be pretty straight forward where strip primitives shift the queue by one, and list primitives clear the queue. Triangle fans keep the first element in the queue.

Fun swizzling formats

A clever idea is that while rendering to 24-bit color or 24-bit depth, there is 8 bits left unused in the MSBs. You can place textures there, because why not. 8H, 4HL, 4HH formats support 8-bit and 4-bit palettes nicely.

Pixel coordinates on PS2 are arranged into “pages”, which are 8 KiB, then subdivided into 32 blocks, and then, the smaller blocks are swizzled into a layout that fits well with a DDA-style renderer. E.g. for 32-bit RGBA, a page is 64×32 pixels, and 32 8×8 blocks are Z-order swizzled into that page.

Framebuffer cache and texture cache

There is a dedicated cache for framebuffer rendering and textures, one page’s worth. Games often abuse this to perform feedback loops, where they render on top of the pixels being sampled from. This is the root cause of extreme pain. N64 avoided this problem by having explicit copies into TMEM (and not really having the bandwidth to do elaborate feedback effects), and other consoles rendered to embedded SRAM (ala a tiler GPU), so these feedbacks aren’t as painful, but the GS is complete YOLO. Dealing with this gracefully is probably the biggest challenge. Combined with the PS2 being a bandwidth monster, developers knew how to take advantage of copious blending and blurring passes …

Texturing

Texturing on the GS is both very familar, and arcane.

On the plus side, the texel-center is at half-pixel, just like modern APIs. It seems like it has 4-bit sub-texel precision instead of 8 however. This is easily solved with some rounding. It also seems to have floor-rounding instead of nearest-rounding for bi-linear.

The bi-linear filter is a normal bi-linear. No weird 3-point N64 filter here.

On the weirder side, there are two special addressing modes.

REGION_CLAMP supports an arbitrary clamp inside a texture atlas (wouldn’t this be nice in normal graphics APIs? :D). It also works with REPEAT, so you can have REPEAT semantics on border, but then clamp slightly into the next “wrap”. This is trivial to emulate.

REGION_REPEAT is … worse. Here we can have custom bit-wise computation per coordinate. So something like u’ = (u & MASK) | FIX. This is done per-coordinate in bi-linear filtering, which is … painful, but solvable. This is another weird PS1 feature that was likely inherited for compatibility. At least on PS1, there was no bi-linear filtering to complicate things 🙂

Mip-mapping is also somewhat curious. Rather than relying on derivatives, the log2 of interpolated Q factor, along with some scaling factors are used to compute the LOD. This is quite clever, but I haven’t really seen any games use it. The down-side is that triangle-setup becomes rather complex if you want to account for correct tri-linear filtering, and it cannot support e.g. anisotropic filtering, but this is 2000, who cares! Not relying on derivatives is a huge boon for the compute implementation.

Formats are always “normalized” to RGBA8_UNORM. 5551 format is expanded to 8888 without bit-replication. There is no RGBA4444 format.

It’s quite feasible to implement the texturing with plain bindless.

CLUT

This is a 1 KiB cache that holds the current palette. There is an explicit copy step from VRAM into that CLUT cache before it can be used. Why hello there, N64 TMEM!

The CLUT is organized such that it can hold one full 256 color palette in 32-bit colors. On the other end, it can hold 32 palettes of 16 colors at 16 bpp.

TEXFLUSH

There is an explicit command that functions like a “sync and invalidate texture cache”. In the beginning I was hoping to rely on this to guide the hazard tracking, but oh how naive I was. In the end, I simply had to ignore TEXFLUSH. Basically, there are two styles of caching we could take with GS.

With “maximal” caching, we can assume that frame buffer caches and texture caches are infinitely large. The only way a hazard needs to be considered is after an explicit flush. This … breaks hard. Either games forget to use TEXFLUSH (because it happened to work on real hardware), or they TEXFLUSH way too much.

With “minimal” caching, we assume there is no caching and hazards are tracked directly. Some edge case handling is considered for feedback loops.

I went with “minimal”, and I believe GSdx did too.

Poking registers with style – GIF

The way to interact with the GS hardware is through the GIF, which is basically a unit that reads data and pokes the correct hardware registers. At the start of a GIF packet, there is a header which configures which registers should be written to, and how many “loops” there are. This maps very well to mesh rendering. We can consider something like one “loop” being:

  • Write RGBA vertex color
  • Write texture coordinate
  • Write position with draw kick

And if we have 300 vertices to render, we’d use 300 loops. State registers can be poked through the Address + Data pair, which just encodes target register + 64-bit payload. It’s possible to render this way too of course, but it’s just inefficient.

Textures are uploaded through the same mechanism. Various state registers are written to set up transfer destinations, formats, etc, and a special register is nudged to transfer 64-bit at a time to VRAM.

Hello Trongle – GS

If you missed the brain-dead simplicity of OpenGL 1.0, this is the API for you! 😀

For testing purposes, I added a tool to generate a .gs dump format that PCSX2 can consume. This is handy for comparing implementation behavior.

First, we program the frame buffer and scissor:

TESTBits test = {};
test.ZTE = TESTBits::ZTE_ENABLED;
test.ZTST = TESTBits::ZTST_GREATER; // Inverse Z, LESS is not supported.
iface.write_register(RegisterAddr::TEST_1, test);

FRAMEBits frame = {};
frame.FBP = 0x0 / PAGE_ALIGNMENT_BYTES;
frame.PSM = PSMCT32;
frame.FBW = 640 / BUFFER_WIDTH_SCALE;
iface.write_register(RegisterAddr::FRAME_1, frame);

ZBUFBits zbuf = {};
zbuf.ZMSK = 0; // Enable Z-write
zbuf.ZBP = 0x118000 / PAGE_ALIGNMENT_BYTES;
iface.write_register(RegisterAddr::ZBUF_1, zbuf);

SCISSORBits scissor = {};
scissor.SCAX0 = 0;
scissor.SCAY0 = 0;
scissor.SCAX1 = 640 - 1;
scissor.SCAY1 = 448 - 1;
iface.write_register(RegisterAddr::SCISSOR_1, scissor);

Then we nudge some registers to draw:

struct Vertex
{
    PackedRGBAQBits rgbaq;
    PackedXYZBits xyz;
} vertices[3] = {};

for (auto &vert : vertices)
{
   vert.rgbaq.A = 0x80;
   vert.xyz.Z = 1;
}

vertices[0].rgbaq.R = 0xff;
vertices[1].rgbaq.G = 0xff;
vertices[2].rgbaq.B = 0xff;

vertices[0].xyz.X = p0.x << SUBPIXEL_BITS;
vertices[0].xyz.Y = p0.y << SUBPIXEL_BITS;
vertices[1].xyz.X = p1.x << SUBPIXEL_BITS;
vertices[1].xyz.Y = p1.y << SUBPIXEL_BITS;
vertices[2].xyz.X = p2.x << SUBPIXEL_BITS;
vertices[2].xyz.Y = p2.y << SUBPIXEL_BITS;

PRIMBits prim = {};
prim.TME = 0; // Turn off texturing.
prim.IIP = 1; // Interpolate RGBA (Gouraud shading)
prim.PRIM = int(PRIMType::TriangleList);

static const GIFAddr addr[] = { GIFAddr::RGBAQ, GIFAddr::XYZ2 };
constexpr uint32_t num_registers = sizeof(addr) / sizeof(addr[0]);
constexpr uint32_t num_loops = sizeof(vertices) / sizeof(vertices[0]);
iface.write_packed(prim, addr, num_registers, num_loops, vertices);

This draws a triangle. We provide coordinates directly in screen-space.

And finally, we need to program the CRTC. Most of this is just copy-pasta from whatever games tend to do.

auto &priv = iface.get_priv_register_state();

priv.pmode.EN1 = 1;
priv.pmode.EN2 = 0;
priv.pmode.CRTMD = 1;
priv.pmode.MMOD = PMODEBits::MMOD_ALPHA_ALP;
priv.smode1.CMOD = SMODE1Bits::CMOD_NTSC;
priv.smode1.LC = SMODE1Bits::LC_ANALOG;
priv.bgcolor.R = 0x0;
priv.bgcolor.G = 0x0;
priv.bgcolor.B = 0x0;
priv.pmode.SLBG = PMODEBits::SLBG_ALPHA_BLEND_BG;
priv.pmode.ALP = 0xff;
priv.smode2.INT = 1;

priv.dispfb1.FBP = 0;
priv.dispfb1.FBW = 640 / BUFFER_WIDTH_SCALE;
priv.dispfb1.PSM = PSMCT32;
priv.dispfb1.DBX = 0;
priv.dispfb1.DBY = 0;
priv.display1.DX = 636; // Magic values that center the screen.
priv.display1.DY = 50; // Magic values that center the screen.
priv.display1.MAGH = 3; // scaling factor = MAGH + 1 = 4 -> 640 px wide.
priv.display1.MAGV = 0;
priv.display1.DW = 640 * 4 - 1;
priv.display1.DH = 448 - 1;

dump.write_vsync(0, iface);
dump.write_vsync(1, iface);

When the GS is dumped, we can load it up in PCSX2 and voila:

And here’s the same .gs dump is played through parallel-gs-replayer with RenderDoc. For debugging, I’ve spent a lot of time making it reasonably convenient. The images are debug storage images where I can store before and after color, depth, debug values for interpolants, depth testing state, etc, etc. It’s super handy to narrow down problem cases. The render pass can be split into 1 or more triangle chunks as needed.

To add some textures, and flex the capabilities of the CRTC a bit, we can try uploading a texture:

int chan;
auto *buf = stbi_load("/tmp/test.png", &w, &h, &chan, 4);
iface.write_image_upload(0x300000, PSMCT32, w, h, buf,
                         w * h * sizeof(uint32_t));
stbi_image_free(buf);

TEX0Bits tex0 = {};
tex0.PSM = PSMCT32;
tex0.TBP0 = 0x300000 / BLOCK_ALIGNMENT_BYTES;
tex0.TBW = (w + BUFFER_WIDTH_SCALE - 1) / BUFFER_WIDTH_SCALE;
tex0.TW = Util::floor_log2(w - 1) + 1;
tex0.TH = Util::floor_log2(h - 1) + 1;
tex0.TFX = COMBINER_DECAL;
tex0.TCC = 1; // Use texture alpha as blend alpha
iface.write_register(RegisterAddr::TEX0_1, tex0);

TEX1Bits tex1 = {};
tex1.MMIN = TEX1Bits::LINEAR;
tex1.MMAG = TEX1Bits::LINEAR;
iface.write_register(RegisterAddr::TEX1_1, tex1);

CLAMPBits clamp = {};
clamp.WMS = CLAMPBits::REGION_CLAMP;
clamp.WMT = CLAMPBits::REGION_CLAMP;
clamp.MINU = 0;
clamp.MAXU = w - 1;
clamp.MINV = 0;
clamp.MAXV = h - 1;
iface.write_register(RegisterAddr::CLAMP_1, clamp);

While PS2 requires POT sizes for textures, REGION_CLAMP is handy for NPOT. Super useful for texture atlases.

struct Vertex
{
    PackedUVBits uv;
    PackedXYZBits xyz;
} vertices[2] = {};

for (auto &vert : vertices)
    vert.xyz.Z = 1;

vertices[0].xyz.X = p0.x << SUBPIXEL_BITS;
vertices[0].xyz.Y = p0.y << SUBPIXEL_BITS;
vertices[1].xyz.X = p1.x << SUBPIXEL_BITS;
vertices[1].xyz.Y = p1.y << SUBPIXEL_BITS;
vertices[1].uv.U = w << SUBPIXEL_BITS;
vertices[1].uv.V = h << SUBPIXEL_BITS;

PRIMBits prim = {};
prim.TME = 1; // Turn on texturing.
prim.IIP = 0;
prim.FST = 1; // Use unnormalized coordinates.
prim.PRIM = int(PRIMType::Sprite);

static const GIFAddr addr[] = { GIFAddr::UV, GIFAddr::XYZ2 };
constexpr uint32_t num_registers = sizeof(addr) / sizeof(addr[0]);
constexpr uint32_t num_loops = sizeof(vertices) / sizeof(vertices[0]);
iface.write_packed(prim, addr, num_registers, num_loops, vertices);

Here we render a sprite with un-normalized coordinates.

Finally, we use the CRTC to do blending against white background.

priv.pmode.EN1 = 1;
priv.pmode.EN2 = 0;
priv.pmode.CRTMD = 1;
priv.pmode.MMOD = PMODEBits::MMOD_ALPHA_CIRCUIT1;
priv.smode1.CMOD = SMODE1Bits::CMOD_NTSC;
priv.smode1.LC = SMODE1Bits::LC_ANALOG;
priv.bgcolor.R = 0xff;
priv.bgcolor.G = 0xff;
priv.bgcolor.B = 0xff;
priv.pmode.SLBG = PMODEBits::SLBG_ALPHA_BLEND_BG;
priv.smode2.INT = 1;

priv.dispfb1.FBP = 0;
priv.dispfb1.FBW = 640 / BUFFER_WIDTH_SCALE;
priv.dispfb1.PSM = PSMCT32;
priv.dispfb1.DBX = 0;
priv.dispfb1.DBY = 0;
priv.display1.DX = 636; // Magic values that center the screen.
priv.display1.DY = 50; // Magic values that center the screen.
priv.display1.MAGH = 3; // scaling factor = MAGH + 1 = 4 -> 640 px wide.
priv.display1.MAGV = 0;
priv.display1.DW = 640 * 4 - 1;
priv.display1.DH = 448 - 1;

Glorious 256×179 logo 😀

Implementation details

The rendering pipeline

Before we get into the page tracker, it’s useful to define a rendering pipeline where synchronization is implied between each stage.

  • Synchronize CPU copy of VRAM to GPU. This is mostly unused, but happens for save state load, or similar
  • Upload data to VRAM (or perform local-to-local copy)
  • Update CLUT cache from VRAM
  • Unswizzle VRAM into VkImages that can be sampled directly, and handle palettes as needed, sampling from CLUT cache
  • Perform rendering
  • Synchronize GPU copy of VRAM back to CPU. This will be useful for readbacks. Then CPU should be able to unswizzle directly from a HOST_CACHED_BIT buffer as needed

This pipeline matches what we expect a game to do over and over:

  • Upload texture to VRAM
  • Upload palette to VRAM
  • Update CLUT cache
  • Draw with texture
    • Trigger unswizzle from VRAM into VkImage if needed
    • Begins building a “render pass”, a batch of primitives

When there are no backwards hazards here, we can happily keep batching and defer any synchronization. This is critical to get any performance out of this style of renderer.

Some common hazards here include:

Copy to VRAM which was already written by copy

This is often a false positive, but we cannot track per-byte. This becomes a simple copy barrier and we move on.

Copy to VRAM where a texture was sampled from, or CLUT cache read from

Since the GS has a tiny 4 MiB VRAM, it’s very common that textures are continuously streamed in, sampled from, and thrown away. When this is detected, we have to submit all vram copy work, all texture unswizzle work and then begin a new batch. Primitive batches are not disrupted.

This means we’ll often see:

  • Copy xN
  • Barrier
  • Unswizzle xN
  • Barrier
  • Copy xN
  • Barrier
  • Unswizzle xN
  • Barrier
  • Rendering

Sample texture that was rendered to

Similar, but here we need to flush out everything. This basically breaks the render pass and we start another one. Too many of these is problematic for performance obviously.

Copy to VRAM where rendering happened

Basically same as sampling textures, this is a full flush.

Other hazards are ignored, since they are implicitly handled by our pipeline.

Page tracker

Arguably, the hardest part of GS emulation is dealing with hazards. VRAM is read and written to with reckless abandon and any potential read-after-write or write-after-write hazard needs to be dealt with. We cannot rely on any game doing this for us, since PS2 GS just deals with sync in most cases, and TEXFLUSH is the only real command games will use (or forget to use).

Tracking per byte is ridiculous, so my solution is to first subdivide the 4 MiB VRAM into pages. A page is the unit for frame buffers and depth buffers, so it is the most meaningful place to start.

PageState

On page granularity, we track:

  • Pending frame buffer write?
  • Pending frame buffer read? (read-only depth)

Textures and VRAM copies have 256 byte alignment, and to avoid a ton of false positives, we need to track on a per-block basis. There are 32 blocks per page, so a u32 bit-mask is okay.

  • VRAM copy writes
  • VRAM copy reads
  • Pending read into CLUT cache or VkImage
  • Blocks which have been clobbered by any write, on next texture cache invalidate, throw away images that overlap

As mentioned earlier, there are also cases where you can render to 24-bit color, while sampling from the upper 8-bits without hazard. We need to optimize for that case too, so there is also:

  • A write mask for framebuffers
  • A read mask for textures

In the example above, FB write mask is 0xffffff and texture cache mask is 0xff000000. No overlap, no invalidate 😀

For host access, there are also timeline semaphore values per page. These values state which sync point to wait for if the host desires mapped read or mapped write access. Mapped write access may require more sync than mapped read if there are pending reads on that page.

Caching textures

Every page contains a list of VkImages which have been associated with it. When a page’s textures has been invalidated, the image is destroyed and has to be unswizzled again from VRAM.

There is a one-to-many relationship with textures and pages. A texture may span more than one page, and it’s enough that only one page is clobbered before the texture is invalidated.

Overall, there are a lot of micro-details here, but the important things to note here is that conservative and simple tracking will not work on PS2 games. Tracking at a 256 byte block level and considering write/read masks is critical.

Special cases

There are various situations where we may have false positives due to how textures work. Since textures are POT sized, it’s fairly common for e.g. a 512×448 texture of a render target to be programmed as a 512×512 texture. The unused region should ideally be clamped out with REGION_CLAMP, but most games don’t. A render target might occupy those unused pages. As long as the game’s UV coordinates don’t extend into the unused red zone, there are no hazards, but this is very painful to track. We would have to analyze every single primitive to detect if it’s sampling into the red zone.

As a workaround, we ignore any potential hazard in that red zone, and just pray that a game isn’t somehow relying on ridiculous spooky-action-at-a-distance hazards to work in the game’s favor.

There are more spicy special cases, especially with texture sampling feedback, but that will be for later.

Updating CLUT in a batched way

Since we want to batch texture uploads, we have to batch CLUT uploads too. To make this work, we have 1024 copies of CLUT, a ring buffer of snapshots.

One workgroup loops through the updates and writes them to an SSBO. I did a similar thing for N64 RDP’s TMEM update, where TMEM was instanced. Fortunately, CLUT update is far simpler than TMEM update.

shared uint tmp_clut[512];

// ...

// Copy from previous instance to allow a
// CLUT entry to be partially overwritten and used later
uint read_index = registers.read_index * CLUT_SIZE_16;
tmp_clut[gl_LocalInvocationIndex] =
    uint(clut16.data[read_index]);
tmp_clut[gl_LocalInvocationIndex + 256u] =
    uint(clut16.data[read_index + 256u]);
barrier();

for (uint i = 0; i < registers.clut_count; i++)
{
  // ...
  if (active_lane)
  {
    // update tmp_clut. If 256 color, all threads participate.
    // 16 color update is a partial update.
  }

  // Flush current CLUT state to SSBO.
  barrier();
  clut16.data[gl_LocalInvocationIndex + clut.instance * CLUT_SIZE_16] =
    uint16_t(tmp_clut[gl_LocalInvocationIndex]);
  clut16.data[gl_LocalInvocationIndex + clut.instance * CLUT_SIZE_16 + 256u] =
    uint16_t(tmp_clut[gl_LocalInvocationIndex + 256u]);
  barrier();
}

One potential optimization is that for 256 color / 32 bpp updates, we can parallelize the CLUT update, since nothing from previous iterations will be preserved, but the CLUT update time is tiny anyway.

Unswizzling textures from VRAM

Since this is Vulkan, we can just allocate a new VkImage, suballocate it from VkDeviceMemory and blast it with a compute shader.

Using Vulkan’s specialization constants, we specialize the texture format and all the swizzling logic becomes straight forward code.

REGION_REPEAT shenanigans is also resolved here, so that the ubershader doesn’t have to consider that case and do manual bilinear filtering.

Even for render targets, we roundtrip through the VRAM SSBO. There is not really a point going to the length of trying to forward render targets into textures. Way too many bugs to squash and edge cases to think about.

Triangle setup and binning

Like paraLLEl-RDP, paraLLEl-GS is a tile-based renderer. Before binning can happen, we need triangle setup. As inputs, we provide attributes in three arrays.

Position
struct VertexPosition
{
  ivec2 pos;
  float z;     // TODO: Should be uint for 32-bit Z.
  int padding; // Free real-estate?
};
Per-Vertex attributes
struct VertexAttribute
{
  vec2 st;
  float q;
  uint rgba; // unpackUnorm4x8
  float fog; // overkill, but would be padding anyway
  u16vec2 uv;
};
Per-primitive attributes
struct PrimitiveAttribute
{
  i16vec4 bb; // Scissor
  // Index into state UBO, as well as misc state bits.
  uint state;
  // Texture state which should be scalarized. Affects code paths.
  // Also holds the texture index (for bindless).
  uint tex;
  // Texture state like lod scaling factors, etc.
  // Does not affect code paths.
  uint tex2;  
  uint alpha; // AFIX / AREF
  uint fbmsk;
  uint fogcol;
};

For rasterization, we have a straight forward barycentric-based rasterizer. It is heavily inspired by https://fgiesen.wordpress.com/2011/07/06/a-trip-through-the-graphics-pipeline-2011-part-6/, which in turn is based on A Parallel Algorithm for Polygon Rasterization (Paneda, 1988) and describes the “standard” way to write a rasterizer with parallel hardware. Of course, the PS2 GS is DDA, i.e. a scanline rasterizer, but in practice, this is just a question of nudging ULPs of precision, and since I’m not aware of a bit-exact description of the GS’s DDA, this is fine. paraLLEl-RDP implements the raw DDA form for example. It’s certainly possible if we have to.

As an extension to a straight-forward triangle rasterizer, I also need to support parallelograms. This is used to implement wide-lines and sprites. Especially wide-line is kinda questionable, but I’m not sure it’s possible to fully solve up-scaling + Bresenham in the general case. At least I haven’t run into a case where this really matters.

Evaluating coverage and barycentric I/J turns into something like this:

bool evaluate_coverage_single(PrimitiveSetup setup,
  bool parallelogram, 
  ivec2 parallelogram_offset,
  ivec2 coord, inout float i, inout float j)
{
  int a = idot3(setup.a, coord);
  int b = idot3(setup.b, coord);
  int c = idot3(setup.c, coord);

  precise float i_result = float(b) * setup.inv_area + setup.error_i;
  precise float j_result = float(c) * setup.inv_area + setup.error_j;
  i = i_result;
  j = j_result;

  if (parallelogram && a.x < 0)
  {
    b += a + parallelogram_offset.x;
    c += a + parallelogram_offset.y;
    a = 0;
  }

  return all(greaterThanEqual(ivec3(a, b, c), ivec3(0)));
}

inv_area is computed in a custom fixed-point RCP, which is ~24.0 bit accurate. Using the standard GPU RCP would be bad since it’s just ~22.5 bit accurate and not consistent across implementations. There is no reason to skimp on reproducibility and accuracy, since we’re not doing work per-pixel.

error_i and error_j terms are caused by the downsampling of the edge equations and tie-break rules. As a side effect of the GS’s [-4k, +4k] pixel range, the range of the cross-product requires 33-bit in signed integers. By downsampling a bit, we can get 32-bit integer math to work just fine with 8 sub-pixel accuracy for super-sampling / multi-sampling. Theoretically, this means our upper up-sampling limit is 8×8, but that’s ridiculous anyway, so we’re good here.

The parallelogram offsets are very small numbers meant to nudge the tie-break rules in our favor as needed. The exact details of the implementation escape me. I wrote that code years ago. It’s not very hard to derive however.

Every primitive gets a struct of transformed attributes as well. This is only read if we actually end up shading a primitive, so it’s important to keep this separate to avoid polluting caches with too much garbage.

struct TransformedAttributes
{
  vec4 stqf0;
  vec4 stqf1;
  vec4 stqf2;
  uint rgba0;
  uint rgba1;
  uint rgba2;
  uint padding;
  vec4 st_bb;
};

Using I/J like this will lead to small inaccuracies when interpolating primitives which expect to land exactly on the top-left corner of a texel with NEAREST filtering. To combat this, a tiny epsilon offset is used when snapping texture coordinates. Very YOLO, but what can you do. As far as I know, hardware behavior is sub-texel floor, not sub-texel round.

precise vec2 uv_1 = uv * scale_1;

// Want a soft-floor here, not round behavior.
const float UV_EPSILON_PRE_SNAP = 1.0 / 16.0;
// We need to bias less than 1 / 512th texel, so that linear filter will RTE to correct subpixel.
// This is a 1 / 1024th pixel bias to counter-act any non-POT inv_scale_1 causing a round-down event.
const float UV_EPSILON_POST_SNAP = 16.0 / 1024.0;

if (sampler_clamp_s)
  uv_1.x = texture_clamp(uv_1.x, region_coords.xz, LOD_1);
if (sampler_clamp_t)
  uv_1.y = texture_clamp(uv_1.y, region_coords.yw, LOD_1);

// Avoid micro-precision issues with UV and flooring + nearest.
// Exact rounding on hardware is somwhat unclear.
// SotC requires exact rounding precision and is hit particularly bad.
// If the epsilon is too high, then FF X save screen is screwed over,
// so ... uh, ye.
// We likely need a more principled approach that is actually HW accurate in fixed point.
uv_1 = (floor(uv_1 * 16.0 + UV_EPSILON_PRE_SNAP) + UV_EPSILON_POST_SNAP) *
       inv_scale_1 * 0.0625;

Binning

This is mostly uninteresting. Every NxN pixel block gets an array of u16 primitive indices to shade. This makes the maximum number of primitives per render pass 64k, but that’s enough for PS2 games. Most games I’ve seen so far tend to be between 10k and 30k primitives for the “main” render pass, but I haven’t tested the real juggernauts of primitive grunt yet, but even so, having to do a little bit of incremental rendering isn’t a big deal.

NxN is usually 32×32, but it can be dynamically changed depending on how heavy the geometry load is. For large resolutions and high primitive counts, the binning and memory cost is unacceptable if the resolution is just 16×16 for example. One subgroup is responsible for iterating through all primitives in a block.

Since binning and triangle is state-less, triangle-setup and binning for back-to-back passes are batched up nicely to avoid lots of silly barriers.

The ubershader

A key difference between N64 and PS2 is fill-rate and per-pixel complexity. For N64, the ideal approach is to specialize the rasterizing shader, write out per-pixel color + depth + coverage + etc, then merge that data in a much simpler ubershader that only needs to consider depth and blend state rather than full texturing state and combiner state. This is very bandwidth intensive on the GPU, but the alternative is the slowest ubershader written by man. We’re saved by the fact that N64 fill-rate is abysmal. Check out this video by Kaze to see how horrible it is.

The GS is a quite different beast. Fill-rate is very high, and per-pixel complexity is fairly low, so a pure ubershader is viable. We can also rely on bindless this time around too, so texturing complexity becomes a fraction of what I had to deal with on N64.

Fine-grained binning

Every tile is 4×4, 4×8 and 8×8 for subgroup sizes 16, 32 and 64 respectively. For super-sampling it’s even smaller (it’s 4×4 / 4×8 / 8×8 in the higher resolution domain instead).

In the outer loop, we pull in up to SubgroupSize’s worth of primitives, and bin them in parallel.

for (int i = 0; i < tile.coarse_primitive_count;
     i += int(gl_SubgroupSize))
{
  int prim_index = i + int(gl_SubgroupInvocationID);
  bool is_last_iteration = i + int(gl_SubgroupSize) >= 
                           tile.coarse_primitive_count;

  // Bin primitives to tile.
  bool binned_to_tile = false;
  uint bin_primitive_index;
  if (prim_index < tile.coarse_primitive_count)
  {
    bin_primitive_index = 
      uint(coarse_primitive_list.data[
           tile.coarse_primitive_list_offset + prim_index]);
    binned_to_tile = primitive_intersects_tile(bin_primitive_index);
  }

  // Iterate per binned primitive, do per pixel work now.
  // Scalar loop.
  uvec4 work_ballot = subgroupBallot(binned_to_tile);

In the inner loop, we can do a scalarized loop which checks coverage per-pixel, one primitive at a time.

// Scalar data
uint bit = subgroupBallotFindLSB(work_ballot);

if (gl_SubgroupSize == 64)
{
  if (bit >= 32)
    work_ballot.y &= work_ballot.y - 1;
  else
    work_ballot.x &= work_ballot.x - 1;
}
else
{
  work_ballot.x &= work_ballot.x - 1;
}

shade_primitive_index = subgroupShuffle(bin_primitive_index, bit);

Early Z

We can take advantage of early-Z testing of course, but we have to be careful if there are rasterized pixels we haven’t resolved yet, and there are Z-writes in flight. In this case we have to defer to late Z to perform test.

// We might have to remove opaque flag.
bool pending_z_write_can_affect_result =
  (pixel.request.z_test || !pixel.request.z_write) &&
  pending_shade_request.z_write;

if (pending_z_write_can_affect_result)
{
  // Demote the pixel to late-Z,
  // it's no longer opaque and we cannot discard earlier pixels.
  // We need to somehow observe the previous results.
  pixel.opaque = false;
}

Deferred on-tile shading

Since we’re an uber-shader, all pixels are “on-chip”, i.e. in registers, so we can take advantage of culling pixels that won’t be visible anyway. The basic idea here is that after rasterization, if a pixel is considered opaque, it will simply replace the shading request that exists for that framebuffer coordinate. It won’t be visible at all anyway.

Lazy pixel shading

We only need to perform shading when we really have to, i.e., we’re shading a pixel that depends on the previous pixel’s results. This can happen for e.g. alpha test (if test fails, we preserve existing data), color write masks, or of course, alpha blending.

If our pixel remains opaque, we can just kill the pending pixel shade request. Very nice indeed. The gain here wasn’t as amazing as I had hoped since PS2 games love blending, but it helps culling out a lot of shading work.

if (pixel.request.coverage > 0)
{
  need_flush = !pixel.opaque && pending_shade_request.coverage > 0;

  // If there is no hazard, we can overwrite the pending pixel.
  // If not, defer the update until we run a loop iteration.
  if (!need_flush)
  {
    set_pending_shade_request(pixel.request, shade_primitive_index);
    pixel.request.coverage = 0;
    pixel.request.z_write = false;
  }
}

If we have flushes that need to happen, we do so if one pixel needs it. It’s just as fast to resolve all pixels anyway.

// Scalar branch
if (subgroupAny(need_flush))
{
  shade_resolve();
  if (has_work && pixel.request.coverage > 0)
    set_pending_shade_request(pixel.request, shade_primitive_index);
}

The resolve is a straight forward waterfall loop that stays in uniform control flow to be well defined on devices without maximal reconvergence support.

while (subgroupAny(has_work))
{
  if (has_work)
  {
    uint state_index =
      subgroupBroadcastFirst(pending_shade_request.state);
    uint tex = subgroupBroadcastFirst(prim_tex);
    if (state_index == pending_shade_request.state && prim_tex == tex)
    {
      has_work = false;
      shade_resolve(pending_primitive_index, state_index, tex);
    }
  }
}

This scalarization ensures that all branches on things like alpha test mode, blend modes, etc, are purely scalar, and GPUs like that. Scalarizing on the texture index is technically not that critical, but it means we end up hitting the same branches for filtering modes, UBOs for scaling factors are loaded uniformly, etc.

When everything is done, the resulting framebuffer color and depth is written out to SSBO. GPU bandwidth is kept to a minimum, just like a normal TBDR renderer.

Super-sampling

Just implementing single sampled rendering isn’t enough for this renderer to be really useful. The software renderer is certainly quite fast, but not fast enough to keep up with intense super-sampling. We can fix that now.

For e.g. 8x SSAA, we keep 10 versions of VRAM on the GPU.

  • 1 copy represents the single-sampled VRAM. It is super-sampled.
  • 1 copy represents the reference value for single-sampled VRAM. This allows us to track when we should discard the super-samples and splat the single sample to all. This can happen if someone copies to VRAM over a render target for whatever reason.
  • 8 copies which each represent the super-samples. Technically, we can reconstruct a higher resolution image from these samples if we really want to, but only the CRTC could easily do that.

When rendering super-sampled, we load the single-sampled VRAM and reference. If they match, we load the super-sampled version. This is important for cases where we’re doing incremental rendering.

On tile completion we use clustered subgroup ops to do multi-sample resolve, then write out the super-samples, and the two single-sampled copies.

uvec4 ballot_color = subgroupBallot(fb_color_dirty);
uvec4 ballot_depth = subgroupBallot(fb_depth_dirty);

// No need to mask, we only care about valid ballot for the
// first sample we write-back.
if (NUM_SAMPLES >= 16)
{
  ballot_color |= ballot_color >> 8u;
  ballot_depth |= ballot_depth >> 8u;
}

if (NUM_SAMPLES >= 8)
{
  ballot_color |= ballot_color >> 4u;
  ballot_depth |= ballot_depth >> 4u;
}

if (NUM_SAMPLES >= 4)
{
  ballot_color |= ballot_color >> 2u;
  ballot_depth |= ballot_depth >> 2u;
}

ballot_color |= ballot_color >> 1u;
ballot_depth |= ballot_depth >> 1u;

// GLSL does not accept cluster reduction as spec constant.
if (NUM_SAMPLES == 16)
  fb_color = packUnorm4x8(subgroupClusteredAdd(
    unpackUnorm4x8(fb_color), 16) / 16.0);
else if (NUM_SAMPLES == 8)
  fb_color = packUnorm4x8(subgroupClusteredAdd(
    unpackUnorm4x8(fb_color), 8) / 8.0);
else if (NUM_SAMPLES == 4)
  fb_color = packUnorm4x8(subgroupClusteredAdd(
    unpackUnorm4x8(fb_color), 4) / 4.0);
else
  fb_color = packUnorm4x8(subgroupClusteredAdd(
    unpackUnorm4x8(fb_color), 2) / 2.0);

fb_color_dirty = subgroupInverseBallot(ballot_color);
fb_depth_dirty = subgroupInverseBallot(ballot_depth);

The main advantage of super-sampling over straight up-scaling is that up-scaling will still have jagged edges, and super-sampling retains a coherent visual look where 3D elements have similar resolution as UI elements. One of my pet peeves is when UI elements have a significantly different resolution from 3D objects and textures. HD texture packs can of course alleviate that, but that’s a very different beast.

Super-sampling also lends itself very well to CRT post-processing shading, which is also a nice bonus.

Dealing with super-sampling artifacts

It’s a fact of life that super-sampling always introduces horrible artifacts if not handled with utmost care. Mitigating this is arguably easier with software renderers over traditional graphics APIs, since we’re not limited by the fixed function interpolators. These tricks won’t make it perfect by any means, but it greatly mitigates jank in my experience, and I already fixed many upscaling bugs that GSdx Vulkan backend does not solve as we shall see later.

Sprite primitives should always render at single-rate

Sprites are always UI elements or similar, and games do not expect us to up-scale them. Doing so either results in artifacts where we sample outside the intended rect, or we risk overblurring the image if bilinear filtering is used.

The trick here is just to force-snap the pixel coordinate we use when rasterizing and interpolating. This is very inefficient of course, but UI shouldn’t take up the entire screen. And if it does (like in a menu), the GPU load is tiny anyway.

const uint SNAP_RASTER_BIT = (1u << STATE_BIT_SNAP_RASTER);
const uint SNAP_ATTR_BIT = (1u << STATE_BIT_SNAP_ATTRIBUTE);

if (SUPER_SAMPLE && (prim_state & SNAP_RASTER_BIT) != 0)
  fb_pixel = tile.fb_pixel_single_rate;

res.request.coverage = evaluate_coverage(
  prim, fb_pixel, i, j,
  res.request.multisample, SAMPLING_RATE_DIM_LOG2);

Flat primitives should interpolate at single-pixel coordinate

Going further, we can demote SSAA interpolation to MSAA center interpolation dynamically. Many UI elements are unfortunately rendered with normal triangles, so we have to be a bit more careful. This snap only affects attribute interpolation, not Z of course.

res.request.st_bb = false;
if (SUPER_SAMPLE &&
    (prim_state & (SNAP_RASTER_BIT | SNAP_ATTR_BIT)) == SNAP_ATTR_BIT)
{
  vec2 snap_ij = evaluate_barycentric_ij(
    prim.b, prim.c, prim.inv_area,
    prim.error_i, prim.error_j, tile.fb_pixel_single_rate,
    SAMPLING_RATE_DIM_LOG2);

  i = snap_ij.x;
  j = snap_ij.y;
  res.request.st_bb = true;
}

Here, we snap interpolation to the top-left pixel. This fixes any artifacts for primitives which align their rendering to a pixel center, but some games are mis-aligned, so this snapping can cause texture coordinates to go outside the expected area. To clean this up, we compute a bounding box of final texture coordinates. Adding bounding boxes can technically cause notorious block-edge artifacts, but that was mostly a thing on PS1 since emulators like to convert nearest sampling to bilinear.

The heuristic for this is fairly simple. If perspective is used, if all vertices in a triangle have exact same Q, we assume it’s a flat UI primitive. The primitive’s Z coordinates must also match. This is done during triangle setup on the GPU. There can of course be false positives here, but it should be rare. In my experience this hack works well enough in the games I tried.

Results

Here’s a good example of up-sampling going awry in PCSX2. This is with Vulkan backend:

Notice the bloom on the glass being mis-aligned and a subtle (?) rectangular pattern being overlaid over the image. This is caused by a post-processing pass rendering in a page-like pattern, presumably to optimize for GS caching behavior.

 

With 8x SSAA in paraLLEl-GS it looks like this instead. There is FSR1 post-upscale in effect here which changes the look a bit, but the usual trappings of bad upscale cannot be observed here. This is another reason to do super-sample; texture mis-alignment has a tendency to fix itself.

Also, if you’re staring at the perf numbers, this is RX 7600 in a low power state :’)

Typical UI issues can be seen in games as well. Here’s native resolution:

and 4x upscale, which … does not look acceptable.

This UI is tricky to render in upscaled mode, since it uses triangles, but the MSAA snap trick above works well and avoids all artifacts. With straight upscale, this is hard to achieve in normal graphics APIs since you’d need interpolateAtOffset beyond 0.5 pixels, which isn’t supported. Perhaps you could do custom interpolation with derivatives or something like that, but either way, this glitch can be avoided. The core message is basically to never upscale UI beyond plain nearest neighbor integer scale. It just looks bad.

There are cases where PCSX2 asks for high blending accuracy. One example is MGS2, and I found a spot where GPU perf is murdered. My desktop GPU cannot keep 60 FPS here at 4x upscale. PCSX2 asks you to turn up blend-accuracy for this game, but …

What happens here is we hit the programmable blending path with barrier between every primitive. Ouch! This wouldn’t be bad for the tiler mobile GPUs, but for a desktop GPU, it is where perf goes to die. The shader in question does subpassLoad and does programmable blending as expected. Barrier, tiny triangle, barrier, tiny triangle, hnnnnnnng.

paraLLEl-GS on the other hand always runs with 100% blend accuracy (assuming no bugs of course). Here’s 16xSSAA (equivalent to 4x upscale). This is just 25 W and 17% GPU utilization on RX 7600. Not bad.

Other difficult cases include texture sampling feedback. One particular case I found was in Valkyrie Profile 2.

This game has a case where it’s sampling it’s own pixel’s alpha as a palette index. Quirky as all hell, and similar to MGS2 there’s a barrier between every pixel.

In paraLLEl-GS, this case is detected, and we emit a magical texture index, which resolved to just looking at in-register framebuffer color instead. Programmable blending go brr. These cases have to be checked per primitive, which is quite rough on CPU time, but it is what it is. If we don’t hit the good path, GPU performance completely tanks.

The trick here is to analyze the effective UV coordinates, and see if UV == framebuffer position. If we fall off this path, we have to go via texture uploads, which is bad.

ivec2 uv0_delta = uv0 - pos[0].pos;
ivec2 uv1_delta = uv1 - pos[1].pos;
ivec2 min_delta = min(uv0_delta, uv1_delta);
ivec2 max_delta = max(uv0_delta, uv1_delta);

if (!quad)
{
  ivec2 uv2_delta = uv2 - pos[2].pos;
  min_delta = min(min_delta, uv2_delta);
  max_delta = max(max_delta, uv2_delta);
}

int min_delta2 = min(min_delta.x, min_delta.y);
int max_delta2 = max(max_delta.x, max_delta.y);

// The UV offset must be in range of [0, 2^SUBPIXEL_BITS - 1].
// This guarantees snapping with NEAREST.
// 8 is ideal. That means pixel centers during interpolation
// will land exactly in the center of the texel.
// In theory we could allow LINEAR if uv delta was
// exactly 8 for all vertices.
if (min_delta2 < 0 || max_delta2 >= (1 << SUBPIXEL_BITS))
  return ColorFeedbackMode::Sliced;

// Perf go brrrrrrr.
return ColorFeedbackMode::Pixel;
if (feedback_mode == ColorFeedbackMode::Pixel)
{
  mark_render_pass_has_texture_feedback(ctx.tex0.desc);
  // Special index indicating on-tile feedback.
  // We could add a different sentinel for depth feedback.
  // 1024k CLUT instances and 32 sub-banks. Fits in 15 bits.
  // Use bit 15 MSB to mark feedback texture.
  return (1u << (TEX_TEXTURE_INDEX_BITS - 1u)) |
         (render_pass.clut_instance * 32 + uint32_t(ctx.tex0.desc.CSA));
}

It’s comfortably full-speed on PCSX2 here, despite the copious number of barriers, but paraLLEl-GS is reasonably close perf-wise, actually. At 8x SSAA.

Overall, we get away with 18 render pass barriers instead of 500+ which was the case without this optimization. You may notice the interlacing artifacts on the swirlies. Silly game has a progressive scan output, but downsamples it on its own to a field before hitting CRTC, hnnnnng 🙁 Redirecting framebuffer locations in CRTC might work as a per-game hack, but either way, I still need to consider a better de-interlacer. Some games actually render explicitly in fields (640×224), which is very annoying.

This scene in the MGS2 intro also exposes some funny edge cases with sampling.

To get the camo effect, it’s sampling its own framebuffer as a texture, with overlapping coordinates, but not pixel aligned, so this raises some serious questions about caching behavior. PCSX2 doesn’t seem to add any barriers here, and I kinda had to do the same thing. It looks fine to me compared to software renderer at least.

if (feedback_mode == ColorFeedbackMode::Sliced)
{
  // If game explicitly clamps the rect to a small region,
  // it's likely doing well-defined feedbacks.
  // E.g. Tales of Abyss main menu ping-pong blurs.
  // This code is quite flawed,
  // and I'm not sure what the correct solution is yet.
  if (desc.clamp.desc.WMS == CLAMPBits::REGION_CLAMP &&
      desc.clamp.desc.WMT == CLAMPBits::REGION_CLAMP)
  {
    ivec4 clamped_uv_bb(
      int(desc.clamp.desc.MINU),
      int(desc.clamp.desc.MINV),
      int(desc.clamp.desc.MAXU),
      int(desc.clamp.desc.MAXV));

    ivec4 hazard_bb(
      std::max<int>(clamped_uv_bb.x, bb.x),
      std::max<int>(clamped_uv_bb.y, bb.y),
      std::min<int>(clamped_uv_bb.z, bb.z),
      std::min<int>(clamped_uv_bb.w, bb.w));

    cache_texture = hazard_bb.x > hazard_bb.z ||
                    hazard_bb.y > hazard_bb.w;
  }
  else
  {
    // Questionable,
    // but it seems almost impossible to do this correctly and fast.
    // Need to emulate the PS2 texture cache exactly,
    // which is just insane.
    // This should be fine.
    cache_texture = false;
  }
}

If we’re in a mode where texture points directly to the frame buffer we should relax the hazard tracking a bit to avoid 2000+ barriers. This is clearly spooky since Tales of Abyss’s bloom effect as shown earlier depends on this to be well behaved, but in that case, at least it uses REGION_CLAMP to explicitly mark the ping-pong behavior. I’m not sure what the proper solution is here.

The only plausible solution to true bit-accuracy with real hardware is to emulate the caches directly, one pixel at a time. You can kiss performance good bye in that case.

One of the worst stress tests I’ve found so far has to be Shadow of the Collosus. Just in the intro, we can make the GPU kneel down to 24 FPS with maximum blend accuracy on PCSX2, at just 2x upscale! Even with normal blending accuracy, it is extremely heavy during the intro cinematic.

At 8x SSAA, perf is still looking pretty good for paraLLEl-GS, but it’s clearly sweating now.

We’re actually still CPU bound on the geometry processing. Optimizing the CPU code hasn’t been a huge priority yet. There’s unfortunately a lot of code that has to run per-primitive, where hazards can happen around every corner that has to be dealt with somehow. I do some obvious optimizations, but it’s obviously not as well-oiled as PCSX2 in that regard.

Deck?

It seems fast enough to comfortably do 4x SSAA. Maybe not in SotC, but … hey. 😀

What now?

For now, the only real way to test this is through GS dumps. There’s a hack-patch for PCSX2 that lets you dump out a raw GS trace, which can be replayed. This works via mkfifo as a crude hack to test in real-time, but some kind of integration into an emulator needs to happen at some point if this is to turn into something that’s useful for end users.

There’s guaranteed to be a million bugs lurking since the PS2 library is ridiculously large and there’s only so much I can be arsed to test myself. At least, paraLLEl-GS has now become my preferred way to play PS2 games, so I can say mission complete.

A potential use case for this is due to its standalone library nature, it may be useful as very old-school rendering API for the old greybeards around that still yearn for the day of PS2 programming for whatever reason :p

Modernizing Granite’s mesh rendering

Granite’s renderer is currently quite old school. It was written with 2017 mobile hardware in mind after all. Little to no indirect drawing, a bindful material system, etc. We’ve all seen that a hundred times before. Granite’s niche ended up exploring esoteric use cases, not high-end rendering, so it was never a big priority for me to fix that.

Now that mesh shading has starting shipping and is somewhat proven in the wild with several games shipping UE5 Nanite, and Alan Wake II – which all rely on mesh shaders to not run horribly slow – it was time to make a more serious push towards rewriting the entire renderer in Granite. This has been a slow burn project that’s been haunting me for almost half a year at this point. I haven’t really had the energy to rewrite a ton of code like this in my spare time, but as always, holidays tend to give me some energy for these things. Video shenanigans have also kept me distracted this fall.

I’m still not done with this rewrite, but enough things have fallen into place, that I think it’s time to write down my findings so far.

Design requirements

Reasonable fallbacks

I had some goals for this new method. Unlike UE5 Nanite and Alan Wake II, I don’t want to hard-require actual VK_EXT_mesh_shader support to run acceptably. Just thinking in terms of meshlets should benefit us in plain multi-draw-indirect (MDI) as well. For various mobile hardware that doesn’t support MDI well (or at all …), I’d also like a fallback path that ends up using regular direct draws. That fallback path is necessary to evaluate performance uplift as well.

What Nanite does to fallback

This is something to avoid. Nanite relies heavily on rendering primitive IDs to a visibility buffer, where attributes are resolved later. In the primary compute software rasterizer, this becomes a 64-bit atomic, and in the mesh shader fallback, a single primitive ID is exported to fragment stage as a per-primitive varying, where fragment shader just does the atomic (no render targets, super fun to debug …). The problem here is that per-primitive varyings don’t exist in the classic vertex -> fragment pipeline. There are two obvious alternatives to work around this:

  • Geometry shaders. Pass-through mode can potentially be used if all the stars align on supported hardware, but using geometry shaders should revoke your graphics programmer’s license.
  • Unroll a meshlet into a non-indexed draw. Duplicate primitive ID into 3 vertices. Use flat shading to pull in the primitive ID.

From my spelunking in various shipped titles, Nanite does the latter, and fallback rendering performance is halved as a result (!). Depending on the game, meshlet fallbacks are either very common or very rare, so real world impact is scene and resolution dependent, but Immortals of Aveum lost 5-15% FPS when I tested it.

The Visibility Buffer: A Cache-Friendly Approach to Deferred Shading suggests rendering out a visibility G-Buffer using InstanceID (fed through some mechanism) and SV_PrimitiveID, which might be worth exploring at some point. I’m not sure why Nanite did not go that route. It seems like it would have avoided the duplicated vertices.

Alan Wake II?

Mesh shaders are basically a hard requirement for this game. It will technically boot without mesh shader support, but the game gives you a stern warning about performance, and they are not kidding. I haven’t dug into what the fallback is doing, but I’ve seen people posting videos demonstrating sub-10 FPS on a 1080 Ti. Given the abysmal performance, I wouldn’t be surprised if they just disabled all culling and draw everything in the fallback.

A compressed runtime meshlet format

While studying https://github.com/zeux/meshoptimizer I found support for compressed meshes, a format that was turned into a glTF EXT. It seems to be designed for decompressing on CPU (completely serial algorithm), which was not all that exciting for me, but this sparked an idea. What if I could decompress meshlets on the GPU instead? There are two ways this can be useful:

  • Would it be fast enough to decompress inline inside the mesh shader? This can potentially save a lot of read bandwidth during rendering and save precious VRAM.
  • Bandwidth amplifier on asset loading time. Only the compressed meshlet format needs to go over PCI-e wire, and we decompress directly into VRAM. Similar idea to GDeflate and other compression formats, except I should be able to come up with something that is way faster than a general purpose algorithm and also give decent compression ratios.

I haven’t seen any ready-to-go implementation of this yet, so I figured this would be my starting point for the renderer. Always nice to have an excuse to write some cursed compute shaders.

Adapting to implementations

One annoying problem with mesh shading is that different vendors have very different fast paths through their hardware. There is no single implementation that fits all. I’ve spent some time testing various parameters and observe what makes NV and AMD go fast w.r.t. mesh shaders, with questionable results. I believe this is the number 1 reason mesh shaders are still considered a niche feature.

Since we’re baking meshlets offline, the format itself must be able to adapt to implementations that prefer 32/64/128/256 primitive meshlets. It must also adapt nicely to MultiDrawIndirect-style rendering.

Random-access

It should be efficient to decode meshlets in parallel, and in complete isolation.

The format

I went through some (read: way too many) design iterations before landing on this design.

256 vert/prim meshlets

Going wide means we get lower culling overhead and emitting larger MDI calls avoids us getting completely bottlenecked on command stream frontend churn. I tried going lower than 256, but performance suffered greatly. 256 seemed like a good compromise. With 256 prim/verts, we can use 8-bit index buffers as well, which saves quite a lot of memory.

Sublets – 8×32 grouping

To consider various hardware implementations, very few will be happy with full, fat 256 primitive meshlets. To remedy this, the encoding is grouped in units of 32 – a “sublet” – where we can shade the 8 groups independently, or have larger workgroups that shade multiple sublets together. Some consideration is key to be performance portable. At runtime we can specialize our shaders to fit whatever hardware we’re targeting.

Using grouping of 32 is core to the format as well, since we can exploit NV warps being 32-wide and force Wave32 on RDNA hardware to get subgroup accelerated mesh shading.

Format header

// Can point to mmap-ed file.
struct MeshView
{
    const FormatHeader *format_header;
    const Bound *bounds;
    const Bound *bounds_256; // Used to cull in units of 256 prims
    const Stream *streams;
    const uint32_t *payload;
    uint32_t total_primitives;
    uint32_t total_vertices;
    uint32_t num_bounds;
    uint32_t num_bounds_256;
};
struct FormatHeader
{
    MeshStyle style;
    uint32_t stream_count;
    uint32_t meshlet_count;
    uint32_t payload_size_words;
};

The style signals type of mesh. This is naturally engine specific.

  • Wireframe: A pure position + index buffer
  • Textured: Adds UV + Normal + Tangent
  • Skinned: Adds bone indices and weights on top

A stream is 32 values encoded in some way.

enum class StreamType
{
    Primitive = 0,
    Position,
    NormalTangentOct8,
    UV,
    BoneIndices,
    BoneWeights,
};

Each meshlet has stream_count number of Stream headers. The indexing is trivial:

streams[RuntimeHeader::stream_offset + int(StreamType)]
// 16 bytes
struct Stream
{
    union
    {
       uint32_t base_value[2];
       struct { uint32_t prim_count; uint32_t vert_count; } counts;
    } u;
    uint32_t bits;
    uint32_t offset_in_words;
};

This is where things get a bit more interesting. I ended up supporting some encoding styles that are tailored for the various attribute formats.

Encoding attributes

There’s two parts to this problem. First is to decide on some N-bit fixed point values, and then find the most efficient way to pull those bits from a buffer. I went through several iterations on the actual bit-stuffing.

Base value + DELTA encoding

A base value is encoded in Stream::base_value, and the decoded bits are an offset from the base. To start approaching speed-of-light decoding, this is about as fancy as we can do it.

I went through various iterations of this model. The first idea had a predictive encoding between neighbor values, where subgroup scan operations were used to complete the decode, but it was too slow in practice, and didn’t really improve bit rates at all.

Index buffer

Since the sublet is just 32-wide, we can encode with 5-bit indices. 15 bit / primitive. There is no real reason to use delta encode here, so instead of storing base values in the stream header, I opted to use those bits to encode vertex/index counts.

Position

This is decoded to 3×16-bit SINT. The shared exponent is stored in top 16 bits of Stream::bits.

vec3 position = ldexp(vec3(i16vec3(decoded)), exponent);

This facilitates arbitrary quantization as well.

UV

Similar idea as position, but 2×16-bit SINT. After decoding similar to position, a simple fixup is made to cater to typical UVs which lie in range of [0, +1], not [-1, +1].

vec2 uv = 0.5 * ldexp(vec2(i16vec2(decoded)), exponent) + 0.5;
Normal/Tangent

Encoded as 4×8-bit SNORM. Normal (XY) and Tangent (ZW) are encoded with Octahedral encoding from meshoptimizer library.

To encode the sign of tangent, Stream::bits stores 2 bits, which signals one of three modes:

  • Uniform W = -1
  • Uniform W = +1
  • LSB of decoded W encodes tangent W. Tangent’s second component loses 1 bit of precision.
Bone index / weight

Basically same as Normal/Tangent, but ignore tangent sign handling.

First (failed?) idea – bitplane encoding

For a long time, I was pursuing bitplane encoding, which is one of the simplest ways to encode variable bitrates. We can encode 1 bit for 32 values by packing them in one u32. To speed up decoding further, I aimed to pack everything into 128-bit aligned loads. This avoids having to wait for tiny, dependent 32-bit loads.

For example, for index buffers:

uint meshlet_decode_index_buffer(
   uint stream_index, uint chunk_index,
   int lane_index)
{
    uint offset_in_b128 =
      meshlet_streams.data[stream_index].offset_in_b128;

    // Fixed 5-bit encoding.
    offset_in_b128 += 4 * chunk_index;

    // Scalar load. 64 bytes in one go.
    uvec4 p0 = payload.data[offset_in_b128 + 0];
    uvec4 p1 = payload.data[offset_in_b128 + 1];
    uvec4 p2 = payload.data[offset_in_b128 + 2];
    uvec4 p3 = payload.data[offset_in_b128 + 3];

    uint indices = 0;

    indices |= bitfieldExtract(p0.x, lane_index, 1) << 0u;
    indices |= bitfieldExtract(p0.y, lane_index, 1) << 1u;
    indices |= bitfieldExtract(p0.z, lane_index, 1) << 2u;
    indices |= bitfieldExtract(p0.w, lane_index, 1) << 3u;

    indices |= bitfieldExtract(p1.x, lane_index, 1) << 8u;
    indices |= bitfieldExtract(p1.y, lane_index, 1) << 9u;
    indices |= bitfieldExtract(p1.z, lane_index, 1) << 10u;
    indices |= bitfieldExtract(p1.w, lane_index, 1) << 11u;

    indices |= bitfieldExtract(p2.x, lane_index, 1) << 16u;
    indices |= bitfieldExtract(p2.y, lane_index, 1) << 17u;
    indices |= bitfieldExtract(p2.z, lane_index, 1) << 18u;
    indices |= bitfieldExtract(p2.w, lane_index, 1) << 19u;

    indices |= bitfieldExtract(p3.x, lane_index, 1) << 4u;
    indices |= bitfieldExtract(p3.y, lane_index, 1) << 12u;
    indices |= bitfieldExtract(p3.z, lane_index, 1) << 20u;

    return indices;
}

On Deck, this ends up looking like

s_buffer_load_dwordx4 x 4
v_bfe_u32 x 15
v_lshl_or_b32 x 15

Thinking about ALU and loads in terms of scalar and vectors can greatly help AMD performance when done right, so this approach felt natural.

For variable bit rates, I’d have code like:

if (bits & 4) { unroll_4_bits_bit_plane(); }
if (bits & 2) { unroll_2_bits_bit_plane(); }
if (bits & 1) { unroll_1_bit_bit_plane(); }

However, I abandoned this idea, since while favoring SMEM so heavily, the VALU with all the bitfield ops wasn’t exactly amazing for perf. I’m still just clocking out one bit per operation here. AMD performance was quite alright compared to what I ended up with in the end, but NVIDIA performance was abysmal, so I went back to the drawing board, and ended up with the absolute simplest solution that would work.

Tightly packed bits

This idea is to just literally pack bits together, clearly a revolutionary idea that noone has ever done before. A VMEM load or two per thread, then some shifts should be all that is needed to move the components into place.

E.g. for index buffers:

uvec3 meshlet_decode_index_buffer(uint stream_index,
   uint chunk_index,
   int lane_index)
{
  uint offset_in_words = 
    meshlet_streams.data[stream_index].offset_in_words;
  return meshlet_decode3(offset_in_words, lane_index, 5);
}

For the actual decode I figured it would be pretty fast if all the shifts could be done in 64-bit. At least AMD has native instructions for that.

uvec3 meshlet_decode3(uint offset_in_words,
   uint index,
   uint bit_count)
{
    const uint num_components = 3;
    uint start_bit = index * bit_count * num_components;
    uint start_word = offset_in_words + start_bit / 32u;
    start_bit &= 31u;
    uint word0 = payload.data[start_word];
    uint word1 = payload.data[start_word + 1u];
    uvec3 v;

    uint64_t word = packUint2x32(uvec2(word0, word1));
    v.x = uint(word >> start_bit);
    start_bit += bit_count;
    v.y = uint(word >> start_bit);
    start_bit += bit_count;
    v.z = uint(word >> start_bit);
    return bitfieldExtract(v, 0, int(bit_count));
}

There is one detail here. For 13, 14 and 15 bit components with uvec3 decode, more than two u32 words may be needed, so in this case, encoder must choose 16 bit. (16-bit works due to alignment.) This only comes up in position encode, and encoder can easily just ensure 12 bit deltas is enough to encode, quantizing a bit more as necessary.

Mapping to MDI

Every 256-wide meshlet can turn into an indexed draw call with VK_INDEX_TYPE_UINT8_EXT, which is nice for saving VRAM. The “task shader” becomes a compute shader that dumps out a big multi-draw indirect buffer. The DrawIndex builtin in Vulkan ends up replacing WorkGroupID in mesh shader for pulling in per-meshlet data.

Performance sanity check

Before going further with mesh shading fun, it’s important to validate performance. I needed at least a ballpark idea of how many primitives could be pumped through the GPU with a good old vkCmdDrawIndexed and the MDI method where one draw call is one meshlet. This was then to be compared against a straight forward mesh shader.

Zeux’s Niagara renderer helpfully has a simple OBJ for us to play with.

When exported to the new meshlet format it looks like:

[INFO]: Stream 0: 54332 bytes. (Index) 15 bits / primitive
[INFO]: Stream 1: 75060 bytes. (Position) ~25 bits / pos
[INFO]: Stream 2: 70668 bytes. (Normal/Tangent) ~23.8 bits / N + T + sign
[INFO]: Total encoded vertices: 23738 // Vertex duplication :(
[INFO]: Average radius 0.037 (908 bounds) // 32-wide meshlet
[INFO]: Average cutoff 0.253 (908 bounds)
[INFO]: Average radius 0.114 (114 bounds) // 256-wide meshlet
[INFO]: Average cutoff 0.697 (114 bounds)
// Backface cone culling isn't amazing for larger meshlets.
[INFO]: Exported meshlet:
[INFO]: 908 meshlets
[INFO]: 200060 payload bytes
[INFO]: 86832 total indices
[INFO]: 14856 total attributes
[INFO]: 703872 uncompressed bytes

One annoying thing about meshlets is attribute duplication when one vertex is reused across meshlets, and using tiny 32-wide meshlets makes this way worse. Add padding on top for encode and the compression ratio isn’t that amazing anymore. The primitive to vertex ratio is ~1.95 here which is really solid, but turning things into meshlets tends to converge to ~1.0.

I tried different sublet sizes, but NVIDIA performance collapsed when I didn’t use 32-wide sublets, and going to 64 primitive / 32 vertex only marginally helped P/V ratios. AMD runtime performance did not like that in my testing (~30% throughput loss), so 32/32 it is!

After writing this section, AMD released a blog post suggesting that the 2N/N structure is actually good, but I couldn’t replicate that in my testing at least and I don’t have the energy anymore to rewrite everything (again) to test that.

Test scene

The classic “instance the same mesh a million times” strategy. This was tested on RTX 3070 (AMD numbers to follow, there are way more permutations to test there …). The mesh is instanced in a 13x13x13 grid. Here we’re throwing 63.59 million triangles at the GPU in one go.

Spam vkCmdDrawIndexed with no culling

5.5 ms

layout(location = 0) in vec3 POS;
layout(location = 1) in mediump vec3 NORMAL;
layout(location = 2) in mediump vec4 TANGENT;
layout(location = 3) in vec2 UV;

layout(location = 0) out mediump vec3 vNormal;
layout(location = 1) out mediump vec4 vTangent;
layout(location = 2) out vec2 vUV;

// The most basic vertex shader.
void main()
{
  vec3 world_pos = (M * vec4(POS, 1.0)).xyz;
  vNormal = mat3(M) * NORMAL;
  vTangent = vec4(mat3(M) * TANGENT.xyz, TANGENT.w);
  vUV = UV;
  gl_Position = VP * vec4(world_pos, 1.0);
}

With per-object frustum culling

This is the most basic thing to do, so for reference.

4.3 ms

One massive MDI

Here we’re just doing basic frustum culling of meshlets as well as back-face cone culling and emitting one draw call per meshlet that passes test.

3.9 ms

Significantly more geometry is rejected now due to back-face cull and tighter frustum cull, but performance isn’t that much better. Once we start considering occlusion culling, this should turn into a major win over normal draw calls. In this path, we have a bit more indirection in the vertex shader, so that probably accounts for some loss as well.

void main()
{
    // Need to index now, but shouldn't be a problem on desktop hardware.
    mat4 M = transforms.data[draw_info.data[gl_DrawIDARB].node_offset];

    vec3 world_pos = (M * vec4(POS, 1.0)).xyz;
    vNormal = mat3(M) * NORMAL;
    vTangent = vec4(mat3(M) * TANGENT.xyz, TANGENT.w);
    vUV = UV;

    // Need to pass down extra data to sample materials, etc.
    // Fragment shader cannot read gl_DrawIDARB.
    vMaterialID = draw_info.data[gl_DrawIDARB].material_index;

    gl_Position = VP * vec4(world_pos, 1.0);
}

Meshlet – Encoded payload

Here, the meshlet will read directly from the encoded payload, and decode inline in the shader. No per-primitive culling is performed.

4.1 ms

Meshlet – Decoded payload

4.0 ms

We’re at the point where we are bound on fixed function throughput. Encoded and Decoded paths are basically both hitting the limit of how much data we can pump to the rasterizer.

Per-primitive culling

To actually make good use of mesh shading, we need to consider per-primitive culling. For this section, I’ll be assuming a subgroup size of 32, and a meshlet size of 32. There are other code paths for larger workgroups, which require some use of groupshared memory, but that’s not very exciting for this discussion.

The gist of this idea was implemented in https://gpuopen.com/geometryfx/. Various AMD drivers adopted the idea as well to perform magic driver culling, but the code here isn’t based on any other code in particular.

Doing back-face culling correctly

This is tricky, but we only need to be conservative, not exact. We can only reject when we know for sure the primitive is not visible.

Perspective divide and clip codes

The first step is to do W divide per vertex and study how that vertex clips against the X, Y, and W planes. We don’t really care about Z. Near-plane clip is covered by negative W tests, and far plane should be covered by simple frustum test, assuming we have a far plane at all.

vec2 c = clip_pos.xy / clip_pos.w;

uint clip_code = clip_pos.w <= 0.0 ? CLIP_CODE_NEGATIVE_W : 0;
if (any(greaterThan(abs(c), vec2(4.0))))
    clip_code |= CLIP_CODE_INACCURATE;
if (c.x <= -1.0)
    clip_code |= CLIP_CODE_NEGATIVE_X;
if (c.y <= -1.0)
    clip_code |= CLIP_CODE_NEGATIVE_Y;
if (c.x >= 1.0)
    clip_code |= CLIP_CODE_POSITIVE_X;
if (c.y >= 1.0)
    clip_code |= CLIP_CODE_POSITIVE_Y;

vec2 window = roundEven(c * viewport.zw + viewport.xy);

There are things to unpack here. The INACCURATE clip code is used to denote a problem where we might start to run into accuracy issues when converting to fixed point, or GPUs might start doing clipping due to guard band exhaustion. I picked the value arbitrarily.

The window coordinate is then computed by simulating the fixed point window coordinate snapping done by real GPUs. Any GPU supporting DirectX will have a very precise way of doing this, so this should be okay in practice. Vulkan also exposes the number of sub-pixel bits in the viewport transform. On all GPUs I know of, this is 8. DirectX mandates exactly 8.

vec4 viewport =
    float(1 << 8 /* shader assumes 8 */) *
        vec4(cmd->get_viewport().x +
               0.5f * cmd->get_viewport().width - 0.5f,
             cmd->get_viewport().y +
               0.5f * cmd->get_viewport().height - 0.5f,
             0.5f * cmd->get_viewport().width,
             0.5f * cmd->get_viewport().height) -
             vec4(1.0f, 1.0f, 0.0f, 0.0f);

This particular way of doing it comes into play later when discussing micro-poly rejection. One thing to note here is that Vulkan clip-to-window coordinate transform does not flip Y-sign. D3D does however, so beware.

Shuffle clip codes and window coordinates

void meshlet_emit_primitive(uvec3 prim, vec4 clip_pos, vec4 viewport)
{
  // ...
  vec2 window = roundEven(c * viewport.zw + viewport.xy);

  // vertex ID maps to gl_SubgroupInvocationID
  // Fall back to groupshared as necessary
  vec2 window_a = subgroupShuffle(window, prim.x);
  vec2 window_b = subgroupShuffle(window, prim.y);
  vec2 window_c = subgroupShuffle(window, prim.z);
  uint code_a = subgroupShuffle(clip_code, prim.x);
  uint code_b = subgroupShuffle(clip_code, prim.y);
  uint code_c = subgroupShuffle(clip_code, prim.z);
}

Early reject or accept

Based on clip codes we can immediately accept or reject primitives.

uint or_code = code_a | code_b | code_c;
uint and_code = code_a & code_b & code_c;
bool culled_planes = (and_code & CLIP_CODE_PLANES) != 0;
bool is_active_prim = false;

if (!culled_planes)
{
    is_active_prim =
        (or_code & (CLIP_CODE_INACCURATE |
                    CLIP_CODE_NEGATIVE_W)) != 0;

    if (!is_active_prim)
        is_active_prim = cull_triangle(window_a,
                                       window_b,
                                       window_c);
}
  • If all three vertices are outside one of the clip planes, reject immediately
  • If any vertex is considered inaccurate, accept immediately
  • If one or two of the vertices have negative W, we have clipping. Our math won’t work, so accept immediately. (If all three vertices have negative W, the first test rejects).
  • Perform actual back-face cull.

Actual back-face cull

bool cull_triangle(vec2 a, vec2 b, vec2 c)
{
  precise vec2 ab = b - a;
  precise vec2 ac = c - a;

  // This is 100% accurate as long as the primitive
  // is no larger than ~4k subpixels, i.e. 16x16 pixels.
  // Normally, we'd be able to do GEQ test, but GE test is conservative,
  // even with FP error in play.

  // Depending on your engine and API conventions, swap these two.
  precise float pos_area = ab.y * ac.x;
  precise float neg_area = ab.x * ac.y;

  // If the pos value is (-2^24, +2^24),
  // the FP math is exact,
  // if not, we have to be conservative.
  // Less-than check is there to ensure that 1.0 delta
  // in neg_area *will* resolve to a different value.
  bool active_primitive;
  if (abs(pos_area) < 16777216.0)
    active_primitive = pos_area > neg_area;
  else
    active_primitive = pos_area >= neg_area;

  return active_primitive;
}

To compute winding, we need a 2D cross product. While noodling with this code, I noticed that we can still do it in FP32 instead of full 64-bit integer math. We’re working with integer-rounded values here, so based on the magnitudes involved we can pick the exact GEQ test. If we risk FP rounding error, we can use GE test. If the results don’t test equal, we know for sure area must be negative, otherwise, it’s possible it could have been positive, but the intermediate values rounded to same value in the end.

3.3 ms

Culling primitives helped as expected. Less pressure on the fixed function units.

Micro-poly rejection

Given how pathologically geometry dense this scene is, we expect that most primitives never trigger the rasterizer at all.

If we can prove that the bounding box of the primitive lands between two pixel grids, we can reject it since it will never have coverage.

if (active_primitive)
{
    // Micropoly test.
    const int SUBPIXEL_BITS = 8;
    vec2 lo = floor(ldexp(min(min(a, b), c), ivec2(-SUBPIXEL_BITS)));
    vec2 hi = floor(ldexp(max(max(a, b), c), ivec2(-SUBPIXEL_BITS)));
    active_primitive = all(notEqual(lo, hi));
}

There is a lot to unpack in this code. If we re-examine the viewport transform:

vec4 viewport = float(1 << 8 /* shader assumes 8 */) *
  vec4(cmd->get_viewport().x +
        0.5f * cmd->get_viewport().width - 0.5f,
      cmd->get_viewport().y +
        0.5f * cmd->get_viewport().height - 0.5f,
      0.5f * cmd->get_viewport().width,
      0.5f * cmd->get_viewport().height) -
      vec4(1.0f, 1.0f, 0.0f, 0.0f);

First, we need to shift by 0.5 pixels. The rasterization test happens at the center of a pixel, and it’s more convenient to sample at integer points. Then, due to top-left rasterization rules on all desktop GPUs (a DirectX requirement), we shift the result by one sub-pixel. This ensures that should a primitive have a bounding box of [1.0, 2.0], we will consider it for rasterization, but [1.0 + 1.0 / 256.0, 2.0] will not. Top-left rules are not technically guaranteed in Vulkan however (it just has to have some rule), so if you’re paranoid, increase the upper bound by one sub-pixel.

1.9 ms

Now we’re only submitting 1.2 M primitives to the rasterizer, which is pretty cool, given that we started with 31 M potential primitives. Of course, this is a contrived example with ridiculous micro-poly issues.

We’re actually at the point here where reporting the invocation stats (one atomic per workgroup) becomes a performance problem, so turning that off:

1.65 ms

With inline decoding there’s some extra overhead, but we’re still well ahead:

2.5 ms

Build active vertex / primitive masks

This is quite straight forward. Once we have the counts, SetMeshOutputCounts is called and we can compute the packed output indices with a mask and popcount.

uint vert_mask = 0u;
if (is_active_prim)
    vert_mask = (1u << prim.x) | (1u << prim.y) | (1u << prim.z);

uvec4 prim_ballot = subgroupBallot(is_active_prim);

shared_active_prim_offset = subgroupBallotExclusiveBitCount(prim_ballot);
shared_active_vert_mask = subgroupOr(vert_mask);

shared_active_prim_count_total = subgroupBallotBitCount(prim_ballot);
shared_active_vert_count_total = bitCount(shared_active_vert_mask);

Special magic NVIDIA optimization

Can we improve things from here? On NVIDIA, yes. NVIDIA seems to under-dimension the shader export buffers in their hardware compared to peak triangle throughput, and their developer documentation on the topic suggests:

  • Replace attributes with barycentrics and allowing the Pixel Shader to fetch and interpolate the attributes

Using VK_KHR_fragment_shader_barycentrics we can write code like:

// Mesh output
layout(location = 0) flat out uint vVertexID[];
layout(location = 1) perprimitiveEXT out uint vTransformIndex[];

// Fragment
layout(location = 0) pervertexEXT in uint vVertexID[];
layout(location = 1) perprimitiveEXT flat in uint vTransformIndex;

// Fetch vertex IDs
uint va = vVertexID[0];
uint vb = vVertexID[1];
uint vc = vVertexID[2];

// Load attributes from memory directly
uint na = attr.data[va].n;
uint nb = attr.data[vb].n;
uint nc = attr.data[vc].n;

// Interpolate by hand
mediump vec3 normal = gl_BaryCoordEXT.x * decode_rgb10a2(na) +
    gl_BaryCoordEXT.y * decode_rgb10a2(nb) +
    gl_BaryCoordEXT.z * decode_rgb10a2(nc);

// Have to transform normals and tangents as necessary.
// Need to pass down some way to load transforms.
normal = mat3(transforms.data[vTransformIndex]) * normal;
normal = normalize(normal);

1.0 ms

Quite the dramatic gain! Nsight Graphics suggests we’re finally SM bound at this point (> 80% utilization), where we used to be ISBE bound (primitive / varying allocation). An alternative that I assume would work just as well is to pass down a primitive ID to a G-buffer similar to Nanite.

There are a lot of caveats with this approach however, and I don’t think I will pursue it:

  • Moves a ton of extra work to fragment stage
    • I’m not aiming for Nanite-style micro-poly hell here, so doing work per-vertex seems better than per-fragment
    • This result isn’t representative of a real scene where fragment shader load would be far more significant
  • Incompatible with encoded meshlet scheme
    • It is possible to decode individual values, but it sure is a lot of dependent memory loads to extract a single value
  • Very awkward to write shader code like this at scale
    • Probably need some kind of meta compiler that can generate code, but that’s a rabbit hole I’m not going down
    • Need fallbacks, barycentrics is a very modern feature
  • Makes skinning even more annoying
    • Loading multiple matrices with fully dynamic index in fragment shader does not scream performance, then combine that with having to compute motion vectors on top …
  • Only seems to help throughput on NVIDIA
  • We’re already way ahead of MDI anyway

Either way, this result was useful to observe.

AMD

Steam Deck

Before running the numbers, we have to consider that the RADV driver already does some mesh shader optimizations for us automatically. The NGG geometry pipeline automatically converts vertex shading workloads into pseudo-meshlets, and RADV also does primitive culling in the driver-generated shader.

To get the raw baseline, we’ll first consider the tests without that path, so we can see how well RADV’s own culling is doing. The legacy vertex path is completely gone on RDNA3 as far as I know, so these tests have to be done on RDNA2.

No culling, plain vkCmdDrawIndexed, RADV_DEBUG=nongg

Even locked to 1600 MHz (peak), GPU is still just consuming 5.5 W. We’re 100% bound on fixed function logic here, the shader cores are sleeping.

44.3 ms

Basic frustum culling

As expected, performance scales as we cull. Still 5.5 W. 27.9 ms

NGG path, no primitive culling, RADV_DEBUG=nonggc

Not too much changed in performance here. We’re still bound on the same fixed function units pumping invisible primitives through. 28.4 ms

Actual RADV path

When we don’t cripple RADV, we get a lot of benefit from driver culling. GPU hits 12.1 W now. 9.6 ms

MDI

Slight win. 8.9 ms

Forcing Wave32 in mesh shaders

Using Vulkan 1.3’s subgroup size control feature, we can force RDNA2 to execute in Wave32 mode. This requires support in

 VkShaderStageFlags requiredSubgroupSizeStages;

The Deck drivers and upstream Mesa ship support for requiredSize task/mesh shaders now which is very handy. AMD’s Windows drivers or AMDVLK/amdgpu-pro do not, however 🙁 It’s possible Wave32 isn’t the best idea for AMD mesh shaders in the first place, it’s just that the format favors Wave32, so I enable it if I can.

Testing various parameters

While NVIDIA really likes 32/32 (anything else I tried fell off the perf cliff), AMD should in theory favor larger workgroups. However, it’s not that easy in practice, as I found.

Decoded meshlet – Wave32 – N/N prim/vert

  • 32/32: 9.3 ms
  • 64/64: 10.5 ms
  • 128/128: 11.2 ms
  • 256/256: 12.8 ms

These results are … surprising.

Encoded meshlet – Wave32 N/N prim/vert

  • 32/32: 10.7 ms
  • 64/64: 11.8 ms
  • 128/128: 12.7 ms
  • 256/256: 14.7 ms

Apparently Deck (or RDNA2 in general) likes small meshlets?

Wave64?

No meaningful difference in performance on Deck.

VertexID passthrough?

No meaningful difference either. This is a very NVIDIA-centric optimization I think.

A note on LocalInvocation output

In Vulkan, there are some properties that AMD sets for mesh shaders.

VkBool32 prefersLocalInvocationVertexOutput;
VkBool32 prefersLocalInvocationPrimitiveOutput;

This means that we should write outputs using LocalInvocationIndex, which corresponds to how RDNA hardware works. Each thread can export one primitive and one vertex and the thread index corresponds to primitive index / vertex index. Due to culling and compaction, we will have to roundtrip through groupshared memory somehow to satisfy this.

For the encoded representation, I found that it’s actually faster to ignore this suggestion, but for the decoded representation, we can just send the vertex IDs through groupshared, and do split vertex / attribute shading. E.g.:

if (meshlet_lane_has_active_vert())
{
    uint out_vert_index = meshlet_compacted_vertex_output();
    uint vert_id = meshlet.vertex_offset + linear_index;

    shared_clip_pos[out_vert_index] = clip_pos;
    shared_attr_index[out_vert_index] = vert_id;
}

barrier();

if (gl_LocalInvocationIndex < shared_active_vert_count_total)
{
    TexturedAttr a =
      attr.data[shared_attr_index[gl_LocalInvocationIndex]];
    mediump vec3 n = unpack_bgr10a2(a.n).xyz;
    mediump vec4 t = unpack_bgr10a2(a.t);
    gl_MeshVerticesEXT[gl_LocalInvocationIndex].gl_Position =
      shared_clip_pos[gl_LocalInvocationIndex];
    vUV[gl_LocalInvocationIndex] = a.uv;
    vNormal[gl_LocalInvocationIndex] = mat3(M) * n;
    vTangent[gl_LocalInvocationIndex] = vec4(mat3(M) * t.xyz, t.w);
}

Only computing visible attributes is a very common optimization in GPUs in general and RADV’s NGG implementation does it roughly like this.

Either way, we’re not actually beating the driver-based meshlet culling on Deck. It’s more or less already doing this work for us. Given how close the results are, it’s possible we’re still bound on something that’s not raw compute. On the positive side, the cost of using encoded representation is very small here, and saving RAM for meshes is always nice.

Already, the permutation hell is starting to become a problem. It’s getting quite obvious why mesh shaders haven’t taken off yet 🙂

RX 7600 numbers

Data dump section incoming …

NGG culling seems obsolete now?

By default RADV disables NGG culling on RDNA3, because apparently it has a much stronger fixed function culling in hardware now. I tried forcing it on with RADV_DEBUG=nggc, but found no uplift in performance for normal vertex shaders. Curious. Here’s with no culling, where the shader is completely export bound.

But, force NGG on, and it still doesn’t help much. Culling path takes as much time as the other, the instruction latencies are just spread around more.

RADV

  • vkCmdDrawIndexed, no frustum culling: 5.9 ms
  • With frustum cull: 3.7 ms
  • MDI: 5.0 ms
Wave32 – Meshlet
  • Encoded – 32/32: 3.3 ms
  • Encoded – 64/64 : 2.5 ms
  • Encoded – 128/128: 2.7 ms
  • Encoded – 256/256: 2.9 ms
  • Decoded – 32/32: 3.3 ms
  • Decoded – 64/64: 2.4 ms
  • Decoded – 128/128: 2.6 ms
  • Decoded – 256/256: 2.7 ms
Wave64 – Meshlet
  • Encoded – 64/64: 2.4 ms
  • Encoded – 128/128: 2.6 ms
  • Encoded – 256/256: 2.7 ms
  • Decoded – 64/64: 2.2 ms
  • Decoded – 128/128: 2.5 ms
  • Decoded – 256/256: 2.7 ms

Wave64 mode is doing quite well here. From what I understand, RADV hasn’t fully taken advantage of the dual-issue instructions in RDNA3 specifically yet, which is important for Wave32 performance, so that might be a plausible explanation.

There was also no meaningful difference in doing VertexID passthrough.

It’s not exactly easy to deduce anything meaningful out of these numbers, other than 32/32 being bad on RDNA3, while good on RDNA2 (Deck)?

AMD doesn’t seem to like the smaller 256 primitive draws on the larger desktop GPUs. I tried 512 and 1024 as a quick test and that improved throughput considerably, still, with finer grained culling in place, it should be a significant win.

amdgpu-pro / proprietary (Linux)

Since we cannot request specific subgroup size, the driver is free to pick Wave32 or Wave64 as it pleases, so I cannot test the difference. It won’t hit the subgroup optimized paths however.

  • vkCmdDrawIndexed, no culling : 6.2 ms
  • With frustum cull: 4.0 ms
  • MDI: 5.3 ms
  • Meshlet – Encoded – 32/32: 2.5 ms
  • Meshlet – Encoded – 64/64 : 2.6 ms
  • Meshlet – Encoded – 128/128: 2.7 ms
  • Meshlet – Encoded – 256/256: 2.6 ms
  • Meshlet – Decoded – 32/32: 2.1 ms
  • Meshlet – Decoded – 64/64: 2.1 ms
  • Meshlet – Decoded – 128/128: 2.1 ms
  • Meshlet – Decoded – 256/256: 2.1 ms

I also did some quick spot checks on AMDVLK, and the numbers are very similar.

The proprietary driver is doing quite well here in mesh shaders. On desktop, we can get significant wins on both RADV and proprietary with mesh shaders, which is nice to see.

It seems like the AMD Windows driver skipped NGG culling on RDNA3 as well. Performance is basically the same.

Task shader woes

The job of task shaders is to generate mesh shader work on the fly. In principle this is nicer than indirect rendering with mesh shaders for two reasons:

  • No need to allocate temporary memory to hold for indirect draw
  • No need to add extra compute passes with barriers

However, it turns out that this shader stage is even more vendor specific when it comes to tuning for performance. So far, no game I know of has actually shipped with task shaders (or the D3D12 equivalent amplification shader), and I think I now understand why.

The basic task unit I settled on was:

struct TaskInfo
{
    uint32_t aabb_instance;  // AABB, for top-level culling
    uint32_t node_instance;  // Affine transform
    uint32_t material_index; // To be eventually forwarded to fragment
    uint32_t mesh_index_count;
    // Encodes count [1, 32] in lower bits.
    // Mesh index is aligned to 32.
    uint32_t occluder_state_offset;
    // For two-phase occlusion state (for later)
};

An array of these is prepared on CPU. Each scene entity translates to one or more TaskInfos. Those are batched up into one big buffer, and off we go.

The logical task shader for me was to have N = 32 threads which tests AABB of N tasks in parallel. For the tasks that pass the test, test 32 meshlets in parallel. This makes it so the task workgroup can emit up to 1024 meshlets.

When I tried this on NVIDIA however …

18.8 ms

10x slowdown … The NVIDIA docs do mention that large outputs are bad, but I didn’t expect it to be this bad:

Avoid large outputs from the amplification shader, as this can incur a significant performance penalty. Generally, we encourage a flexible implementation that allows for fine-tuning. With that in mind, there are a number of generic factors that impact performance:

  • Size of the payloads. The AS payload should preferably stay below 108 bytes, but if that is not possible, then keep it at least under 236 bytes.

If we remove all support for hierarchical culling, the task shader runs alright again. 1 thread emits 0 or 1 meshlet. However, this means a lot of threads dedicated to culling, but it’s similar in performance to plain indirect mesh shading.

AMD however, is a completely different story. Task shaders are implemented by essentially emitting a bunch of tiny indirect mesh shader dispatches anyway, so the usefulness of task shaders on AMD is questionable from a performance point of view. While writing this blog, AMD released a new blog on the topic, how convenient!

When I tried NV-style task shader on AMD, performance suffered quite a lot.

However, the only thing that gets us to max perf on both AMD and NV is to forget about task shaders and go with vkCmdDrawMeshTasksIndirectCountEXT instead. While the optimal task shader path for each vendor gets close to indirect mesh shading, having a universal fast path is good for my sanity. The task shader loss was about 10% for me even in ideal situations on both vendors, which isn’t great. As rejection ratios increase, this loss grows even more. This kind of occupancy looks way better 🙂

The reason for using multi-indirect-count is to deal with the limitation that we can only submit about 64k workgroups in any dimension, similar to compute. This makes 1D atomic increments awkward, since we’ll easily blow past the 64k limit. One alternative is to take another tiny compute pass that prepares a multi-indirect draw, but that’s not really needed. Compute shader code like this works too:

// global_offset = atomicAdd() in thread 0

if (gl_LocalInvocationIndex == 0 && draw_count != 0)
{
  uint max_global_offset = global_offset + draw_count - 1;
  // Meshlet style.
  // Only guaranteed to get 0xffff meshlets,
  // so use 32k as cutoff for easy math.
  // Allocate the 2D draws in-place, avoiding an extra barrier.
  uint multi_draw_index = max_global_offset / 0x8000u;
  uint local_draw_index = max_global_offset & 0x7fffu;
  const int INC_OFFSET = NUM_CHUNK_WORKGROUPS == 1 ? 0 : 1;
  atomicMax(output_draws.count[1], multi_draw_index + 1);
  atomicMax(output_draws.count[
    2 + 3 * multi_draw_index + INC_OFFSET],
    local_draw_index + 1);

  if (local_draw_index <= draw_count)
  {
    // This is the thread that takes us over the threshold.
    output_draws.count[
      2 + 3 * multi_draw_index + 1 - INC_OFFSET] =
      NUM_CHUNK_WORKGROUPS;
    output_draws.count[2 + 3 * multi_draw_index + 2] = 1;
  }

  // Wrapped around, make sure last bucket sees 32k meshlets.
  if (multi_draw_index != 0 && local_draw_index < draw_count)
  {
    atomicMax(output_draws.count[
      2 + 3 * (multi_draw_index - 1) +
      INC_OFFSET], 0x8000u);
  }
}

This prepares a bunch of (8, 32k, 1) dispatches that are processed in one go. No chance to observe a bunch of dead dispatches back-to-back like task shaders can cause. In the mesh shader, we can use DrawIndex to offset the WorkGroupID by the appropriate amount (yay, Vulkan). A dispatchX count of 8 is to shade the full 256-wide meshlet through 8x 32-wide workgroups. As the workgroup size increases to handle more sublets per group, dispatchX count decreases similarly.

Occlusion culling

To complete the meshlet renderer, we need to consider occlusion culling. The go-to technique for this these days is two-phase occlusion culling with HiZ depth buffer. Some references:

Basic gist is to keep track of which meshlets are considered visible. This requires persistent storage of 1 bit per unit of visibility. Each pass in the renderer needs to keep track of its own bit-array. E.g. shadow passes have different visibility compared to main scene rendering.

For Granite, I went with an approach where 1 TaskInfo points to one uint32_t bitmask. Each of the 32 meshlets within the TaskInfo gets 1 bit. This makes the hierarchical culling nice too, since we can just test for visibility != 0 on the entire word. Nifty!

First phase

Here we render all objects which were considered visible last frame. It’s extremely likely that whatever was visible last frame is visible this frame, unless there was a full camera cut or similar. It’s important that we’re actually rendering to the framebuffer now. In theory, we’d be done rendering now if there were no changes to camera or objects in the scene.

HiZ pass

Based on the objects we drew in phase 1, build a HiZ depth map. This topic is actually kinda tricky. Building the mip-chain in one pass is great for performance, but causes some problems. With NPOT textures and single pass, there is no obvious way to create a functional HiZ, and the go-to shader for this, FidelityFX SPD, doesn’t support that use case.

The problem is that the size of mip-maps round down, so if we have a 7×7 texture, LOD 1 is 3×3 and LOD 2 is 1×1. In LOD2, we will be able to query a 4×4 depth region, but the edge pixels are forgotten.

The “obvious” workaround is to pad the texture to POT, but that is a horrible waste of VRAM. The solution I went with instead was to fold in the neighbors as the mips are reduced. This makes it so that the edge pixels in each LOD also remembers depth information for pixels which were truncated away due to NPOT rounding.

I rolled a custom HiZ shader similar to SPD with some extra subgroup shenanigans because why not (SubgroupShuffleXor with 4 and 8).

Second phase

In this pass we submit for rendering any object which became visible this frame, i.e. the visibility bit was not set, but it passed occlusion test now. Again, if camera did not change, and objects did not move, then nothing should be rendered here.

However, we still have to test every object, in order to update the visibility buffer for next frame. We don’t want visibility to remain sticky, unless we have dedicated proxy geometry to serve as occluders (might still be a thing if game needs to handle camera cuts without large jumps in rendering time).

In this pass we can cull meshlet bounds against the HiZ.

Because I cannot be arsed to make a fancy SVG for this, the math to compute a tight AABB bound for a sphere is straight forward once the geometry is understood.

The gist is to figure out the angle, then rotate the (X, W) vector with positive and negative angles. X / W becomes the projected lower or upper bound. Y bounds are computed separately.

vec2 project_sphere_flat(float view_xy, float view_z, float radius)
{
    float len = length(vec2(view_xy, view_z));
    float sin_xy = radius / len;

    float cos_xy = sqrt(1.0 - sin_xy * sin_xy);
    vec2 rot_lo = mat2(cos_xy, sin_xy, -sin_xy, cos_xy) *
      vec2(view_xy, view_z);
    vec2 rot_hi = mat2(cos_xy, -sin_xy, +sin_xy, cos_xy) *
      vec2(view_xy, view_z);

    return vec2(rot_lo.x / rot_lo.y, rot_hi.x / rot_hi.y);
}

The math is done in view space where the sphere is still a sphere, which is then scaled to window coordinates afterwards. To make the math easier to work with, I use a modified view space in this code where +Y is down and +Z is in view direction.

bool hiz_cull(vec2 view_range_x, vec2 view_range_y, float closest_z)
// view_range_x: .x -> lower bound, .y -> upper bound
// view_range_y: same
// closest_z: linear depth. ViewZ - Radius for a sphere

First, convert to integer coordinates.

// Viewport scale first applies any projection scale in X/Y
// (without Y flip).
// The scale also does viewport size / 2 and then
// offsets into integer window coordinates.

vec2 range_x = view_range_x *
  frustum.viewport_scale_bias.x +
  frustum.viewport_scale_bias.z;
vec2 range_y = view_range_y *
  frustum.viewport_scale_bias.y +
  frustum.viewport_scale_bias.w;

ivec2 ix = ivec2(range_x);
ivec2 iy = ivec2(range_y);

ix.x = clamp(ix.x, 0, frustum.hiz_resolution.x - 1);
ix.y = clamp(ix.y, ix.x, frustum.hiz_resolution.x - 1);
iy.x = clamp(iy.x, 0, frustum.hiz_resolution.y - 1);
iy.y = clamp(iy.y, iy.x, frustum.hiz_resolution.y - 1);

Figure out a LOD where we only have to sample a 2×2 footprint. findMSB to the rescue.

int max_delta = max(ix.y - ix.x, iy.y - iy.x);
int lod = min(findMSB(max_delta - 1) + 1, frustum.hiz_max_lod);
ivec2 lod_max_coord = max(frustum.hiz_resolution >> lod, ivec2(1)) - 1;

// Clamp to size of the actual LOD.
ix = min(ix >> lod, lod_max_coord.xx);
iy = min(iy >> lod, lod_max_coord.yy);

And finally, sample:

ivec2 hiz_coord = ivec2(ix.x, iy.x);

float d = texelFetch(uHiZDepth, hiz_coord, lod).x;
bool nx = ix.y != ix.x;
bool ny = iy.y != iy.x;

if (nx)
    d = max(d, texelFetchOffset(uHiZDepth,
      hiz_coord, lod,
      ivec2(1, 0)).x);

if (ny)
    d = max(d, texelFetchOffset(uHiZDepth,
      hiz_coord, lod,
      ivec2(0, 1)).x);

if (nx && ny)
    d = max(d, texelFetchOffset(uHiZDepth,
      hiz_coord, lod,
      ivec2(1, 1)).x);

return closest_z < d;

Trying to get up-close, it’s quite effective.

Without culling:

With two-phase:

As the culling becomes more extreme, GPU go brrrrr. Mostly just bound on HiZ pass and culling passes now which can probably be tuned a lot more.

Conclusion

I’ve spent way too much time on this now, and I just need to stop endlessly tuning various parameters. This is the true curse of mesh shaders, there’s always something to tweak. Given the performance I’m getting, I can call this a success, even if there might be some wins left on the table by tweaking some more. Now I just need to take a long break from mesh shaders before I actually rewrite the renderer to use this new code … And maybe one day I can even think about how to deal with LODs, then I would truly have Nanite at home!

The “compression” format ended up being something that can barely be called a compression format. To chase decode performance of tens of billions of primitives per second through, I suppose that’s just how it is.

Vulkan video shenanigans – FFmpeg + RADV integration experiments

Vulkan video is finally here and it’s a fierce battle to get things working fully. The leaders of the pack right now with the full release is RADV (Dave Airlie) and FFmpeg (Lynne).

In Granite, I’ve been wanting a solid GPU video decoding solution and I figured I’d work on a Vulkan video implementation over the holidays to try helping iron out any kinks with real-world application integration. The goal was achieving everything a 3D engine could potentially want out of video decode.

  • Hardware accelerated
  • GPU decode to RGB without round-trip through system memory (with optional mip generation when placed in a 3D world)
  • Audio decode
  • A/V synchronization

This blog is mostly here to demonstrate the progress in FFmpeg + RADV. I made a neat little sample app that fully uses Vulkan video to do a simple Sponza cinema. It supports A/V sync and seeking, which covers most of what a real media player would need. Ideally, this can be used as a test bench.

Place a video feed as a 3D object inside Sponza, why not?

Introduction blog post – read this first

This blog post by Lynne summarizes the state of Vulkan video at the time it was written. Note that none of this is merged upstream as of writing and APIs are changing rapidly.

Building FFmpeg + RADV + Granite

FFmpeg

Make sure to install the very latest Vulkan headers. On Arch Linux, install vulkan-headers-git from AUR for example.

Check out the branch in the blog and build. Make sure to install it in some throwaway prefix, e.g.

./configure --disable-doc --disable-shared --enable-static --disable-ffplay --disable-ffprobe --enable-vulkan --prefix=$HOME/ffmpeg-vulkan

Mesa

Check out https://gitlab.freedesktop.org/airlied/mesa/-/commits/radv-vulkan-video-decode. Then build with:

mkdir build
cd build
meson setup .. -Dvideo-codecs=h264dec,h265dec --buildtype release
ninja

Granite

git clone https://github.com/Themaister/Granite
cd Granite
git submodule update --init
mkdir build
cd build
cmake .. -DCMAKE_BUILD_TYPE=Release -DGRANITE_FFMPEG=ON -DGRANITE_AUDIO=ON -DGRANITE_FFMPEG_VULKAN=ON -G Ninja -DCMAKE_PREFIX_PATH=$HOME/ffmpeg-vulkan
ninja video-player

Running test app

Basic operation, a weird video player where the image is a flat 3D object floating in space. For fun the video is also mip-mapped and the plane is anisotropically filtered, because why not.

RADV_PERFTEST=video_decode GRANITE_FFMPEG_VULKAN=1 ./tests/video-player /tmp/test.mkv

Controls
  • WASD: move camera
  • Arrow keys: rotate camera
  • Space: Toggle pause
  • HJKL: Vim style for seeking

If you have https://github.com/KhronosGroup/glTF-Sample-Models checked out you can add a glTF scene as well for fun. I hacked it together with Sponza in mind, so:

RADV_PERFTEST=video_decode GRANITE_FFMPEG_VULKAN=1 ./tests/video-player $HOME/git/glTF-Sample-Models/2.0/Sponza/glTF/Sponza.gltf /tmp/test.mkv

and then you get the screenshot above with whatever video you’re using 🙂

Integration API

The Granite implementation can be found in https://github.com/Themaister/Granite/blob/master/video/ffmpeg_decode.cpp. It will probably be different in the final upstreamed version, so beware. I’m not an FFmpeg developer either FWIW, so take this implementation with a few grains of salt.

To integrate with Vulkan video, there are some steps we need to take. This assumes some familiarity with FFmpeg APIs. This is mostly interesting for non-FFmpeg developers. I had to figure this out with help from Lynne, spelunking in mpv and looking over the hardware decode samples in FFmpeg upstream.

Creating shared device

Before opening the decode context with:

avcodec_open2(ctx, codec, nullptr)

we will provide libavcodec with a hardware device context.

avcodec_get_hw_config(codec, index)

to scan through until you find a Vulkan configuration.

AVBufferRef *hw_dev = av_hwdevice_ctx_alloc(config->device_type);
auto *hwctx = reinterpret_cast<AVHWDeviceContext *>(hw_dev->data);
auto *vk = static_cast<AVVulkanDeviceContext *>(hwctx->hwctx);

hwctx->user_opaque = this; // For callbacks later.

To interoperate with FFmpeg, we have to provide it our own Vulkan device and lots of information about how we created the device.

vk->get_proc_addr = Vulkan::Context::get_instance_proc_addr();
vk->inst = device->get_instance();
vk->act_dev = device->get_device();
vk->phys_dev = device->get_physical_device();
vk->device_features = *device->get_device_features().pdf2;
vk->enabled_inst_extensions =
  device->get_device_features().instance_extensions;
vk->nb_enabled_inst_extensions =
  int(device->get_device_features().num_instance_extensions);
vk->enabled_dev_extensions =
  device->get_device_features().device_extensions;
vk->nb_enabled_dev_extensions =
  int(device->get_device_features().num_device_extensions);

Fortunately, I had most of this query scaffolding in place for Fossilize integration already. Vulkan 1.3 core is required here as well, so I had to bump that too when Vulkan video is enabled.

auto &q = device->get_queue_info();

vk->queue_family_index =
  int(q.family_indices[Vulkan::QUEUE_INDEX_GRAPHICS]);
vk->queue_family_comp_index =
  int(q.family_indices[Vulkan::QUEUE_INDEX_COMPUTE]);
vk->queue_family_tx_index =
  int(q.family_indices[Vulkan::QUEUE_INDEX_TRANSFER]);
vk->queue_family_decode_index =
  int(q.family_indices[Vulkan::QUEUE_INDEX_VIDEO_DECODE]);

vk->nb_graphics_queues = int(q.counts[Vulkan::QUEUE_INDEX_GRAPHICS]);
vk->nb_comp_queues = int(q.counts[Vulkan::QUEUE_INDEX_COMPUTE]);
vk->nb_tx_queues = int(q.counts[Vulkan::QUEUE_INDEX_TRANSFER]);
vk->nb_decode_queues = int(q.counts[Vulkan::QUEUE_INDEX_VIDEO_DECODE]);

vk->queue_family_encode_index = -1;
vk->nb_encode_queues = 0;

We need to let FFmpeg know about how it can query queues. Close match with Granite, but I had to add some extra APIs to make this work.

We also need a way to lock Vulkan queues:

vk->lock_queue = [](AVHWDeviceContext *ctx, int, int) {
   auto *self = static_cast<Impl *>(ctx->user_opaque);
   self->device->external_queue_lock();
};

vk->unlock_queue = [](AVHWDeviceContext *ctx, int, int) {
   auto *self = static_cast<Impl *>(ctx->user_opaque);
   self->device->external_queue_unlock();
};

For integration purposes, not making vkQueueSubmit internally synchronized in Vulkan was a mistake I think, oh well.

Once we’ve created a hardware context, we can let the codec context borrow it:

hw.device = av_hwdevice_ctx_init(hw_dev); // Unref later.
ctx->hw_device_ctx = av_buffer_ref(hw.device);

We also have to override get_format() and return the hardware pixel format.

ctx->opaque = this;
ctx->get_format = [](
    AVCodecContext *ctx,
    const enum AVPixelFormat *pix_fmts) -> AVPixelFormat {
  auto *self = static_cast<Impl *>(ctx->opaque);
  while (*pix_fmts != AV_PIX_FMT_NONE)
  {
    if (*pix_fmts == self->hw.config->pix_fmt)
      return *pix_fmts;
    pix_fmts++;
  }

  return AV_PIX_FMT_NONE;
};

This will work, but we’re also supposed to create a frames context before returning from get_format(). This also lets us configure how Vulkan images are created.

int ret = avcodec_get_hw_frames_parameters(
      ctx, ctx->hw_device_ctx,
      AV_PIX_FMT_VULKAN, &ctx->hw_frames_ctx);
// Check error.

auto *frames =
  reinterpret_cast<AVHWFramesContext *>(ctx->hw_frames_ctx->data);
auto *vk = static_cast<AVVulkanFramesContext *>(frames->hwctx);

vk->img_flags |= VK_IMAGE_CREATE_MUTABLE_FORMAT_BIT;

ret = av_hwframe_ctx_init(ctx->hw_frames_ctx);
// Check error.

The primary motivation for overriding image creation was that I wanted to do YCbCr to RGB conversion in a more unified way, i.e. using individual planes. That would be compatible with non-Vulkan video as well, but taking plane views of an image requires VK_IMAGE_CREATE_MUTABLE_FORMAT_BIT.

Using per-plane views is important, as we’ll see later. YCbCr samplers fall flat when dealing with practical video use cases.

Processing AVFrames

In FFmpeg, decoding works by sending AVPackets to a codec and it spits out AVFrame objects. If these frames are emitted by a software codec, we just poke at AVFrame::data[] directly, but with hardware decoders, AVFrame::pix_fmt is an opaque type.

There are two ways we can deal with this. For non-Vulkan hardware decoders, just read-back and upload planes to a VkBuffer staging buffer later, ewwww.

AVFrame *sw_frame = av_frame_alloc();

if (av_hwframe_transfer_data(sw_frame, av_frame, 0) < 0)
{
   LOGE("Failed to transfer HW frame.\n");
   av_frame_free(&sw_frame);
   av_frame_free(&av_frame);
}
else
{
   sw_frame->pts = av_frame->pts;
   av_frame_free(&av_frame);
   av_frame = sw_frame;
}

Each hardware pixel format lets you reinterpret AVFrame::data[] in a “magical” way if you’re willing to poke into low-level data structures. For VAAPI, VDPAU and APIs like that there are ways to use buffer sharing somehow, but the details are extremely hairy and is best left to experts. For Vulkan, we don’t even need external memory!

First, we need to extract the decode format:

auto *frames =
  reinterpret_cast<AVHWFramesContext *>(ctx->hw_frames_ctx->data);
active_upload_pix_fmt = frames->sw_format;

Then we can query the VkFormat if we want to stay multi-plane.

auto *hwdev =
  reinterpret_cast<AVHWDeviceContext *>(hw.device->data);
const VkFormat *fmts = nullptr;
VkImageAspectFlags aspects;
VkImageUsageFlags usage;
int nb_images;

int ret = av_vkfmt_from_pixfmt2(hwdev, active_upload_pix_fmt,
                                VK_IMAGE_USAGE_SAMPLED_BIT, &fmts,
                                &nb_images, &aspects, &usage);

However, this has some pitfalls in practice. Video frames tend to be aligned to a macro-block size or similar, meaning that the VkImage dimension might not be equal to the actual size we’re supposed to display. Even 1080p falls in this category for example since 1080 does not cleanly divide into 16×16 macro blocks. The only way to resolve this without extra copies is to view planes separately with VK_IMAGE_ASPECT_PLANE_n_BIT and do texture coordinate clamping manually. This way we avoid sampling garbage when converting to RGB. av_vkfmt_from_pixfmt can help here to deduce the per-plane Vulkan formats, but I just did it manually either way.

// Real output size.
ubo.resolution = uvec2(video.av_ctx->width, video.av_ctx->height);

if (video.av_ctx->hw_frames_ctx && hw.config &&
    hw.config->device_type == AV_HWDEVICE_TYPE_VULKAN)
{
   // Frames (VkImages) may be padded.
   auto *frames = reinterpret_cast<AVHWFramesContext *>(
       video.av_ctx->hw_frames_ctx->data);
   ubo.inv_resolution = vec2(
       1.0f / float(frames->width),
       1.0f / float(frames->height));
}
else
{
   ubo.inv_resolution = vec2(1.0f / float(video.av_ctx->width),
                             1.0f / float(video.av_ctx->height));
}

// Have to emulate CLAMP_TO_EDGE to avoid filtering against garbage.
ubo.chroma_clamp =
  (vec2(ubo.resolution) - 0.5f * float(1u << plane_subsample_log2[1])) *
  ubo.inv_resolution;

Processing the frame itself starts with magic casts:

auto *frames =
  reinterpret_cast<AVHWFramesContext *>(ctx->hw_frames_ctx->data);
auto *vk = static_cast<AVVulkanFramesContext *>(frames->hwctx);
auto *vk_frame = reinterpret_cast<AVVkFrame *>(av_frame->data[0]);

We have to lock the frame while accessing it, FFmpeg is threaded.

vk->lock_frame(frames, vk_frame);
// Do stuff
vk->unlock_frame(frames, vk_frame);

Now, we have to wait on the timeline semaphore (note that Vulkan 1.3 is required, so this is guaranteed to be supported).

// Acquire the image from FFmpeg.
if (vk_frame->sem[0] != VK_NULL_HANDLE && vk_frame->sem_value[0])
{
   // vkQueueSubmit(wait = sem[0], value = sem_value[0])
}

Create a VkImageView from the provided image. Based on av_vkfmt_from_pixfmt2 or per-plane formats from earlier, we know the appropriate Vulkan format to use when creating a view.

Queue family ownership transfer is not needed. FFmpeg uses CONCURRENT for sake of our sanity.

Transition the layout:

cmd->image_barrier(
    *wrapped_image,
    vk_frame->layout[0],
    VK_IMAGE_LAYOUT_SHADER_READ_ONLY_OPTIMAL,
    VK_PIPELINE_STAGE_COMPUTE_SHADER_BIT /* sem wait stage */, 0,
    VK_PIPELINE_STAGE_COMPUTE_SHADER_BIT,
    VK_ACCESS_2_SHADER_SAMPLED_READ_BIT);

vk_frame->layout[0] = VK_IMAGE_LAYOUT_SHADER_READ_ONLY_OPTIMAL;

Now, we can convert this to RGB as we desire. I went with an async compute formulation. If this were a pure video player we could probably blit this directly to screen with some fancy scaling filters.

When we’re done, we have to “release” the image back to FFmpeg.

// Release the image back to FFmpeg.
if (vk_frame->sem[0] != VK_NULL_HANDLE)
{
   vk_frame->sem_value[0] += 1;
   // vkQueueSubmit(signal = sem[0], value = sem_value[0]);
}

And that’s it!

Test results

I tried various codec configurations to see state of things.

RADV
  • H.264 – 8bit: Works
  • H.264 – 10bit: Not supported by hardware
  • H.265 – 8bit: Works
  • H.265 – 10bit: Works
nvidia
  • H.264: Broken
  • H.265: Seems to work
ANV

There’s a preliminary branch by Airlie again, but it doesn’t seem to have been updated for final spec yet.

Conclusion

Exciting times for Vulkan video. The API is ridiculously low level and way too complicated for mere graphics programming mortals, which is why having first class support in FFmpeg and friends will be so important to make the API usable.

Compressed GPU texture formats – a review and compute shader decoders – part 1

Compressed texture formats is one of the esoteric aspects of graphics programming almost no one cares all that much about. Neither did I, however, I’ve recently taken an academic interest in the zoo of compressed texture formats.

During development in Granite, I occasionally find it useful to test scenes which target mobile on desktop and vice versa, and in Vulkan, where there are no fallback paths for unsupported compression formats, we gotta roll our own decompression.

While it really isn’t all that useful to write a decoder for these formats, my goal is to create a suite of reasonably understandable compute shader kernels which can decode all of the standard formats I care about. Of course, I could just use a Frankenstein decoder which merges together a lot of C reference decoders and call it a day, but that’s not aesthetically pleasing or interesting to me. By implementing these formats straight from the Khronos Data Format specification, I learned a lot of things I would not otherwise know about these formats.

There are several major families of formats we can consider multi-vendor and standardized. Each of them fill their own niche. Unfortunately, desktop and mobile each have their own timelines with different texture compression standards, which is not fully resolved to this day in GPU hardware. (Basis Universal is something I will need to study eventually as well as it aims to solve this problem in software.)

By implementing all these formats, I got to see the evolution of block compression formats, see the major differences and design decisions that went into each format.

The major format families

First, it is useful to summarize all the families of texture compression I’ve looked at.

S3TC / DXT

The simplest family of formats. These formats are also known as the “BC” formats in Vulkan, or rather, BC 1, 2 and 3. This is the granddad of texture compression, similar to how I view MPEG1 in the video compression world.

These formats are firmly rooted in desktop GPUs. They are basically non-existent on mobile GPUs, probably for historical patent reasons.

RGTC

A very close relative of S3TC. These formats are very simple formats which specialize in encoding 1 and 2 uncorrelated channels, perfect for normal maps, metallic-roughness maps, etc. It is somewhat questionable to call these a separate family of formats (the Data Format specification separates them), since the basic format is basically exactly equal to the alpha format of S3TC, except that it extends the format to also support SNORM (-1, 1 range) alongside UNORM. These formats represent BC4 and BC5 in Vulkan.

These formats are firmly rooted in desktop GPUs. They are basically non-existent on mobile GPUs.

ETC

The ETC family of formats is very similarly laid out to S3TC in how different texture types are supported, but the implementation detail is quite different (and ETC2 is quite the interesting format). To support encoding full depth alpha and 1/2-component textures, there is the EAC format, which mirrors the RGTC formats.

These formats are firmly rooted in mobile GPUs. ETC1 was originally the only mandated format for OpenGLES 2.0 implementations, and ETC2 was mandated for OpenGLES 3.0 GPUs. It has almost no support on desktop GPUs. Intel iGPU is an exception here.

BPTC

This is where complexity starts to explode and where things get interesting. BC6 and BC7 are designed to compress high quality color images at 8bpp. BC6 adds support for HDR, which is to this day, one of only two ways to compress HDR images.

On desktop, BPTC is the state of the art in texture compression and was introduced around 2010.

ASTC

ASTC is the final boss of texture compression, and is the current state of the art in texture compression. Its complexity is staggering and aims to dominate the world with 128 bits. Mere mortals are not supposed to understand this format.

ASTC’s roots are on Mali GPUs, but it was always a Khronos standard, and is widely supported now on mobile Vulkan implementation (and Intel iGPU :3), at least the LDR profile. What you say, profiles in a texture compression format? Yes … yes, this is just the beginning of the madness that is ASTC.

PVRTC?

PVRTC is a PowerVR-exclusive format that has had some staying power due to iOS and I will likely ignore it in this series. However, it seems like a very different kind of format to all the others and studying it might be interesting. However, there is zero reason to use this format in Granite, and I don’t want to chew over too much.

What is a texture compression format anyway?

In a texture compression format, the specification describes a process for taking random bits given to it, and how to decode the bit-soup into texels. There are fundamental constraints in texture compression which is unique to this problem domain, and these restrictions heavily influence the design of the formats in question.

Fixed block size

To be able to randomly access any texel in a texture, there must be an O(1) mapping from texture coordinate to memory address. The only reasonable way to do this is to have a fixed block size. In all formats, 4×4 is the most common one. (As you can guess, ASTC can do odd-ball block sizes like 6×5).

Similarly, for reasons of random access, the number of bits spent per block must be constant. The typical block sizes are 64-bits and 128-bits, which is 4bpp and 8bpp respectively at 4×4 block size.

Image and video compression has none of these restrictions. That is a major reason why image and video compression is so much more efficient.

A set of coding tools

Each format has certain things it can do. The more complex the operations the format can do, the more expensive the decoding hardware becomes (and complex a software decoder becomes), so there’s always a challenge to balance complexity with quality per bit when standardizing a format. The most typical way to add coding tools is to be able to select between different modes of operation based on the content of the block, where each mode is suited to certain patterns of input. Use the right tool for the job! As we will see in this study, the number of coding tools will increase exponentially, and it starts to become impossible to make good use of all the tools given to you by the format.

Encoding becomes an optimization task where we aim to figure out the best coding tools to use among the ones given to us. In simpler formats, there are very few things to try, and approaching the optimal solution becomes straight forward, but as we get into the more esoteric formats, the real challenge is to prune dead ends early, since brute forcing our way through a near-infinite configuration space is not practical (but maybe it is with GPU encode? :3)

Commonalities across formats

Image compression and video compression uses the Discrete Cosine Transform (DCT) even to this day. This fundamental compression technique has been with us since the 80s and refuses to die. All the new compression formats just keep piling on complexity on top of more complexity, but in the center of it all, we find the DCT, or some equivalent of it.

Very similarly, texture compression achieves its compression through interpolation between two color values. Somehow, the formats all let us specify two endpoints which are constant over the 4×4 block and interpolation weights, which vary per pixel. Most of the innovation in the formats all comes down to how complicated and esoteric we can make the process of generating the endpoints and weights.

The weight values are typically expressed with very few bits of precision per texel (usually 2 or 3), and this is the main way we will keep bits spent per pixel down. This snippet is the core coding tool in all the formats I have studied:

decoded_texel = mix(endpoint0, endpoint1, weight_between_0_and_1);

To correlate, or not to correlate?

The endpoint model blends all components in lock-step. Typically the endpoint will be an RGB value. We call this correlated, because this interpolation will only work well if chrominance remains fairly constant with luminance being the only component which varies significantly. In uncorrelated input, say, RGB with an alpha mask, many formats let us express decorrelated inputs with two sets of weights.

decoded_rgb = mix(endpoint0_rgb, endpoint1_rgb, rgb_weight);
decoded_alpha = mix(endpoint0_alpha, endpoint1_alpha, alpha_weight);

This costs a lot more bits to encode since alpha_weight is very different from rgb_weight, but it should be worth it.

Many formats let us express if there is correlation or not. Correlation should always be exploited.

Working around the horrible endpoint interpolation artifacts

Almost all formats beyond the most trivial ones try really hard to come up with ways to work around the fact that endpoint interpolation leads to horrible results in all but the simplest input. The most common approach here is to split the block into partitions, where each partition has its own endpoints.

S3TC – The basics

A compute shader decoder:

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/s3tc.comp

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/rgtc.h

BC1 – 4×4 – 64 bits

The BC1 format is extremely simple and a good starting point. 32 bits is used to encode two RGB endpoints in RGB565 format. The other 32 bits encode 16 weights, with 2 bits allocated to each texel.

This lets us represent interpolation weights of 0, 1/3, 2/3 and 1.

Since there is a symmetry in this design, i.e.:

mix(a, b, l) == mix(b, a, 1.0 - l)

there would be two ways to specify the same block, where we swap endpoints and invert the weights to compensate. This is an extra bit of information we can exploit. Based on the integer representation of the two endpoints, we can check if one of greater than the other, and use a different decoding mode based on that information. This exploitation of symmetry will pop up again in many formats later! In the secondary mode, we add support for 1-bit alpha, called a punch-through in most formats. In this mode, the interpolation weights become 0, 1/2, 1 and BLACK. This lets us represent fully transparent pixels. However, color becomes BLACK, so this will only work with pre-multiplied alpha schemes, otherwise there will be black rings around textures. I don’t think this mode is used all that much these days, but it is an option.

That is it for this format, it really is that simple.

One thing to note is that the specification is defined in terms of floating point with under-specified requirements for precision, and thus there is no bit-exact representation of the decoded values. Almost all hardware decoders of this format will give slightly different results, which is unfortunate. MPEG1 and MPEG2 also made the same mistake back in the day, where the DCT is specified in terms of floating point.

BC2 – 4×4 – 128 bits

BC2 is a format which adds alpha support by splicing together two blocks. A BC1 block describes color, and a second block adds an alpha plane with 4-bit UNORM. This format is quite obscure since the next format, BC3, generally does a much better job at compressing alpha. A curious side effect of BC2 and 3 is that the punch-through mode in the BC1 block no longer exists, i.e. the symmetry exists, so we lose 1 bit of information. I wonder why that information bit was not used to toggle between BC2 alpha decode (noisy alpha) and BC3 alpha decode (smooth alpha) …

BC3 – 4×4 – 128 bits

The alpha encoding of BC3 is very similar to how BC1 works. It also forms the basis of RGTC. The 64 bits of the alpha block spends 16 bits to encode 2 8-bit endpoints for alpha. We now get 3 bits as interpolation weights instead of 2 since there’s 48 bits left for this purpose. Similar to BC1, we can exploit symmetry and get two different modes based on the endpoints. In the first mode, the interpolation weights are as expected, [0, 1, 2, …, 7] / 7. In the second mode, we have 6 interpolated values [0, 1, 2, …, 5] / 5, and two values are reserved to represent 0.0 and 1.0. This can be very useful. This is essentially a very early attempt to introduce partitions into the mix, as we can essentially split up a block into 3 partitions on demand: (Fully opaque texels, fully transparent texels, the in-betweens). This can let us specify a tighter range for the endpoints as there is never a need to use a full [0, 0xff] endpoint range.

Summary

The S3TC formats are very simple, but there are certainly things to note about them. Alpha support is just bolted on, it is not integrated into the format. This means that even though the block is 128 bits, there is no way to spend more than 64 bits on color, even if the alpha plane has a completely flat value.

RGTC

RGTC (red-green) is basically BC3’s alpha block format turned into its own thing. Their main use is with non-color textures, e.g. normal maps, metallic-roughness maps, luminance maps, etc. It is very simple, and quality is quite good.

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/rgtc.comp

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/rgtc.h

BC4 – 4×4 – 64 bits

This is BC3’s alpha block format as its own format, which returns one component. The only real difference from BC3 alpha is that it also supports an SNORM variant, which is very useful for normal maps, although I only bothered with UNORM, since my shaders need to assume input can be from any format.

BC5 – 4×4 – 128 bits

RGTC assumes uncorrelated channels, and thus the only sensible choice was to just slap together two BC4 blocks side-by-side, and voila, we can encode 2 channels instead of 1.

Summary

RGTC is simple and nice. It only needs to consider single channels of data, and writing encoders for it is very easy, and it is probably the simplest format out there. For what they do, I really like these formats.

Like S3TC, there is no bit-exact decoding, which is rather unfortunate.

ETC – Refining S3TC

ETC, or Ericsson Texture Compression is a family with multiple generations. ETC2 is backwards compatible with ETC1 in that all valid ETC1 blocks will decode the same way in ETC2, but ETC2 exploits some undefined behavior of ETC1 to extend the format into something more interesting.

ETC has a bit-exact decode, which makes verification very easy. 😀

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/etc2.comp

ETC1 – 4×4 – 64 bits

In many ways, ETC1 is quite similar to BC1, but there are some key differences. Just like BC1, 32 bits are spent to encode endpoints, and 32 bits are spent to give 16 texels 2 bits of weight information each. The main difference between ETC1 and BC1 is how endpoints work.

Sub-blocks

As a very crude form of partitioning, ETC1 allows you to split a block into either 2×4 or 4×2 sub-blocks, where each sub-block has its own endpoints. To do this, endpoints are expressed in a more compact way than BC1. Rather than specifying two RGB values, ETC in general likes to express endpoints as RGB +/- delta-intensity, where delta-intensity is described by a table. This makes things far more compact since we enforce constant chrominance. By saving so many bits, we can express 4 endpoints in total, 2 for each sub-block.

Uncorrelated or correlated endpoints?

Since we have to specify two sets of endpoints, the format gives us a way to specify if the two endpoints have completely different colors, or if the endpoints should be specified in base + offset form. This is controlled with a single bit, which changes the encoding from ep0 = RGB444, ep1 = RGB444 to ep0 = RGB555, ep1 = ep0 + sign_extend(RGB333). These values are not allowed to overflow in any way, which is something ETC2 exploits to great effect later.

NOTE: I found it more instructional to call it uncorrelated and correlated endpoint modes, but the specification calls it “individual” and “differential” modes.

Summary

ETC1 is somewhat different than BC1, but overall, it’s quite similar. It turns out, if you add a few restrictions on top of ETC1, you get ETC1S, which can trivially be transcoded to BC1. Basically, enforce the correlated endpoint mode, set the RGB333 bits to 0, enforce that delta-luma is the same for both sub-blocks.

There is also no way to express alpha with ETC1, which is unfortunate, but ETC1 is completely obsolete for my use cases either way.

ETC2 RGB – 4×4 – 64 bits

As mentioned earlier, ETC1 is a sub-set of ETC2. ETC2 adds a bunch of new and curious modes which gives us a small glimpse into more flexible ways to express endpoints, and even adds a mode which I have never seen in any other formats ever since.

Exploiting undefined behavior

When you select the correlated (differential) endpoint mode, there were some restrictions on overflow. We can exploit this fact in order to gain 3 new modes of operation for the ETC2 color codec!

First, we check if R + sign_extend(dR) is outside [0, 31] range. If so, we activate the so-called “T” mode. In this mode, we essentially add a partitioning scheme to the codec. We now remove the concept of two sub-blocks and let all texels access all available endpoints. We encode two RGB444 values (A, B), and a delta value (d). We form a T-shape by specifying 4 possible color values as A, B, B + d, B – d. This can be useful if the block is smooth, except for some weird outliers. A would be the outlier color, and B represents the middle of the smooth colors.

If G + sign_extend(dG) overflows, we enter a very similar “H” mode. In this mode, we do the exact same thing, except that the 4 possible colors become A + d, A – d, B + d, B – d.

If B + sign_extend(dB) overflows, we enter a very interesting mode, which I have never seen again in future formats. I’m not sure why, since it seems very useful for expressing smooth gradients. Essentially, in this mode we don’t encode weights per texel, but rather express RGB at texel (0, 0), at texel (4, 0) and texel(0, 4), and just bilinearly interpolate across the block to obtain the actual color. This is very different from the other endpoint interpolation we’ve seen earlier, because that flattens everything into a single line in the color space, but now we can access an entire 2D plane in the space instead.

Punch-through alpha

Like BC1’s 1-bit alpha scheme, ETC2 is very similar. When this format is enabled, we remove the capability for uncorrelated endpoints (individual mode), and replace the bit with a selection to select if all texels are Opaque, or potentially transparent. This idea is the exact same as BC1. In the transparent mode, code == 2 marks the texel as being transparent black. It does not work in planar mode though, this bit is ignored there.

Alpha support

Very similar to BC3, ETC2 also supports full 8-bit alpha by slapping together a separate block alongside the color block. The way this works is very similar to how RGTC works, but instead of two endpoints, ETC2 encodes a center point, and then uses tables to expand that range into 8 possible values using a table selector and a multiplier. These 8 possible values for the block are then selected with 3bpp indices. We lose the capability to cleanly represent 0.0 and 1.0 though, which is somewhat curious.

EAC – 4×4 – 64 bits

https://github.com/Themaister/Granite/blob/master/assets/shaders/decode/eac.comp

EAC is ETC2’s version of RGTC, it is designed as a way to encode 1 and 2 de-correlated channels, basically the exact same approach as RGTC where the alpha block format is reused for the 1/2-channel formats. EAC is a bit different in that the internal precision is 11 bits (for some reason).

Unfortunately, EAC is kinda awkward since it’s technically bit-exact in the final fixed point value between [0, 2047], but it specifies many different ways how this can be converted to a floating point value.

Next up, BPTC

S3TC, RGTC and ETC represents the simpler formats, and hopefully I’ve summarized what these formats can do, next up, I’ll go through the BPTC formats, which significantly increases complexity.

Clustered shading evolution in Granite

Shading many lights in a 3D engine is kinda hard once you step outside the bubble of classic deferred shading. Granite supports a fair amount of different ways to do lighting, mostly because I like to experiment with different rendering structures.

I presented some work on this topic at SIGGRAPH 2018 in the Moving Mobile Graphics course when I still worked at Arm: https://community.arm.com/cfs-file/__key/communityserver-blogs-components-weblogfiles/00-00-00-20-66/5127.siggraph_2D00_2018_2D00_mmg_2D00_3_2D00_rendering_2D00_hanskristian.pdf. Unfortunately, some slides have videos embedded, and they didn’t seem to have made the transition unscathed to PDF.

Since then I’ve looked into VK_EXT_descriptor_indexing in order to remove some critical limitations with my older implementation and I ended up with some uncommon implementation details.

Classic deferred

Just to summarize what I mean by this, this is the good old method of rendering light volumes on top of a G-buffer and light is accumulated per-light by blending on top of the frame buffer. This method is considered completely obsolete on desktop these days, but it’s still quite viable on mobile with on-chip G-buffers.

Where classic deferred breaks down

  • Very high bandwidth on desktop due to fill-rate/bandwidth
  • No forward shading support (well, duh)
  • No transparent objects
  • No volumetric lighting

Thus, I explored some alternatives …

Why I don’t like tile deferred shading

This is a somewhat overloaded term, but this is the method where you send the G-buffer to a compute shader. The G-buffer is split into tiles, assigned to a workgroup, depth range (or multiple ranges) is found for the tile, which then forms a frustum. Lights are then culled against this frustum, and shaded in one go. This technique was a really great fit for the PS3 SPUs, as shown by Battlefield 3 back in the day.

It still doesn’t really solve the underlying issues of classic deferred except that it is far more bandwidth efficient on desktop-class GPUs.

Forward shading isn’t feasible unless you split the algorithm into multiple stages with Z-prepass -> build light list per tile -> resubmit geometry and shade, but then it’s probably called something entirely different … (Is this what’s known as Forward+ perhaps?)

Transparency isn’t doable unless you render all transparent geometry in a separate pass to find min/max depth per pixel or something, and forget about volumetrics.

On mobile TBDR architectures you lose all bandwidth savings from staying on-tile, and doing FRAGMENT -> COMPUTE barriers on a tile-based architecture is usually a terrible idea for pipelining. The exception here is tile shaders on recent iPhone hardware which seems almost designed to do this algorithm in hardware.

Clustered shading – the old implementation

Clustered shading is really nice in that it is completely agnostic to a depth buffer, so all the problems mentioned earlier just go away. Lights are assigned in 3D-space rather than screen-space. The original paper on the subject is from around the same time tile deferred was getting popular.

Abandoning the frustum structure

In Granite, I chose early on to abandon the “frustum” layout. Culling spot lights and point lights against elongated frustums analytically is very hard and computationally expensive, and getting the Z slices correct is also very fiddly.

A common workaround for culling is using conservative rasterization, but that is a feature I cannot rely on. I figured that using a more grid-like structure I could get away with much simpler culling math. Since all elements in the grid-based cluster are near-perfect cubes, I could get away with treating the elements as spheres, and https://bartwronski.com/2017/04/13/cull-that-cone/ makes spot light-to-sphere culling very cheap and quite tight. Here’s a visualization of the structure. This structure is stored in a 3D texture. Each “mip-level” is packed in the Z dimension as the resolution of each level is the same.

Bitscan loops instead of light lists

Light lists is the approach where each element in the cluster contains a (start, count) and all lights found in that range are shaded. Computing this list on GPU is rather messy. The memory footprint for a single element is unknown in CPU timeline, and we cannot deal properly with worst-case scenarios. This is easier when we cluster on the CPU instead, but that’s boring!

I really wanted to cluster lights on the GPU, so I landed on a bitmask approach instead. The worst case storage is just 1 bit per light per element rather than 32 bit.

The main limitation of this technique is still the number of lights we could feasibly support. With a bitmask structure we need to allocate for worst-case and it can get out of hand when we consider worst case with 1000+ lights. I only had modest ideas in the beginning, so I supported 32 spot lights and 32 point lights, which were encoded in RG32UI per element in a 3D texture. At a resolution of the cluster at 64x32x16x9, culling on the GPU is very fast, even on mobile. We can set the ceiling higher of course if we expand to RGBA32UI or use more texels per element.

Bitscan loops are great for scalarization

A thing I realized quickly when doing clustered forward shading is the importance of keeping VGPRs down on AMD hardware. The trick to move VGPRs to more plentiful SGPRs is to ensure that values are uniform across a subgroup. E.g. instead of doing this:

// VGPR
int light_bitmask = fetch_bitmask_for_world_coord(coord);

vec3 color = vec3(0.0);
while (light_bitmask != 0)
{
    int lsb = findLSB(light_bitmask);

    // All light data must be loaded in VGRPs since lsb is a VGPR.
    color += shade_light(lsb);
    light_bitmask &= ~(1 << lsb);
}

we can do a simple trick with subgroup operations:

// VGPR
int light_bitmask = fetch_bitmask_for_world_coord(coord);

// OR over all active threads.
// As this is the same value for all threads, compiler promoted to SGPR.
light_bitmask = subgroupOr(light_bitmask);

vec3 color = vec3(0.0);
while (light_bitmask != 0)
{
    int lsb = findLSB(light_bitmask);

    // All light data can be loaded into SGPRs instead.
    // Far better occupancy, much amaze, wow!
    color += shade_light(lsb);
    light_bitmask &= ~(1 << lsb);
}

Uniformly loading light data from buffers is excellent. I’ve observed up to 15% uplift on AMD by doing this. The light list approach mentioned earlier has a much harder time employing this kind of optimization. We would have to scalarize on the cluster element itself, which could lead to very bad worst-case performance.

No bindless – ugly atlasing

Another problem with clustered shading (and tile deferred for that matter) is that we need to shade a lot of lights in one go, and those lights can have shadow maps. Without bindless, all shadow maps for spot lights must fit into one texture, and point lights must fit into one texture. Atlassing is the classic solution here, but it is a little too messy for my taste. As the number of lights was rather low, I just had a plain 2D texture for spot lights, and a cube array for point lights. Implementing variable resolution with an atlas is also rather annoying, and for point lights, I would be forced to flatten the cube down to 6 2D rects and do manual cube lookup instead without proper seam filtering, ugh.

Scaling to “arbitrary” number of lights

While performance for reasonable number of lights was quite excellent compared to alternative techniques, I couldn’t really scale number of lights arbitrarily, and it has been nagging me a bit. Memory becomes a concern, and while the “list of lights” approach is likely less memory hungry in the average case, it has even worse worst-case memory requirements, and it’s not very friendly for GPU culling.

Kinda clustered shading? New bindless hotness

The talk on Call of Duty’s renderer in Advances in Real-Time Rendering 2017 presents a very fresh idea on how to do shading, and it hits all the right buttons. Culling on GPU, bitscan loops, scalarization, scales to a lot of lights, a lot of things to like here.

I spent some time this holiday season implementing a new path for clustered shading based on this technique, and I ended up deviating in a few places. I’ll go through my implementation of it, but it will make more sense once you study the presentation first.

Decoupling XY culling from Z

The key feature of the Call of Duty implementation is how it partitions space. Rather than a full 3D cluster, we decouple XY and Z dimensions, so rather than O(X * Y * Z) we get O(X * Y + Z).

Z is also binned linearly, and we can have several thousand bins in Z. This makes everything much nicer to deal with later. Culling here is trivial, since we compute min/max Z in view space for each light, which is very simple.

For each Z-slice, we just need to figure out the minimum light index and maximum light index which hits our Z-slice. Of course, to make the ranges as tight as possible, we sort lights by Z distance.

Data structures

  • Each tile in XY needs a bitmask array, u32[ceil(num_lights / 32)] for each tile. This can be tightly packed in a single buffer.
  • A buffer containing Z-slices as described above.
  • Per-light information: position/radius/cone/light type/shadow matrix/etc

Going back to frustum culling

Now that we cluster in XY and Z separately, I went back to frustum partitioning, and now we need a way to do frustum culling against spot and point lights … Conservative rasterization really is the perfect extension to use here, just a shame it’s not widely available yet on all relevant hardware.

The presentation has an alternative for conservative rasterization as current consoles do not support this feature, which is mostly to render light volumes a-la classic deferred (not completely dead yet!) at full-resolution and splatting out bits with atomics as fragment threads are spawned. If you have depth information you can eliminate coverage using classic deferred techniques. However, I went in a completely different direction without using depth at all.

Compute shader conservative rasterization

It felt natural to do all the culling work in compute shaders. This is where most of the “fun” is. This is also a rather esoteric way of doing it, but I like doing esoteric stuff with compute shaders, see https://themaister.net/blog/2019/10/12/emulating-a-fake-retro-gpu-in-vulkan-compute/ as proof of that.

Conservative sphere rasterization

To solve this problem, we’re going to tackle it geometrically instead. First, we solve the screen-space bounding box problem. Fortunately, this problem is separable, so we can compute screen-space bounds in X and Y separately.

(Behold glorious Inkscape skillz ._.)

What we want to do is to figure out where P_lo and P_hi intersect the near plane. P can be rotated in 2D by the half-angle in two directions. This way we find tangent points on the circle. sin(theta) is conveniently equal to r / length(L), so building a 2×2 rotation matrix is very easy. After rotating P, we can project P_lo and P_hi, and now we have clip space bounds in one dimension. Compute separately for XZ and YZ dimensions and we have screen-space boundaries.

Projecting a sphere to clip space creates an ellipsis, and to compute the ellipsis, we need to rotate view space such that the sphere center lies perfectly on the X or Y axis. For simplicity, we orient it on the +X axis. We can then perform the range test, and an ellipsis is formed. We can now test any point directly against this ellipsis. If the sphere intersects the near plane in any way, we can fall back to screen-space bounding box.

Here’s a ShaderToy demonstrating the math. Of course, Inigo Quilez did it first :p

In the real implementation, we only need to compute setup data once for each point light, to rasterize a pixel, we apply a transform, and perform a conservative ellipsis test, which is rather straight forward.

Conservative spot light rasterization

I tried an analytical approach, but I gave up. Spot-Frustum culling is really hard if you want tight culling, so I went with the simpler approach of just straight up rasterizing 6 triangles which forms a spot light. We can rasterize in floating point since we don’t care about water-tight rasterization rules, and it’s conservative anyways. It’s not the prettiest thing in the world to do primitive clipping inside a shader, but you gotta do what you gotta do …

The … peculiar shader can be found here.

Binning shader

Once we have setup data for point lights and spot lights, we do the classic culling optimization where we bin N lights in parallel over a tile, broadcast the results to all threads, which then computes the relevant lights per pixel. Subgroup ballots is a nice trick here, which replaces the old shared memory approach. Each workgroup preferably works on 32 lights at a time to compute a 32-bit bitmask.

Shader: https://github.com/Themaister/Granite/blob/master/assets/shaders/lights/clusterer_bindless_binning.comp

The shading loop

The final loop to shade becomes something like:

uint cluster_mask_range(uint mask, uvec2 range, uint start_index)
{
	range.x = clamp(range.x, start_index, start_index + 32u);
	range.y = clamp(range.y + 1u, range.x, start_index + 32u);

	uint num_bits = range.y - range.x;
	uint range_mask = num_bits == 32 ?
		0xffffffffu :
		((1u << num_bits) - 1u) << (range.x - start_index);
	return mask & uint(range_mask);
}

vec3 shade_clustered(Material material, vec3 world_pos)
{
    ivec2 cluster_coord = compute_clustered_coord(gl_FragCoord.xy);
    int linear_cluster_coord = linearize_coord(cluster_coord);
    int z = compute_z_slice(dot(world_pos - camera_pos, camera_front));

    uvec2 z_range = cluster_range[z];

    // Find min/max light we need to consider when shading slice Z.
    // Make this uniform across subgroup. SGPR.
    int z_start = int(subgroupMin(z_range.x) >> 5u);
    int z_end = int(subgroupMax(z_range.y) >> 5u);

    for (int i = z_start; i <= z_end; i++)
    {
        // SGPR
        uint mask =
            cluster_bitmask[linear_cluster_coord *
                            num_lights_div_32 + i];

        // Restrict to lights within our Z-range. VGPR now.
        mask = cluster_mask_range(mask, z_range, 32u * i);
        // SGPR again.
        mask = subgroupOr(mask);

        // SGPR
        uint type_mask = cluster_transforms.type_mask[i];

        // Good old scalarized loop <3
        while (mask != 0u)
        {
            int bit_index = findLSB(mask);
            int light_index = 32 * i + bit_index;
            if ((type_mask & (1 << bit_index)) != 0)
            {
                result += compute_point_light(light_index,
                                              material,
                                              world_pos);
            }
            else
            {
                result += compute_spot_light(light_index,
                                             material,
                                             world_pos);
            }
        }
    }
}

Shader: https://github.com/Themaister/Granite/blob/master/assets/shaders/lights/clusterer_bindless.h

Potential performance problems

By decoupling XY and Z in culling there’s a lot of potential for false positives where large lights might dominate how Z-ranges are computed and trigger a lot of over-shading. I haven’t done much testing here though, but this is probably the only real weakness I can think of with this technique. Regular tile-deferred has similar issues.

Culling tightness

I’m using smaller lights here to demonstrate. Red or green light signifies that a light was computed for that pixel:

Placing point light in volumetric fog for good measure. The weird red/green “artifacting” around the edges is caused by the forward shading and subgroupOr logic when shading to ensure subgroup uniform behavior.

Potential improvements

Clipping Z-range calculation against view frustum might help a fair bit, since Z-range can be way too conservative for large positional lights like these. Especially spot lights which point to the side like the image above. Classic deferred has a very similar problem case unless ancient stencil culling techniques are used to get double sided depth tests.

Conclusion

I’m happy with the implementation. Performance seems very good, but I haven’t dug deep in analysis there yet. I was mostly concerned getting it to work. Just waiting for mobile GPU vendors to support bindless, so I can test there as well I guess.

The weird world of shader divergence and LOD

Mip-mapping is hard – importance of keeping your quads full

Sampling textures with mip-mapping is ancient, but it’s still hard apparently. Implicit LOD calculation is the first instance where we poke a hole into the “single threaded” abstraction of high level shading languages like GLSL and HLSL and dive into the maddening world of warps, waves, quads and everything in-between. For fragment shading, at least a group of 2×2 threads (a quad) need to run side by side so we can have gradient information over the screen. On any modern GPU, these 2×2 threads are actually running in lock-step as we’ll see when looking at GPU ISA later …

The Vulkan/GL/GLES ecosystem has always specified that implicit LOD instructions must happen in dynamically uniform control flow. Dynamically uniform just means that either all threads have to execute a texture() instruction, or no threads do. This ensures that there is always 4 valid texture coordinates from which to compute derivatives. The easiest way to ensure that this guarantee holds is simply to never sample in control flow, but that’s not really practical in more interesting shaders.

If you’re sampling in control flow you better make sure you uphold the guarantees of the spec.

Having to be dynamically uniform over an entire draw call is a bit silly, so the Vulkan specification recently tightened the scope such that if you have subgroupSize >= 4, you only need to be dynamically uniform on a per-quad granularity. This makes sense. We only need correct derivatives in the quad, we shouldn’t have to care if some unrelated quad or even triangle is diverging.

An interesting case came up recently where apparently developers expect that you actually can sample with implicit LOD in diverging control flow. Apparently HLSL “defines” this control flow to be valid code.

vec2 uv = from_somewhere();
if (weight > 0.0)
    sum += weight * texture(Texture, uv);

The idea is that we shouldn’t have to sample the texture unless we’re going to use it, but it’s still nice to provide UV for LOD purposes. Unfortunately, there is no obvious way to express this optimization in high level languages. UV is well defined in the outer scope which is dynamically uniform, so that’s something … Intuitively, this code makes sense, but it gets really murky once we dig deeper.

With subgroup ops, we can probably get a good approximation on the HLL side.

bool quadAny(bool value)
{
    // Perhaps this can be translated into s_wqm on AMD
    // if compiler checks this pattern?
    return subgroupClusteredOr(int(value), 4) != 0;
}

vec2 uv = from_somewhere();
// Hoist texture sampling out of branch and force quad uniformity.
vec4 tex;
if (quadAny(weight > 0.0))
    tex = texture(Texture, uv);
if (weight > 0.0)
    sum += weight * tex;

Querying gradients and then sampling with that in the branch is fine as well, but it is slow, and not really a fix, at best a workaround.

HLSL seems a bit murky about if this kind of code is legal, it’s all “that one app did this thing that one time and now we’re screwed”. From my understanding compilers can do some heroics here to work around this in applications.

I wanted to try this kind of code on all Vulkan devices I have available to see what happens. We’re in undefined territory as far as LOD goes, so anything can happen. There’s three outcomes I’m looking for which seem like plausible HW behavior:

  • It just happens to work. This is kinda scary, since it’ll probably break in 5 years anyways.
  • The LOD computed is garbage.
  • The LOD is forced to some value on divergence.

Here’s the concrete shader I’m using, from https://github.com/Themaister/Granite/blob/master/tests/assets/shaders/divergent_lod.frag. A test to run the shader is https://github.com/Themaister/Granite/blob/master/tests/divergent_lod_test.cpp

#version 450

layout(location = 0) in vec2 vUV;
layout(location = 0) out vec4 FragColor;
layout(set = 0, binding = 1) uniform sampler2D uSampler;
layout(set = 0, binding = 0, std140) uniform Weights
{
    vec4 weights[4];
};

void main()
{
    vec3 tex = vec3(0.0);
    float lod = -10.0;
    vec2 uv = vUV;
    if (weights[int(gl_FragCoord.x) + 2 * int(gl_FragCoord.y)].x > 0.0)
    {
        tex = texture(uSampler, uv).rgb;
        lod = textureQueryLod(uSampler, uv).y;
    }

    FragColor = vec4(tex, lod);
}

I render this on a 2×2 frame buffer with a full-screen “expanded triangle” to not get any helper lane shenanigans. Let’s try to run this across a wide range of hardware and see what happens. NOTE: any result here is equally valid in Vulkan, this is intentionally going out of spec.

AMD

I tested this on a Navi card. RDNA ISA seems similar enough to GCN … We effectively have 4 driver stacks for AMD cards now to test.

RADV (LLVM 10)

Garbage LOD

main:
BB16_0:
	s_mov_b64 s[0:1], exec    
	s_wqm_b64 exec, exec
	v_cvt_i32_f32_e32 v3, v3
	s_mov_b32 s6, s3
	s_movk_i32 s7, 0x8000          
	v_cvt_i32_f32_e32 v2, v2     
	v_mov_b32_e32 v5, 0xc1200000    
	s_load_dwordx4 s[8:11], s[6:7], 0x0    
	v_lshlrev_b32_e32 v3, 1, v3  
	v_add_lshl_u32 v2, v3, v2, 4   
	s_waitcnt lgkmcnt(0)    
	buffer_load_dword v4, v2, s[8:11], 0 offen    
	v_mov_b32_e32 v2, 0  
	v_mov_b32_e32 v3, v2      
	s_waitcnt vmcnt(0)        
	v_cmp_lt_f32_e32 vcc, 0, v4 against 0       
	v_mov_b32_e32 v4, v2     
	s_and_saveexec_b64 s[8:9], vcc
	s_cbranch_execz BB16_2
BB16_1:
	s_mov_b32 s3, s7      
	s_add_i32 s6, s2, 0x50           
	s_mov_b32 m0, s4           
	s_load_dwordx8 s[12:19], s[2:3], 0x0       
	s_load_dwordx4 s[20:23], s[6:7], 0x0          
	v_interp_p1_f32_e32 v7, v0, attr0.x
	v_interp_p1_f32_e32 v8, v0, attr0.y      
	v_interp_p2_f32_e32 v7, v1, attr0.x          
	v_interp_p2_f32_e32 v8, v1, attr0.y   
	s_waitcnt lgkmcnt(0)     
	image_sample v[2:4], v[7:8], s[12:19], s[20:23] dmask:0x7 dim:SQ_RSRC_IMG_2D
	image_get_lod v5, v[7:8], s[12:19], s[20:23] dmask:0x2 dim:SQ_RSRC_IMG_2D
BB16_2:
	v_nop                 
	s_or_b64 exec, exec, s[8:9]              
	s_and_b64 exec, exec, s[0:1]         
	s_waitcnt vmcnt(0)                
	exp mrt0 v2, v3, v4, v5 done vm          
	s_endpgm               

We see that v7 and v8 hold the UV coordinates, but they are actually only computed inside the branch (v_interp). The optimizer is allowed to place UV computation inside the branch here. If there is divergence in a quad, the disabled lanes won’t get correct values for v7 and v8 (since execution is masked), and LOD becomes garbage.

RADV (ACO)

Coupled with Navi cards, this is probably the most bleeding edge setup you can run. It’s a completely new compiler backend for AMD cards, not based on LLVM.

Just happens to work

BB0:
	s_wqm_b64 exec, exec 
	s_mov_b32 s0, s3    
	s_movk_i32 s1, 0x8000  
	s_load_dwordx4 s[8:11], s[0:1], 0x0 
	s_mov_b32 m0, s4 
	v_interp_p1_f32_e32 v4, v0, attr0.y 
	v_cvt_i32_f32_e32 v2, v2   
	v_cvt_i32_f32_e32 v3, v3   
	v_lshl_add_u32 v2, v3, 1, v2 
	v_lshlrev_b32_e32 v2, 4, v2 
	s_waitcnt lgkmcnt(0)   
	buffer_load_dword v2, v2, s[8:11], 0 offen 
	v_interp_p2_f32_e32 v4, v1, attr0.y 
	v_interp_p1_f32_e32 v0, v0, attr0.x  
	v_interp_p2_f32_e32 v0, v1, attr0.x 
	v_mov_b32_e32 v1, v4  
	s_waitcnt vmcnt(0)   
	v_cmp_lt_f32_e32 vcc, 0, v2  
	s_and_saveexec_b64 s[0:1], vcc     
	s_cbranch_execz BB3  
BB1:
	s_movk_i32 s3, 0x8000 
	s_load_dwordx8 s[4:11], s[2:3], 0x0  
	s_load_dwordx4 s[12:15], s[2:3], 0x50 
	s_waitcnt lgkmcnt(0)   
	image_sample v[2:4], v[0:1], s[4:11], s[12:15] dmask:0x7 dim:SQ_RSRC_IMG_2D
	image_get_lod v0, v[0:1], s[4:11], s[12:15] dmask:0x2 dim:SQ_RSRC_IMG_2D
BB3:
	s_andn2_b64 exec, s[0:1], exec 
	s_cbranch_execz BB6  
BB4:
	v_mov_b32_e32 v0, 0xc1200000 
	v_mov_b32_e32 v2, 0  
	v_mov_b32_e32 v3, 0   
	v_mov_b32_e32 v4, 0  
BB6:
	s_mov_b64 exec, s[0:1] 
	s_waitcnt vmcnt(0)   
	exp mrt0 v2, v3, v4, v0 done vm 
	s_endpgm

This time, UV is interpolated outside the branch, so sampling in divergent control flow ends up working after all. The registers are well defined as they enter the branch. For AMD, it seems like it just comes down to whether or not the lanes have correct values placed in them and not having them be clobbered by the time we get around to sampling. There doesn’t seem to be any hardware level checks for divergence.

AMDVLK

Garbage LOD

AMDVLK uses the same LLVM stack that RADV LLVM uses, and no surprise, same result, and basically same exact ISA is generated.

Windows

Also just happens to work

I guess it’s the exact same case as the ACO compiler here. No need to paste disassembly.

Intel

Tested on UHD 620 (8th gen mobile CPU I think).

Anvil (Mesa)

The Mesa compiler can spit out assembly, which is nice.

Just happens to work

ISA (a little too wide to embed): https://gist.github.com/Themaister/7c5b011cde3c7585459b089f80f897e2

From what I can make out of the ISA, the UV is interpolated outside control flow, and then only the sampling takes place in control flow. It seems like Intel has similar behavior as AMD here, in that just as long as the registers are valid, divergent sampling “works”.

Windows

Just happens to work

Doesn’t seem to be a way to get ISA from Windows driver, but I suppose it’s same as ANV here.

Nvidia

Tested on a Turing GPU on Linux. Didn’t bother testing on Windows as well considering the driver stack is basically the same.

LOD is clamped to 0, textureQueryLod returns -32.0.

Apparently, now we start seeing interesting behavior. Unfortunately, there is no public ISA to look at. The -32.0 LOD might look weird, but this is kind of expected. This is apparently the smallest possible representable LOD on this GPU. LOD is usually represented in some kind of fixed point, log2(0) = -inf after all.

I confirmed it worked as expected when using non-divergent execution as a sanity check.

Arm

Tested on Mali-G72.

LOD is clamped to 0, textureQueryLod returns -128.0.

Very similar behavior to Nvidia here, except the LOD is -128.0 rather than -32.0. I confirmed it worked as expected when using non-divergent execution as a sanity check.

QCOM

Tested on Adreno 506.

Garbage LOD

Again, no ISA to look at. I confirmed it worked as expected when using non-divergent execution as a sanity check.

Conclusion

Never ever rely on LOD behavior with divergent quads (EDIT: at least the way it’s specced out and implemented on Vulkan drivers right now) . You’d be contributing to the pain and suffering of compiler engineers the world over. Staying quad-uniform is fine though.