raymarching.sttf

{"links":[{"end":"RenderOutput","filter":"Linear","slot":0,"start":"Image","wrapMode":"Repeat"}],"metadata":{"Author":"iq","Description":"A set of raw primitives. All except the ellipsoid are exact euclidean distances. More info here: [url=https://iquilezles.org/articles/distfunctions]https://iquilezles.org/articles/distfunctions[/url]","Name":"Raymarching - Primitives","ShaderToyURL":"https://www.shadertoy.com/view/Xds3zN"},"nodes":[{"class":"RenderOutput","name":"RenderOutput"},{"class":"GLSLShader","name":"Image","source":"// Copyright Inigo Quilez, 2016 - https://iquilezles.org/\n// I am the sole copyright owner of this Work.\n// You cannot host, display, distribute or share this Work in any form,\n// including physical and digital. You cannot use this Work in any\n// commercial or non-commercial product, website or project. You cannot\n// sell this Work and you cannot mint an NFTs of it.\n// I share this Work for educational purposes, and you can link to it,\n// through an URL, proper attribution and unmodified screenshot, as part\n// of your educational material. If these conditions are too restrictive\n// please contact me and we'll definitely work it out.\n\n// A list of useful distance function to simple primitives. All\n// these functions (except for ellipsoid) return an exact\n// euclidean distance, meaning they produce a better SDF than\n// what you'd get if you were constructing them from boolean\n// operations (such as cutting an infinite cylinder with two planes).\n\n// List of other 3D SDFs:\n//    https://www.shadertoy.com/playlist/43cXRl\n// and\n//    https://iquilezles.org/articles/distfunctions\n\n#if HW_PERFORMANCE==0\n#define AA 1\n#else\n#define AA 2   // make this 2 or 3 for antialiasing\n#endif\n\n//------------------------------------------------------------------\nfloat dot2( in vec2 v ) { return dot(v,v); }\nfloat dot2( in vec3 v ) { return dot(v,v); }\nfloat ndot( in vec2 a, in vec2 b ) { return a.x*b.x - a.y*b.y; }\n\nfloat sdPlane( vec3 p )\n{\n\treturn p.y;\n}\n\nfloat sdSphere( vec3 p, float s )\n{\n    return length(p)-s;\n}\n\nfloat sdBox( vec3 p, vec3 b )\n{\n    vec3 d = abs(p) - b;\n    return min(max(d.x,max(d.y,d.z)),0.0) + length(max(d,0.0));\n}\n\nfloat sdBoxFrame( vec3 p, vec3 b, float e )\n{\n       p = abs(p  )-b;\n  vec3 q = abs(p+e)-e;\n\n  return min(min(\n      length(max(vec3(p.x,q.y,q.z),0.0))+min(max(p.x,max(q.y,q.z)),0.0),\n      length(max(vec3(q.x,p.y,q.z),0.0))+min(max(q.x,max(p.y,q.z)),0.0)),\n      length(max(vec3(q.x,q.y,p.z),0.0))+min(max(q.x,max(q.y,p.z)),0.0));\n}\nfloat sdEllipsoid( in vec3 p, in vec3 r ) // approximated\n{\n    float k0 = length(p/r);\n    float k1 = length(p/(r*r));\n    return k0*(k0-1.0)/k1;\n}\n\nfloat sdTorus( vec3 p, vec2 t )\n{\n    return length( vec2(length(p.xz)-t.x,p.y) )-t.y;\n}\n\nfloat sdCappedTorus(in vec3 p, in vec2 sc, in float ra, in float rb)\n{\n    p.x = abs(p.x);\n    float k = (sc.y*p.x>sc.x*p.y) ? dot(p.xy,sc) : length(p.xy);\n    return sqrt( dot(p,p) + ra*ra - 2.0*ra*k ) - rb;\n}\n\nfloat sdHexPrism( vec3 p, vec2 h )\n{\n    vec3 q = abs(p);\n\n    const vec3 k = vec3(-0.8660254, 0.5, 0.57735);\n    p = abs(p);\n    p.xy -= 2.0*min(dot(k.xy, p.xy), 0.0)*k.xy;\n    vec2 d = vec2(\n       length(p.xy - vec2(clamp(p.x, -k.z*h.x, k.z*h.x), h.x))*sign(p.y - h.x),\n       p.z-h.y );\n    return min(max(d.x,d.y),0.0) + length(max(d,0.0));\n}\n\nfloat sdOctogonPrism( in vec3 p, in float r, float h )\n{\n  const vec3 k = vec3(-0.9238795325,   // sqrt(2+sqrt(2))/2 \n                       0.3826834323,   // sqrt(2-sqrt(2))/2\n                       0.4142135623 ); // sqrt(2)-1 \n  // reflections\n  p = abs(p);\n  p.xy -= 2.0*min(dot(vec2( k.x,k.y),p.xy),0.0)*vec2( k.x,k.y);\n  p.xy -= 2.0*min(dot(vec2(-k.x,k.y),p.xy),0.0)*vec2(-k.x,k.y);\n  // polygon side\n  p.xy -= vec2(clamp(p.x, -k.z*r, k.z*r), r);\n  vec2 d = vec2( length(p.xy)*sign(p.y), p.z-h );\n  return min(max(d.x,d.y),0.0) + length(max(d,0.0));\n}\n\nfloat sdCapsule( vec3 p, vec3 a, vec3 b, float r )\n{\n\tvec3 pa = p-a, ba = b-a;\n\tfloat h = clamp( dot(pa,ba)/dot(ba,ba), 0.0, 1.0 );\n\treturn length( pa - ba*h ) - r;\n}\n\nfloat sdRoundCone( in vec3 p, in float r1, float r2, float h )\n{\n    vec2 q = vec2( length(p.xz), p.y );\n    \n    float b = (r1-r2)/h;\n    float a = sqrt(1.0-b*b);\n    float k = dot(q,vec2(-b,a));\n    \n    if( k < 0.0 ) return length(q) - r1;\n    if( k > a*h ) return length(q-vec2(0.0,h)) - r2;\n        \n    return dot(q, vec2(a,b) ) - r1;\n}\n\nfloat sdRoundCone(vec3 p, vec3 a, vec3 b, float r1, float r2)\n{\n    // sampling independent computations (only depend on shape)\n    vec3  ba = b - a;\n    float l2 = dot(ba,ba);\n    float rr = r1 - r2;\n    float a2 = l2 - rr*rr;\n    float il2 = 1.0/l2;\n    \n    // sampling dependant computations\n    vec3 pa = p - a;\n    float y = dot(pa,ba);\n    float z = y - l2;\n    float x2 = dot2( pa*l2 - ba*y );\n    float y2 = y*y*l2;\n    float z2 = z*z*l2;\n\n    // single square root!\n    float k = sign(rr)*rr*rr*x2;\n    if( sign(z)*a2*z2 > k ) return  sqrt(x2 + z2)        *il2 - r2;\n    if( sign(y)*a2*y2 < k ) return  sqrt(x2 + y2)        *il2 - r1;\n                            return (sqrt(x2*a2*il2)+y*rr)*il2 - r1;\n}\n\nfloat sdTriPrism( vec3 p, vec2 h )\n{\n    const float k = sqrt(3.0);\n    h.x *= 0.5*k;\n    p.xy /= h.x;\n    p.x = abs(p.x) - 1.0;\n    p.y = p.y + 1.0/k;\n    if( p.x+k*p.y>0.0 ) p.xy=vec2(p.x-k*p.y,-k*p.x-p.y)/2.0;\n    p.x -= clamp( p.x, -2.0, 0.0 );\n    float d1 = length(p.xy)*sign(-p.y)*h.x;\n    float d2 = abs(p.z)-h.y;\n    return length(max(vec2(d1,d2),0.0)) + min(max(d1,d2), 0.);\n}\n\n// vertical\nfloat sdCylinder( vec3 p, vec2 h )\n{\n    vec2 d = abs(vec2(length(p.xz),p.y)) - h;\n    return min(max(d.x,d.y),0.0) + length(max(d,0.0));\n}\n\n// arbitrary orientation\nfloat sdCylinder(vec3 p, vec3 a, vec3 b, float r)\n{\n    vec3 pa = p - a;\n    vec3 ba = b - a;\n    float baba = dot(ba,ba);\n    float paba = dot(pa,ba);\n\n    float x = length(pa*baba-ba*paba) - r*baba;\n    float y = abs(paba-baba*0.5)-baba*0.5;\n    float x2 = x*x;\n    float y2 = y*y*baba;\n    float d = (max(x,y)<0.0)?-min(x2,y2):(((x>0.0)?x2:0.0)+((y>0.0)?y2:0.0));\n    return sign(d)*sqrt(abs(d))/baba;\n}\n\n// vertical\nfloat sdCone( in vec3 p, in vec2 c, float h )\n{\n    vec2 q = h*vec2(c.x,-c.y)/c.y;\n    vec2 w = vec2( length(p.xz), p.y );\n    \n\tvec2 a = w - q*clamp( dot(w,q)/dot(q,q), 0.0, 1.0 );\n    vec2 b = w - q*vec2( clamp( w.x/q.x, 0.0, 1.0 ), 1.0 );\n    float k = sign( q.y );\n    float d = min(dot( a, a ),dot(b, b));\n    float s = max( k*(w.x*q.y-w.y*q.x),k*(w.y-q.y)  );\n\treturn sqrt(d)*sign(s);\n}\n\nfloat sdCappedCone( in vec3 p, in float h, in float r1, in float r2 )\n{\n    vec2 q = vec2( length(p.xz), p.y );\n    \n    vec2 k1 = vec2(r2,h);\n    vec2 k2 = vec2(r2-r1,2.0*h);\n    vec2 ca = vec2(q.x-min(q.x,(q.y < 0.0)?r1:r2), abs(q.y)-h);\n    vec2 cb = q - k1 + k2*clamp( dot(k1-q,k2)/dot2(k2), 0.0, 1.0 );\n    float s = (cb.x < 0.0 && ca.y < 0.0) ? -1.0 : 1.0;\n    return s*sqrt( min(dot2(ca),dot2(cb)) );\n}\n\nfloat sdCappedCone(vec3 p, vec3 a, vec3 b, float ra, float rb)\n{\n    float rba  = rb-ra;\n    float baba = dot(b-a,b-a);\n    float papa = dot(p-a,p-a);\n    float paba = dot(p-a,b-a)/baba;\n\n    float x = sqrt( papa - paba*paba*baba );\n\n    float cax = max(0.0,x-((paba<0.5)?ra:rb));\n    float cay = abs(paba-0.5)-0.5;\n\n    float k = rba*rba + baba;\n    float f = clamp( (rba*(x-ra)+paba*baba)/k, 0.0, 1.0 );\n\n    float cbx = x-ra - f*rba;\n    float cby = paba - f;\n    \n    float s = (cbx < 0.0 && cay < 0.0) ? -1.0 : 1.0;\n    \n    return s*sqrt( min(cax*cax + cay*cay*baba,\n                       cbx*cbx + cby*cby*baba) );\n}\n\n// c is the sin/cos of the desired cone angle\nfloat sdSolidAngle(vec3 pos, vec2 c, float ra)\n{\n    vec2 p = vec2( length(pos.xz), pos.y );\n    float l = length(p) - ra;\n\tfloat m = length(p - c*clamp(dot(p,c),0.0,ra) );\n    return max(l,m*sign(c.y*p.x-c.x*p.y));\n}\n\nfloat sdOctahedron(vec3 p, float s)\n{\n    p = abs(p);\n    float m = p.x + p.y + p.z - s;\n\n    // exact distance\n    #if 0\n    vec3 o = min(3.0*p - m, 0.0);\n    o = max(6.0*p - m*2.0 - o*3.0 + (o.x+o.y+o.z), 0.0);\n    return length(p - s*o/(o.x+o.y+o.z));\n    #endif\n    \n    // exact distance\n    #if 1\n \tvec3 q;\n         if( 3.0*p.x < m ) q = p.xyz;\n    else if( 3.0*p.y < m ) q = p.yzx;\n    else if( 3.0*p.z < m ) q = p.zxy;\n    else return m*0.57735027;\n    float k = clamp(0.5*(q.z-q.y+s),0.0,s); \n    return length(vec3(q.x,q.y-s+k,q.z-k)); \n    #endif\n    \n    // bound, not exact\n    #if 0\n\treturn m*0.57735027;\n    #endif\n}\n\nfloat sdPyramid( in vec3 p, in float h )\n{\n    float m2 = h*h + 0.25;\n    \n    // symmetry\n    p.xz = abs(p.xz);\n    p.xz = (p.z>p.x) ? p.zx : p.xz;\n    p.xz -= 0.5;\n\t\n    // project into face plane (2D)\n    vec3 q = vec3( p.z, h*p.y - 0.5*p.x, h*p.x + 0.5*p.y);\n   \n    float s = max(-q.x,0.0);\n    float t = clamp( (q.y-0.5*p.z)/(m2+0.25), 0.0, 1.0 );\n    \n    float a = m2*(q.x+s)*(q.x+s) + q.y*q.y;\n\tfloat b = m2*(q.x+0.5*t)*(q.x+0.5*t) + (q.y-m2*t)*(q.y-m2*t);\n    \n    float d2 = min(q.y,-q.x*m2-q.y*0.5) > 0.0 ? 0.0 : min(a,b);\n    \n    // recover 3D and scale, and add sign\n    return sqrt( (d2+q.z*q.z)/m2 ) * sign(max(q.z,-p.y));;\n}\n\n// la,lb=semi axis, h=height, ra=corner\nfloat sdRhombus(vec3 p, float la, float lb, float h, float ra)\n{\n    p = abs(p);\n    vec2 b = vec2(la,lb);\n    float f = clamp( (ndot(b,b-2.0*p.xz))/dot(b,b), -1.0, 1.0 );\n\tvec2 q = vec2(length(p.xz-0.5*b*vec2(1.0-f,1.0+f))*sign(p.x*b.y+p.z*b.x-b.x*b.y)-ra, p.y-h);\n    return min(max(q.x,q.y),0.0) + length(max(q,0.0));\n}\n\nfloat sdHorseshoe( in vec3 p, in vec2 c, in float r, in float le, vec2 w )\n{\n    p.x = abs(p.x);\n    float l = length(p.xy);\n    p.xy = mat2(-c.x, c.y, \n              c.y, c.x)*p.xy;\n    p.xy = vec2((p.y>0.0 || p.x>0.0)?p.x:l*sign(-c.x),\n                (p.x>0.0)?p.y:l );\n    p.xy = vec2(p.x,abs(p.y-r))-vec2(le,0.0);\n    \n    vec2 q = vec2(length(max(p.xy,0.0)) + min(0.0,max(p.x,p.y)),p.z);\n    vec2 d = abs(q) - w;\n    return min(max(d.x,d.y),0.0) + length(max(d,0.0));\n}\n\nfloat sdU( in vec3 p, in float r, in float le, vec2 w )\n{\n    p.x = (p.y>0.0) ? abs(p.x) : length(p.xy);\n    p.x = abs(p.x-r);\n    p.y = p.y - le;\n    float k = max(p.x,p.y);\n    vec2 q = vec2( (k<0.0) ? -k : length(max(p.xy,0.0)), abs(p.z) ) - w;\n    return length(max(q,0.0)) + min(max(q.x,q.y),0.0);\n}\n\n//------------------------------------------------------------------\n\nvec2 opU( vec2 d1, vec2 d2 )\n{\n\treturn (d1.x<d2.x) ? d1 : d2;\n}\n\n//------------------------------------------------------------------\n\n#define ZERO (min(iFrame,0))\n\n//------------------------------------------------------------------\n\nvec2 map( in vec3 pos )\n{\n    vec2 res = vec2( pos.y, 0.0 );\n\n    // bounding box\n    if( sdBox( pos-vec3(-2.0,0.3,0.25),vec3(0.3,0.3,1.0) )<res.x )\n    {\n      res = opU( res, vec2( sdSphere(    pos-vec3(-2.0,0.25, 0.0), 0.25 ), 26.9 ) );\n\t  res = opU( res, vec2( sdRhombus(  (pos-vec3(-2.0,0.25, 1.0)).xzy, 0.15, 0.25, 0.04, 0.08 ),17.0 ) );\n    }\n\n    // bounding box\n    if( sdBox( pos-vec3(0.0,0.3,-1.0),vec3(0.35,0.3,2.5) )<res.x )\n    {\n\tres = opU( res, vec2( sdCappedTorus((pos-vec3( 0.0,0.30, 1.0))*vec3(1,-1,1), vec2(0.866025,-0.5), 0.25, 0.05), 25.0) );\n    res = opU( res, vec2( sdBoxFrame(    pos-vec3( 0.0,0.25, 0.0), vec3(0.3,0.25,0.2), 0.025 ), 16.9 ) );\n\tres = opU( res, vec2( sdCone(        pos-vec3( 0.0,0.45,-1.0), vec2(0.6,0.8),0.45 ), 55.0 ) );\n    res = opU( res, vec2( sdCappedCone(  pos-vec3( 0.0,0.25,-2.0), 0.25, 0.25, 0.1 ), 13.67 ) );\n    res = opU( res, vec2( sdSolidAngle(  pos-vec3( 0.0,0.00,-3.0), vec2(3,4)/5.0, 0.4 ), 49.13 ) );\n    }\n\n    // bounding box\n    if( sdBox( pos-vec3(1.0,0.3,-1.0),vec3(0.35,0.3,2.5) )<res.x )\n    {\n\tres = opU( res, vec2( sdTorus(      (pos-vec3( 1.0,0.30, 1.0)).xzy, vec2(0.25,0.05) ), 7.1 ) );\n    res = opU( res, vec2( sdBox(         pos-vec3( 1.0,0.25, 0.0), vec3(0.3,0.25,0.1) ), 3.0 ) );\n    res = opU( res, vec2( sdCapsule(     pos-vec3( 1.0,0.00,-1.0),vec3(-0.1,0.1,-0.1), vec3(0.2,0.4,0.2), 0.1  ), 31.9 ) );\n\tres = opU( res, vec2( sdCylinder(    pos-vec3( 1.0,0.25,-2.0), vec2(0.15,0.25) ), 8.0 ) );\n    res = opU( res, vec2( sdHexPrism(    pos-vec3( 1.0,0.2,-3.0), vec2(0.2,0.05) ), 18.4 ) );\n    }\n\n    // bounding box\n    if( sdBox( pos-vec3(-1.0,0.35,-1.0),vec3(0.35,0.35,2.5))<res.x )\n    {\n\tres = opU( res, vec2( sdPyramid(    pos-vec3(-1.0,-0.6,-3.0), 1.0 ), 13.56 ) );\n\tres = opU( res, vec2( sdOctahedron( pos-vec3(-1.0,0.15,-2.0), 0.35 ), 23.56 ) );\n    res = opU( res, vec2( sdTriPrism(   pos-vec3(-1.0,0.15,-1.0), vec2(0.3,0.05) ),43.5 ) );\n    res = opU( res, vec2( sdEllipsoid(  pos-vec3(-1.0,0.25, 0.0), vec3(0.2, 0.25, 0.05) ), 43.17 ) );\n    res = opU( res, vec2( sdHorseshoe(  pos-vec3(-1.0,0.25, 1.0), vec2(cos(1.3),sin(1.3)), 0.2, 0.3, vec2(0.03,0.08) ), 11.5 ) );\n    }\n\n    // bounding box\n    if( sdBox( pos-vec3(2.0,0.3,-1.0),vec3(0.35,0.3,2.5) )<res.x )\n    {\n    res = opU( res, vec2( sdOctogonPrism(pos-vec3( 2.0,0.2,-3.0), 0.2, 0.05), 51.8 ) );\n    res = opU( res, vec2( sdCylinder(    pos-vec3( 2.0,0.14,-2.0), vec3(0.1,-0.1,0.0), vec3(-0.2,0.35,0.1), 0.08), 31.2 ) );\n\tres = opU( res, vec2( sdCappedCone(  pos-vec3( 2.0,0.09,-1.0), vec3(0.1,0.0,0.0), vec3(-0.2,0.40,0.1), 0.15, 0.05), 46.1 ) );\n    res = opU( res, vec2( sdRoundCone(   pos-vec3( 2.0,0.15, 0.0), vec3(0.1,0.0,0.0), vec3(-0.1,0.35,0.1), 0.15, 0.05), 51.7 ) );\n    res = opU( res, vec2( sdRoundCone(   pos-vec3( 2.0,0.20, 1.0), 0.2, 0.1, 0.3 ), 37.0 ) );\n    }\n    \n    return res;\n}\n\n// https://iquilezles.org/articles/boxfunctions\nvec2 iBox( in vec3 ro, in vec3 rd, in vec3 rad ) \n{\n    vec3 m = 1.0/rd;\n    vec3 n = m*ro;\n    vec3 k = abs(m)*rad;\n    vec3 t1 = -n - k;\n    vec3 t2 = -n + k;\n\treturn vec2( max( max( t1.x, t1.y ), t1.z ),\n\t             min( min( t2.x, t2.y ), t2.z ) );\n}\n\nvec2 raycast( in vec3 ro, in vec3 rd )\n{\n    vec2 res = vec2(-1.0,-1.0);\n\n    float tmin = 1.0;\n    float tmax = 20.0;\n\n    // raytrace floor plane\n    float tp1 = (0.0-ro.y)/rd.y;\n    if( tp1>0.0 )\n    {\n        tmax = min( tmax, tp1 );\n        res = vec2( tp1, 1.0 );\n    }\n    //else return res;\n    \n    // raymarch primitives   \n    vec2 tb = iBox( ro-vec3(0.0,0.4,-0.5), rd, vec3(2.5,0.41,3.0) );\n    if( tb.x<tb.y && tb.y>0.0 && tb.x<tmax)\n    {\n        //return vec2(tb.x,2.0);\n        tmin = max(tb.x,tmin);\n        tmax = min(tb.y,tmax);\n\n        float t = tmin;\n        for( int i=0; i<70 && t<tmax; i++ )\n        {\n            vec2 h = map( ro+rd*t );\n            if( abs(h.x)<(0.0001*t) )\n            { \n                res = vec2(t,h.y); \n                break;\n            }\n            t += h.x;\n        }\n    }\n    \n    return res;\n}\n\n// https://iquilezles.org/articles/rmshadows\nfloat calcSoftshadow( in vec3 ro, in vec3 rd, in float mint, in float tmax )\n{\n    // bounding volume\n    float tp = (0.8-ro.y)/rd.y; if( tp>0.0 ) tmax = min( tmax, tp );\n\n    float res = 1.0;\n    float t = mint;\n    for( int i=ZERO; i<24; i++ )\n    {\n\t\tfloat h = map( ro + rd*t ).x;\n        float s = clamp(8.0*h/t,0.0,1.0);\n        res = min( res, s );\n        t += clamp( h, 0.01, 0.2 );\n        if( res<0.004 || t>tmax ) break;\n    }\n    res = clamp( res, 0.0, 1.0 );\n    return res*res*(3.0-2.0*res);\n}\n\n// https://iquilezles.org/articles/normalsSDF\nvec3 calcNormal( in vec3 pos )\n{\n#if 0\n    vec2 e = vec2(1.0,-1.0)*0.5773*0.0005;\n    return normalize( e.xyy*map( pos + e.xyy ).x + \n\t\t\t\t\t  e.yyx*map( pos + e.yyx ).x + \n\t\t\t\t\t  e.yxy*map( pos + e.yxy ).x + \n\t\t\t\t\t  e.xxx*map( pos + e.xxx ).x );\n#else\n    // inspired by tdhooper and klems - a way to prevent the compiler from inlining map() 4 times\n    vec3 n = vec3(0.0);\n    for( int i=ZERO; i<4; i++ )\n    {\n        vec3 e = 0.5773*(2.0*vec3((((i+3)>>1)&1),((i>>1)&1),(i&1))-1.0);\n        n += e*map(pos+0.0005*e).x;\n      //if( n.x+n.y+n.z>100.0 ) break;\n    }\n    return normalize(n);\n#endif    \n}\n\n// https://iquilezles.org/articles/nvscene2008/rwwtt.pdf\nfloat calcAO( in vec3 pos, in vec3 nor )\n{\n\tfloat occ = 0.0;\n    float sca = 1.0;\n    for( int i=ZERO; i<5; i++ )\n    {\n        float h = 0.01 + 0.12*float(i)/4.0;\n        float d = map( pos + h*nor ).x;\n        occ += (h-d)*sca;\n        sca *= 0.95;\n        if( occ>0.35 ) break;\n    }\n    return clamp( 1.0 - 3.0*occ, 0.0, 1.0 ) * (0.5+0.5*nor.y);\n}\n\n// https://iquilezles.org/articles/checkerfiltering\nfloat checkersGradBox( in vec2 p, in vec2 dpdx, in vec2 dpdy )\n{\n    // filter kernel\n    vec2 w = abs(dpdx)+abs(dpdy) + 0.001;\n    // analytical integral (box filter)\n    vec2 i = 2.0*(abs(fract((p-0.5*w)*0.5)-0.5)-abs(fract((p+0.5*w)*0.5)-0.5))/w;\n    // xor pattern\n    return 0.5 - 0.5*i.x*i.y;                  \n}\n\nvec3 render( in vec3 ro, in vec3 rd, in vec3 rdx, in vec3 rdy )\n{ \n    // background\n    vec3 col = vec3(0.7, 0.7, 0.9) - max(rd.y,0.0)*0.3;\n    \n    // raycast scene\n    vec2 res = raycast(ro,rd);\n    float t = res.x;\n\tfloat m = res.y;\n    if( m>-0.5 )\n    {\n        vec3 pos = ro + t*rd;\n        vec3 nor = (m<1.5) ? vec3(0.0,1.0,0.0) : calcNormal( pos );\n        vec3 ref = reflect( rd, nor );\n        \n        // material        \n        col = 0.2 + 0.2*sin( m*2.0 + vec3(0.0,1.0,2.0) );\n        float ks = 1.0;\n        \n        if( m<1.5 )\n        {\n            // project pixel footprint into the plane\n            vec3 dpdx = ro.y*(rd/rd.y-rdx/rdx.y);\n            vec3 dpdy = ro.y*(rd/rd.y-rdy/rdy.y);\n\n            float f = checkersGradBox( 3.0*pos.xz, 3.0*dpdx.xz, 3.0*dpdy.xz );\n            col = 0.15 + f*vec3(0.05);\n            ks = 0.4;\n        }\n\n        // lighting\n        float occ = calcAO( pos, nor );\n        \n\t\tvec3 lin = vec3(0.0);\n\n        // sun\n        {\n            vec3  lig = normalize( vec3(-0.5, 0.4, -0.6) );\n            vec3  hal = normalize( lig-rd );\n            float dif = clamp( dot( nor, lig ), 0.0, 1.0 );\n          //if( dif>0.0001 )\n        \t      dif *= calcSoftshadow( pos, lig, 0.02, 2.5 );\n\t\t\tfloat spe = pow( clamp( dot( nor, hal ), 0.0, 1.0 ),16.0);\n                  spe *= dif;\n                  spe *= 0.04+0.96*pow(clamp(1.0-dot(hal,lig),0.0,1.0),5.0);\n                //spe *= 0.04+0.96*pow(clamp(1.0-sqrt(0.5*(1.0-dot(rd,lig))),0.0,1.0),5.0);\n            lin += col*2.20*dif*vec3(1.30,1.00,0.70);\n            lin +=     5.00*spe*vec3(1.30,1.00,0.70)*ks;\n        }\n        // sky\n        {\n            float dif = sqrt(clamp( 0.5+0.5*nor.y, 0.0, 1.0 ));\n                  dif *= occ;\n            float spe = smoothstep( -0.2, 0.2, ref.y );\n                  spe *= dif;\n                  spe *= 0.04+0.96*pow(clamp(1.0+dot(nor,rd),0.0,1.0), 5.0 );\n          //if( spe>0.001 )\n                  spe *= calcSoftshadow( pos, ref, 0.02, 2.5 );\n            lin += col*0.60*dif*vec3(0.40,0.60,1.15);\n            lin +=     2.00*spe*vec3(0.40,0.60,1.30)*ks;\n        }\n        // back\n        {\n        \tfloat dif = clamp( dot( nor, normalize(vec3(0.5,0.0,0.6))), 0.0, 1.0 )*clamp( 1.0-pos.y,0.0,1.0);\n                  dif *= occ;\n        \tlin += col*0.55*dif*vec3(0.25,0.25,0.25);\n        }\n        // sss\n        {\n            float dif = pow(clamp(1.0+dot(nor,rd),0.0,1.0),2.0);\n                  dif *= occ;\n        \tlin += col*0.25*dif*vec3(1.00,1.00,1.00);\n        }\n        \n\t\tcol = lin;\n\n        col = mix( col, vec3(0.7,0.7,0.9), 1.0-exp( -0.0001*t*t*t ) );\n    }\n\n\treturn vec3( clamp(col,0.0,1.0) );\n}\n\nmat3 setCamera( in vec3 ro, in vec3 ta, float cr )\n{\n\tvec3 cw = normalize(ta-ro);\n\tvec3 cp = vec3(sin(cr), cos(cr),0.0);\n\tvec3 cu = normalize( cross(cw,cp) );\n\tvec3 cv =          ( cross(cu,cw) );\n    return mat3( cu, cv, cw );\n}\n\nvoid mainImage( out vec4 fragColor, in vec2 fragCoord )\n{\n    vec2 mo = iMouse.xy/iResolution.xy;\n\tfloat time = 32.0 + iTime*1.5;\n\n    // camera\t\n    vec3 ta = vec3( 0.25, -0.75, -0.75 );\n    vec3 ro = ta + vec3( 4.5*cos(0.1*time + 7.0*mo.x), 2.2, 4.5*sin(0.1*time + 7.0*mo.x) );\n    // camera-to-world transformation\n    mat3 ca = setCamera( ro, ta, 0.0 );\n\n    vec3 tot = vec3(0.0);\n#if AA>1\n    for( int m=ZERO; m<AA; m++ )\n    for( int n=ZERO; n<AA; n++ )\n    {\n        // pixel coordinates\n        vec2 o = vec2(float(m),float(n)) / float(AA) - 0.5;\n        vec2 p = (2.0*(fragCoord+o)-iResolution.xy)/iResolution.y;\n#else    \n        vec2 p = (2.0*fragCoord-iResolution.xy)/iResolution.y;\n#endif\n\n        // focal length\n        const float fl = 2.5;\n        \n        // ray direction\n        vec3 rd = ca * normalize( vec3(p,fl) );\n\n         // ray differentials\n        vec2 px = (2.0*(fragCoord+vec2(1.0,0.0))-iResolution.xy)/iResolution.y;\n        vec2 py = (2.0*(fragCoord+vec2(0.0,1.0))-iResolution.xy)/iResolution.y;\n        vec3 rdx = ca * normalize( vec3(px,fl) );\n        vec3 rdy = ca * normalize( vec3(py,fl) );\n        \n        // render\t\n        vec3 col = render( ro, rd, rdx, rdy );\n\n        // gain\n        // col = col*3.0/(2.5+col);\n        \n\t\t// gamma\n        col = pow( col, vec3(0.4545) );\n\n        tot += col;\n#if AA>1\n    }\n    tot /= float(AA*AA);\n#endif\n    \n    fragColor = vec4( tot, 1.0 );\n}\n","type":"Image"}]}