     return fract(dot(rnd/m, vec4(1.0, -1.0, 1.0, -1.0)));
 }
 
-// find first hit over specified segment
-bool traceDistance(in vec3 start, in vec3 dir, float maxDist,
-                   inout vec3 hit, inout int material)
-{
-    float minMarch = minLengthScale;
-    const float HUGE_VAL = 1.0e20;
-    float sdf = HUGE_VAL;
-#ifdef HAS_SURFACE
-    float sdf_surfa = abs(SDF_SURFACE(start));    sdf = min(sdf, sdf_surfa);
-#endif
-#ifdef HAS_METAL
-    float sdf_metal = abs(SDF_METAL(start));      sdf = min(sdf, sdf_metal);
-#endif
-#ifdef HAS_DIELECTRIC
-    float sdf_diele = abs(SDF_DIELECTRIC(start)); sdf = min(sdf, sdf_diele);
-#endif
-#ifdef HAS_VOLUME
-    float sdf_volum = abs(SDF_VOLUME(start));     sdf = min(sdf, sdf_volum);
-#endif
-
-    float InitialSign = sign(sdf);
-    float t = InitialSign * sdf; // (always take the first step along the ray direction)
-    if (t>=maxDist) return false;
-    for (int n=0; n<__MAX_MARCH_STEPS__; n++)
-    {
-        vec3 pW = start + t*dir;
-        sdf = HUGE_VAL;
-#ifdef HAS_SURFACE
-        sdf_surfa = abs(SDF_SURFACE(pW));    if (sdf_surfa<minMarch) { material = MAT_SURFA; hit = start + t*dir; return true; } sdf = min(sdf, sdf_surfa);
-#endif
-#ifdef HAS_METAL
-        sdf_metal = abs(SDF_METAL(pW));      if (sdf_metal<minMarch) { material = MAT_METAL; hit = start + t*dir; return true; } sdf = min(sdf, sdf_metal);
-#endif
-#ifdef HAS_DIELECTRIC
-        sdf_diele = abs(SDF_DIELECTRIC(pW)); if (sdf_diele<minMarch) { material = MAT_DIELE; hit = start + t*dir; return true; } sdf = min(sdf, sdf_diele);
-#endif
-#ifdef HAS_VOLUME
-        sdf_volum = abs(SDF_VOLUME(pW));     if (sdf_volum<minMarch) { material = MAT_VOLUM; hit = start + t*dir; return true; } sdf = min(sdf, sdf_volum);
-#endif
-        // With this formula, the ray advances whether sdf is initially negative or positive --
-        // but on crossing the zero isosurface, sdf flips allowing bracketing of the root.
-        t += InitialSign * sdf;
-        if (t>=maxDist) return false;
-    }
-    return !maxStepsIsMiss;
-}
-
-// find first hit along infinite ray
-bool traceRay(in vec3 start, in vec3 dir,
-              inout vec3 hit, inout int material, float maxMarchDist)
-{
-    material = MAT_INVAL;
-    return traceDistance(start, dir, maxMarchDist, hit, material);
-}
-
-vec3 normal(in vec3 pW, int material)
-{
-    // Compute normal as gradient of SDF
-    float normalEpsilon = 2.0*minLengthScale;
-    vec3 e = vec3(normalEpsilon, 0.0, 0.0);
-    vec3 xyyp = pW+e.xyy; vec3 xyyn = pW-e.xyy;
-    vec3 yxyp = pW+e.yxy; vec3 yxyn = pW-e.yxy;
-    vec3 yyxp = pW+e.yyx; vec3 yyxn = pW-e.yyx;
-    vec3 N;
-#ifdef HAS_SURFACE
-    if (material==MAT_SURFA) { N = vec3(   SDF_SURFACE(xyyp) -    SDF_SURFACE(xyyn),    SDF_SURFACE(yxyp) -    SDF_SURFACE(yxyn),    SDF_SURFACE(yyxp) -    SDF_SURFACE(yyxn)); return normalize(N); }
-#endif
-#ifdef HAS_METAL
-    if (material==MAT_METAL) { N = vec3(     SDF_METAL(xyyp) -      SDF_METAL(xyyn),      SDF_METAL(yxyp) -      SDF_METAL(yxyn),      SDF_METAL(yyxp) -      SDF_METAL(yyxn)); return normalize(N); }
-#endif
-#ifdef HAS_DIELECTRIC
-    if (material==MAT_DIELE) { N = vec3(SDF_DIELECTRIC(xyyp) - SDF_DIELECTRIC(xyyn), SDF_DIELECTRIC(yxyp) - SDF_DIELECTRIC(yxyn), SDF_DIELECTRIC(yyxp) - SDF_DIELECTRIC(yyxn)); return normalize(N); }
-#endif
-#ifdef HAS_VOLUME
-    if (material==MAT_VOLUM) { N = vec3(    SDF_VOLUME(xyyp) -     SDF_VOLUME(xyyn),     SDF_VOLUME(yxyp) -     SDF_VOLUME(yxyn),     SDF_VOLUME(yyxp) -     SDF_VOLUME(yyxn)); return normalize(N); }
-#endif
-}
-
-// MC-estimates the amount of light transmitted along an infinite ray
-float Visibility(in vec3 pW, in vec3 rayDir, inout vec4 rnd, bool inVolume)
-{
-    vec3 pW_next;
-    int hitMaterial;
-    bool hit = traceRay(pW, rayDir, pW_next, hitMaterial, maxLengthScale);
-#ifndef HAS_VOLUME
-    return float(!hit);
-#else
-    // Deal with all volume interactions until the ray either hits a surface or escapes
-    const int maxVolHits = 16;
-    int volHits = 0;
-    float sigma_t_max = VOLUME_EXTINCTION_MAX();
-    float inv_sigma_t_max = 1.0/sigma_t_max;
-    float eps = 3.0*minLengthScale;
-    const int maxWoodcockSteps = 128;
-    while ((inVolume || hitMaterial==MAT_VOLUM) && volHits<maxVolHits)
-    {
-        // If ray lies in the volume interior
-        if (inVolume)
-        {
-            // Do Woodcock tracking to check for a possible volume scattering event along the in-volume segment pW -> pW_next
-            // (except if the current ray lies inside a dielectric, in which case there is no volumetric scattering, i.e. the dielectric displaces the volume)
-            vec3 pScatter;
-            float sigma_t;
-            bool interacted = false;
-            float segmentLength = length(pW_next - pW);
-            float distanceMarched = -log(rand(rnd)) * inv_sigma_t_max;
-            int steps = 0;
-            while (distanceMarched<segmentLength-eps && steps<maxWoodcockSteps)
-            {
-                pScatter = pW + distanceMarched*rayDir;
-                sigma_t = VOLUME_EXTINCTION(pScatter);
-                interacted = bool(rand(rnd) < sigma_t * inv_sigma_t_max);
-                if (interacted) break;
-                distanceMarched += -log(rand(rnd)) * inv_sigma_t_max;
-                steps++;
-            }
-
-            if (interacted) return 0.0; // on absorption/scattering, return zero visibility
-            else
-            {
-                pW = pW_next;
-                if (hitMaterial==MAT_VOLUM) // if hitMaterial is volume boundary, displace pW into volume exterior
-                {
-                    vec3 nW = normal(pW, MAT_VOLUM);
-                    pW += nW * eps;
-                    inVolume = false;
-                }
-                else // if hitMaterial is a surface, return zero visibility
-                    return 0.0;
-            }
-        }
-
-        // else if ray passed through the volume exterior and hit the volume boundary
-        else // !inVolume
-        {
-            pW = pW_next; // advance to the volume boundary
-            vec3 nW = normal(pW, MAT_VOLUM); // displace pW into volume interior
-            pW -= nW * eps;
-            inVolume = true;
-        }
-
-        // Re-trace from new pW
-        hit = traceRay(pW, rayDir, pW_next, hitMaterial, maxLengthScale);
-        volHits++;
-    }
-    return float(!hit || volHits==maxVolHits); // if no hit or scattering, return visibility 1.0 (transparent)
-#endif // HAS_VOLUME
-}
+//////////////////////////////////////////////////////////////
+// Dynamically injected code
+//////////////////////////////////////////////////////////////
+
+__DEFINES__
+
+__SHADER__
+
+__IOR_FUNC__
 
 
 /////////////////////////////////////////////////////////////////////////


vec3 samplePhaseFunction(in vec3 rayDir, inout vec4 rnd)
{
    float g = 0.0; // @todo: parameter
    float costheta;
    if (abs(g)<1.0e-3)
        costheta = 1.0 - 2.0*rand(rnd);

                if (interacted)
                {
                    vec3 wiW = samplePhaseFunction(rayDir, rnd); // sample scattered dir (NB, PF is the angle PDF so the MC PDF denom. cancels it in the L estimator)
                    float emission = VOLUME_EMISSION_EVAL(pScatter, XYZ);
                    float albedo = VOLUME_ALBEDO_EVAL(pScatter, XYZ);
                    absorb = (rand(rnd) > albedo);
                    if (absorb)
                    {
