Code

float sdPlane( in vec3 p, in vec4 s )
{
    return p.y + s.w;
}

float sdSphere( in vec3 p, in float s )
{
    return sqrt(p.x*p.x+p.y*p.y+p.z*p.z)-s;
}

float sdBox( vec3 p, vec3 b )
{
  vec3 d = abs(p) - b;
  return min(max(d.x,max(d.y,d.z)),0.0) + length(max(d,0.0));
}


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

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

vec2 map( in vec3 pos )
{
    vec2 res = vec2( sdPlane(pos,vec4(0.0,1.0,0.0,0.0)), 13.0 );

    res = opU( res, vec2( sdSphere( pos-vec3(-2.0,0.40,-2.0), 0.40 ), 0.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-2.0,0.25,-1.0), 0.25 ), 0.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-2.0,0.40, 0.0), 0.40 ), 1.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-2.0,0.25, 1.0), 0.25 ), 1.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-2.0,0.40, 2.0), 0.40 ), 2.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-1.0,0.25,-2.0), 0.25 ), 2.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-1.0,0.40,-1.0), 0.40 ), 3.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-1.0,0.25, 0.0), 0.25 ), 3.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-1.0,0.40, 1.0), 0.40 ), 4.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3(-1.0,0.25, 2.0), 0.25 ), 4.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 0.0,0.40,-2.0), 0.40 ), 5.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 0.0,0.25,-1.0), 0.25 ), 5.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 0.0,0.40, 0.0), 0.40 ), 6.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 0.0,0.25, 1.0), 0.25 ), 6.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 0.0,0.40, 2.0), 0.40 ), 7.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 1.0,0.25,-2.0), 0.25 ), 7.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 1.0,0.40,-1.0), 0.40 ), 8.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 1.0,0.25, 0.0), 0.25 ), 8.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 1.0,0.40, 1.0), 0.40 ), 9.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 1.0,0.25, 2.0), 0.25 ), 9.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 2.0,0.40,-2.0), 0.40 ), 10.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 2.0,0.25,-1.0), 0.25 ), 10.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 2.0,0.40, 0.0), 0.40 ), 11.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 2.0,0.25, 1.0), 0.25 ), 11.0 ) );
    res = opU( res, vec2( sdSphere( pos-vec3( 2.0,0.40, 2.0), 0.40 ), 12.0 ) );

    return res;
}

vec2 castRay( in vec3 ro, in vec3 rd, in float tmin, in float tmax )
{
    float t = tmin;
    float m = 0.0;
    for( int i=0; i<=63; i++ )
    {
        vec2 res = map( ro+rd*t );
        if( res.x<0.0 || t>tmax ) break;
        t += res.x;
        m = res.y;
    }

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


float shadow( in vec3 ro, in vec3 rd, in float mint, in float tmax )
{
    float t = mint;
    for( int i=0; i<64; i++ )
    {
        vec2 res = map( ro+rd*t );
        if( res.x<0.0 || t>tmax ) break;
        t += res.x;
    }

    if( t>tmax ) return 1.0;
    return 0.0;
}

float calcOcclusionArlo( in vec3 pos, in vec3 nor )
{
    float occ = 0.0;
    for( int i=0; i<8; i++ )
    {
        float h = 0.005 + 0.25*float(i)/7.0;
        vec3 dir = normalize( sin( float(i)*73.4 + vec3(0.0,2.1,4.2) ));
        dir = normalize( nor + dir );
        occ += (h-map( pos + h*dir ).x);
    }
    return clamp( 1.0 - 9.0*occ/8.0, 0.0, 1.0 );
}

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);
}

