//Based on: https://www.shadertoy.com/view/sdy3RW

#define m(p) (length(mod(p,4.)-2.)-2.)*z

void mainImage(out vec4 O, vec2 I)
{
    vec3 a = iResolution, b=iMouse.xyz,
    d = vec3(I+I-a.xy,a),
    p = vec3(I-I,iTime),c = p;

    d.xz *= mat2(sin(vec4(0,11,33,0)+3.*(b/a-.5*sign(b)).x));
    d /= length(d);

    float e = I.x/a.x*4.,
    z = exp2(ceil(++e));
    for(int i = 0; i++<999 && m(a=ceil(c*z+d*step(a,b))/z-.5/z)>.8; )
        a = sign(d),
        b = fract(-c*z*a),
        a *= (b+step(b,-b))/d,
        b = min(a.xxx,min(a.y,a.z)),
        c += d*b/z;

    O.rgb = mod(a*z*2.+.5,3.)/(9.+dot(p-c,p-c))*m(c);
    O *= 3.2-2.*pow(abs(mod(e+e,2.)-1.),3.);
}

/*
///Original 527 chars:

#define V vec3
#define m(p) (length(mod(p,4.)-2.)-2.)

void mainImage(out vec4 O, vec2 I)
{
    V a = iResolution, b=iMouse.xyz,
    c = 3.*(b/a-.5)*sign(b),
    d = V(a.y/.5,I+I-a.xy),
    p = d-d; p.x=iTime;
    d.xy *= mat2(sin(vec4(0,11,33,0)+c.x));
    d /= length(d);

    c = p;
    float e = I.x/a.x*4.,
    L = float(4<<int(e)),
    z = L;
    for(int i = 0; i++<999 && z<=L;)
    {
        a = sign(d);
        b = fract(-c*z*a);
        b += step(b,-b);
        a *= b/d;
        b = min(a.xxx,min(a.y,a.z));
        c += d*b/z;
        z *= m(ceil(c*z+d*step(a,b))/z-.5/z)*z<.8?2.:1.;
    }

    O.rgb = (mod(ceil(c*z/2.),3.)+.5)/V(8.+dot(p-c,p-c))*(m(c)*z);
    O *= 1.2-.8*pow(abs(mod(e+e,2.)-1.),3.);
}
*/