/* Creative Commons Licence Attribution-NonCommercial-ShareAlike 
   phreax 2022
*/

#define PI 3.141592
#define TAU (2.*PI)
#define SIN(x) (sin(x)*.5+.5)
#define PHI 1.618033988749895
#define hue(v) ( .6 + .6 * cos( 6.3*(v) + vec3(0,23,21) ) )
//#define KALEIDOSCOPIC 2
//#define POLARKALEI 
#define THICKNESS .1
#define AA 0
#define BUMP_EPS 0.004
//#define KALEI


float tt, g_mat;

mat2 rot(float a) { return mat2(cos(a), -sin(a), sin(a), cos(a)); }

float fadeInOut(float t1, float t2, float fadeTime) {
    return smoothstep(t1, t1+fadeTime, iTime)-smoothstep(t2-fadeTime, t2, iTime);
}

// by Nusan
float curve(float t, float d) {
  t/=d;
  return mix(floor(t), floor(t)+1., pow(smoothstep(0.,1.,fract(t)), 10.));
}

////////------------ RAINBOW TEMPLE ------------///////////

vec3 transformTemple(vec3 p) {
    vec2 uv = p.xy;
   
    uv *= 2.;
       // uv*= 1.5 + sin(.4*iTime);

        
    #ifdef KALEI
    float N = 6.0;
    float t = atan(uv.x,-uv.y);
    t = mod(t+PI/N,2.0*PI/N)-PI/N;
    uv = length(uv)*vec2(cos(t),sin(t));
    uv = abs(uv)-0.25;
    uv = abs(uv)-0.24;
    #endif

  p.xy = uv;
 /*
  if (iMouse.x > 0.0) {
    // Full range of rotation across the screen.
    float phi = (2.0*iMouse.x-iResolution.x)/iResolution.x*PI;
    float theta = (2.0*iMouse.y-iResolution.y)/iResolution.y*PI;
    p.yz = rotate(p.yz,-theta);
    p.zx = rotate(p.zx,phi);
  }
  */
   //p.yz *= rot(PI*curve(tt, 8.));
   //p.xz *= rot(PI+PI*curve(tt, 8.));
  return p;
}

vec3 temple(vec3 ro, vec3 rd) {
  vec3 d = transformTemple(rd);
  vec3 p = ro;
  p.z += 0.2*iTime;
  vec3 q, p0 = p;
  vec3 lightdir = vec3(1,1,1);
  
  //lightdir.xz = rotate(lightdir.xz,0.123*iTime);
  for (int i = 0; i < 80; i++) {
  
    float lfactor = 0.11;
    float overlap = 0.65;
    //overlap += 0.1*sin(0.1*iTime);
    float floorheight = mix(0.2, .2, SIN(.3*tt));
    
    float s = overlap-length(fract(p-.5)-0.5);
    
    
   // s = min(s, p.y +floorheight);
    p += lfactor*d*s;
 
    if (i == 40) {
     
      // We must have hit the surface by now.
      q = p; // Remember position
      p -= d*.01; // Back up a little to stop self-shadowing
      d = abs(lightdir-.5); // Light direction (for shadows)
    }
 
  }
  
  float z = length(q-p0);
  ivec3 u = ivec3(abs(q)*5e2); // u in [0..255]
  u = u.yzx^u.zxy;             // Mix up for chequer effect
  u &= 255;
  vec3 col = vec3(u)/255.0 - vec3(0.082,0.118,0.200);
  col *= .6*(length(p-q) + .7);
  col += 0.05*(length(p-p0));
  
  p = 8.*(abs(q-.5));
  
  p.z *= 2.;
  
  p = abs(mod(p, 0.125)-0.0125);
    
  float x = min(p.z,min(p.x, p.y))/0.03125;
  

  col += .2*(1.-smoothstep(.1, 0.1+fwidth(x), x)) * smoothstep(1.3, 0.5, z);;

  col = clamp(col, vec3(0), vec3(1));
  col *= mix(1., .8, (1.5-pow(dot(d*2., d*2.), .5)));
  col = pow(col, vec3(1.5));
  return col;
}

// zucconis spectral palette https://www.alanzucconi.com/2017/07/15/improving-the-rainbow-2/
vec3 bump3y (vec3 x, vec3 yoffset)
{
    vec3 y = 1. - x * x;
    y = clamp((y-yoffset), vec3(0), vec3(1));
    return y;
}


