/*original https://www.shadertoy.com/view/lslyRn https://www.shadertoy.com/view/wly3Dy https://www.shadertoy.com/view/tllfRX */
#define iterations 17
#define formuparam 0.53

#define volsteps 20
#define stepsize 0.1

#define zoom   0.800
#define tile   0.850
#define speed  0.000

#define brightness 0.0015
#define darkmatter 0.300
#define distfading 0.730
#define saturation 0.850

#define NUM_LAYERS 8.
#define TAU 6.28318
#define PI 3.141592
#define Velocity .025 //modified value to increse or decrease speed, negative value travel backwards
#define StarGlow 0.025
#define StarSize 02.
#define CanvasView 20.


float Star(vec2 uv, float flare){
    float d = length(uv);
  	float m = sin(StarGlow*1.2)/d;
    float rays = max(0., .5-abs(uv.x*uv.y*1000.));
    m += (rays*flare)*2.;
    m *= smoothstep(1., .1, d);
    return m;
}

float Hash21(vec2 p){
    p = fract(p*vec2(123.34, 456.21));
    p += dot(p, p+45.32);
    return fract(p.x*p.y);
}
#define LAYERS 6.0

float hash21( vec2 id ){
   id = fract( id * vec2( 123.34 , 456.21 ) );
   id += dot( id , id + 425.321);
   return fract( id.x * id.y );
}

mat2 rotate( float a ){
    // Takes an angle in Radians
    float s = sin(a);
    float c = cos(a);

    // Returns a 2d Rotation Matrix
    return mat2(c, -s , s , c);
}
float star( vec2 uv , float flare , float rand ){
    // Simulating a light like shape
    float d = length(uv);
    float starLight = 0.52/d;

    vec2 rotUv = uv * rotate( 3.1415 / 4.0 );
    float rotRays = max( 1.0 , 1.0 - abs( rotUv.x * rotUv.y * 1000.0 ) ) * 1.05 * flare;

    float allLight =  starLight + rotRays;

    allLight *= smoothstep( 0.3 , 0.2 , d );
    return max( allLight , 0.0 );
}
vec3 StarLayer(vec2 uv){
     vec2 id = floor(uv);
    vec2 gv = fract(uv);

    vec3 col = vec3( 0 );

    for( int y=-1; y<=1; y++ ) {
        for( int x=-1; x<=1; x++ ) {
            vec2 offset = vec2(x , y);
            vec2 thisId = offset + id;

            float starRand = Hash21( thisId );
            float size = fract(starRand *405.32 );
            vec2 starPos = gv - offset - vec2( sin(starRand * 0.3423) * 0.5 + 0.5 , fract( starRand * 1.123 ) ) + 0.5;

            float thisStar = star( starPos , smoothstep( .9 , 1.0, 1.0 ) , starRand  );

            vec3 color = sin( vec3(0.02 , 0.01 , 0.09) * fract( starRand * 10.2 ) * 123.2 ) * 0.5 + 0.5;
            color = max( color, vec3( 0.04 , 1.0 , 0.4 ) );

            color *= vec3(0.1 , 0.6 , 1.0 + size );

            col += thisStar *0.01 * color;
        }
    }

    return col;
}


#define BURST
#define NUM_LAYERS2 10.

mat2 Rot(float a) {
    float s=sin(a), c=cos(a);
    return mat2(c,-s,s,c);
}

float Star2(vec2 uv, float a, float sparkle) {
    vec2 av1 = abs(uv);
 	vec2 av2 = abs(uv*Rot(a));
    vec2 av = min(av1, av2);

    vec3 col = vec3(0);
    float d = length(uv);
    float star = av1.x*av1.y;
    star = max(av1.x*av1.y, av2.x*av2.y);
    star = max(0., 1.-star*1e3);

    float m = min(5., 1e-2/d);

    return m+pow(star, 4.)*sparkle;
}

float Hash212(vec2 p) {
    p = fract(p*vec2(123.34,145.54));
    p += dot(p, p+45.23);
    return fract(p.x*p.y);
}

vec3 StarLayer2(vec2 uv, float t, float sparkle) {
    vec2 gv = fract(uv)-.5;
    vec2 id = floor(uv);
	vec3 col = vec3(0);

    #ifndef BURST
    t = 0.;
    #endif

    for(int y=-1; y<=1; y++) {
        for(int x=-1; x<=1; x++) {
            vec2 offs = vec2(x, y);
            float n = Hash21(id-offs);
			vec3 N = fract(n*vec3(10,100,1000));
            vec2 p = (N.xy-.5)*.7;

            float brightness2 = Star2(gv-p+offs, n*6.2831+t, sparkle);
            vec3 star = brightness2*vec3(.6+p.x, .4, .6+p.y)*N.z*N.z;



            star *= 1.+sin((t+n)*20.)*smoothstep(sin(t)*.5+.5, 1., fract(10.*n));

            float d = length(gv+offs);

            col += star*smoothstep(1.5, .8, d);
        }
    }
    return col;
}


