// uncomment to use square pattern instead of circles
//#define SQUARES

#define TAU 6.283185307179586

//
// GLSL textureless classic 3D noise "cnoise",
// with an RSL-style periodic variant "pnoise".
// Author:  Stefan Gustavson (stefan.gustavson@liu.se)
// Version: 2011-10-11
//
// Many thanks to Ian McEwan of Ashima Arts for the
// ideas for permutation and gradient selection.
//
// Copyright (c) 2011 Stefan Gustavson. All rights reserved.
// Distributed under the MIT license. See LICENSE file.
// https://github.com/stegu/webgl-noise
//

vec3 mod289(vec3 x)
{
  return x - floor(x * (1.0 / 289.0)) * 289.0;
}

vec4 mod289(vec4 x)
{
  return x - floor(x * (1.0 / 289.0)) * 289.0;
}

vec4 permute(vec4 x)
{
  return mod289(((x*34.0)+1.0)*x);
}

vec4 taylorInvSqrt(vec4 r)
{
  return 1.79284291400159 - 0.85373472095314 * r;
}

vec3 fade(vec3 t) {
  return t*t*t*(t*(t*6.0-15.0)+10.0);
}

// Classic Perlin noise
float cnoise(vec3 P)
{
  vec3 Pi0 = floor(P); // Integer part for indexing
  vec3 Pi1 = Pi0 + vec3(1.0); // Integer part + 1
  Pi0 = mod289(Pi0);
  Pi1 = mod289(Pi1);
  vec3 Pf0 = fract(P); // Fractional part for interpolation
  vec3 Pf1 = Pf0 - vec3(1.0); // Fractional part - 1.0
  vec4 ix = vec4(Pi0.x, Pi1.x, Pi0.x, Pi1.x);
  vec4 iy = vec4(Pi0.yy, Pi1.yy);
  vec4 iz0 = Pi0.zzzz;
  vec4 iz1 = Pi1.zzzz;

  vec4 ixy = permute(permute(ix) + iy);
  vec4 ixy0 = permute(ixy + iz0);
  vec4 ixy1 = permute(ixy + iz1);

  vec4 gx0 = ixy0 * (1.0 / 7.0);
  vec4 gy0 = fract(floor(gx0) * (1.0 / 7.0)) - 0.5;
  gx0 = fract(gx0);
  vec4 gz0 = vec4(0.5) - abs(gx0) - abs(gy0);
  vec4 sz0 = step(gz0, vec4(0.0));
  gx0 -= sz0 * (step(0.0, gx0) - 0.5);
  gy0 -= sz0 * (step(0.0, gy0) - 0.5);

  vec4 gx1 = ixy1 * (1.0 / 7.0);
  vec4 gy1 = fract(floor(gx1) * (1.0 / 7.0)) - 0.5;
  gx1 = fract(gx1);
  vec4 gz1 = vec4(0.5) - abs(gx1) - abs(gy1);
  vec4 sz1 = step(gz1, vec4(0.0));
  gx1 -= sz1 * (step(0.0, gx1) - 0.5);
  gy1 -= sz1 * (step(0.0, gy1) - 0.5);

  vec3 g000 = vec3(gx0.x,gy0.x,gz0.x);
  vec3 g100 = vec3(gx0.y,gy0.y,gz0.y);
  vec3 g010 = vec3(gx0.z,gy0.z,gz0.z);
  vec3 g110 = vec3(gx0.w,gy0.w,gz0.w);
  vec3 g001 = vec3(gx1.x,gy1.x,gz1.x);
  vec3 g101 = vec3(gx1.y,gy1.y,gz1.y);
  vec3 g011 = vec3(gx1.z,gy1.z,gz1.z);
  vec3 g111 = vec3(gx1.w,gy1.w,gz1.w);

  vec4 norm0 = taylorInvSqrt(vec4(dot(g000, g000), dot(g010, g010), dot(g100, g100), dot(g110, g110)));
  g000 *= norm0.x;
  g010 *= norm0.y;
  g100 *= norm0.z;
  g110 *= norm0.w;
  vec4 norm1 = taylorInvSqrt(vec4(dot(g001, g001), dot(g011, g011), dot(g101, g101), dot(g111, g111)));
  g001 *= norm1.x;
  g011 *= norm1.y;
  g101 *= norm1.z;
  g111 *= norm1.w;

  float n000 = dot(g000, Pf0);
  float n100 = dot(g100, vec3(Pf1.x, Pf0.yz));
  float n010 = dot(g010, vec3(Pf0.x, Pf1.y, Pf0.z));
  float n110 = dot(g110, vec3(Pf1.xy, Pf0.z));
  float n001 = dot(g001, vec3(Pf0.xy, Pf1.z));
  float n101 = dot(g101, vec3(Pf1.x, Pf0.y, Pf1.z));
  float n011 = dot(g011, vec3(Pf0.x, Pf1.yz));
  float n111 = dot(g111, Pf1);

  vec3 fade_xyz = fade(Pf0);
  vec4 n_z = mix(vec4(n000, n100, n010, n110), vec4(n001, n101, n011, n111), fade_xyz.z);
  vec2 n_yz = mix(n_z.xy, n_z.zw, fade_xyz.y);
  float n_xyz = mix(n_yz.x, n_yz.y, fade_xyz.x); 
  return 2.2 * n_xyz;
}

