const float E = 0.0001;
const float FAR = 40.0;
const int STEPS = 64;

vec3 glow = vec3(0.0);

float speed = 2.0;
float tunnelrot = 1.4;
float tunneltwst = .6;
float camyaw = 0.4;
float morphspd = 0.25;


float noise(vec3 p)
{
	vec3 ip=floor(p);
    p-=ip; 
    vec3 s=vec3(41,157,353);
    vec4 h=vec4(0.,s.yz,s.y+s.z)+dot(ip,s);
    p=p*p*(3.-2.*p); 
    h=mix(fract(sin(h)*43758.5),fract(sin(h+s.x)*43758.5),p.x);
    h.xy=mix(h.xz,h.yw,p.y);
    return mix(h.x,h.y,p.z); 
}

float noise2(vec3 p)
{
	vec3 ip=floor(p);
    p-=ip; 
    vec3 s=vec3(11,313,701);
    vec4 h=vec4(0.,s.yz,s.y+s.z)+dot(ip,s);
    p=p*p*(3.-2.*p); 
    h=mix(fract(sin(h)*43758.5),fract(sin(h+s.x)*43758.5),p.x);
    h.xy=mix(h.xz,h.yw,p.y);
    return mix(h.x,h.y,p.z); 
}

void rot(inout vec2 p, float a) {
    p = cos(a)*p + sin(a)*vec2(p.y, -p.x);
}

float box (vec3 p, vec3 r){
    vec3 d = abs(p) - r;
    return length(max(d,0.0) + min(max(d.x, max(d.y, d.z)), 0.0));
}

float oct( vec3 p, float s)
{
  p = abs(p);
  return (p.x+p.y+p.z-s)*0.57735027;
}

vec3 twistz(vec3 p, float k) {
    float c = cos(k*p.z);
    float s = sin(k*p.z);
    mat2  m = mat2(c,-s,s,c);
    return vec3(m*p.xy,p.z);
}

float tun2(vec3 p){
    vec3 pp = p;
    vec3 t = vec3(1.) - abs(vec3(length(pp.xz),length(p.xy),1.0));
    return max(t.x,t.y)+.1;
    return min(max(t.x,t.y),0.0);
}

float tun(vec3 p){
    vec3 t = vec3(1.) - abs(twistz(p,tunneltwst));
    return min(t.x,t.y);
}

float op(vec3 p){
    float d1 = tun2(p);
    float d2 = tun(p);
    float d3 = mix(noise(p),noise2(p),sin(iTime));
    float t = clamp(iTime*0.1,0.0,1.);
    return mix(d2,d1,t)+d3;
}

float scene(vec3 p)
{
    vec3 pp = p;
    pp -= vec3(0.,0.,iTime*speed+1.5);
    for (int i = 0; i < 2; ++i)
    {
        rot(pp.xz,iTime);
        rot(pp.xy,iTime*2.);
        pp -= float(i)*vec3(0.1,0.1,0.0);
    }
    float a = oct(pp,.5);
    float c = box(pp,vec3(.2));
    float d =  mix(c,a,cos(3.+iTime*morphspd)*1.2);
    rot(p.yx,iTime*tunnelrot);
    float b = op(p);
    d = max(d,0.01);
    
    if (mod(p.z,noise2(p)) <= 1.) {
            glow += vec3(.0,.0,.1)*0.001/(0.01+abs(b));
        }
    glow += vec3(0.6,0.3,0.0)* 0.01 / (0.01+abs(d));
    float r = min(d,b);
    return r;
}

float march(in vec3 ro, in vec3 rd, out vec3 p)
{
    p = ro;
    float t = E;
    for (int i = 0; i < STEPS; ++i) {
        float d = scene(p);
        t += d;        
        if (d < E || t > FAR) {
            break;
        }
        p -= rd*d;
    }
    return t;
}


void mainImage( out vec4 fragColor, in vec2 fragCoord )
{

    vec3 cp = vec3(sin(iTime)*camyaw,.0,iTime*speed);
    vec3 ct = vec3(0,0,-1);
    
    vec2 uv = fragCoord.xy / iResolution.xy;
    vec2 q = -1.0+2.0*uv;
    q.x *= iResolution.x/iResolution.y;

    vec3 cf = normalize(ct-cp);
    vec3 cr = normalize(cross(vec3(0.0,1.0,0.0),cf));
    vec3 cu = normalize(cross(cf,cr));
    
    vec3 rd = normalize(mat3(cr,cu,cf)*vec3(q,radians(60.0)));
   
    vec3 col = vec3(1.);
    vec3 p = vec3(0);
    float t = march(cp,rd,p);
    if (t < FAR) {
            col = mix(vec3(.2,.2,.6),vec3(.3,.3,.8),step(0.9,fract(p.z)))*exp((cp.z-p.z)*0.3);       
    }
    col += glow;
    col = smoothstep(0.0,1.0,col);
    col *= (1.0-vec3(t/FAR));

    // Output to screen
    fragColor = vec4(col,1.0);
}