vec3 spectral_zucconi6(float x) {
    x = fract(x);
    const vec3 c1 = vec3(3.54585104, 2.93225262, 2.41593945);
    const vec3 x1 = vec3(0.69549072, 0.49228336, 0.27699880);
    const vec3 y1 = vec3(0.02312639, 0.15225084, 0.52607955);
    const vec3 c2 = vec3(3.90307140, 3.21182957, 3.96587128);
    const vec3 x2 = vec3(0.11748627, 0.86755042, 0.66077860);
    const vec3 y2 = vec3(0.84897130, 0.88445281, 0.73949448);
    return bump3y(c1 * (x - x1), y1) + bump3y(c2 * (x - x2), y2) ;
}


float rect( vec2 p, vec2 b, float r ) {
    vec2 d = abs(p) - (b - r);
    return length(max(d, 0.)) + min(max(d.x, d.y), 0.) - r;
}

vec3 kalei(vec3 p) {
    float iter = 3.;
    p.yz *= rot(tt);
    p.xy = abs(p.xy)+.9;
    
    #ifdef POLARKALEI
    p.xy = vec2(length(p.xy)*.8, atan(p.y, p.x));
    #endif
   
    float s= 1.;
    for(float i=0.; i< iter; i++) {
        p = abs(p);  
        p.xz *= rot(i/iter*PI+.2*tt);
        p.xz += .2 - .1*i + sin(p.y*.2) + .2*exp(.01+p.z*.2);
        p.z -= .2*sin(p.z*.5*i/iter*2.*PI);
        p -= .2;
         
    }
    
    return clamp(p, -1e5, 1e5);

}

// fold space for toroid
vec3 transform(vec3 p) {
   
    p.xy *= 1.2 + 1.2*sin(.5*tt);

    // n-fold symmetry by mla
    float N = 6.0;
    float t = atan(p.x,-p.y);
    t = mod(t+PI/N,2.0*PI/N)-PI/N;
    p.xy = length(p.xy)*vec2(cos(t),sin(t));
    p.xy = abs(p.xy)-0.25;
    p.xy = abs(p.xy)-0.24;
    
    p.yz *= rot(PI*.5);
    p.xz *= rot(.2*tt);
 
    #ifdef KALEIDOSCOPIC
    p = kalei(p);
    #endif

    float r = 1.9;
    vec2 cp = vec2(length(p.xz)-r, p.y);
    

    float rev = 2.5;
    float a = atan(p.z, p.x);
    
 
    cp *= rot(rev*a+.3*tt);
    cp= abs(cp) - mix(.1, 1., SIN(tt));
    cp *= rot(.5*tt);


    return vec3(cp, p.z);
}


float map(vec3 p) {   
  
    vec3 bp = p;
    float edge = 0.01;
   
    vec2 cp = transform(p).xy;
    float dr = rect(cp.xy, vec2(.8, THICKNESS), edge);
    
    // g_mat = db < dr ? 1. : 0.;
    float d = dr;
    return .2*d;
}


vec3 getNormal(vec3 p) {

    vec2 eps = vec2(0.001, 0.0);
    return normalize(vec3(map(p + eps.xyy) - map(p - eps.xyy),
                          map(p + eps.yxy) - map(p - eps.yxy),
                          map(p + eps.yyx) - map(p - eps.yyx)
                         )
                     );
}

float bumpSurf3D( in vec3 p){

    p.z += .3*tt;
    p = abs(mod(p*2., 2.*0.125)-0.0125);
    
    float x = min(p.x,min(p.z, p.y))/0.03125;

    return clamp(x, 0., 1.);


}

// Standard function-based bump mapping function (from Shane)
vec3 doBumpMap(in vec3 p, in vec3 nor, float bumpfactor){
    
    const float eps = BUMP_EPS;
    float ref = bumpSurf3D(p);                 
    vec3 grad = vec3( bumpSurf3D(vec3(p.x-eps, p.y, p.z))-ref,
                      bumpSurf3D(vec3(p.x, p.y-eps, p.z))-ref,
                      bumpSurf3D(vec3(p.x, p.y, p.z-eps))-ref )/eps;                     
          
    grad -= nor*dot(nor, grad);          
                      
    return normalize( nor + bumpfactor*grad );
	
}

