// "Active Flowers" by dr2 - 2017
// License: Creative Commons Attribution-NonCommercial-ShareAlike 3.0 Unported License

vec3 HsvToRgb (vec3 c);
vec2 Rot2D (vec2 q, float a);
float SmoothBump (float lo, float hi, float w, float x);
float Hashfv2 (vec2 p);
float Fbm2 (vec2 p);
vec3 VaryNf (vec3 p, vec3 n, float f);

vec3 sunDir, qHit;
vec2 gId;
float tCur, dstFar, gSize, szFac, rotAng, rotAngF, dTwist;
int idObj;
const float pi = 3.14159;

#define SQRT3 1.7320508

vec2 PixToHex (vec2 p)
{
  vec3 c, r, dr;
  c.xz = vec2 ((1./SQRT3) * p.x - (1./3.) * p.y, (2./3.) * p.y);
  c.y = - c.x - c.z;
  r = floor (c + 0.5);
  dr = abs (r - c);
  r -= step (dr.yzx, dr) * step (dr.zxy, dr) * dot (r, vec3 (1.));
  return r.xz;
}

vec2 HexToPix (vec2 h)
{
  return vec2 (SQRT3 * (h.x + 0.5 * h.y), (3./2.) * h.y);
}

float ObjDf (vec3 p)
{
  vec3 q, qq;
  float dMin, d, a, r, h;
  dMin = dstFar / szFac;
  p.xz -= HexToPix (gId);
  p.xz = Rot2D (p.xz, rotAngF);
  p.z -= 0.1;
  p /= szFac;
  q = p;
  h = 1.5;
  d = 0.8 * max (length (Rot2D (q.xz - vec2 (0., 0.07), 4. * pi * dTwist * q.y / h) -
     vec2 (0., 0.03)) - (0.06 - 0.015 * q.y / h + 0.005 * sin (4. * pi * q.y)), abs (q.y - h) - h);
  if (d < dMin) { dMin = d;  idObj = 1;  qHit = q; }
  q = p;
  r = length (q.xz);
  h = 0.2 * (0.8 - r);
  d = 0.8 * max (max (r - 0.8, abs (q.y - h) - h), - q.y);
  if (d < dMin) { dMin = d;  idObj = 2;  qHit = q; }
  p.xz = Rot2D (p.xz, rotAng);
  qq = p;
  qq.y -= 3.1;
  qq.yz = Rot2D (qq.yz, 0.3 * pi + 0.1 * rotAng);
  q = qq;
  r = length (q.xz);
  d = max (r - 0.25, abs (q.y) - 0.1 * (1. - 1.5 * r));
  if (d < dMin) { dMin = d;  idObj = 3;  qHit = q; }
  d = length (q - vec3 (0., clamp (q.y, -0.2, 0.), 0.)) - 0.15;
  if (d < dMin) { dMin = d;  idObj = 4;  qHit = q; }
  for (float k = 0.; k < 2.; k ++) {
    q = qq;
    if (k > 0.) q.xz = Rot2D (q.xz, pi / 6.);
    a = atan (q.z, - q.x);
    q.xz = Rot2D (q.xz, 2. * pi * (floor (6. * a / (2. * pi) + 0.5) / 6.));
    q.x += 0.2;
    q.y -= 0.01 * (1. - 2. * k) + 0.4 * dot (q.xz, q.xz);
    q.z *= 2.5;
    d = 0.5 * max (length (q.xz + vec2 (0.4, 0.)) - 0.5, abs (q.y) - 0.03 * (1.2 - length (q.xz)));
    if (d < dMin) { dMin = d;  idObj = 5;  qHit = q; }
  }
  return szFac * dMin;
}

void SetFlwConf ()
{
  float s;
  s = length (vec3 (gId.xy, gId.x + gId.y)) / gSize;
  szFac = 0.7 * (1. - 0.2 * s * s + 0.2 * (Hashfv2 (17. * gId + 99.) - 0.5));
  s = Hashfv2 (23. * gId + 99.) - 0.5;
  dTwist = sign (s);
  s = max (0.1, abs (s)) * dTwist;
  rotAngF = 0.2 * pi * s;
  rotAng = 0.3 * pi * (sin (pi * s * (tCur + 10.)) + s) - rotAngF;
}

