float waterline = 1.0;

vec4 hash4( vec4 n ) { return fract(abs(sin(n))*1399763.5453123); }
vec3 hash3( vec3 n ) { return fract(abs(sin(n))*1399763.5453123); }
vec2 NC0=vec2(1.0,136.0);
vec3 NC1=vec3(1.0,136.0,31.0);
float seed = 0.0;

vec2 getBL(vec3 pos) { return vec2(atan(pos.x,pos.z),atan(length(pos.xz),pos.y))/3.1415; }
vec3 setBL(vec2 bl) { bl*=3.1415; return vec3(cos(bl.y)*cos(bl.x),sin(bl.y),cos(bl.y)*sin(bl.x)); }

// Simple 2d interpolated noise for mercator projection, using local seed.
float sphereNoise2Ds(vec2 bl,float scale,float s)
{
    vec2 x = bl;
    x.x = clamp(x.x,-1.0,0.9999);
    vec2 xc = vec2(floor(x.y*scale)) + vec2(1.0,0.0);
    xc = floor(sin(3.1415 * xc/scale)* scale);
    vec4 p = vec4(xc.xy*x.x,x.y*scale,x.y*scale+1.0);
    vec4 f = fract(p);
    f=f*f*(3.0-2.0*f);
    p = floor(p);
    p.zw +=  vec2(scale);
    vec4 pp = vec4(p.x+1.0>=xc.x? -xc.x:p.x+1.0,p.y+1.0>=xc.y? -xc.y:p.y+1.0,p.zw);

    vec4 n= hash4(vec4(dot(p.yz,NC0),dot(pp.yz,NC0),dot(p.xw,NC0),dot(pp.xw,NC0))+s);
    n.xy = mix(n.xz,n.yw,f.yx);
    return mix(n.x,n.y,f.z);
}

vec2 sphereNoise2D2s(vec2 bl,float scale,vec2 s)
{
    vec2 x = bl;
    x.x = clamp(x.x,-1.0,0.9999);
    vec2 xc = vec2(floor(x.y*scale)) + vec2(1.0,0.0);
    xc = floor(sin(3.1415 * xc/scale)* scale);
    vec4 p = vec4(xc.xy*x.x,x.y*scale,x.y*scale+1.0);
    vec4 f = fract(p);
    f=f*f*(3.0-2.0*f);
    p = floor(p);
    p.zw +=  vec2(scale);
    vec4 pp = vec4(p.x+1.0>=xc.x? -xc.x:p.x+1.0,p.y+1.0>=xc.y? -xc.y:p.y+1.0,p.zw);

    vec4 n1= hash4(vec4(dot(p.yz,NC0),dot(pp.yz,NC0),dot(p.xw,NC0),dot(pp.xw,NC0))+s.x);
    vec4 n2= hash4(vec4(dot(p.yz,NC0),dot(pp.yz,NC0),dot(p.xw,NC0),dot(pp.xw,NC0))+s.y);
    n1 = vec4(mix(n1.xz,n1.yw,f.yx),mix(n2.xz,n2.yw,f.yx));

    return mix(n1.xz,n1.yw,f.zz);
}

// Simple 2d interpolated noise for mercator projection, using global seed.
float sphereNoise2D(vec2 bl,float scale) { return sphereNoise2Ds(bl,scale,seed); }

// Simple 2d interpolated noise for mercator projection, using global seed animated.
float sphereNoise2D(vec3 blt,float scale)
{
    float z = floor(blt.z)*54.0;
    float f = fract(blt.z);
    //f=f*f*(3.0-2.0*f);
    return mix(sphereNoise2Ds(blt.xy,scale,seed+z),
    sphereNoise2Ds(blt.xy,scale,seed+z+54.0),f);
}

// Simple 2d interpolated distortion noise for mercator projection, using global seed.
float sphereNoise2Do(vec2 bl,float scale,float force)
{
    bl = bl+vec2(sphereNoise2Ds(bl,scale*2.0f,3.0),sphereNoise2Ds(bl,scale*2.0f,20.0))*force/scale;
    bl.x = bl.x>=1.0 ? bl.x - 2.0 : bl.x;
    bl.y = clamp(bl.y,0.0,1.0);
    return sphereNoise2D(bl,scale);
}

