  // animated noise
vec4 NC0=vec4(0.0,157.0,113.0,270.0);
vec4 NC1=vec4(1.0,158.0,114.0,271.0);
vec4 WS=vec4(0.25,0.25,0.25,0.25);
// mix noise for alive animation, full source
vec4 hash4( vec4 n ) { return fract(sin(n)*1399763.5453123); }
float noise4q(vec4 x)
{
  vec4 n3 = vec4(0,0.25,0.5,0.75);
  vec4 p2 = floor(x.wwww+n3);
  vec4 b = floor(x.xxxx+n3) + floor(x.yyyy+n3)*157.0 + floor(x.zzzz +n3)*113.0;
  vec4 p1 = b + fract(p2*0.00390625)*vec4(164352.0, -164352.0, 163840.0, -163840.0);
  p2 = b + fract((p2+1.0)*0.00390625)*vec4(164352.0, -164352.0, 163840.0, -163840.0);
  vec4 f1 = fract(x.xxxx+n3);
  vec4 f2 = fract(x.yyyy+n3);
  f1=f1*f1*(3.0-2.0*f1);
  f2=f2*f2*(3.0-2.0*f2);
  vec4 n1 = vec4(0,1.0,157.0,158.0);
  vec4 n2 = vec4(113.0,114.0,270.0,271.0);
  vec4 vs1 = mix(hash4(p1), hash4(n1.yyyy+p1), f1);
  vec4 vs2 = mix(hash4(n1.zzzz+p1), hash4(n1.wwww+p1), f1);
  vec4 vs3 = mix(hash4(p2), hash4(n1.yyyy+p2), f1);
  vec4 vs4 = mix(hash4(n1.zzzz+p2), hash4(n1.wwww+p2), f1);
  vs1 = mix(vs1, vs2, f2);
  vs3 = mix(vs3, vs4, f2);
  vs2 = mix(hash4(n2.xxxx+p1), hash4(n2.yyyy+p1), f1);
  vs4 = mix(hash4(n2.zzzz+p1), hash4(n2.wwww+p1), f1);
  vs2 = mix(vs2, vs4, f2);
  vs4 = mix(hash4(n2.xxxx+p2), hash4(n2.yyyy+p2), f1);
  vec4 vs5 = mix(hash4(n2.zzzz+p2), hash4(n2.wwww+p2), f1);
  vs4 = mix(vs4, vs5, f2);
  f1 = fract(x.zzzz+n3);
  f2 = fract(x.wwww+n3);
  f1=f1*f1*(3.0-2.0*f1);
  f2=f2*f2*(3.0-2.0*f2);
  vs1 = mix(vs1, vs2, f1);
  vs3 = mix(vs3, vs4, f1);
  vs1 = mix(vs1, vs3, f2);
  float r=dot(vs1,vec4(0.25));
  return r*r*(3.0-2.0*r);
}


// rays of a star
float ringRayNoise(vec3 ray,vec3 pos,float r,float size,float anim)
{
    float b = dot(ray,pos);
    vec3 pr=ray*b-pos;
    float c=length(pr);
    pr=normalize(pr);
    float s=max(0.0,(1.0-size*abs(r-c)));
    float nd=noise4q(vec4(pr*1.0,-anim+c))*2.0;
    nd=pow(nd,2.0);
    float n=0.4;
    float ns=1.0;
    if (c>r) {
        n=noise4q(vec4(pr*10.0,-anim+c));
        ns=noise4q(vec4(pr*50.0,-anim*2.5+c*2.0))*2.0;
    }
    n=n*n*nd*ns;

    return pow(s,4.0)+s*s*n;
}



float sdSphere( vec3 p, float s )
{
    return length(p)-s;
}



//----------------------------------------------------------------------


vec2 opU( vec2 d1, vec2 d2 )
{
        return (d1.x<d2.x) ? d1 : d2;
}

const float S=2.5;
const float au=2.0;
const float t1=1./0.241,t2=1./0.615,t3=1.,t4=1./1.88,t5=1./11.86; // planets speed
const float r1=.38*au,r2=.72*au,r3=1.*au,r4=1.52*au,r5=5.20*au; // planets orbit radius
const float s1=S*.024,s2=S*.061,s3=S*.063,s4=S*.035,s5=S*0.713; // planets radius


//----------------------------------------------------------------------
const float rsun = 0.695*.3;

vec2 map( in vec3 pos )
{
    float a=iGlobalTime+15.;
    vec2 res = vec2( sdSphere( pos , rsun ), 100.0 );
         res = opU( res, vec2( sdSphere(    pos-vec3( r1*cos(t1*a),0., r1*sin(t1*a)), s1 ), 47.1 ) );
         res = opU( res, vec2( sdSphere(    pos-vec3( r2*cos(t2*a),0., r2*sin(t2*a)), s2 ), 10.3 ) );
         res = opU( res, vec2( sdSphere(    pos-vec3( r3*cos(t3*a),0., r3*sin(t3*a)), s3 ), 80.2 ) );
         res = opU( res, vec2( sdSphere(    pos-vec3( r4*cos(t4*a),0., r4*sin(t4*a)), s4 ), 57.5 ) );
//         res = opU( res, vec2( sdSphere(    pos-vec3( r5*cos(t5*a),0., r5*sin(t5*a)), s5 ), 60. ) );


    return res;
}

