Corrs cause CORS

## Thursday, 5 November 2015

### sRGB Integer Conversions

A comment on a previous post on this subject got me thinking; 'Anonymous' mentioned using a lookup table (LUT) for converting sRGB to and from linear values. This would be less than ideal on a GPU, but a good compromise on a CPU, particularly if integer values were expected.

But first, I tried to wring some more mileage out of the polynomial approximation of the sRGB gamma curve, which, as you'll remember, is:

if (srgb <= 0.04045)
lin = srgb / 12.92;
else
lin = pow((srgb + 0.055) / 1.055, 2.4);

It is common that the input 'srgb' values are bytes (0..255 for each of red, green and blue channels) and so are the output 'lin' values. My initial thought was to use a quartic integer polynomial of the form:

uint32 x = srgb;
uint32 y = a*x*x*x*x + b*x*x*x + c*x*x + d*x + e;
lin = (uint8)(y >> 24);

That is, pick the constants 'a', 'b', 'c', 'd' and 'e' such that, after multiplying through as 32-bit unsigned integers, the result is conveniently in the uppermost eight bits of 'y'. I dutifully prepared a table of values and fed them into an on-line tool such as the excellent Polynomial Regression Data Fit page by Paul Lutus. It came up with the following polynomial for 'y' in terms of 'x':

y = -1.296·10-1 x4 + b x3 + c x2 + d x + e

This was a little disappointing: the value 'a' turned out to be about -0.1296. That's definitely not an integer coefficient. So I thought I'd be clever and approximate the first term as

-(x * x * x * x / 8)

Serendipitously, the division can be performed using a right-shift. This approximation is all very well, but the remaining coefficients are no longer correct. So I rearranged the table so that the residuals could be approximated as a cubic polynomial:

y + x*x*x*x/8 = b*x*x*x + c*x*x + d*x + e

The next data fit suggested a polynomial with a value of 'b' of approximately 142. I continued the process of fitting, rounding to an integer coefficient and computing the residuals until I got all the values:

a = -1/8
b = 142
c = 35322
d = 602000
e = 11918016

Plugging these numbers back into the formula and testing all 256 input possibilities, I got exact results for 238 inputs and a maximum error of ±1 for the remaining 18. That's not too bad, but hardly stellar, especially as a lookup table of 256 bytes would be more accurate and almost certainly faster!

However, if the input/output integers were 16 bits, instead of 8, a polynomial solution might be worthwhile as the size of a 65536-element array would be cache-unfriendly. This is not as unlikely as it seems: high quality graphics processing is expected to deal with 16-bits-per-channel sRGB inputs.

So, going through the whole process again, but this time using 16- and 64-bit integers, I came up with:

uint64 x1 = srgb;
uint64 x2 = x1 * x1;
uint64 x4 = x2 * x2;
uint64 y = (x1 + 3441) * (x1 + 72046) * x1 * 32627 - x4 / 10;
lin = (uint16)(y >> 48);

If you forgive the unfortunate integer division by 10, this boils down to nine arithmetic operations and one shift. The error is never more than about 1 part in a 1000, but just about all 65536 conversions are slightly off. Still, not too bad.

Finally, I bit the bullet and went for the somewhat more prosaic piecewise linear approximation. This is an undervalued technique that allows you to balance the speed of table-based lookups with the immediacy of full computation. Here's the algorithm:

const uint32 lut[256][2];
uint8 i = (uint8)(srgb >> 8);
uint32 y = lut[i][0] * srgb + lut[i][1];
lin = (uint16)(y >> 16);

One kilobyte of constant lookup table data isn't too cache-unfriendly and the code has a certain beautiful simplicity.

I divided the 65536-element sRGB curve table into 256 equal parts. Each part was fitted to the appropriate portion of the curve using simple linear regression, being careful to handle rounding correctly. The gradient coefficients are in 'lut[i][0]' and the y-intercepts are in 'lut[i][1]'. They are biased so that the required 16-bit result is in the upper half of 'y'. For completeness, the full table is listed at the end of this post.

So, how does it perform? Actually, very well. The error is never more than ±1 in 65535, with 61009 of the 65536 mappings being exact results. As an aside, if you replaced the integer entries in 'lut' with appropriate floating-point values, you could also construct a fairly accurate integer-sRGB-to-float-linear conversion.

