// Based on Shanes' Fiery Spikeball https://www.shadertoy.com/view/4lBXzy (I think that his implementation is more understandable than the original :) )
// Relief come from Siggraph workshop by Beautypi/2015 https://www.shadertoy.com/view/MtsSRf
// License Creative Commons Attribution-NonCommercial-ShareAlike 3.0

uniform mat3 geo2loc;
uniform float blend;
uniform vec4 star; // xyz- star origin(m)  , w - time
uniform sampler2D tex0;
out vec4 fragColor;
in  vec3 dir;
in  vec2 texc;

const float internal_rad1=8.;
const float internal_rad2=256.;


// IQ's noise
float pn( in vec3 p )
{
    vec3 ip = floor(p);
    p = fract(p);
    p *= p*(3.0-2.0*p);
    vec2 uv = (ip.xy+vec2(37.0,17.0)*ip.z) + p.xy;
    uv = texture( tex0, (uv+ 0.5)/256.0).yx;
    return mix( uv.x, uv.y, p.z );
}

// FBM
float fpn(vec3 p) {
    return pn(p*.06125)*.57 + pn(p*.125)*.28 + pn(p*.25)*.15;
}

float rand(vec2 co){// implementation found at: lumina.sourceforge.net/Tutorials/Noise.html
	return fract(sin(dot(co*0.123,vec2(12.9898,78.233))) * 43758.5453);
}

float cosNoise( in vec2 p )
{
    return 0.5*( sin(p.x) + sin(p.y) );
}

const mat2 m2 = mat2(1.6,-1.2,
                     1.2, 1.6);

float sdTorus( vec3 p, vec2 t )
{
  return length( vec2(length(p.xz)-t.x*1.2,p.y) )-t.y;
}

float smin( float a, float b)
{
	float h = clamp( 0.5 + 0.5*(b-a), 0.0, 1.0 );
	return mix( b, a, h ) - h*(1.0-h);
}

float SunSurface( in vec3 pos )
{
    float h = 0.0;
    vec2 q = pos.xy*0.25;

    float s = 0.5;


    for( int i=0; i<6; i++ )
    {
        h += s*cosNoise( q );
        q = m2*q*0.85;
        q += vec2(2.41,8.13);
        s *= 0.48 + 0.2*h;
    }
    h *= 2.0;

    float d1 = pos.z - h;

    // rings
    vec3 r1 = mod(2.3+pos+1.0,10.0)-5.0;
    r1.z = pos.z-0.1 - 0.7*h + 0.5*sin( 3.0*star.w+pos.x + 3.0*pos.y);
    r1.xy=cos(pos.x)*r1.xy+vec2(r1.y, -r1.x);
    float d2 = sdTorus( r1.xyz, vec2(clamp(abs(pos.x/pos.y),0.7,2.5), 0.20) );
    return smin( d1, d2);
}


float map(vec3 p,bool linear)
{
   float n=fpn(p*50.+star.w*5.) * 0.45;
   float b;
   if (linear) b=SunSurface(geo2loc*p+vec3(10.0,10.0,-internal_rad2));
   else          b=length(p)-internal_rad1;
   return b+n;
}


// See "Combustible Voronoi"
// https://www.shadertoy.com/view/4tlSzl
vec3 firePalette(float i)
{
    float T = 3400. +1300.*i-blend*2000.; // Temperature range (in Kelvin).
    vec3 L = vec3(7.4, 5.6, 4.4); // Red, green, blue wavelengths (in hundreds of nanometers).
    L = pow(L,vec3(5.0)) * (exp(1.43876719683e5/(T*L))-1.0);
    return 1.0-exp(-5e8/L); // Exposure level. Set to "50." For "70," change the "5" to a "7," etc.
}

//----------------------------------
void main( )
{
   const float target=1.1;
   float internal_rad;
   float a=1.0-length(texc); if (a<=0.01) discard;
   float ds = length(star.xyz);
   bool linear=ds<target;
   float flash = 1.0-smoothstep(0.01,0.02,abs(target-ds));

   if (linear) internal_rad=internal_rad2;
   else        internal_rad=internal_rad1;
   vec3 rd  =  normalize(dir);   // rd: direction of the ray
   vec3 ro  = -star.xyz*internal_rad;  // ro: position on the ray
   float offset = max(1.0,(ds-2.)*internal_rad);

   // ld, td: local, total density
   // w: weighting factor
   float ld=0., td=0., w=0.;

   // t: length of the ray
   // d: distance function
   float d=1., t=offset;

   // Distance threshold.
   const float h = .1;

   // total color
   float tc = 0.;

   vec2 seed = rd.xy + fract(star.w);

   // rm loop
   for (int i=0; i<56; i++) {

      // Loop break conditions. Seems to work, but let me
      // know if I've overlooked something.
      if (t>10.*internal_rad+offset||td>(1.-1./200.)) break;

      // evaluate distance function
      d = map (ro+t*rd,linear);

      // check whether we are close enough (step)
      // compute local density and weighting factor
      //const float h = .1;
      ld = (h - d) * step(d, h);
      w = (1. - td) * ld;

      // accumulate color and density
      tc += w*w + 1./50.;  // Different weight distribution.
      td += w + 1./200.;

	  // dithering implementation come from Eiffies' https://www.shadertoy.com/view/MsBGRh
      // if (linear)
      // {
      //  d = max(d, 0.04);
      //  d=abs(d)*(1.+0.28*rand(seed*vec2(i)));
      // }
      // else
      // {
       // add in noise to reduce banding and create fuzz
       d=abs(d)*(.8+0.28*rand(seed*vec2(i)));
       d = max(d, 0.03);
      // }
      // step forward
      t += d*0.5;

   }

   // Fire palette.
   vec3 color = firePalette(tc);

   // if (linear)
   //  color *= 1.05 / exp( ld * 2.82 ) ;

   fragColor.rgb = pow( color+flash, vec3(2.2) );
   fragColor.a   = smoothstep(0.01,0.55,a);
}