vec2 castRay( in vec3 ro, in vec3 rd)
{
    float tmin = 1.0;
    float tmax = 20.0;

    float precis = 0.0002;
    float t = tmin;
    float m = -1.0;
    for( int i=0; i<100; i++ )
    {
        vec3 vx  = ro+rd*t;
        vec2 res = map( vx );
        if( res.x<precis || t>tmax ) break;
        t += res.x;
        m = res.y;
    }

    if( t>tmax ) m=-1.0;
    return vec2( t, m );
}



vec3 calcNormal( in vec3 pos )
{
        vec3 eps = vec3( 0.001, 0.0, 0.0 );
        vec3 nor = vec3(
        map(pos+eps.xyy).x - map(pos-eps.xyy).x,
        map(pos+eps.yxy).x - map(pos-eps.yxy).x,
        map(pos+eps.yyx).x - map(pos-eps.yyx).x );
        return normalize(nor);
}





vec4 render( in vec3 ro, in vec3 rd )
{
    const vec4 colx = vec4(0.01,0.01,0.01,0.2);
    vec4 col = vec4(0.0);
    vec2 res = castRay(ro,rd);
    float t = res.x;
    float m = res.y;
    float t0 = -ro.y/rd.y;
    float l0 = length(ro + t0*rd);
    float S =.5 * t0*t0/iResolution.x;

    float d0 = 1.0-smoothstep(0.0,S,abs(l0-r1));
    d0      += 1.0-smoothstep(0.0,S,abs(l0-r2));
    d0      += 1.0-smoothstep(0.0,S,abs(l0-r3));
    d0      += 1.0-smoothstep(0.0,S,abs(l0-r4));
    if (l0<r4) col = colx;
    if (t0<t&&d0>0.0)
    {
     col.rgb = vec3(max(colx.r,d0*.1));
     col.a   = max(colx.a,d0);
    }
    else
     if( m>99. ) col = vec4(1., 1.,.75,1.0);
     else
      if( m>0. )
      {
        vec3 pos = ro + t*rd;
        vec3 nor = calcNormal( pos );
        vec3 ref = reflect( rd, nor );

        // material
        col.rgb = 0.45 + 0.3*sin( vec3(0.05,0.08,0.10)*(m-1.0) );
        col.a   = 1.0;

        // lighitng
        vec3  lig = normalize( -pos );
        float dif = clamp( dot( nor, lig ), 0.0, 1.0 );
        float spe = pow(clamp( dot( ref, lig ), 0.0, 1.0 ),3.0);


        vec3 lin = vec3(0.0);
        lin += 1.20*dif*vec3(1.00,0.85,0.55);
        lin += 1.20*spe*vec3(1.00,0.85,0.55)*dif;
        col.rgb = col.rgb*lin;
        col.rgb = mix( col.rgb, vec3(0.8,0.9,1.0), 1.0-exp( -0.002*t*t ) );
      }
    return clamp(col,0.0,1.0);
}

mat3 setCamera( in vec3 ro, in vec3 ta)
{
        vec3 cw = normalize(ta-ro);
        vec3 cu = normalize( cross(cw,vec3(0., 1.,0.)) );
        vec3 cv = normalize( cross(cu,cw) );
    return mat3( cu, cv, cw );
}

vec4 doit(in vec2 q)
{
    vec2 p = -1.0+2.0*q;
    p.x *= iResolution.x/iResolution.y;
    vec2 mo = iMouse.xy/iResolution.xy;
        // camera
    vec3 ro = vec3( .0, 3., 4. );
    vec3 ta = vec3( 0., 0., 1.2 );

        // camera-to-world transformation
    mat3 ca = setCamera( ro, ta);

    // ray direction
    vec3 rd = ca * normalize( vec3(p.xy,2.0) );

    // render
    vec4 col = render( ro, rd );

    float time=iGlobalTime*2.;
    float s3=ringRayNoise(rd,ro,rsun,1.0,time);
    col+=s3*vec4(1.0,1.0,.5,1.0);
    return col;
}

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    // render this with four sampels per pixel
    vec4 col0 = doit( (fragCoord.xy+vec2(0.0,0.0) )/iResolution.xy );
    vec4 col1 = doit( (fragCoord.xy+vec2(0.5,0.0) )/iResolution.xy );
    vec4 col2 = doit( (fragCoord.xy+vec2(0.0,0.5) )/iResolution.xy );
    vec4 col3 = doit( (fragCoord.xy+vec2(0.5,0.5) )/iResolution.xy );
    vec4 col = 0.25*(col0 + col1 + col2 + col3);

    col.rgb = pow( col.rgb, vec3(0.4545) );
    fragColor = col;
}