vec2 polar(vec2 p_rect) {
    return vec2(atan(p_rect.y, p_rect.x), length(p_rect));
}

vec2 rect(vec2 p_polar) {
    return vec2(p_polar.y * cos(p_polar.x), p_polar.y * sin(p_polar.x));
}

vec2 rotate(float angle, vec2 p) {
    return rect(polar(p) + vec2(angle, 0.0));
}

float map(float l0, float r0, float l1, float r1, float x) {
  return (x - l0) / (r0 - l0) * (r1 - l1) + l1;
}

vec2 map(vec2 l0, vec2 r0, vec2 l1, vec2 r1, vec2 p) {
  return vec2(map(l0.x, r0.x, l1.x, r1.x, p.x), map(l0.y, r0.y, l1.y, r1.y, p.y));
}

vec2 map(float l0, float r0, float l1, float r1, vec2 p) {
  return map(vec2(l0), vec2(r0), vec2(l1), vec2(r1), p);
}

vec3 rgb2hsl(vec3 rgb) {
    float r = rgb.r;
    float g = rgb.g;
    float b = rgb.b;
    float v, m, vm, r2, g2, b2;
    float h = 0.0;
    float s = 0.0;
    float l = 0.0;
    v = max(max(r, g), b);
    m = min(min(r, g), b);
    l = (m + v) / 2.0;
    if(l > 0.0) {
        vm = v - m;
        s = vm;
        if(s > 0.0) {
            s /= (l <= 0.5) ? (v + m) : (2.0 - v - m);
            r2 = (v - r) / vm;
            g2 = (v - g) / vm;
            b2 = (v - b) / vm;
            if(r == v) {
                h = (g == m ? 5.0 + b2 : 1.0 - g2);
            } else if(g == v) {
                h = (b == m ? 1.0 + r2 : 3.0 - b2);
            } else {
                h = (r == m ? 3.0 + g2 : 5.0 - r2);
            }
        }
    }
    h /= 6.0;
    return vec3(h, s, l);
}

vec3 hsl2rgb(vec3 hsl) {
    float h = hsl.x;
    float s = hsl.y;
    float l = hsl.z;
    float r = l;
    float g = l;
    float b = l;
    float v = (l <= 0.5) ? (l * (1.0 + s)) : (l + s - l*s);
    if(v > 0.0) {
        float m, sv;
        int sextant;
        float fract, vsf, mid1, mid2;
        m = l + l - v;
        sv = (v - m) / v;
        h *= 6.0;
        sextant = int(h);
        fract = h - float(sextant);
        vsf = v * sv * fract;
        mid1 = m + vsf;
        mid2 = v - vsf;
        if(sextant == 0) {
            r = v;
            g = mid1;
            b = m;
        } else if(sextant == 1) {
            r = mid2;
            g = v;
            b = m;
        } else if(sextant == 2) {
            r = m;
            g = v;
            b = mid1;
        } else if(sextant == 3) {
            r = m;
            g = mid2;
            b = v;
        } else if(sextant == 4) {
            r = mid1;
            g = m;
            b = v;
        } else if(sextant == 5) {
            r = v;
            g = m;
            b = mid2;
        }
    }
    return vec3(r, g, b);
}

vec3 hueshift(float dh, vec3 color) {
  vec3 hsl = rgb2hsl(color);
  hsl.x = fract(hsl.x + 1.0 + dh);
  return hsl2rgb(hsl);
}

vec2 aspect_correct(float aspect, vec2 p) {
  return aspect >= 1.0
    ? vec2(p.x, p.y/aspect)
    : vec2(p.x*aspect, p.y);
}

