precision highp float;
uniform vec2 resolution;
uniform vec2 mouse;
uniform float time;
uniform sampler2D backbuffer;
out vec4 outColor;

//#define DEBUG

#ifdef DEBUG
  #define AA 0
#else
//-----------------GPUの処理能力に応じて変更-----------------
  #define AA 3
//-----------------------------------------------------------
#endif

// 乱数
#define hash(x) fract(sin(x) * 43758.5453123)

// オブジェクトの種類
#define SPHERE 0
#define WALL 1
#define STAIR_UP 2
#define STAIR_SIDE 3

// 文字のサイズ
#define CHAR_SIZE ivec2(4, 8)

const int maxDepth = 6;               // レイの経路の深さの最大値
const float PI = acos(-1.);           // 円周率
const float PI2 = PI * 2.;            // τ
const float EPS = 0.0001;             // パストレーシングなどに使う微小量
const float FAR = 1e3;                // 最大描画距離
const float stairWidth = 2.;          // 階段の横幅
const float tileGrooveWidth = 0.005;  // タイルの溝の幅
const float tileRaCurvature = 0.03;   // タイルの縁の曲率半径
const float lightSize = 0.4;          // 光源（球体）の半径
const vec3 lightColor = vec3(30.);    // 光源（球体）の発光色
vec3 lightPos = vec3(0);              // 光源（球体）の中心の座標
float pathSeed = 0.;                  // パストレーシングで使う乱数のシード

// ピクセルフォント
// 参考：美咲フォント "misaki_gothic_2nd_4x8"
// https://littlelimit.net/misaki.htm
const uint ch_2 = 0x4A2448E0u;
const uint ch_0 = 0x4AAEAA40u;
//const uint ch_4 = 0x26AAE220u;
const uint ch_4 = 0x26AAF220u;
const uint ch_C = 0x4A888A40u;
const uint ch_E = 0xE88C88E0u;
const uint ch_I = 0xE44444E0u;
const uint ch_L = 0x888888E0u;
const uint ch_M = 0xAEEEAAA0u;
const uint ch_N = 0xCAAAAAA0u;
const uint ch_O = 0x4AAAAA40u;
const uint ch_S = 0x4A842A40u;
const uint ch_T = 0xE4444440u;
const uint ch_W = 0xAAAEEEA0u;

// nの各ビットを出力
bool extract_bit(uint n, int b) {
  return bool(1U & (n >> b));
}

// 単一の文字を出力
bool sprite(uint spr, ivec2 I) {
  int bit = (CHAR_SIZE.x - 1 - I.x) + I.y * CHAR_SIZE.x;
  bool bounds = 0 <= I.x && I.x < CHAR_SIZE.x && 0 <= I.y && I.y < CHAR_SIZE.y;
  return bounds && extract_bit(spr, bit);
}

// 1Dの乱数（シードを更新）
float random() {
  return hash(pathSeed++);
}

// 3Dの乱数
float hash13(vec3 p) {
  return hash(dot(p, vec3(127.1, 311.7, 74.7)));
}

// HSVからRGBへの変換
vec3 hsv(float h, float s, float v) {
  vec3 res = fract(h + vec3(0, 2, 1) / 3.);
  res = clamp(abs(res * 6. - 3.) - 1., 0., 1.);
  res = (res - 1.) * s + 1.;
  return res * v;
}

// 画面上の座標からレイの方向を算出
vec3 rayDir(vec2 uv, inout vec3 dir, float fov) {
  dir = normalize(dir);
  vec3 u = abs(dir.y) < 0.999 ? vec3(0, 1, 0) : vec3(0, 0, 1);
  vec3 side = normalize(cross(dir, u));
  vec3 up = cross(side, dir);
  return normalize(uv.x * side + uv.y * up + dir / tan(fov / 360. * PI));
}

// ベクトルvを回転軸axで角度aだけ回転
vec3 rot3D(vec3 v, float a, vec3 ax) {
  ax = normalize(ax);
  return mix(dot(ax, v) * ax, v, cos(a)) - sin(a) * cross(ax, v);
}