void mainVR( out vec4 fragColor, in vec2 fragCoord, in vec3 ro, in vec3 rd )
{

	vec3 dir=rd;
	vec3 from=ro;

	//volumetric rendering
	float s=0.1,fade=1.;
	vec3 v=vec3(0.);
	for (int r=0; r<volsteps; r++) {
		vec3 p=from+s*dir*.5;
		p = abs(vec3(tile)-mod(p,vec3(tile*2.))); // tiling fold
		float pa,a=pa=0.;
		for (int i=0; i<iterations; i++) {
			p=abs(p)/dot(p,p)-formuparam;
            p.xy*=mat2(cos(iTime*0.1),sin(iTime*0.1),-sin(iTime*0.1),cos(iTime*0.1));// the magic formula
			a+=abs(length(p)-pa); // absolute sum of average change
			pa=length(p);
		}
		float dm=max(0.,darkmatter-a*a*.001); //dark matter
		a*=a*a; // add contrast
		if (r>6) fade*=1.-dm; // dark matter, don't render near
		//v+=vec3(dm,dm*.5,0.);
		v+=fade;
		v+=vec3(s,s*s,s*s*s*s)*a*brightness*fade; // coloring based on distance
		fade*=distfading; // distance fading
		s+=stepsize;
	}
	v=mix(vec3(length(v)),v,saturation); //color adjust
	fragColor = vec4(v*.03,1.);
}

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
	//get coords and direction
	vec2 uv=fragCoord.xy/iResolution.xy-.5;
    	vec2 uv2=fragCoord.xy/iResolution.xy-.5;
	uv.y*=iResolution.y/iResolution.x;
	vec3 dir=vec3(uv*zoom,1.);
	float time=iTime*speed+.25;




	float t = iTime*.3;

    float twirl = sin(t*.1);
    twirl *= twirl*twirl*sin(dot(uv,uv));
    uv2 *= Rot(-t*.2);

    uv2 *= 20.+sin(t*.05);

    vec3 col = vec3(0);
    float speed2 = -.2;

    speed2 = .1;
    float bla = sin(t+sin(t+sin(t)*.5))*.5+.5;
    float d = dot(uv2,uv2);

    float a = atan(uv2.x, uv2.y);

    float burst = sin(-iTime*.05);


    float stp = 1./NUM_LAYERS;

    for(float i=0.; i<1.; i+=stp) {
    	float lt = fract(t*speed+i);
        float scale = mix(10., .25, lt);
        float fade = smoothstep(0., .4, lt)*smoothstep(1., .95, lt);
        vec2 sv = uv2*scale+i*134.53+iTime;

        col += StarLayer2(sv, t, fade)*fade;
    }



    float burstFade = smoothstep(0., .02, abs(burst));
    float size = .9*sin(t)+1.;
    size = max(size, sqrt(size));
    float fade = size/d;
    col *= mix(1., fade, burstFade);
    col += fade*.2*vec3(1., .5, .1)*bla*burstFade;

    t*=1.5;

    a -= -iTime*.1;
    float rays = sin(a*5.+t*3.)-cos(a*7.-t);
    rays *= sin(a+t+sin(a*4.)*10.)*.5+.5;
    col += rays*bla*.1*burstFade;
    col += 1.-burstFade;

    col *= 4.;



    float t2 = iTime*Velocity;
    vec3 col2 = vec3(0);
    for(float i=0.; i<1.; i+=1./NUM_LAYERS){
        float depth = fract(i*t);
        float scale = mix(CanvasView, .5, depth);

        col += StarLayer(uv*scale+i*453.2-iTime*.05);}
	float a1=.5+iMouse.x/iResolution.x*2.;
	float a2=.8+iMouse.y/iResolution.y*2.;
	mat2 rot1=mat2(cos(a1),sin(a1),-sin(a1),cos(a1));
	mat2 rot2=mat2(cos(a2),sin(a2),-sin(a2),cos(a2));
	dir.xz*=rot1;
	dir.xy*=rot2;
	vec3 from=vec3(1.,.5,0.5);
	from+=vec3(time*2.,time,-2.);
	from.xz*=rot1;
	from.xy*=rot2;

	mainVR(fragColor, fragCoord, from, dir);
    fragColor+=vec4(col,1.);
}