uniform sampler2D SCREEN_TEXTURE;

// Number of frequencies to combine, can't be a parameter/uniform else it causes problems in compatibility mode
// I have no idea why
const highp int N = 32;

uniform highp float spatial_scale; // spatial scaling of modes, higher = fine turbulence, lower = coarse turbulence
uniform highp float strength_scale; // distortion strength
uniform highp float speed_scale; // scaling factor on the speed of the animation
// Matrix to convert screen coordinates into grid coordinates
// This is to "pin" the effect to the grid, so that it does not shimmer as you move
uniform highp mat3 grid_ent_from_viewport_local;

const highp float TWO_PI = 6.28318530718;
 // This is just the default target values so that the external parameters can be normalized to 1
const highp float strength_factor = 0.0005;
const highp float spatial_factor = 22.0;

// 1D pseudo-random function
highp float random_1d(highp float n) {
    return fract(sin(n * 12.9898) * 43758.5453);
}

// Kolmogorov amplitude, power spectrum goes as k^(–11/6)
highp float kolAmp(highp float k) {
    return pow(k, -11.0 / 6.0);
}

void fragment() {

    highp vec2 ps = vec2(1.0/SCREEN_PIXEL_SIZE.x, 1.0/SCREEN_PIXEL_SIZE.y);
    highp float aspectratio = ps.x / ps.y;

    // scale the scale factor with the number of modes just cuz it works reasonably
    highp float s_scale = spatial_scale * spatial_factor / sqrt(float(N));

    // Coordinates to use to calculate the effects, convert to grid coordinates
    highp vec2 uvW = (grid_ent_from_viewport_local * vec3(UV.x, UV.y, 1.0)).xy;
    // Scale the coordinates
    uvW *= s_scale;

    // accumulate phase gradienta
    highp vec2 grad = vec2(0.0);

    for (lowp int i = 0; i < N; i++) {
        // float cast of the index
        highp float fi = float(i);

        // Pick a random direction
        highp float ang = random_1d(fi + 1.0) * TWO_PI;
        highp vec2 dir = vec2(cos(ang), sin(ang));

        // Pick a random spatial frequency from 0.5 to 30
        highp float k = mix(0.5, 30.0, random_1d(fi + 17.0));

        // Pick a random speed from 0.05 to 0.20
        highp float speed = mix(3., 8., random_1d(fi + 33.0));

        // Pick a random phase offset
        highp float phi_0 = random_1d(fi + 49.0) * TWO_PI;

        // phase argument
        highp float t = dot(dir, uvW) * k + TIME * speed * speed_scale + phi_0;

        // analytical gradient: ∇[sin(t)] = cos(t) * ∇t
        // ∇t = k * dir * scale (scale is factored out)
        grad += kolAmp(k) * cos(t) * k * dir;
    }
    // Spatial scaling (coarse or fine turbulence)
    grad *= s_scale;

    // The texture should have been blurred using a previous operation
    // We use the alpha channel to cut off the blur that bleeds outside the tile, then we rescale
    // the mask back up to 0.0 to 1.0
    highp float mask = clamp((zTexture(UV).a - 0.5)*2.0, 0.00, 1.0);

    // Calculate warped UV using the turbulence gradient
    // The strength of the turbulence is encoded into the red channel of TEXTURE
    // Give it a little polynomial boost: https://www.wolframalpha.com/input?i=-x%5E2+%2B2x+from+0+to+1
    highp float heatStrength = zTexture(UV).r*1.0;
    heatStrength = clamp(-heatStrength*heatStrength + 2.0*heatStrength, 0.0, 1.0);
    highp vec2 uvDist = UV + (strength_scale * strength_factor * heatStrength * mask) * grad;

    // Apply to the texture
    COLOR = texture2D(SCREEN_TEXTURE, uvDist);

    // Uncomment the following two lines to view the strength buffer directly
    // COLOR.rgb = vec3(heatStrength * mask);
    // COLOR.a = mask;
}