🕗 💤 👨🏿‍🔬 Descriptografia do ray tracer do tamanho de um cartão postal 😒 🦇 🧒🏼

“Ele fez de novo!” - foi o que me ocorreu quando olhei para o verso do folheto da Pixar ^[1] , completamente preenchido com código. Um conjunto de construções e expressões foi assinado no canto inferior direito por ninguém menos que Andrew Kensler. Para quem não o conhece, direi: Andrew é um programador que inventou um ray tracer de 1337 bytes em 2009.

Desta vez, Andrew apresentou algo mais volumoso, mas com um resultado visual muito mais interessante. Desde que terminei de escrever meus livros negros do Game Engine sobre Wolf3D e DOOM , tive tempo de aprender o interior de seu código enigmático. E quase imediatamente, fiquei literalmente fascinado pelas técnicas descobertas nele. Eles eram muito diferentes do trabalho anterior de Andrew, baseado em um traçador de raios "padrão". Eu estava interessado em aprender sobre marchar com raios, recursos de geometria volumétrica construtiva, renderização / rastreamento de caminhos de Monte Carlo, bem como muitos outros truques que ele usou para espremer código em um pedaço de papel tão pequeno.

Código fonte

A frente do folheto é um anúncio para o departamento de recrutamento da Pixar. No verso, são impressos 2.037 bytes de código C ++, ofuscados para ocupar a menor superfície possível.

#include <stdlib.h> // card > pixar.ppm #include <stdio.h> #include <math.h> #define R return #define O operator typedef float F;typedef int I;struct V{F x,y,z;V(F v=0){x=y=z=v;}V(F a,F b,F c=0){x=a;y=b;z=c;}V O+(V r){RV(x+rx,y+ry,z+rz);}VO*(V r){RV(x*rx,y*r. y,z*rz);}FO%(V r){R x*r.x+y*r.y+z*rz;}VO!(){R*this*(1/sqrtf(*this%*this) );}};FL(F l,F r){R l<r?l:r;}FU(){R(F)rand()/RAND_MAX;}FB(V p,V l,V h){l=p +l*-1;h=h+p*-1;RL(L(L(lx,hx),L(ly,hy)),L(lz,hz));}FS(V p,I&m){F d=1\ e9;V f=p;fz=0;char l[]="5O5_5W9W5_9_COC_AOEOA_E_IOQ_I_QOUOY_Y_]OWW[WaOa_aW\ eWa_e_cWiO";for(I i=0;i<60;i+=4){V b=V(l[i]-79,l[i+1]-79)*.5,e=V(l[i+2]-79,l [i+3]-79)*.5+b*-1,o=f+(b+e*L(-L((b+f*-1)%e/(e%e),0),1))*-1;d=L(d,o%o);}d=sq\ rtf(d);V a[]={V(-11,6),V(11,6)};for(I i=2;i--;){V o=f+a[i]*-1;d=L(d,ox>0?f\ absf(sqrtf(o%o)-2):(o.y+=oy>0?-2:2,sqrtf(o%o)));}d=powf(powf(d,8)+powf(pz, 8),.125)-.5;m=1;F r=L(-L(B(p,V(-30,-.5,-30),V(30,18,30)),B(p,V(-25,17,-25),V (25,20,25))),B(V(fmodf(fabsf(px),8),py,pz),V(1.5,18.5,-25),V(6.5,20,25))) ;if(r<d)d=r,m=2;F s=19.9-py;if(s<d)d=s,m=3;R d;}IM(V o,V d,V&h,V&n){I m,s= 0;F t=0,c;for(;t<100;t+=c)if((c=S(h=o+d*t,m))<.01||++s>99)R n=!V(S(h+V(.01,0 ),s)-c,S(h+V(0,.01),s)-c,S(h+V(0,0,.01),s)-c),m;R 0;}VT(V o,V d){V h,n,r,t= 1,l(!V(.6,.6,1));for(I b=3;b--;){I m=M(o,d,h,n);if(!m)break;if(m==1){d=d+n*( n%d*-2);o=h+d*.1;t=t*.2;}if(m==2){F i=n%l,p=6.283185*U(),c=U(),s=sqrtf(1-c), g=nz<0?-1:1,u=-1/(g+nz),v=nx*ny*u;d=V(v,g+ny*ny*u,-ny)*(cosf(p)*s)+V( 1+g*nx*nx*u,g*v,-g*nx)*(sinf(p)*s)+n*sqrtf(c);o=h+d*.1;t=t*.2;if(i>0&&M(h +n*.1,l,h,n)==3)r=r+t*V(500,400,100)*i;}if(m==3){r=r+t*V(50,80,100);break;}} R r;}I main(){I w=960,h=540,s=16;V e(-22,5,25),g=!(V(-3,4,0)+e*-1),l=!V(gz, 0,-gx)*(1./w),u(gy*lz-gz*ly,gz*lx-gx*lz,gx*ly-gy*lx);printf("P\ 6 %d %d 255 ",w,h);for(I y=h;y--;)for(I x=w;x--;){V c;for(I p=s;p--;)c=c+T(e ,!(g+l*(xw/2+U())+u*(yh/2+U())));c=c*(1./s)+14./241;V o=c+1;c=V(cx/ox,c. y/oy,cz/oz)*255;printf("%c%c%c",(I)cx,(I)cy,(I)cz);}}// Andrew Kensler