vec2 raymarch(vec3 ro, vec3 rd, float steps) {

    float mat = 0.,
          t   = 0.,
          d   = 0.;
    vec3 p = ro;
    for(float i=.0; i<steps; i++) {
    
        d = map(p);
        mat = g_mat;  // save global material
        
        if(abs(d) < 0.0001 || t > 100.) break;
        
        t += d;
        p += rd*d;
    }
    
    return vec2(t, mat);
}

float n21(vec2 p) {
      return fract(sin(dot(p, vec2(524.423,123.34)))*3228324.345);
}

float noise(vec2 n) {
    const vec2 d = vec2(0., 1.0);
    vec2 b = floor(n);
    vec2 f = mix(vec2(0.0), vec2(1.0), fract(n));
    return mix(mix(n21(b), n21(b + d.yx), f.x), mix(n21(b + d.xy), n21(b + d.yy), f.x), f.y);
}

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

    vec3 tot = vec3(0);

   	vec2 uv = (fragCoord - .5*iResolution.xy)/iResolution.y;

    tt = iTime + 17.;
    vec3 ro = vec3(uv*6.,-4.),
          rd = vec3(0,0,.6),
          lp = vec3(2., 2., -10),
          lp2 = vec3(-2., -3., -5);

    vec3 ro2 = vec3(0);
    vec3 rd2 = vec3(uv, 1);
    vec3 col;
       
    float mat = 0.,
          t   = 0.,
          d   = 0.;
 

    vec2 e = vec2(0.0035, -0.0035);
     
    // background color
    vec3 c1 = vec3(0.000,1.000,0.835);
    vec3 c2 = vec3(1.000,0.494,0.078);
    
    // light color
    vec3 lc1 = vec3(0.961,0.816,0.737);
    vec3 lc2 = vec3(0.588,0.992,0.945);
    
    vec3 bg = temple(ro2, rd2);
    
    vec3 al = mix(vec3(0.639,1.000,0.843), vec3(0.714,0.522,1.000), length(uv));
    
    
    // currently only one pass
    for(float i = 0.; i < 1.; i++) {
        float steps = i > 0. ? 50. : 250.;
        vec2 rm = raymarch(ro, rd, steps);
        mat = rm.y;
        
        
        vec3 p = ro + rm.x*rd;
        
     
        
        vec3 n = normalize( e.xyy*map(p+e.xyy) + e.yyx*map(p+e.yyx) +
                                e.yxy*map(p+e.yxy) + e.xxx*map(p+e.xxx));

        vec3 pt = transform(p);
        vec3 np = transform(n);
        n = doBumpMap(pt, n, .004);
    
     
        if(rm.x < 50.) {
        
            vec3 l = normalize(lp-p);
            vec3 l2 = normalize(lp2-p);
            float dif = max(dot(n, l), .0);
            float dif2 = max(dot(n, l), .0);
            float spe = pow(max(dot(reflect(-rd, n), -l), .0),40.);

            float sss = smoothstep(0., 1., map(p + l * .4)) / .4;
            vec3 n2 = n;
            n2.xy += noise(p.xy) * .5 - .025;
            n2 = normalize(n2);
            float height = atan(n2.y, n2.x);

            vec3 iri = spectral_zucconi6(height*.8)*smoothstep(.8, .2, abs(n2.z));
       
            col += al*dif*lc1*.6 + al*.8*dif2*lc2 + .7*iri + .3*al;
            col *= .9;
           
            if(mat == 0.) {
                n += .5*texture(iChannel1, n.xy*100.).rgb;
                rd = reflect(rd, n);
         
                vec3 refl = texture(iChannel0, rd).rgb;
                refl *= .5*spectral_zucconi6(n.x*n.y);
                
                vec3 refl2 = .7*temple(p, rd);
                refl = mix(refl, refl2, .5);
                ; // reflect rainbows too
               
                col = mix(col, refl, .6);
                
                col = mix(col, bg, (1.3*SIN(.5*tt+PI)-.3)*(1.9*length(uv)));
               // col *= 4.;
              
            } 

        } else {
            col =  .8*bg;
               
        } 
    
    }
    
    //col = bg;
    
    col *= mix(.2, 1., (1.5-pow(dot(uv, uv), .5))); // vignette
    col = pow(col, vec3(.6)); // gamma
  
   
    fragColor = vec4(col, 1.0);
}