// 被写界深度
void dof(inout vec3 ro, inout vec3 rd, vec3 dir, float L, float factor) {
  float r = sqrt(hash(random()));
  float a = hash(random()) * PI2;
  vec3 v = vec3(r * vec2(cos(a), sin(a)) * factor, 0);
  vec3 diro = vec3(0, 0, -1);
  float d = dot(dir, diro);
  float va = acos(d);
  vec3 vax = abs(d) < 0.999 ? cross(dir, diro) : dir.yzx;
  v = rot3D(v, va, vax);
  ro += v;
  rd = normalize(rd * L - v);
}

// 物体の色を算出
vec3 objColor(vec3 p, int type, bool isGroove, vec3 n) {
  if(isGroove) {
    return vec3(0.8);
  }
  
  if(type == SPHERE) {
    return vec3(0);
  }
  
  vec3 ID = floor(p * float(CHAR_SIZE.y) + n * 0.5);
  float h = hash13(mod(ID, 1e3) + float(type) * 0.1328); // 長時間経過時の精度対策でmod()
  if(ID.y + ID.z > h * 40. + 8.) {
    return vec3(0.9, 0.7, 0.5);
  }
  
  vec3 col = hsv(0.4 - h * 0.1, 0.9, 0.1 + h * 0.5);
  if(type != STAIR_SIDE) {
    return col;
  }
  
  ivec2 I = ivec2(ID.x, mod(ID.y, float(CHAR_SIZE.y * 5)));
  
  const uint[] str = uint[](ch_W, ch_E, ch_L, ch_C, ch_O, ch_M, ch_E, // len = 7
                            ch_T, ch_O, // len = 2
                            ch_S, ch_E, ch_S, ch_S, ch_I, ch_O, ch_N, ch_S, // len = 8
                            ch_2, ch_0, ch_2, ch_4); // len = 4
  const int numLine = 4;
  const int[numLine] lenArr = int[](7, 2, 8, 4);
  const int[numLine] posArr = int[](0, lenArr[0], lenArr[0] + lenArr[1], lenArr[0] + lenArr[1] + lenArr[2]);
  
  int line = numLine - I.y / CHAR_SIZE.y;
  if(line == numLine) {
    return col;
  }
  int len = lenArr[line];
  I.x += len * CHAR_SIZE.x / 2 - 1;
  int index = I.x / CHAR_SIZE.x;
  if(index < 0 || index >= len) {
    return col;
  }
  
  I.y = I.y % CHAR_SIZE.y;
  I.x -= index * CHAR_SIZE.x;
  index += posArr[line];
  
  return sprite(str[index], I) ? vec3(0.9) : col;
}

// エミッションの色
vec3 emission(int type) {
  return type == SPHERE ? lightColor : vec3(0);
}

// 3次元空間上の点ceを通り、normalを法線ベクトルとする平面とレイの交差判定関数
// （normalは正規化されていなくてもよい）
float plaIntersect(vec3 ro, vec3 rd, vec3 ce, vec3 normal) {
  float front1 = -dot(rd, normal);
  vec3 oc = ro - ce;
  float front2 = dot(oc, normal);
  return front2 / front1;
}

// 3次元空間上の点ceを通り、normalを法線ベクトルとする平面とレイの交差判定関数
// （normalは正規化されていなくてもよい）
float plaIntersect2(vec3 ro, vec3 rd, vec3 ce, vec3 normal) {
  float front1 = -dot(rd, normal);
  vec3 oc = ro - ce;
  float front2 = dot(oc, normal);
  if(front1 <= 0. || front2 <= 0.) {
    return -1.; // 平面の裏側が見えている場合、または平面がカメラの後方にある場合は負の値を返す
  }
  return front2 / front1;
}

// 球面とレイの交差判定関数
float sphIntersect2(vec3 ro, vec3 rd, vec3 ce, float ra) {
  vec3 oc = ro - ce;
  float c = dot(oc, oc) - ra * ra;
  if(c <= 0.) { // 球面の裏側は描画しない
    return -1.;
  }
  float b = dot(oc, rd);
  float h = b * b - c;
  if(h < 0.) {
    return -1.;
  }
  return -b - sqrt(h);
}