vec3 render( in vec3 ro, in vec3 rd )
{
    vec3 col = vec3(0.7, 0.9, 1.0) + rd.y*vec3(0.8,0.8,0.8);
    vec2 res = castRay(ro,rd, 0.0, 20.0);
    float t = res.x;
    float m = res.y;
    if( m>-0.5 )
    {
        vec3  pos = ro + t*rd;
        vec3  nor = calcNormal( pos );
        vec3  ref = reflect( rd, nor );
        float ang = 5.0*atan( pos.x, pos.z );
        float coa = cos(ang);

             if( m< 0.5 ) col = vec3(1.00, 0.00, 0.00);
        else if( m< 1.5 ) col = vec3(0.92, 0.03, 0.06);
        else if( m< 2.5 ) col = vec3(0.84, 0.06, 0.12);
        else if( m< 3.5 ) col = vec3(0.76, 0.09, 0.18);
        else if( m< 4.5 ) col = vec3(0.68, 0.12, 0.24);
        else if( m< 5.5 ) col = vec3(0.60, 0.15, 0.30);
        else if( m< 6.5 ) col = vec3(0.52, 0.18, 0.36);
        else if( m< 7.5 ) col = vec3(0.44, 0.21, 0.42);
        else if( m< 8.5 ) col = vec3(0.36, 0.24, 0.48);
        else if( m< 9.5 ) col = vec3(0.28, 0.27, 0.54);
        else if( m<10.5 ) col = vec3(0.20, 0.30, 0.60);
        else if( m<11.5 ) col = vec3(0.12, 0.33, 0.66);
        else if( m<12.5 ) col = vec3(0.04, 0.36, 0.72);
        else
        {
            if( mod( mod(floor(pos.x),2.0) + mod(floor(pos.z),2.0), 2.0 ) < 0.5 )
            {
                col = vec3(0.4, 0.4, 0.4);
            }
            else
            {
                col = vec3(0.5, 0.5, 0.5);
                vec2 uv = 256.0*0.1*pos.xz + vec2(37.0,17.0);
                col *= texture( iChannel0, (floor(uv)+ 0.5)/256.0, -100.0 ).y;
            }

            if( max( abs(pos.x), abs(pos.z)) < 1.5 )
            {
                vec2 q = mod(pos.xz+0.5,1.0)-0.5;
                col *= smoothstep( 0.0, 0.5, sqrt(q.x*q.x+q.y*q.y) );
            }

            col *= pow(abs(tan(ang)*coa),2.0) + pow(abs(coa),2.0);
        }


        // lighitng
        vec3  lig = normalize( vec3(-0.6, 0.7, -0.5) );
        float amb = 0.5+0.5*nor.y; if( amb<0.0 ) amb=0.0; else if( amb>1.0 ) amb=1.0;
        float dif = dot( normalize(nor), lig ); if( dif<0.0 ) dif=0.0; else if( dif>1.0 ) dif=1.0;
        float spe = pow(clamp( dot(normalize(ref), lig ), 0.0, 1.0 ),16.0);

        dif *= shadow( pos, lig, 0.01, 2.5 );

        vec3 lin = vec3(0.0, 0.0, 0.0);
        lin += 1.20*dif*vec3(1.00,0.85,0.55);
        lin += 0.50*amb*vec3(0.50,0.70,1.00);
        col = col*lin;
        col += 1.0*spe*dif*vec3(1.0,1.0,1.0);
    }

    return col;
}

vec3 cross_product( in vec3 a, in vec3 b )
{
    return vec3( a.y*b.z - a.z*b.y,
                -a.x*b.z + a.z*b.x,
                 a.x*b.y - a.y*b.x );
}

mat3 setCamera( in vec3 ro, in vec3 ta, float cr )
{
    vec3 cw = normalize(ta-ro);
    vec3 cp = vec3(sin(cr), cos(cr),0.0);
    vec3 cu = normalize( cross_product(cw,cp) );
    vec3 cv = normalize( cross_product(cu,cw) );
    return mat3( cu, cv, cw );
}

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

    float time = 50.0;

    // camera
    vec3 ro = vec3( 3.0*cos(0.1*time), 2.0, 3.0*sin(0.1*time) );
    vec3 ta = vec3( 0.0, -0.2, 0.0 );

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

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

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

    col = pow( col, vec3(0.5) );

    fragColor=vec4( col, 1.0 );
}