float sphereNoise2Do2(vec2 bl,float scale,float force, float s)
{
    bl = bl+vec2(sphereNoise2Ds(bl,scale*s,3.0),sphereNoise2Ds(bl,scale*s,20.0))*force/scale;
    bl.x = bl.x>=1.0 ? bl.x - 2.0 : bl.x;
    bl.y = clamp(bl.y,0.0,1.0);
    return sphereNoise2D(bl,scale);
}

vec2 sphereIntersect(vec3 ray,vec3 pos,float r)
{
    float r2=r*r;
    float cd = dot(pos,ray);
    vec3 p = ray * cd-pos;
    float t2 = dot(p,p);
    if (t2>=r2) return vec2(-1.0, -1.0);
    float d = sqrt(r2-t2);
    return cd+vec2(- d,d);
}

float sphereIntersectIn(vec3 ray,vec3 pos,float r)
{
    float r2=r*r;
    float cd = dot(pos,ray);
    vec3 p = ray * cd-pos;
    float t2 = dot(p,p);
    if (t2>=r2) return -1.0;
    return cd + sqrt(r2-t2);
}

vec3 sphereNormal(vec3 pos, vec3 surface)
{
    return normalize(surface-pos);
}

float cloud(vec2 bl,float scale)
{
    float f = 0.0;
    float m = 1.0;
    float s = 0.32;
    for (int i=0;i<20;i++) {
        f += sphereNoise2D(bl,m*scale)*s;
        m = m * 2.0;
        s = s * 0.7;
    }
    f = f;
    return clamp(f,0.,1.);
}

float spherePerlin(vec2 bl,float scale)
{
    return (sphereNoise2D(bl,512.0*scale)*0.8+sphereNoise2D(bl,256.0*scale)*0.8
     +sphereNoise2D(bl,128.0*scale)+sphereNoise2D(bl,64.0*scale)
     +sphereNoise2D(bl,32.0*scale)+sphereNoise2D(bl,16.0*scale)
     +sphereNoise2D(bl,8.0*scale)+sphereNoise2D(bl,4.0*scale)
     +sphereNoise2D(bl,2.0*scale)*1.2+sphereNoise2D(bl,scale)*1.2)*0.1;
}

float spherePerlinA(vec3 blt,float scale)
{
    return (sphereNoise2D(blt,16.0*scale)
     +sphereNoise2D(blt,8.0*scale)+sphereNoise2D(blt,4.0*scale)
     +sphereNoise2D(blt,2.0*scale)*1.2+sphereNoise2D(blt,scale)*1.2)*0.2;
}

float spherePerlinHalf(vec2 bl,float scale)
{
    return (sphereNoise2D(bl,16.0*scale)+sphereNoise2D(bl,8.0*scale)+sphereNoise2D(bl,4.0*scale)
     +sphereNoise2D(bl,2.0*scale)+sphereNoise2D(bl,scale))*0.2;
}

float spherePerlinHalfs(vec2 bl,float scale, float seed)
{
    return (sphereNoise2Ds(bl,16.0*scale,seed)+sphereNoise2Ds(bl,8.0*scale,seed)+sphereNoise2Ds(bl,4.0*scale,seed)
     +sphereNoise2Ds(bl,2.0*scale,seed)+sphereNoise2Ds(bl,scale,seed))*0.2;
}

vec2 spherePerlinHalf2s(vec2 bl,float scale, vec2 seed)
{
    return (sphereNoise2D2s(bl,16.0*scale,seed)+sphereNoise2D2s(bl,8.0*scale,seed)+sphereNoise2D2s(bl,4.0*scale,seed)
     +sphereNoise2D2s(bl,2.0*scale,seed)+sphereNoise2D2s(bl,scale,seed))*0.2;
}