float ObjRay (vec3 ro, vec3 rd)
{
  vec3 vri, vf, hv, p;
  vec2 edN[3], pM, gIdP;
  float dHit, d, s;
  edN[0] = vec2 (1., 0.);
  edN[1] = 0.5 * vec2 (1., SQRT3);
  edN[2] = 0.5 * vec2 (1., - SQRT3);
  for (int k = 0; k < 3; k ++) edN[k] *= sign (dot (edN[k], rd.xz));
  vri = 1. / vec3 (dot (rd.xz, edN[0]), dot (rd.xz, edN[1]), dot (rd.xz, edN[2]));
  vf = (SQRT3/2.) - vec3 (dot (ro.xz, edN[0]), dot (ro.xz, edN[1]),
     dot (ro.xz, edN[2]));
  pM = HexToPix (PixToHex (ro.xz));
  gIdP = vec2 (-99.);
  dHit = 0.;
  for (int j = 0; j < 200; j ++) {
    hv = (vf + vec3 (dot (pM, edN[0]), dot (pM, edN[1]), dot (pM, edN[2]))) * vri;
    s = min (hv.x, min (hv.y, hv.z));
    p = ro + dHit * rd;
    gId = PixToHex (p.xz);
    if (length (vec3 (gId.xy, gId.x + gId.y)) <= gSize) {
      if (gId.x != gIdP.x || gId.y != gIdP.y) {
        gIdP = gId;
        SetFlwConf ();
      }
      d = ObjDf (p);
    } else d = dstFar;
    if (dHit + d < s) {
      dHit += d;
    } else {
      dHit = s + 0.005;
      pM += SQRT3 * ((s == hv.x) ? edN[0] : ((s == hv.y) ? edN[1] : edN[2]));
    }
    if (d < 0.0005 || dHit > dstFar) break;
  }
  return dHit;
}

vec3 ObjNf (vec3 p)
{
  vec4 v;
  vec3 e = vec3 (0.002, -0.002, 0.);
  v = vec4 (ObjDf (p + e.xxx), ObjDf (p + e.xyy), ObjDf (p + e.yxy), ObjDf (p + e.yyx));
  return normalize (vec3 (v.x - v.y - v.z - v.w) + 2. * v.yzw);
}

float ObjSShadow (vec3 ro, vec3 rd)
{
  vec3 p;
  vec2 gIdP;
  float sh, d, h;
  sh = 1.;
  gIdP = vec2 (-99.);
  d = 0.01;
  for (int j = 0; j < 30; j ++) {
    p = ro + rd * d;
    gId = PixToHex (p.xz);
    if (length (vec3 (gId.xy, gId.x + gId.y)) <= gSize) {
      if (gId.x != gIdP.x || gId.y != gIdP.y) {
        gIdP = gId;
        SetFlwConf ();
      }
      h = ObjDf (p);
      sh = min (sh, smoothstep (0., 0.05 * d, h));
      d += clamp (h, 0.05, 0.5);
    } else d += 0.2;
    if (sh < 0.05) break;
  }
  return 0.5 + 0.5 * sh;
}

vec3 ShowScene (vec3 ro, vec3 rd)
{
  vec3 col, vn;
  vec2 vf;
  float dstObj, a, sh, hue, spec, glit;
  bool isBg;
  gSize = 12.;
  isBg = true;
  spec = 0.1;
  glit = 0.;
  dstObj = ObjRay (ro, rd);
  if (dstObj < dstFar) {
    isBg = false;
    ro += dstObj * rd;
    gId = PixToHex (ro.xz);
    vn = ObjNf (ro);
    hue = Hashfv2 (gId * vec2 (17., 27.) + 0.5);
    if (idObj >= 3) a = atan (qHit.z, - qHit.x) / (2. * pi);
    if (idObj == 1) {
      col = vec3 (0.7, 0.4, 0.);
      vf = vec2 (64., 2.);
    } else if (idObj == 2) {
      col = mix (vec3 (0.4, 0.1, 0.), vec3 (0.1, 0.4, 0.), Fbm2 (16. * ro.xz));
      vf = vec2 (16., 10.);
    } else if (idObj == 3) {
      col = HsvToRgb (vec3 (mod (hue + 0.6, 1.), 1., 0.7)) *
         (0.7 + 0.3 * SmoothBump (0.05, 0.95, 0.01, mod (16. * a, 1.)));
      vf = vec2 (32., 1.);
    } else if (idObj == 4) {
      col = HsvToRgb (vec3 (mod (hue + 0.1, 1.), 0.7, 1.)) *
         (0.7 + 0.3 * SmoothBump (0.05, 0.95, 0.01, mod (8. * a, 1.)));
      vf = vec2 (64., 1.);
    } else if (idObj == 5) {
      col = HsvToRgb (vec3 (mod (hue, 1.), 0.7 + 0.3 * mod (7. * hue, 1.), 1.));
      col *= 1. - 0.1 * SmoothBump (0.05, 0.95, 0.01, mod (16. * a, 1.));
      spec = 0.05;
      glit = 500. * step (0.01, max (0., dot (vn, sunDir))) *
         pow (max (0., dot (sunDir, reflect (rd, vn))), 16.) *
         pow (1. - 0.6 * abs (dot (normalize (sunDir - rd), VaryNf (1000. * ro, vn, 2.))), 8.);
      vf = vec2 (128., 0.1);
    }
    vn = VaryNf (vf.x * qHit, vn, vf.y);
    sh = ObjSShadow (ro, sunDir);
  } else if (rd.y < 0.) {
    isBg = false;
    dstObj = - ro.y / rd.y;
    ro += dstObj * rd;
    vn = VaryNf (16. * ro, vec3 (0., 1., 0.), 5. * (1. - smoothstep (0.5, 1., dstObj / dstFar)));
    sh = (dstObj < dstFar) ? ObjSShadow (ro, sunDir) : 1.;
    col = mix (vec3 (0.4, 0.1, 0.), vec3 (0.1, 0.4, 0.), Fbm2 (16. * ro.xz));
  }
  if (isBg) col = vec3 (0., 0., 0.2);
  else {
    col = col * (0.2 + 0.8 * sh * max (dot (vn, sunDir), 0.)) + vec3 (1., 1., 0.5) * sh *
       (glit + spec * pow (max (dot (normalize (sunDir - rd), vn), 0.), 64.));
  }
  return pow (clamp (col, 0., 1.), vec3 (0.8));
}