// カメラ（ro）から物体表面まで飛ばしたレイの長さ
float castRay(vec3 ro, vec3 rd) {
  float t = FAR;
  
  float w = 0.5 - step(0., rd.x);
  float tWall = plaIntersect2(ro, rd, vec3(-w * 2. * stairWidth, 0, 0), vec3(w, 0, 0));
  if(tWall > 0.) {
    t = min(t, tWall);
  }
  
  float tLight = sphIntersect2(ro, rd, lightPos, lightSize);
  if(tLight > 0.) {
    t = min(t, tLight);
  }
  
  float tBound = 0.;
  if(ro.y > -ro.z) {
    vec2 rp = ro.yz + t * rd.yz;
    if(rp.x > -rp.y) {
      return t;
    }
    tBound = plaIntersect2(ro, rd, vec3(0), vec3(0, 1, 1));
  }
  
  vec2 rp = ro.yz + tBound * rd.yz;
  float i = floor((rp.x - rp.y) * 0.5);
  vec3 ce = vec3(0, i, -1. - i);
  float tStairUp = plaIntersect2(ro, rd, ce, vec3(0, 1, 0));
  float tStairSide = plaIntersect2(ro, rd, ce, vec3(0, 0, 1));
  float tStair = (tStairUp > 0. && tStairSide > 0.) ? min(tStairUp, tStairSide) : max(tStairUp, tStairSide);
  
  float temp = t;
  if(tStair > 0.) {
    t = min(t, tStair);
  }
  rp = ro.yz + t * rd.yz;
  if(rp.x > -rp.y) {
    return temp;
  }
  
  return t;
}

// 被写界深度の焦点まで伸ばしたレイの長さ
float focusIntersect(vec3 ro, vec3 rd) {
  float t = FAR;
  // ro.yとro.zがともに整数のとき、この関数内でplaIntersectの代わりにplaIntersect2を使うとなぜかバグる
  float w = 0.5 - step(0., rd.x);
  float tWall = plaIntersect(ro, rd, vec3(-w * 2. * stairWidth, 0, 0), vec3(w, 0, 0));
  if(tWall > 0.) {
    t = min(t, tWall);
  }
  
  const float height = -0.5;
  float tBound = plaIntersect(ro, rd, vec3(0, height, 0), vec3(0, 1, 1));
  if(tBound > 0.) {
    t = min(t, tBound);
  }
  
  return t;
}

// 物体の種類
int objType(vec3 p) {
  if(abs(p.x) > stairWidth - EPS) {
    return WALL;
  }
  
  vec3 q = p - lightPos;
  const float r = lightSize + EPS;
  if(dot(q, q) < r * r) {
    return SPHERE;
  }
  
  p.yz = abs(fract(p.yz) - 0.5);
  if(p.y < p.z) {
    return STAIR_SIDE;
  }
  
  return STAIR_UP;
}

// 法線ベクトルを計算する（タイルの凹凸は含まない）
vec3 calcNormal(vec3 p, int type) {
  if(type == WALL) {
    return vec3(-sign(p.x), 0, 0);
  }
  
  if(type == STAIR_SIDE) {
    return vec3(0, 0, 1);
  }
  
  if(type == STAIR_UP) {
    return vec3(0, 1, 0);
  }
  
  return (p - lightPos) / lightSize;
}

// タイルの凹凸
float bumpFunc(vec3 p, vec3 n) {
  p = abs(p);
  vec2 uv = n.y > 0.5 ? p.zx : n.z > 0.5 ? p.xy : p.yz;
  if(uv.y > uv.x) {
    uv = uv.yx;
  }
  uv.y -= tileGrooveWidth + tileRaCurvature;
  const float ra2 = tileRaCurvature * tileRaCurvature;
  if(uv.y > 0.) {
    return EPS + ra2;
  }
  
  if(uv.y < -tileRaCurvature) {
    return EPS;
  }
  
  float h = -uv.y * uv.y;
  return EPS + ra2 + h;
}