// Perlin distortion noise
float spherePerlinHalfd(vec2 bl,float scale,float force)
{
    return (sphereNoise2Do(bl,16.0*scale,force)+sphereNoise2Do(bl,8.0*scale,force)
     +sphereNoise2Do(bl,4.0*scale,force)
     +sphereNoise2Do(bl,2.0*scale,force)+sphereNoise2Do(bl,scale,force))*0.2;
}

// Perlin distortion noise with increasing coefficient
float spherePerlinHalfds(vec2 bl,float scale,float force)
{
    return sphereNoise2Do(bl,16.0*scale,force)*0.05
    +sphereNoise2Do(bl,8.0*scale,force)*0.15+sphereNoise2Do(bl,4.0*scale,force)*0.2
    +sphereNoise2Do(bl,2.0*scale,force)*0.25+sphereNoise2Do(bl,scale,force)*0.35;
}

float spherePerlinHalfds2(vec2 bl,float scale,float force,float s)
{
    return sphereNoise2Do2(bl,16.0*scale,force,16.*s)*0.05
    +sphereNoise2Do2(bl,8.0*scale,force,8.*s)*0.15+sphereNoise2Do2(bl,4.0*scale,force,4.*s)*0.2
    +sphereNoise2Do2(bl,2.0*scale,force,2.*s)*0.25+sphereNoise2Do2(bl,scale,force,s)*0.35;
}

vec3 wind3D(vec3 pos,float rings,float skip)
{
    vec3 xv = normalize(cross(pos,vec3(0.,1.,0.)));
    vec3 yv = normalize(cross(pos,xv));
    vec3 f = fract(pos)-0.5;
    vec3 px = floor(pos);
    float r = length(pos);
    vec3 wind =vec3(0.);

    for (int z = -0;z<2;z++)
    for (int y = -0;y<2;y++)
    for (int x = -0;x<2;x++) {
        vec3 off = vec3(float(x),float(y),float(z));
        vec3 xyz = px + off;
        float p = dot(xyz,NC1);
        vec3 add = hash3(vec3(p,p+30.0,p+23.0));
        if (add.z<skip) continue;
        vec3 fxyz = f + add*0.8-off;
        float ml = clamp(1.0-length(fxyz)*1.5,0.,1.);
        vec3 wp = normalize(pos - fxyz) * r;
        fxyz = wp-pos;
        float l1 = length(fxyz);
        float l2 = 1.-l1;
        vec3 w = vec3(0.);
        if (l1>0.) {
            w = normalize(cross(fxyz,pos));
        }
        wind +=w*ml*clamp(0.5-l1,0.,1.)*clamp(1.0-l2*l2*l2,0.,1.);
    }
    return wind;
}

float atmosphere(vec3 pos, vec3 atmosSurface, vec3 ray, vec3 lnorm, float planetRadius, float atmosRadius, float atmosDepth)
{
    float wall = sqrt(atmosRadius*atmosRadius - planetRadius*planetRadius);
    float acc = clamp(atmosDepth/wall,0.0,1.0);
    acc = acc*clamp(dot(lnorm,normalize(atmosSurface-pos)),0.0,1.0);
    acc = acc*acc;
    acc = acc*acc;

    return acc;
}