Ele trabalha mesmo?

Com o código, há uma instrução para seu lançamento. A idéia é redirecionar a saída padrão para um arquivo. Por extensão, podemos assumir que o formato de saída é um formato de imagem de texto chamado NetPBM ^[2] .

  $ clang -o card2 -O3 raytracer.cpp
 $ time ./card> pixar.ppm

 2m58.524s reais
 usuário 2m57.567s
 sys 0m0.415s

Após dois minutos e cinquenta e oito segundos ^[3] , a seguinte imagem é gerada. É incrível como pouco código é necessário para isso.

Você pode extrair muito da imagem acima. Grit é um sinal óbvio de um "rastreador de caminho". Esse tipo de renderizador difere do rastreamento de raios, pois os raios não são rastreados de volta às fontes de luz. Nesse método, milhares de raios por pixel são emitidos das fontes e o programa os monitora, esperando que eles encontrem a fonte de luz. Essa é uma técnica interessante que, muito melhor do que o traçado de raios, pode lidar com a renderização de oclusão ambiental, sombras suaves, cáusticas e radiosidade.

Vamos dividir o código em partes

Passar entrada para o CLion formata o código (veja a saída aqui ) e o divide em partes / tarefas menores.

  #include <stdlib.h> // cartão> pixar.ppm 
  #include <stdio.h> 
  #include <math.h>

  #define R return 
  #define O operator 
  typedef float F; typedef int I;

  struct V {F x, y, z; V (F v = 0) {x = y = z = v;} V (F a, F b, F 
  c = 0) {x = a; y = b; z = c;} V O + (V r) {RV (x + rx, y + ry, z + rz);} VO * (V r) {RV ( x * rx, y * r. 
  y, z * rz);} FO% (V r) {R x * r.x + y * r.y + z * rz;} VO! () {R * this * (1 / sqrtf (* this% * isto) 
  );

  FL (F l, F r) {R l <r? L: r;} FU () {R (F) rand () / RAND_MAX;} FB (V p, V l, V h) {l = p 
  + l * -1; h = h + p * -1; RL (L (L (L (lx, hx), L (ly, hy)), L (lz, hz));}

  FS (Vp, I & m) {F d = 1 \ 
  e9; V f = p; fz = 0; char l [] = "5O5_5W9W5_9_COC_AOEOA_E_IOQ_I_QOUOY_Y_] OWW [WaOa_aW \ 
  eWa_e_cWiO "; para (I i = 0; i <60; i + = 4) {V b = V (l [i] -79, l [i + 1] -79) *. 5, e = V (l [ i + 2] -79, l 
  [i + 3] -79) *. 5 + b * -1, o = f + (b + e * L (-L ((b + f * -1)% e / (e% e), 0), 1)) * - 1; d = L (d, o% o);} d = sq \ 
  rtf (d); V a [] = {V (-11,6), V (11,6)}; para (I i = 2; i -;) {V o = f + a [i] * -1; d = L (d, ox> 0? F \ 
  absf (sqrtf (o% o) -2) :( o.y + = oy> 0? -2: 2, sqrtf (o% o)));} d = powf (powf (d, 8) + powf (pz , 
  8), 125) -. 5; m = 1; F r = L (-L (B (p, V (-30, -. 5, -30), V (30,18,30)), B (p, V (-25,17, -25), V 
  (25,20,25))), B (V (fmodf (fabsf (px), 8), py, pz), V (1,5,18,5, -25), V (6,5,20,25))) 
  ; if (r <d) d = r, m = 2; F s = 19,9-py; if (s <d) d = s, m = 3; R d;}

  IM (V o, V d, V e h, V e n) {I m, s = 
  0; F t = 0, c; para (; t <100; t + = c) se ((c = S (h = o + d * t, m)) <. 01 || ++ s> 99) R n =! V (S (h + V (0,01,0 
  ), s) -c, S (h + V (0, 01), s) -c, S (h + V (0,0, 0,01), s) -c), m; R 0;}

  VT (V o, V d) {V h, n, r, t = 
  1, l (! V (.6, .6,1)); para (I b = 3; b -;) {I m = M (o, d, h, n); se (! M) quebrar ; se (m == 1) {d = d + n * ( 
  n% d * -2); o = h + d * .1; t = t * .2;} se (m == 2) {F i = n% l, p = 6,283185 * U (), c = U (), s = sqrtf (1-c), 
  g = nz <0? -1: 1, u = -1 / (g + nz), v = nx * ny * u; d = V (v, g + ny * ny * u, -ny) * (cosf (p) * s) + V ( 
  1 + g * nx * nx * u, g * v, -g * nx) * (sinf (p) * s) + n * sqrtf (c); o = h + d * .1; t = t *. 2; se (i> 0 && M (h 
  + n * .1, l, h, n) == 3) r = r + t * V (500,400,100) * i;} se (m == 3) {r = r + t * V (50,80,100) ; break;}} 
  R r;

  I main () {I w = 960, h = 540, s = 16; V e (-22,5,25), g =! (V (-3,4,0) + e * -1), l =! V (gz, 
  0, -gx) * (1./w), u (gy * lz-gz * ly, gz * lx-gx * lz, gx * ly-gy * lx); printf ("P \ 
  6% d% d 255 ", w, h); para (I y = h; y -;) para (I x = w; x -;) {V c; para (I p = s; p- -;) c = c + T (e 
  , (g + l * (xw / 2 + U ()) + u * (yh / 2 + U ()))); c = c * (1./s) + 14. / 241; V o = c + 1; c = V (cx / ox, c. 
  y / oy, cz / oz) * 255; printf ("% c% c% c", (I) cx, (I) cy, (I) cz);}}

  // Andrew Kensler

Cada uma das seções é descrita em detalhes no restante do artigo:
■ - truques comuns, ■ - classe vetorial, ■ - código auxiliar, ■ - banco de dados, ■ - marcha de raios, ■ - amostragem, ■ - código principal.

Truques comuns com #define e typedef

Truques comuns estão usando #define e typedef para reduzir significativamente a quantidade de código. Aqui denotamos F = float, I = int, R = return e O = operador. A engenharia reversa é trivial.

Classe v

Em seguida, vem a classe V, que renomeei para Vec (embora, como veremos abaixo, ela também seja usada para armazenar canais RGB no formato float).

 struct Vec { float x, y, z; Vec(float v = 0) { x = y = z = v; } Vec(float a, float b, float c = 0) { x = a; y = b; z = c;} Vec operator+(Vec r) { return Vec(x + rx, y + ry, z + rz); } Vec operator*(Vec r) { return Vec(x * rx, y * ry, z * rz); } // dot product float operator%(Vec r) { return x * rx + y * ry + z * rz; } // inverse square root Vec operator!() {return *this * (1 / sqrtf(*this % *this) );} };

Observe que não há operador de subtração (-); portanto, em vez de escrever "X = A - B", "X = A + B * -1" é usado. A raiz quadrada inversa é útil mais tarde para normalizar os vetores.

Função principal

main () é o único caractere que não pode ser ofuscado porque é chamado pela função _start da biblioteca libc. Geralmente, vale a pena começar com isso, porque será mais fácil trabalhar dessa maneira. Demorei um pouco para descobrir o significado das primeiras letras, mas ainda assim consegui criar algo legível.

 int main() { int w = 960, h = 540, samplesCount = 16; Vec position(-22, 5, 25); Vec goal = !(Vec(-3, 4, 0) + position * -1); Vec left = !Vec(goal.z, 0, -goal.x) * (1. / w); // Cross-product to get the up vector Vec up(goal.y * left.z - goal.z * left.y, goal.z * left.x - goal.x * left.z, goal.x * left.y - goal.y * left.x); printf("P6 %d %d 255 ", w, h); for (int y = h; y--;) for (int x = w; x--;) { Vec color; for (int p = samplesCount; p--;) color = color + Trace(position, !(goal + left * (x - w / 2 + randomVal())+ up * (y - h / 2 + randomVal()))); // Reinhard tone mapping color = color * (1. / samplesCount) + 14. / 241; Vec o = color + 1; color = Vec(color.x / ox, color.y / oy, color.z / oz) * 255; printf("%c%c%c", (int) color.x, (int) color.y, (int) color.z); } }

Observe que literais flutuantes não contêm a letra "f" e a parte fracionária é descartada para economizar espaço. O mesmo truque é usado abaixo, onde a parte inteira é descartada (float x = 0,5). Também incomum é a construção "for" com uma expressão de iteração inserida dentro da condição de interrupção.

Esta é uma função principal bastante padrão para um rastreador de raio / caminho. Os vetores da câmera são definidos aqui e os raios são emitidos para cada pixel. A diferença entre o traçador de raios e o traçador de caminhos é que vários raios são emitidos por pixel no TP, que são levemente deslocados aleatoriamente. Então a cor obtida para cada raio em um pixel é acumulada em três canais de flutuação R, B, G. No final, é realizada a correção tonal do resultado do método Reinhardt.

A parte mais importante é sampleCount, que teoricamente pode ser definido como 1 para acelerar a renderização e a iteração. Aqui estão representações de amostra com valores de amostra de 1 a 2048.

Título de spoiler

32.

128

256

512

1024

2048

Código auxiliar

Outro trecho simples de código são as funções auxiliares. Nesse caso, temos uma função trivial min (), um gerador de valor aleatório no intervalo [0,1] e um boxTest muito mais interessante, que faz parte do sistema Construct Solid Geometry (CSG) usado para cortar o mundo. O CSG é discutido na próxima seção.

 float min(float l, float r) { return l < r ? l : r; } float randomVal() { return (float) rand() / RAND_MAX; } // Rectangle CSG equation. Returns minimum signed distance from // space carved by lowerLeft vertex and opposite rectangle // vertex upperRight. float BoxTest(Vec position, Vec lowerLeft, Vec upperRight) { lowerLeft = position + lowerLeft * -1; upperRight = upperRight + position * -1; return -min( min( min(lowerLeft.x, upperRight.x), min(lowerLeft.y, upperRight.y) ), min(lowerLeft.z, upperRight.z)); }

Funções da geometria volumétrica construtiva

Não há vértices no código. Tudo é feito usando as funções CSG. Se você não estiver familiarizado com elas, basta dizer que essas são funções que descrevem se a coordenada está dentro ou fora do objeto. Se a função retornar uma distância positiva, o ponto estará dentro do objeto. Uma distância negativa indica que o ponto está fora do objeto. Existem muitas funções para descrever objetos diferentes, mas, para simplificar, tomemos, por exemplo, uma esfera e dois pontos, A e B.

 // Signed distance point(p) to sphere(c,r) float testSphere(Vec p, Vec c, float r) { Vec delta = c - p; float distance = sqrtf(delta%delta); return radius - distance; } Vec A {4, 6}; Vec B {3, 2}; Vec C {4, 2}; float r = 2.; testSphere(A, C, r); // == -1 (outside) testSphere(B, C, r); // == 1 (inside)

A função testSphere () retorna -1 para o ponto A (ou seja, está fora) e 1 para B (ou seja, está dentro). Sinais a distâncias são apenas um truque, permitindo que você obtenha duas informações em vez de uma no caso de um único valor. Um tipo semelhante de função pode ser escrito para descrever um paralelogramo (é exatamente isso que é executado na função BoxTest).

  // Signed distance point(p) to Box(c1,c2) float testRectangle(Vec p, Vec c1, Vec c2) { c1 = p + c1 * -1; c2 = c2 + position * -1; return min( min( min(c1.x, c2.x), min(c1.y, c2.y)), min(c1.z, c2.z)); } Vec A {3, 3}; Vec B {4, 6}; Vec C1 {2, 2}; Vec C2 {5, 4}; testRectangle(A, C1, C2); // 1.41 (inside) testRectangle(B, C1, C2); // -2.23 (outside)

Agora vamos ver o que acontece se você virar o sinal do valor de retorno.

  // Signed distance point(p) to carved box(c1,c2) float testCarveBox(Vec p, Vec c1, Vec c2) { c1 = p + c1 * -1; c2 = c2 + position * -1; return -min( min( min(c1.x, c2.x), min(c1.y, c2.y)), min(c1.z, c2.z)); } Vec A {3, 3}; Vec B {4, 6}; Vec C1 {2, 2}; Vec C2 {5, 4}; testCarveBox(A, C1, C2); // == -1.41 (outside) testCarveBox(B, C1, C2); // == 2.23 (inside)

Agora não descrevemos um objeto sólido, mas declaramos o mundo inteiro sólido e cortamos um espaço vazio nele. As funções podem ser usadas como tijolos de construção, que quando combinados podem descrever formas mais complexas. Usando o operador de adição lógica (função min), podemos cortar um par de retângulos um acima do outro e o resultado será semelhante a este.

  // Signed distance point to room float testRoom(Vec p) { Vec C1 {2, 4}; Vec C2 {5, 2}; // Lower room Vec C3 {3, 5}; Vec C4 {4, 4}; // Upper room // min() is the union of the two carved volumes. return min(testCarvedBox(p, C1, C2), testCarvedBox(p, C3, C4)); } Vec A {3, 3}; Vec B {4, 6}; testRoom(A, C1, C2); // == -1.41 (outside) testRoom(B, C1, C2); // == 1.00 (inside)

Se você pensar bem, parece a sala que estamos estudando, porque a sala inferior é expressa exatamente dessa maneira - com a ajuda de dois paralelogramos cortados.

Agora, tendo dominado o poderoso conhecimento do CSG, podemos retornar ao código e considerar a função do banco de dados, que é a mais difícil de lidar.

 #define HIT_NONE 0 #define HIT_LETTER 1 #define HIT_WALL 2 #define HIT_SUN 3 // Sample the world using Signed Distance Fields. float QueryDatabase(Vec position, int &hitType) { float distance = 1e9; Vec f = position; // Flattened position (z=0) fz = 0; char letters[15*4+1] = // 15 two points lines "5O5_" "5W9W" "5_9_" // P (without curve) "AOEO" "COC_" "A_E_" // I "IOQ_" "I_QO" // X "UOY_" "Y_]O" "WW[W" // A "aOa_" "aWeW" "a_e_" "cWiO"; // R (without curve) for (int i = 0; i < sizeof(letters); i += 4) { Vec begin = Vec(letters[i] - 79, letters[i + 1] - 79) * .5; Vec e = Vec(letters[i + 2] - 79, letters[i + 3] - 79) * .5 + begin * -1; Vec o = f + (begin + e * min(-min((begin + f * -1) % e / (e % e), 0), 1) ) * -1; distance = min(distance, o % o); // compare squared distance. } distance = sqrtf(distance); // Get real distance, not square distance. // Two curves (for P and R in PixaR) with hard-coded locations. Vec curves[] = {Vec(-11, 6), Vec(11, 6)}; for (int i = 2; i--;) { Vec o = f + curves[i] * -1; distance = min(distance, ox > 0 ? fabsf(sqrtf(o % o) - 2) : (oy += oy > 0 ? -2 : 2, sqrtf(o % o)) ); } distance = powf(powf(distance, 8) + powf(position.z, 8), .125) - .5; hitType = HIT_LETTER; float roomDist ; roomDist = min(// min(A,B) = Union with Constructive solid geometry //-min carves an empty space -min(// Lower room BoxTest(position, Vec(-30, -.5, -30), Vec(30, 18, 30)), // Upper room BoxTest(position, Vec(-25, 17, -25), Vec(25, 20, 25)) ), BoxTest( // Ceiling "planks" spaced 8 units apart. Vec(fmodf(fabsf(position.x), 8), position.y, position.z), Vec(1.5, 18.5, -25), Vec(6.5, 20, 25) ) ); if (roomDist < distance) distance = roomDist, hitType = HIT_WALL; float sun = 19.9 - position.y ; // Everything above 19.9 is light source. if (sun < distance)distance = sun, hitType = HIT_SUN; return distance; }

Você pode ver aqui a função de "cortar" o paralelogramo, no qual apenas dois retângulos são usados para construir toda a sala (nosso cérebro faz o resto, representa paredes). A escada horizontal é uma função CSG um pouco mais complexa usando a divisão restante. E, finalmente, as letras da palavra PIXAR são compostas por 15 linhas com um par "origem / delta" e dois casos especiais para curvas nas letras P e R.

Marcha de raio

Tendo um banco de dados de funções CSG descrevendo o mundo, é suficiente ignorarmos todos os raios emitidos na função main (). A marcha de raio usa a função de distância. Isso significa que a posição de amostragem se move para frente uma distância até o obstáculo mais próximo.

 // Perform signed sphere marching // Returns hitType 0, 1, 2, or 3 and update hit position/normal int RayMarching(Vec origin, Vec direction, Vec &hitPos, Vec &hitNorm) { int hitType = HIT_NONE; int noHitCount = 0; float d; // distance from closest object in world. // Signed distance marching for (float total_d=0; total_d < 100; total_d += d) if ((d = QueryDatabase(hitPos = origin + direction * total_d, hitType)) < .01 || ++noHitCount > 99) return hitNorm = !Vec(QueryDatabase(hitPos + Vec(.01, 0), noHitCount) - d, QueryDatabase(hitPos + Vec(0, .01), noHitCount) - d, QueryDatabase(hitPos + Vec(0, 0, .01), noHitCount) - d) , hitType; // Weird return statement where a variable is also updated. return 0; }

A idéia de marchar com raios com base na distância é avançar uma distância até o objeto mais próximo. No final, o feixe se aproximará tanto da superfície que pode ser considerado um ponto de incidência.

Observe que a marcha de raios não retorna uma interseção verdadeira com a superfície, mas uma aproximação. É por isso que a marcha para no código quando d <0,01f.

Juntando tudo: amostragem

A investigação do traçador de caminhos está quase completa. Está faltando uma ponte que conecte a função main () ao ray marcher. Esta última parte, que renomei de “Trace”, é o “cérebro” no qual os raios saltam ou param, dependendo do que encontram.

 Vec Trace(Vec origin, Vec direction) { Vec sampledPosition, normal, color, attenuation = 1; Vec lightDirection(!Vec(.6, .6, 1)); // Directional light for (int bounceCount = 3; bounceCount--;) { int hitType = RayMarching(origin, direction, sampledPosition, normal); if (hitType == HIT_NONE) break; // No hit. This is over, return color. if (hitType == HIT_LETTER) { // Specular bounce on a letter. No color acc. direction = direction + normal * ( normal % direction * -2); origin = sampledPosition + direction * 0.1; attenuation = attenuation * 0.2; // Attenuation via distance traveled. } if (hitType == HIT_WALL) { // Wall hit uses color yellow? float incidence = normal % lightDirection; float p = 6.283185 * randomVal(); float c = randomVal(); float s = sqrtf(1 - c); float g = normal.z < 0 ? -1 : 1; float u = -1 / (g + normal.z); float v = normal.x * normal.y * u; direction = Vec(v, g + normal.y * normal.y * u, -normal.y) * (cosf(p) * s) + Vec(1 + g * normal.x * normal.x * u, g * v, -g * normal.x) * (sinf(p) * s) + normal * sqrtf(c); origin = sampledPosition + direction * .1; attenuation = attenuation * 0.2; if (incidence > 0 && RayMarching(sampledPosition + normal * .1, lightDirection, sampledPosition, normal) == HIT_SUN) color = color + attenuation * Vec(500, 400, 100) * incidence; } if (hitType == HIT_SUN) { // color = color + attenuation * Vec(50, 80, 100); break; // Sun Color } } return color; }

Eu experimentei um pouco com essa função para alterar o número máximo de reflexões de feixe permitidas. O valor "2" confere às letras uma cor Vantablack lacada surpreendentemente bela ^[4] .