vec2 aspect_decorrect(float aspect, vec2 p) {
  return aspect >= 1.0
    ? vec2(p.x, p.y*aspect)
    : vec2(p.x/aspect, p.y);
}

bool checkers(vec2 d, vec2 uv) {
    bool px = mod(uv.x, 2.0*d.x) < d.x;
    bool py = mod(uv.y, 2.0*d.y) < d.y;
    return px != py;
}

bool circlers(vec2 d, float kr, vec2 uv) {
    float x = fract(uv.x/d.x);
    float y = fract(uv.y/d.y);

    return distance(vec2(x, y), vec2(0.5)) < 0.5*kr;
}

float quant(float n, float x) {
    return floor(n * x)/n;
}

float cool_pow(float b, float e) {
  return (b < 0.0 ? -1.0 : 1.0) * pow(abs(b), e);
}

vec3 shatter(float time, vec2 xy) {
    float aspect = iResolution.x / iResolution.y;
    vec2 uv = aspect_correct(aspect, map(vec2(0.0), iResolution.xy, vec2(-1.0), vec2(1.0), xy));
    vec2 uv_p = polar(uv);
    uv_p.x += uv_p.x < 0.0 ? TAU : 0.0;
    
    float turn = map(0.0, TAU, 0.0, 1.0, uv_p.x);
    
    float wob1 = sin(0.47*time) * sin(7.0*uv_p.x + 1.47*time);
    float wob2 = cos(-0.66*time)  * cos(5.0*uv_p.x - 1.21*time);
    vec2 uv1_p = uv_p * vec2(1.0, 1.0 + 0.3 * wob1);
    vec2 uv2_p = uv_p * vec2(1.0, 1.0 + 0.28 * wob2);
    vec2 uv1 = rect(uv1_p);
    vec2 uv2 = rect(uv2_p);

    vec2 uv1a = rotate(0.05*time, uv1);
    vec2 uv1b = rotate(-0.036*time, uv2);

    float hhh = 90.0;
    float ggg = 4.0;
    float n1 = quant(ggg, map(-1.0, 1.0, 0.0, 1.0, cnoise(vec3(uv1a * hhh, 0.4*time))));
    float n2 = quant(ggg, map(-1.0, 1.0, 0.0, 1.0, cnoise(vec3(-0.53*time, uv1b * hhh))));
    vec3 c1 = vec3(abs(uv1a.x), 1.0-abs(uv1a.x), 0.0) * n1;
    vec3 c2 = vec3(0.0, 0.0, abs(uv1b.y)) * n2;
    float kr = 0.5;

    #ifdef SQUARES
    bool p1 = checkers(vec2(0.012), uv1a + time * vec2(0.005, 0.01));
    bool p2 = checkers(vec2(0.016), uv1b + time * vec2(-0.004, 0.0065));
    #else
    bool p1 = circlers(vec2(0.012), kr, uv1a + time * vec2(0.005, 0.01));
    bool p2 = circlers(vec2(0.02), kr, uv1b + time * vec2(-0.004, 0.0065));    
    #endif
    
    vec3 color = (p1 != p2) ? c1 : c2;
    float sh = map(-1.0, 1.0, 0.0, 1.0, wob1) * map(-1.0, 1.0, 0.0, 1.0, wob2);
    color = hueshift(0.1*time + sh*0.3 + turn, color) * pow(sh, map(-1.0, 1.0, 0.01, 0.15, sin(0.12*time)));

    vec3 color_hsl = rgb2hsl(color);
    color_hsl.y *= map(0.0, 1.0, 0.25, 1.0, pow(sh, 1.0));
    color_hsl.z *= map(0.0, 1.0, 0.25, 1.0, pow(sh, 1.0));

    color = hsl2rgb(color_hsl);

    return color;
}

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    vec2 xy = fragCoord.xy;

    vec2 xy_p = polar(xy - 0.5*iResolution.xy);
    xy_p.x += 0.00008*xy_p.y * cool_pow(sin(-0.13*iTime), 1.3) * sin(0.89*pow(xy_p.y, 0.5) - 4.0*iTime);
    xy_p.x += 0.0005*sin(0.23*iTime)*xy_p.y;
    vec2 xy1 = rect(xy_p) + 0.5*iResolution.xy;
    
    vec3 color = shatter(iTime, xy1);
    fragColor = vec4(color*3.0, 0.0);
}