vec3 sphereTexture(vec3 realSurface,vec3 ray, vec3 snorm,vec3 lnorm,float time)
{
    seed = 1.0;
    vec2 bl = getBL(realSurface);

    seed = 2.0;
    float d = spherePerlinHalfds(bl,8.0,0.5);
    float c = spherePerlin(bl,2.5);
    float can = clamp(1.0-abs(d-0.5)*10.0,0.0,1.0) * clamp(c*5.0-1.6,0.0,1.0);
    float ll = clamp(1.0-abs(d-0.5)*4.0,0.0,1.0);
    float can2 = can*can;//*0.34;
    float can3 = can2*can;

    seed = 10.0;
    float e = spherePerlinHalfds(bl,6.0,0.5);
    float mc = clamp(e*6.0-2.4, 0.0, 1.0);

    seed = 5.0;
    float fa = (spherePerlinA(vec3(bl,time),200.0)+spherePerlinA(vec3(bl,time+2.5),200.0))*0.5;

    seed = 5.0;
    float f = spherePerlin(bl,4.0);
    f = mix(1.0-pow(1.0-f*1.7,2.0),f*1.7,0.3);
    f = f-can3*0.34;

    float canl = (can3+fa-0.5);

    float light = clamp(dot(snorm,lnorm),0.0,1.0);
    vec3 col = clamp(mix(vec3(2.0,0.2,0.0),vec3(0.02,0.0,0.0),f)*f,vec3(0.,0.,0.),vec3(1.,1.,1.));
    col = mix(col,mix(vec3(1.0,0.5,0.0),vec3(1.0,0.8,0.33),mc),0.2);
    vec3 lav = mix(vec3(1.,0.0,0.),mix(vec3(1.,1.0,0.),vec3(1.,1.,1.),clamp(abs(0.5-canl)*3.0,0.,1.)),clamp(canl,0.,1.));

    return clamp(mix(col*light,lav*2.0,can3),vec3(0.,0.,0.),vec3(1.,1.,1.));//mix(col*light,mix(lav,col*light,clamp((1.0-ll*ll),0.,1.)),ll);

    //vec3 col = clamp(mix(vec3(0.94,0.63,0.33),vec3(1.0,0.85,0.57),f)*f,vec3(0.,0.,0.),vec3(1.,1.,1.));
    //return mix(col,mix(vec3(1.0,0.5,0.0),vec3(1.0,0.8,0.33),mc),0.2);
    //return vec3(0.5);
}

vec2 cloudTexture(vec3 realSurface, float time)
{
   seed = 14.0;
   vec2 bl = getBL(realSurface);

   vec3 w = wind3D(realSurface*1.,1.5,0.);
   w += wind3D(realSurface*2.,1.5,0.)*0.2;
   w += wind3D(realSurface*6.,1.5,0.)*0.1;
   vec2 t = time*vec2(0.2,0.5);
   vec2 bld1 = getBL(normalize(realSurface-w*(1.4+0.4*fract(t.x))));
   vec2 bld2 = getBL(normalize(realSurface-w*(1.4+0.4*fract(t.x+0.5))));
   vec2 blc1 = getBL(normalize(realSurface-w*(1.4+0.4*fract(t.y))));
   vec2 blc2 = getBL(normalize(realSurface-w*(1.4+0.4*fract(t.y+0.5))));

   seed = 14.0+floor(t.x)*0.1;
   float d1 = spherePerlinHalf(bld1, 6.0);
   seed = 14.05+floor(t.x+0.5)*0.1;
   float d2 = spherePerlinHalf(bld2, 6.0);

   seed = 14.0+floor(t.y)*2.0;
   float c1 = spherePerlinHalf(blc1, 192.0);
   seed = 15.0+floor(t.y+0.5)*2.0;
   float c2 = spherePerlinHalf(blc2, 192.0);
   vec2 fdc = abs(0.5-fract(t))*2.0;
   float d = mix(d1,d2,fdc.x);
   float c = mix(c1,c2,fdc.y);

   float f = clamp((d*d)*4.0-c-0.5,0.0,1.0);
   f = 1.-f;
   return vec2(1.-f*f,0);
   //return vec2(0.);
}