void mainImage (out vec4 fragColor, vec2 fragCoord)
{
  mat3 vuMat;
  vec4 mPtr;
  vec3 ro, rd;
  vec2 canvas, uv, ori, ca, sa;
  float el, az, zmFac;
  canvas = iResolution.xy;
  uv = 2. * fragCoord.xy / canvas - 1.;
  uv.x *= canvas.x / canvas.y;
  tCur = iTime;
  mPtr = iMouse;
  mPtr.xy = mPtr.xy / canvas - 0.5;
  az = 0.1 * pi;
  el = -0.13 * pi;
  if (mPtr.z > 0.) {
    az += 2. * pi * mPtr.x;
    el += pi * mPtr.y;
  } else {
    az +=  pi * (SmoothBump (0.25, 0.75, 0.25, mod (0.01 * tCur, 1.)) - 0.5);
  }
  el = clamp (el, - 0.4 * pi, -0.1 * pi);
  ori = vec2 (el, az);
  ca = cos (ori);
  sa = sin (ori);
  vuMat = mat3 (ca.y, 0., - sa.y, 0., 1., 0., sa.y, 0., ca.y) *
          mat3 (1., 0., 0., 0., ca.x, - sa.x, 0., sa.x, ca.x);
  ro = vuMat * vec3 (0., 0., -45.);
  zmFac = 15. - 10. * SmoothBump (0.25, 0.75, 0.15, mod (az / pi, 1.));
  rd = vuMat * normalize (vec3 (uv, zmFac));
  dstFar = 100.;
  sunDir = vuMat * normalize (vec3 (0.5, 1., -1.));
  fragColor = vec4 (ShowScene (ro, rd), 1.);
}

vec3 HsvToRgb (vec3 c)
{
  vec3 p;
  p = abs (fract (c.xxx + vec3 (1., 2./3., 1./3.)) * 6. - 3.);
  return c.z * mix (vec3 (1.), clamp (p - 1., 0., 1.), c.y);
}

vec2 Rot2D (vec2 q, float a)
{
  return q * cos (a) + q.yx * sin (a) * vec2 (-1., 1.);
}

float SmoothBump (float lo, float hi, float w, float x)
{
  return (1. - smoothstep (hi - w, hi + w, x)) * smoothstep (lo - w, lo + w, x);
}

const float cHashM = 43758.54;

float Hashfv2 (vec2 p)
{
  return fract (sin (dot (p, vec2 (37., 39.))) * cHashM);
}
vec2 Hashv2v2 (vec2 p)
{
  vec2 cHashVA2 = vec2 (37., 39.);
  return fract (sin (vec2 (dot (p, cHashVA2), dot (p + vec2 (1., 0.), cHashVA2))) * cHashM);
}

float Noisefv2 (vec2 p)
{
  vec2 t, ip, fp;
  ip = floor (p);  
  fp = fract (p);
  fp = fp * fp * (3. - 2. * fp);
  t = mix (Hashv2v2 (ip), Hashv2v2 (ip + vec2 (0., 1.)), fp.y);
  return mix (t.x, t.y, fp.x);
}

float Fbm2 (vec2 p)
{
  float f, a;
  f = 0.;
  a = 1.;
  for (int i = 0; i < 5; i ++) {
    f += a * Noisefv2 (p);
    a *= 0.5;
    p *= 2.;
  }
  return f * (1. / 1.9375);
}

float Fbmn (vec3 p, vec3 n)
{
  vec3 s;
  float a;
  s = vec3 (0.);
  a = 1.;
  for (int i = 0; i < 3; i ++) {
    s += a * vec3 (Noisefv2 (p.yz), Noisefv2 (p.zx), Noisefv2 (p.xy));
    a *= 0.5;
    p *= 2.;
  }
  return dot (s, abs (n));
}

vec3 VaryNf (vec3 p, vec3 n, float f)
{
  vec3 e = vec3 (0.1, 0., 0.);
  float s = Fbmn (p, n);
  vec3 g = vec3 (Fbmn (p + e.xyy, n) - s,
     Fbmn (p + e.yxy, n) - s, Fbmn (p + e.yyx, n) - s);
  return normalize (n + f * (g - n * dot (n, g)));
}