// 参考: "Maze Lattice" by Shane
// https://www.shadertoy.com/view/llGGzh
// タイルの凹凸ベクトル
// コードは改変済み
vec3 bumpMap(vec3 p, vec3 n, out bool isGroove) {
  const float bumpFactor = 0.5;
  const vec2 e = vec2(EPS, 0);
  
  const float interval = 1. / 8.;
  const float h = interval * 0.5;
  p = mod(p - h, interval) - h;
  float ref = bumpFunc(p, n);
  isGroove = (ref == EPS);
  vec3 grad = (vec3(bumpFunc(p - e.xyy, n),
                    bumpFunc(p - e.yxy, n),
                    bumpFunc(p - e.yyx, n)) - ref) / e.x;
  grad /= 2. * sqrt(ref); // sqrt(x) の微分
  grad -= n * dot(n, grad);
  return normalize(n + grad * bumpFactor);
}

// 三角波 周期:2, 振幅:1
float triWave(float x) {
  x -= 0.5;
  x *= 0.5;
  float res = abs(fract(x) - 0.5) - 0.25;
  return res * 4.;
}

// 範囲[-1, 1.]の間で等間隔にn個の値を取る滑らかな階段状のノイズ
float stepNoise(float x, float n) {
  const float factor = 0.3;
  float i = floor(x);
  float f = x - i; // fract(x)
  float u = smoothstep(0.5 - factor, 0.5 + factor, f);
  float res = mix(floor(hash(i) * n), floor(hash(i + 1.) * n), u);
  res /= (n - 1.) * 0.5;
  return res - 1.;
}

// 光源（球体）の位置を動かす
void moveLightPos() {
  float T = time - 0.7829;
  lightPos.y += abs(sin(T * PI2)) * 1.3;
  T += triWave(T / 4.) * 2.; // 時間停止のために三角波を足す
  lightPos += vec3(0, 1, -1) * T;
  lightPos.x += cos(T * PI) * 0.5;
}

// カメラの位置・向きを動かす
void moveCamera(inout vec3 ro, inout vec3 dir) {
  ro.yz += vec2(1, -1) * time;
  ro.y += abs(sin(ro.z * PI)) * 0.13;
  ro += sin(vec3(5, 7, 9) * time) * 0.01;
  dir.x += stepNoise(time * 0.4, 3.) * 0.8;
  dir.y += stepNoise(time * 0.4 - 500., 3.) * 0.4;
  dir += sin(vec3(7, 9, 5) * time - 500.) * 0.01;
  dir = normalize(dir);
}

// vを天頂として極座標phi, thetaだけ回転させたベクトル
vec3 jitter(vec3 v, float phi, float sinTheta, float cosTheta) {
  vec3 zAxis = normalize(v);
  vec3 xAxis = normalize(cross(zAxis.yzx, zAxis));
  vec3 yAxis = cross(zAxis, xAxis);
  return (xAxis * cos(phi) + yAxis * sin(phi)) * sinTheta + zAxis * cosTheta;
}