I'll leave the reverse transformation (linear to sRGB) as an exercise for the reader.

Here's that table...

const uint32 lut[256][2] = {
{   5076, 0x00007CEA },
{   5068, 0x000086E6 },
{   5071, 0x00008368 },
{   5069, 0x00008C66 },
{   5078, 0x000067EB },
{   5073, 0x00007B68 },
{   5071, 0x00008A68 },
{   5070, 0x000090E7 },
{   5065, 0x0000C064 },
{   5078, 0x00004BEB },
{   5192, 0xFFFBB824 },
{   5497, 0xFFEEA33C },
{   5798, 0xFFE08653 },
{   6103, 0xFFD1136C },
{   6415, 0xFFC00C08 },
{   7059, 0xFF99184A },
{   7391, 0xFF830BF0 },
{   7719, 0xFF6BF894 },
{   8054, 0xFF531BBB },
{   8390, 0xFF38D863 },
{   8734, 0xFF1C9F0F },
{   9079, 0xFEFEF63C },
{   9428, 0xFEDF956A },
{   9778, 0xFEBEC519 },
{  10138, 0xFE9B9CCD },
{  10500, 0xFE76D482 },
{  10870, 0xFE4FC93B },
{  11227, 0xFE28B46E },
{  11600, 0xFDFE71A8 },
{  11972, 0xFDD2D862 },
{  12347, 0xFDA5709E },
{  12726, 0xFD7612DB },
{  13114, 0xFD440F9D },
{  13505, 0xFD1021E0 },
{  13891, 0xFCDB5BA2 },
{  14287, 0xFCA3AB68 },
{  14673, 0xFC6BE828 },
{  15082, 0xFC2F3775 },
{  15489, 0xFBF134C0 },
{  15894, 0xFBB1EC0B },
{  16273, 0xFB753D48 },
{  16699, 0xFB2F7A1E },
{  17126, 0xFAE7C973 },
{  17549, 0xFA9F13C6 },
{  17963, 0xFA565196 },
{  18384, 0xFA0AAFE8 },
{  18823, 0xF9BA1644 },
{  19246, 0xF96AC597 },
{  19679, 0xF917E4F0 },
{  20107, 0xF8C44FC6 },
{  20545, 0xF86D16A0 },
{  20989, 0xF812EA7E },
{  21418, 0xF7BA20D5 },
{  21845, 0xF7603D2A },
{  22327, 0xF6F8AE1C },
{  22784, 0xF694B180 },
{  23228, 0xF631D25E },
{  23678, 0xF5CBE33F },
{  24141, 0xF56135A6 },
{  24608, 0xF4F3C210 },
{  25070, 0xF485A8F7 },
{  25541, 0xF4139362 },
{  25986, 0xF3A611C1 },
{  26477, 0xF32B5836 },
{  26945, 0xF2B48320 },
{  27423, 0xF2394610 },
{  27896, 0xF1BD7A7C },
{  28387, 0xF13B0CF2 },
{  28864, 0xF0BA7360 },
{  29347, 0xF0365ED2 },
{  29830, 0xEFB06A43 },
{  30318, 0xEF272B37 },
{  30815, 0xEE9972B0 },
{  31308, 0xEE0AEF26 },
{  31800, 0xED7ACB1C },
{  32770, 0xEC592501 },
{  33305, 0xEBB6168C },
{  33798, 0xEB1DF803 },
{  34309, 0xEA7E5182 },
{  34817, 0xE9DD9A80 },
{  35333, 0xE9385682 },
{  35844, 0xE892AB02 },
{  36361, 0xE7E90B84 },
{  36886, 0xE73AB80B },
{  37403, 0xE68D058E },
{  37919, 0xE5DDA990 },
{  38447, 0xE5282998 },
{  38970, 0xE472551D },
{  39493, 0xE3BA7BA2 },
{  40021, 0xE2FED4AA },
{  40556, 0xE23E9636 },
{  41099, 0xE1794FC6 },
{  41633, 0xE0B53CD0 },
{  42167, 0xDFEF0FDC },
{  42710, 0xDF236F6B },
{  43255, 0xDE54E7FC },
{  43812, 0xDD7F9892 },
{  44339, 0xDCB3CD1A },
{  44892, 0xDBDBC4AE },
{  45439, 0xDB03F040 },
{  45983, 0xDA2B3150 },
{  46545, 0xD9491468 },
{  47082, 0xD86EEDF5 },
{  47654, 0xD7845613 },
{  48211, 0xD69DB2AA },
{  48764, 0xD5B6923E },
{  49891, 0xD3D8F7F2 },
{  50459, 0xD2E4EE0E },
{  51033, 0xD1EC112C },
{  51599, 0xD0F47248 },
{  52156, 0xCFFE9BDE },
{  52745, 0xCEF85C84 },
{  53324, 0xCDF44726 },
{  53896, 0xCCF11344 },
{  54487, 0xCBE2FBEC },
{  55641, 0xC9CCC12C },
{  56217, 0xC8BEC54C },
{  56799, 0xC7ABB970 },
{  57394, 0xC6903619 },
{  57990, 0xC571DA43 },
{  58585, 0xC4519FEC },
{  59169, 0xC3347710 },
{  59770, 0xC20CA4BD },
{  60357, 0xC0E97262 },
{  60947, 0xBFC2798A },
{  61554, 0xBE90A439 },
{  62151, 0xBD6180E4 },
{  62758, 0xBC2AEA93 },
{  63359, 0xBAF50C40 },
{  63976, 0xB9B482F4 },
{  64565, 0xB880359A },
{  65313, 0xB6F5B510 },
{  65594, 0xB6609D1D },
{  66392, 0xB4B5DDAC },
{  67048, 0xB35433F4 },
{  67654, 0xB20B1C23 },
{  68259, 0xB0C03BD2 },
{  68881, 0xAF69A808 },
{  69496, 0xAE1485BC },
{  70120, 0xACB7FAF4 },
{  70740, 0xAB5B3F2A },
{  71367, 0xA9F822E4 },
{  71987, 0xA896921A },
{  72614, 0xA72E93D3 },
{  73252, 0xA5BDC612 },
{  73882, 0xA44F214D },
{  74524, 0xA2D6F48E },
{  75163, 0xA15E054E },
{  75799, 0x9FE4628C },
{  76435, 0x9E683ECA },
{  77063, 0x9CEE7B04 },
{  77702, 0x9B6B9DC3 },
{  78338, 0x99E81901 },
{  79002, 0x9850E64D },
{  79643, 0x96C5460E },
{  80292, 0x95322DD2 },
{  80938, 0x939E6E15 },
{  81602, 0x91FCD061 },
{  82225, 0x90728F18 },
{  82889, 0x8ECBCE64 },
{  83541, 0x8D2A20AA },
{  84197, 0x8B8352F2 },
{  84856, 0x89D806BC },
{  85512, 0x882C1B04 },
{  86176, 0x86785C50 },
{  86832, 0x84C75298 },
{  87509, 0x8305BC6A },
{  88164, 0x81503732 },
{  88831, 0x7F901500 },
{  90164, 0x7C08B21A },
{  90840, 0x7A3A906C },
{  91493, 0x7879A832 },
{  92187, 0x7699D98E },
{  92859, 0x74C695DE },
{  93538, 0x72EBD0B1 },
{  94212, 0x7111E102 },
{  94887, 0x6F34A1D4 },
{  95566, 0x6D51ECA7 },
{  96251, 0x6B683D7E },
{  96941, 0x69784BD6 },
{  97615, 0x67913128 },
{  98304, 0x659CA500 },
{  98993, 0x63A547D8 },
{  99670, 0x61B421AB },
{ 100369, 0x5FB01288 },
{ 101055, 0x5DB2ECE0 },
{ 101744, 0x5BB0E1B8 },
{ 102443, 0x59A4A496 },
{ 103835, 0x55887B4E },
{ 104527, 0x53795EA8 },
{ 105232, 0x515D9D88 },
{ 105935, 0x4F40A268 },
{ 106634, 0x4D23F945 },
{ 107327, 0x4B094C20 },
{ 108041, 0x48DB7F84 },
{ 108731, 0x46BDC1DE },
{ 109469, 0x4477684E },
{ 110164, 0x4250472A },
{ 110872, 0x401C138C },
{ 111573, 0x3DEABC6A },
{ 112282, 0x3BB03E4D },
{ 113000, 0x396BB7B4 },
{ 113717, 0x3725289A },
{ 114417, 0x34E9B4F8 },
{ 115149, 0x32914B66 },
{ 115865, 0x30432DCC },
{ 116593, 0x2DE85338 },
{ 117317, 0x2B8DEAA2 },
{ 118038, 0x2933358B },
{ 118760, 0x26D4D1F4 },
{ 119487, 0x246F66E0 },
{ 120214, 0x220721CB },
{ 120938, 0x1F9E9935 },
{ 121667, 0x1D2EF422 },
{ 122401, 0x1AB82990 },
{ 123141, 0x18394902 },
{ 123858, 0x15CB7DE9 },
{ 124586, 0x13515B55 },
{ 125319, 0x10D00244 },
{ 126063, 0x0E421EB8 },
{ 126796, 0x0BBB07A6 },
{ 127535, 0x092BC398 },
{ 128289, 0x068C3910 },
{ 129043, 0x03E9AC8A },
{ 129778, 0x01553C79 },
{ 130452, 0xFEF517CA },
{ 131275, 0xFC0B4AE6 },
{ 131988, 0xF9825BCA },
{ 132739, 0xF6D3F2C2 },
{ 133496, 0xF41D06BC },
{ 134231, 0xF17774AC },
{ 134991, 0xEEB7E328 },
{ 135745, 0xEBFAE3A0 },
{ 136491, 0xE9426F16 },
{ 137238, 0xE686280B },
{ 137998, 0xE3BAB987 },
{ 138758, 0xE0EC5103 },
{ 139518, 0xDE1AE87F },
{ 140265, 0xDB52F474 },
{ 141036, 0xD8711D76 },
{ 141788, 0xD59E83EE },
{ 142558, 0xD2B7A06F },
{ 143327, 0xCFCEAA70 },
{ 144092, 0xCCE6966E },
{ 144853, 0xC9FF696A },
{ 145619, 0xC71066EA },
{ 146384, 0xC41F5CE8 },
{ 147190, 0xC102CA7B },
{ 147924, 0xBE2A77EA },
{ 148707, 0xBB1E88F2 }
};