void mainImage( out vec4 fragColor, in vec2 fragCoord )
{
    // resolution
    float res = 1.0 / iResolution.y;
    vec2 p = (-iResolution.xy + 2.0*fragCoord.xy) *res;

    // bace const
    float time=iTime*1.0;
    vec3 ray = normalize(vec3(p,2.0));
    vec3 col = vec3(0.0);
    vec2 mouse = (iMouse.xy*2.0-iResolution.xy)*res;

    // segments
    vec2 scroll = iMouse.z>0.0 ? mouse : vec2(floor(p.x*2.5)*0.4+0.2,0.0);
    scroll = clamp(scroll*0.5+0.5,vec2(0.),vec2(1.));
    // test scroll
    //const vec3 pal[6] = vec3[6](vec3(1.0,0.0,0.0),vec3(1.0,1.0,0.0),vec3(0.0,1.0,0.0),vec3(0.0,1.0,1.0),vec3(0.0,0.0,1.0),vec3(1.0,0.0,1.0));
    //col = mix(pal[int(floor(scroll.x*5.0))],pal[int(floor(scroll.x*5.0))+1],fract(scroll.x*5.0));

    // transforms
    vec2 rotate = iMouse.z>0.0 ? mouse*2. : vec2(0.,-0.6);
    //vec2 rotate = vec2(0.,-0.6);
    float protate = time*0.025;
    vec4 mcs=vec4(sin(rotate),cos(rotate));
    vec2 pcs=vec2(cos(protate),sin(protate)); // planet rotation
    mat3 mr=mat3(vec3(mcs.z,0.0,mcs.x),vec3(0.0,1.0,0.0),vec3(-mcs.x,0.0,mcs.z));
    mr=mat3(vec3(1.0,0.0,0.0),vec3(0.0,mcs.w,mcs.y),vec3(0.0,-mcs.y,mcs.w))*mr;
    mat3 pmr=mat3(vec3(pcs.x,0.0,pcs.y),vec3(0.0,1.0,0.0),vec3(-pcs.y,0.0,pcs.x));

    // object data
    float dist=2.3;
    float rplanet = 1.0;
    float ratmos = 1.15;
    float rcloud = 1.005;
    vec3 spos = vec3(0.0,0.0,dist);
    vec3 lpos = vec3(4.5,4.5,-9.0+dist);
    lpos = mr*(lpos-spos) + spos;
    vec3 lnorm = normalize(lpos - spos);

    // render calc
    float dplanet = sphereIntersect(ray,spos,rplanet).x;
    vec2 dcloud = sphereIntersect(ray,spos,rcloud);
    float datmos = sphereIntersect(ray,spos,ratmos).x;
    if (dplanet>0.0) {
        vec3 ssurface = ray*dplanet;
        vec3 snorm = sphereNormal(spos,ssurface);
        vec3 lsnorm = normalize(lpos - ssurface);
        float dcs = sphereIntersectIn(lnorm,spos-ssurface,rcloud);
        vec3 scloud = ssurface+lnorm*dcs;
        float cs = dot(snorm,lnorm)>0. ? clamp(1.0-cloudTexture((scloud-spos)*mr*pmr,time*0.2).x,0.0,1.0) : 1.0;

        col = sphereTexture((ssurface-spos)*mr*pmr,ray,snorm,lsnorm,time);
        col = mix(col,col*cs,0.7);
    }
    if (datmos>0.0) {
        float adepth = dplanet>0.0 ? (dplanet-datmos) : abs(dot(ray,spos)-datmos);
        vec3 ssurface = ray*datmos;
        float ap = atmosphere(spos, ssurface, ray, lnorm, rplanet, ratmos, adepth);
        col = mix(col,mix(vec3(1.0,0.3,0.0),vec3(1.0,0.99,0.94),ap),ap);
    }
    if (dcloud.x>0.0) {
        vec3 sback = ray*dcloud.y;
        /*if (dplanet<0.0) {
            vec3 bnorm = sphereNormal(spos,sback);
            col = mix(col, vec3(0.0), cloudTexture((sback-spos)*mr*pmr,time*0.2) * clamp(dot(bnorm,lnorm),0.0,1.0));
        }*/
        vec3 ssurface = ray*dcloud.x;
        vec3 snorm = sphereNormal(spos,ssurface);
        vec3 lsnorm = normalize(lpos - ssurface);
        vec2 al = cloudTexture((ssurface-spos+lsnorm*0.01)*mr*pmr,time*0.2);
        vec2 a = cloudTexture((ssurface-spos)*mr*pmr,time*0.2);
        float light = clamp(dot(snorm,lsnorm),0.0,1.0);
        col = mix(col,mix(vec3(0.0,0.0,0.0),vec3(0.0,0.0,0.0),a.y)*light,clamp(a.x,0.,1.))
            + vec3(clamp((a.x-al.x)*0.7,0.,1.)+clamp((a.y-al.y)*1.*a.x+a.x*0.3*light*light,0.,1.))*light;
    }
    // Output to screen
    fragColor = vec4(col,1.0);
}