// パストレーシングをする
// 参考: "パストレーシング - Computer Graphics - memoRANDOM" by Shocker_0x15
// https://rayspace.xyz/CG/contents/path_tracing/
// 参考: "GLSL smallpt" by Zavie
// https://www.shadertoy.com/view/4sfGDB
// コードは改変済み
vec3 pathTrace(vec3 ro, vec3 rd) {
  vec3 acc = vec3(0);
  vec3 mask = vec3(1);
  
  for(int depth = 0; depth < maxDepth; depth++) {
    float t = castRay(ro, rd);
    if(t >= FAR) {
      break;
    }
    ro += t * rd;
    int type = objType(ro);
    vec3 n0 = calcNormal(ro, type);
    vec3 n = n0;
    bool isGroove = false;
    if(type != SPHERE) {
      n = bumpMap(ro, n0, isGroove);
    }
    
    ro += n0 * EPS;
    vec3 objC = objColor(ro, type, isGroove, n0);
    vec3 objE = emission(type);
    
    vec3 e = vec3(0);
    vec3 l0 = lightPos - ro;
    float cosA_max = sqrt(1. - lightSize * lightSize / dot(l0, l0));
    float cosA = mix(cosA_max, 1., random());
    vec3 l = jitter(l0, random() * PI2, sqrt(1. - cosA * cosA), cosA);
    
    float tl = castRay(ro, l);
    type = objType(ro + tl * l);
    
    if(type == SPHERE) {
      vec3 em = emission(type);
      float omega = PI2 * (1. - cosA_max);
      e += (em * max(dot(l, n), 0.) * omega) / PI;
    }
    
    const float rough = 0.01;
    vec3 nl = reflect(rd, n);
    //if(isGroove) {rough = 0.9; nl = n;}
    
    const float c = 1.; // GlossyなマテリアルのNext Event Estimationについて要検証
    //float c = rough;
    acc += mask * (objE + objC * e * c);
    mask *= objC;
    
    float ur = random() * rough; // 要検証
    rd = jitter(nl, random() * PI2, sqrt(ur), sqrt(1. - ur));
  }
  
  return acc;
}

void main(){
  vec3 col = vec3(0); // 画面上の色
  
  lightPos = vec3(0, lightSize, 0.);
  moveLightPos(); // 光源（球体）を動かす
  
  vec3 ro = vec3(0, 0, 5); // カメラの位置（レイの原点）
  vec3 dir = vec3(0, 0, -1); // カメラの向き
  moveCamera(ro, dir); // カメラの位置と向きを動かす
  float fov = 60.; // FOV（視野角）範囲：(0., 180.)
  
  float L = focusIntersect(ro, dir); // 被写界深度の焦点まで伸ばしたレイの長さ
  
  for(int m = 0; m < AA; m++) {
    for(int n = 0; n < AA; n++) {
      // アンチエイリアシング
      vec2 of = vec2(m, n) / float(AA) - 0.5; // オフセット
      vec2 uv = ((gl_FragCoord.xy + of) * 2. - resolution) / min(resolution.x, resolution.y); // 画面上の座標を正規化
      pathSeed += hash13(vec3(uv + of, fract(time * 0.1)));
      
      vec3 rd = rayDir(uv, dir, fov); // レイの向き
      
      // 被写界深度
      vec3 ros = ro;
      vec3 rds = rd;
      dof(ros, rds, dir, L, 0.15);
      
      col += pathTrace(ros, rds); // パストレーシング
    }
  }
  col /= float(AA * AA); // 色の平均をとる
  
  //if(col.r > 1. || col.g > 1. || col.b > 1.) col = vec3(1, 0, 0); // saturationチェック
  col = clamp(col, 0., 1.);
  
  col = pow(col, vec3(1. / 2.2)); // ガンマ補正
  
  // 口径食（vignetting）
  vec2 p = gl_FragCoord.xy / resolution;
  col *= 0.5 + 0.5 * pow(16.0 * p.x * p.y * (1.0 - p.x) * (1.0 - p.y), 0.5);
  
  // パストレーシングで得られた色の分散を低減するために前のフレームの色を合成する
  //col = mix(col, texture(backbuffer, gl_FragCoord.xy / resolution).rgb, 0.2);
  
#ifdef DEBUG
  col = vec3(0);
  vec2 uv = (gl_FragCoord.xy * 2. - resolution) / min(resolution.x, resolution.y); // 画面上の座標を正規化
  vec3 rd = rayDir(uv, dir, fov); // レイの向き
  float t = castRay(ro, rd);
  ro += t * rd;
  int type = objType(ro);
  vec3 n0 = calcNormal(ro, type);
  bool isGroove = false;
  vec3 n = n0;
  if(type != SPHERE) {
    n = bumpMap(ro, n0, isGroove);
  }
  col += exp(-t * t * 0.03) * (n * 0.5 + 0.5);
#endif
  
  outColor = vec4(col, 1);
}