## Sunday, 1 November 2015

### Cocktail Chord Diagram

I've had the third and probably final instalment of my investigations into HTML 5 sitting on my hard disk for quite some time, so I've decided to bite the bullet, tidy it up and upload.

This experiment uses my JSON cocktail database to display a cocktail chord diagram. Chord diagrams have been very trendy in the last couple of years, to the point where I've seen them turn up as art exhibits in galleries and museums. I've always been very sceptical about how much information they actually convey, so I thought I'd experiment with the cocktail dataset.

I've grouped cocktail ingredients together in a fairly arbitrary manner and then pruned any group that is used in fewer than ten of the 289 cocktails. The width of a chord between two ingredient groups is proportional to the number of cocktails in the database that uses both groups.

It turns out that these diagrams are quite tricky to make aesthetically pleasing, especially the order and form of the chords (see the JavaScript if you're curious).

The result is mildly gratifying, but I'm still not convinced about these diagrams' utility, so I didn't spend a huge amount of time tweaking it.

## Saturday, 3 October 2015

### RGB/HSV in HLSL 7

I got a very kind email from David Schaeffer yesterday with some bug fixes for the shader code I wrote for converting RGB to the HCY colour space.

My original implementation was

float3 RGBtoHCY(in float3 RGB)
{
float3 HCV = RGBtoHCV(RGB);
float Y = dot(RGB, HCYwts);
if (HCV.y != 0)
{
float Z = dot(HUEtoRGB(HCV.x), HCYwts);
if (Y > Z)
{
Y = 1 - Y;
Z = 1 - Z;
}
HCV.y *= Z / Y;
}
return float3(HCV.x, HCV.y, Y);
}

As he points out:
It seems to me that you have an error in your RGB to HCY code, though.  I was having some pixels always saturate to white, and after looking at Kuzma Shapran's original code I don't think the value of Y should actually be changed, only using the adjusted value of Y to calculate HCV.y.  Additionally I think the check for 0 in the original code was meant to prevent divide by 0 errors, so should be using Y instead of HCV.y.
I'll have to bow to David's superior knowledge on HCY usage as I've never actually used it in anger. So here's his improved version:

float3 RGBtoHCY(in float3 RGB)
{
float3 HCV = RGBtoHCV(RGB);
float Y = dot(RGB, HCYwts);
float Z = dot(HUEtoRGB(HCV.x), HCYwts);
if (Y < Z)
{
HCV.y *= Z / (Epsilon + Y);
}
else
{
HCV.y *= (1 - Z) / (Epsilon + 1 - Y);
}
return float3(HCV.x, HCV.y, Y);
}

I've updated the snippets page accordingly.

## Thursday, 1 October 2015

When I was at primary school, they used to run a book club where you could buy books with your pocket money. One of the purchases I made was Denis Gifford's "Monsters of the Movies" for the princely sum of 45p.

The book fascinated me and probably inspired my love of horror films. However, re-reading the book recently there are a few wrinkles I hadn't noticed.

The book was published in 1977 and had 46 entries:

 The Alligator People in The Alligator People (1959) Certificate X The Ape Man in The Ape Man (1942) Certificate X (now PG) Barnabas Collins in House of Dark Shadows (1970) Certificate X (now 18) The Beast in Beauty and the Beast (1946) Certificate A (now PG) The Blood Beast Terror in The Blood Beast Terror (1967) Certificate X (now 12) The Bride of Frankenstein in The Bride of Frankenstein (1935) Certificate A (now PG) Carmilla in The Vampire Lovers (1969) Certificate X (now 15) The Cat in The Cat and the Canary (1939) Certificate A (now PG) Count Yorga in Count Yorga, Vampire (1970) Certificate X Countess Dracula in Countess Dracula (1971) Certificate X (now 18) The Creature in The Creature from the Black Lagoon (1954) Certificate A (now PG) The Creeper in The Brute Man (1946) Certificate A The Demon in The Curse of the Demon (1958) Certificate X (now PG) Dr Caligari in The Cabinet of Dr Caligari (1919) Certificate A (now U) Dr Jekyll and Mr Hyde in Dr Jekyll and Mr Hyde (1931) Certificate A (now 12A) Dr Moreau in The Island of Lost Souls (1932) Certificate - (now PG) Dr Phibes in The Abominable Dr Phibes (1971) Certificate X (now 15) Dr X in Doctor X (1932) Certificate A Dracula in Dracula (1931) Certificate A (now PG) The Electric Man in Man-Made Monster (1941) Certificate A The Fly in The Fly (1958) Certificate X (now PG) Frankenstein's Monster in Frankenstein (1931) Certificate A (now PG) Fu Manchu in The Mask of Fu Manchu (1932) Certificate A The Ghoul in The Ghoul (1933) Certificate A (now PG) The Golem in The Golem (1920) Certificate - (now PG) The Gorgon in The Gorgon (1964) Certificate X (now 12) The Hunchback of Notre Dame in The Hunchback of Notre Dame (1923) Certificate A (now PG) The Invisible Man in The Invisible Man (1933) Certificate A (now PG) King Kong in King Kong (1933) Certificate A (now PG) The Mad Monster in The Mad Monster (1941) Certificate - The Manster in The Split (1962) Certificate X The Mummy in The Mummy's Hand (1941) Certificate A (now PG) The Munsters in Munster, Go Home! (1966) Certificate U (now U) Nosferatu in Nosferatu (1922) Certificate A (now PG) Orlac in The Hands of Orlac (1935) Certificate A The Phantom of the Opera in The Phantom of the Opera (1925) Certificate U (now PG) The Raven in The Raven (1935) Certificate A (now 15) The Reptile in The Reptile (1966) Certificate X (now 15) Teenage Werewolf in I Was a Teenage Werewolf (1957) Certificate X (now 15) The Walking Dead in The Walking Dead (1936) Certificate A The Werewolf of London in Werewolf of London (1934) Certificate A (now PG) White Zombie in White Zombie (1932) Certificate A (now PG) The Wild Woman in The Jungle Captive (1944) Certificate H The Wolf Man in The Wolf Man (1941) Certificate A (now PG) Zaroff in The Most Dangerous Game (1932) Certificate A (now 12) The Zombie in The Plague of the Zombies (1966) Certificate X (now 12A)

Here are my thoughts:

• The latest films listed were released in 1971, so perhaps Mr Gifford wrote the book a few years before it was finally published.
• The fact that's there 46 entries, as opposed to some nice "round" number like 50 is probably due to the limitations of the cheap printing: there are 48 leaves between the soft covers.
• Some of the films are a bit suspect for a book obviously aimed at children: "The Vampire Lovers" and "The Abominable Dr Phibes" being prime examples.
• "The Island of Lost Souls" was rejected by the BBFC and only given an X certificate in 1958.
• "The Golem" was never release in the UK.
• "The Mad Monster" was initially rejected by the BBFC and only given an X certificate in 1952.
• "Countess Dracula" is currently listed as a 18 certificate, but this is presumably because it hasn't been re-submitted. It's the equivalent of a PG in the USA.
• All the zombies listed are "voodoo zombies"; "Night of the Living Dead" (1968) isn't listed.
• Some of the films aren't even monster films and a few are actually comedies.
• I've managed to see just over half of the films on the list.
• The most obscure is surely "The Split" which I'm dying to see.
If I were going to "round up" the list to 50 by adding four more, which monsters would I choose? Funnily enough, they're all 1950s films (which seems to be under-represented in the original list):
1. "The Thing from Another World" (1951)
2. "The Beast from 20,000 Fathoms" (1953) not because it's a good film, but because of fond memories staying up late to watch it on TV during the school holidays.
3. "Godzilla" (1954) - how did this not get on the list?
4. "The Blob" (1958)
Two films that just missed the cut were "Forbidden Planet" (1956) and "The Day of the Triffids" (1962).

## Tuesday, 29 September 2015

### Anglo-Americanisms

Was it George Bernard Shaw who said that England and America are two countries divided by a common language?

Some words have very different meanings on either side of the pond. So while an Englishman puts his trunks in the boot of his car before driving off on holiday, an American would put his boots in the trunk.

Similarly, an Englishman puts his pants on before his shorts.

But my all-time favourite [or should that be favorite?] only really works when spoken aloud:
In England, we pay our bills with a cheque. In America, they pay their checks with a bill.

## Sunday, 27 September 2015

### Dominoes, Triangles, Squares and Tetris

My maternal grandfather was a coal miner from the North West of England. Unsurprisingly, for a working class man of that era, he played dominoes. Dominoes is often thought of as a children's game; but that's very unfair. Indeed there are championships played throughout the world. By subtle play, an experienced player can discern which tiles are in their opponents' hands. For a double-six deck, this is supposedly easy; for a double-nine deck, it is somewhat more difficult.

Dominoes are a favourite of mathematicians. Consider the number of tiles in a deck. If we arrange the tiles into rows where the value of the greater "side" is zero, one, two, three, four, five and six in turn, we get:

00
01 11
02 12 22
03 13 23 33
04 14 24 34 44
05 15 25 35 45 55
06 16 26 36 46 56 66

This obviously forms a triangle: the number of tiles in a complete double-n deck is a triangular number:

tiles in a double-n deck, D(n) = T(n+1) = (n + 1) × (n + 2) ÷ 2

For a double-six deck, D(6) = 28. For a double-nine deck, D(9) = 55.

In terms of construction, dominoes are two squares bolted together with a common edge:

There's no real choice about how we join the two squares, but when we go to three squares, there are two arrangements (ignoring rotations):

These are known as trominoes, not to be confused with triominoes.

When we get to four squares, we get the tetrominoes, The one-sided variants (allowing reflection) also being known as the Tetris pieces:

If we take the seven Tetris pieces and add the two trominoes and the single domino, we get ten pieces made up of a total of 36 squares (7×4 + 2×3 + 1×2):

Thirty-six is interesting: it is the first number after one that is both a square and a triangular number. Therefore, we can arrange the ten pieces above into both a square:

And into a triangle:

You can also arrange them into 2-by-18, 3-by-12 and 4-by-9 rectangles: