#ifndef ASTRAL_COLOSSUS_PLANET_GLSL
#define ASTRAL_COLOSSUS_PLANET_GLSL

uniform sampler2D astralColossusAtmosphere;

/* Astral Colossus
 *
 * A true analytic sphere placed behind the real sun. The image sampler is
 * only atmospheric material: silhouette, depth, surface normal, lighting,
 * differential jet motion, storm evolution, time palette, limb and corona
 * are all rebuilt in shader space.
 *
 * RGB is premultiplied radiance and alpha is displaced sky energy:
 *     sky = sky * (1.0 - planet.a) + planet.rgb;
 */

float colossus_hash12(vec2 p) {
    vec3 p3 = fract(vec3(p.xyx) * 0.1031);
    p3 += dot(p3, p3.yzx + 33.33);
    return fract((p3.x + p3.y) * p3.z);
}

float colossus_mirror(float value) {
    float folded = fract(value * 0.5) * 2.0;
    return 1.0 - abs(folded - 1.0);
}

vec3 colossus_rotate_z(vec3 p, float angle) {
    float c = cos(angle);
    float s = sin(angle);
    return vec3(c * p.x - s * p.y,
        s * p.x + c * p.y, p.z);
}

vec3 colossus_rotate_y(vec3 p, float angle) {
    float c = cos(angle);
    float s = sin(angle);
    return vec3(c * p.x + s * p.z,
        p.y, -s * p.x + c * p.z);
}

vec2 colossus_atmosphere_uv(
        vec3 normalDirection, float longitudeFlow,
        float latitudeFlow) {
    vec3 n = normalize(normalDirection);
    float longitude = atan(n.z, n.x) * 0.15915494309189535
        + 0.5 + longitudeFlow;
    float latitude = asin(clamp(n.y, -1.0, 1.0))
        * 0.3183098861837907 + 0.5 + latitudeFlow;
    return vec2(
        colossus_mirror(longitude),
        clamp(latitude, 0.002, 0.998));
}

float colossus_visibility(vec3 sunDirection) {
    // Retain the planet through golden hour but hand the sky completely to
    // the rift at night. Any established rain system closes the daytime
    // event completely; it must never remain visible inside a storm deck.
    float daylight = day_blend_float(0.38, 1.0, 0.0);
    daylight = smoothstep(0.025, 0.66, daylight);
    float sunHeight = smoothstep(-0.16, 0.075, sunDirection.y);
    float rainVeil = 1.0 - smoothstep(
        0.008, 0.085, clamp(rainStrength, 0.0, 1.0));
    return daylight * sunHeight * rainVeil;
}

void colossus_time_palette(
        out vec3 deepColor,
        out vec3 midColor,
        out vec3 highColor,
        out vec3 energyColor,
        out vec3 atmosphereColor) {
    float time = clamp(day_moment, 0.0, 1.0);
    float dawnToNoon = smoothstep(0.015, 0.245, time);
    float noonToDusk = smoothstep(0.285, 0.505, time);

    const vec3 dawnDeep = vec3(0.095, 0.025, 0.105);
    const vec3 dawnMid = vec3(0.39, 0.135, 0.18);
    const vec3 dawnHigh = vec3(0.94, 0.48, 0.21);
    const vec3 dawnEnergy = vec3(1.12, 0.40, 0.16);
    const vec3 dawnAtmosphere = vec3(0.92, 0.36, 0.30);

    const vec3 noonDeep = vec3(0.018, 0.075, 0.19);
    const vec3 noonMid = vec3(0.055, 0.33, 0.46);
    const vec3 noonHigh = vec3(0.58, 0.84, 0.90);
    const vec3 noonEnergy = vec3(0.18, 0.78, 1.24);
    const vec3 noonAtmosphere = vec3(0.25, 0.70, 1.08);

    const vec3 duskDeep = vec3(0.105, 0.018, 0.145);
    const vec3 duskMid = vec3(0.43, 0.075, 0.26);
    const vec3 duskHigh = vec3(0.92, 0.32, 0.20);
    const vec3 duskEnergy = vec3(1.16, 0.25, 0.34);
    const vec3 duskAtmosphere = vec3(0.72, 0.23, 0.72);

    deepColor = mix(dawnDeep, noonDeep, dawnToNoon);
    midColor = mix(dawnMid, noonMid, dawnToNoon);
    highColor = mix(dawnHigh, noonHigh, dawnToNoon);
    energyColor = mix(dawnEnergy, noonEnergy, dawnToNoon);
    atmosphereColor = mix(
        dawnAtmosphere, noonAtmosphere, dawnToNoon);

    deepColor = mix(deepColor, duskDeep, noonToDusk);
    midColor = mix(midColor, duskMid, noonToDusk);
    highColor = mix(highColor, duskHigh, noonToDusk);
    energyColor = mix(energyColor, duskEnergy, noonToDusk);
    atmosphereColor = mix(
        atmosphereColor, duskAtmosphere, noonToDusk);
}

vec4 getAstralColossus(
        vec3 viewDirection,
        vec3 sunWorldDirection,
        vec3 baseSkyColor) {
    #ifndef ASTRAL_COLOSSUS_PLANET
        return vec4(0.0);
    #else
        vec3 viewDir = normalize(viewDirection);
        vec3 sunDir = normalize(sunWorldDirection);
        float visibility = colossus_visibility(sunDir);
        if (visibility < 0.002 || viewDir.y < -0.055)
            return vec4(0.0);

        float clock = frameTimeCounter * ASTRAL_COLOSSUS_MOTION;
        float scale = max(ASTRAL_COLOSSUS_SCALE, 0.35);

        // The sun remains in front of the lower hemisphere. Only a very slow
        // positional parallax occurs; all prominent animation belongs to the
        // atmosphere, storms, energy fronts and corona rather than a sprite
        // orbiting around the screen.
        vec3 frameUp = abs(sunDir.y) > 0.94
            ? vec3(1.0, 0.0, 0.0)
            : vec3(0.0, 1.0, 0.0);
        vec3 solarX = normalize(cross(frameUp, sunDir));
        vec3 solarY = normalize(cross(sunDir, solarX));
        vec2 centerOffset = vec2(
            sin(clock * 0.0053) * 0.035,
            0.185 + cos(clock * 0.0037) * 0.022);
        vec3 planetDirection = normalize(
            sunDir + solarX * centerOffset.x
                + solarY * centerOffset.y);

        vec3 planetUp = abs(planetDirection.y) > 0.94
            ? solarX : vec3(0.0, 1.0, 0.0);
        vec3 planetX = normalize(cross(planetUp, planetDirection));
        vec3 planetY = normalize(cross(planetDirection, planetX));
        vec2 planetPlane = vec2(
            dot(viewDir, planetX), dot(viewDir, planetY));
        float planetDot = clamp(
            dot(viewDir, planetDirection), -1.0, 1.0);
        float radial = sqrt(max(2.0 - 2.0 * planetDot, 0.0));
        float planetRadius = 0.34 * scale;
        float atmosphereReach = 0.155 * scale;

        // Reject most of the sky before either texture sample. This keeps the
        // giant event affordable even though its visible detail is animated.
        if (radial > planetRadius + atmosphereReach)
            return vec4(0.0);

        float edgeAA = max(fwidth(radial), 0.0012);
        float sphereMask = 1.0 - smoothstep(
            planetRadius - edgeAA * 1.25,
            planetRadius + edgeAA * 1.25, radial);

        float angularRadiusSin = planetRadius * sqrt(max(
            1.0 - planetRadius * planetRadius * 0.25, 0.02));
        vec2 discCoord = planetPlane / max(angularRadiusSin, 0.001);
        float discLength = length(discCoord);
        vec2 surfaceCoord = discCoord / max(discLength, 1.0);
        float surfaceZ = sqrt(max(
            1.0 - dot(surfaceCoord, surfaceCoord), 0.0));
        vec3 surfaceNormal = normalize(vec3(
            surfaceCoord.x, surfaceCoord.y, surfaceZ));

        // A fixed axial tilt creates a recognizable planetary orientation.
        // Differential longitude shear then drives every latitude at a
        // different speed, so the atmosphere cannot read as one rotating map.
        vec3 tiltedNormal = colossus_rotate_z(
            surfaceNormal, -0.34);
        float latitude = tiltedNormal.y;
        float latitudeShearA = sin(
            latitude * 15.0 + clock * 0.115) * 0.020
            + sin(latitude * 31.0 - clock * 0.071) * 0.007;
        float longitudeFlowA = clock
            * (0.010 + abs(latitude) * 0.006)
            + latitudeShearA;
        float latitudeFlowA = sin(
            clock * 0.043 + latitude * 8.0) * 0.007;
        vec3 layerNormalA = colossus_rotate_y(
            tiltedNormal, clock * 0.012 + latitudeShearA * 4.0);
        vec2 atmosphereUvA = colossus_atmosphere_uv(
            layerNormalA, longitudeFlowA, latitudeFlowA);
        vec3 atmosphereA = texture2D(
            astralColossusAtmosphere, atmosphereUvA).rgb;

        #if PROFILE_QUALITY == 1
            vec3 atmosphereField = atmosphereA;
            float upperCloudSource = dot(
                atmosphereA, vec3(0.2126, 0.7152, 0.0722));
        #else
            // The high atmosphere moves against the deep jets, with local
            // displacement borrowed from layer A. This supplies parallax and
            // real morphing rather than a second copy rotating uniformly.
            float atmosphereALuma = dot(
                atmosphereA, vec3(0.2126, 0.7152, 0.0722));
            float latitudeShearB = sin(
                latitude * 21.0 - clock * 0.132
                    + atmosphereALuma * 3.8) * 0.024;
            vec3 layerNormalB = colossus_rotate_y(
                colossus_rotate_z(surfaceNormal, -0.29),
                -clock * 0.008 + latitudeShearB * 3.2);
            vec2 atmosphereUvB = colossus_atmosphere_uv(
                layerNormalB,
                -clock * (0.007 + abs(latitude) * 0.004)
                    + latitudeShearB
                    + (atmosphereA.b - atmosphereA.g) * 0.024,
                cos(clock * 0.037 + latitude * 11.0) * 0.009
                    + (atmosphereA.r - atmosphereA.b) * 0.010);
            vec3 atmosphereB = texture2D(
                astralColossusAtmosphere, atmosphereUvB).rgb;
            float layerExchange = 0.48
                + sin(clock * 0.067 + latitude * 7.0
                    + atmosphereALuma * 4.3) * 0.15;
            vec3 atmosphereField = mix(
                atmosphereA, atmosphereB, layerExchange);
            float upperCloudSource = dot(
                atmosphereB, vec3(0.2126, 0.7152, 0.0722));
        #endif

        float atmosphereLuma = dot(
            atmosphereField, vec3(0.2126, 0.7152, 0.0722));
        float atmosphereMax = max(atmosphereField.r,
            max(atmosphereField.g, atmosphereField.b));
        float atmosphereMin = min(atmosphereField.r,
            min(atmosphereField.g, atmosphereField.b));
        float chroma = atmosphereMax - atmosphereMin;
        float highClouds = smoothstep(
            0.42, 0.79, upperCloudSource + chroma * 0.18);
        float pressureValleys = 1.0 - smoothstep(
            0.10, 0.43, atmosphereLuma);

        float longitude = atan(
            tiltedNormal.z, tiltedNormal.x);

        // Cinematic lower-hemisphere reveal. The upper half dissolves into
        // the sky through a broad atmospheric transition instead of ending
        // at a mechanical straight cut. Slow pressure ripples and the live
        // atmosphere field gently disturb the dissolve boundary.
        float cinematicDissolve = (atmosphereLuma - 0.5) * 0.075
            + sin(longitude * 3.0 + clock * 0.037) * 0.018
            + sin(longitude * 7.0 - clock * 0.052) * 0.009;
        float cinematicLowerReveal = 1.0 - smoothstep(
            -0.105, 0.145, discCoord.y + cinematicDissolve);
        float travelingFront = 0.5 + 0.5 * sin(
            longitude * 3.2 + latitude * 12.0
                - clock * 0.29 + atmosphereLuma * 5.7);
        travelingFront = smoothstep(0.43, 0.87, travelingFront);
        float counterFront = 0.5 + 0.5 * sin(
            -longitude * 5.1 + latitude * 18.0
                + clock * 0.21 + upperCloudSource * 4.1);
        counterFront = smoothstep(0.57, 0.92, counterFront);

        // A narrow energy front sweeps across the full sphere. Unlike map
        // rotation it changes local light, pressure and cloud height as it
        // passes, producing a large readable animation within a few seconds.
        float energyFront = smoothstep(0.56, 0.92,
            travelingFront * 0.64 + counterFront * 0.36);
        float veinSignal = max(
            atmosphereField.b - atmosphereField.r * 0.46,
            atmosphereField.r - atmosphereField.b * 0.58);
        float energyVeins = smoothstep(
            0.075, 0.31, veinSignal + chroma * 0.24)
            * smoothstep(0.20, 0.72, atmosphereLuma + chroma);

        float eventClock = clock * 0.115;
        float eventLocal = fract(eventClock);
        float eventEnvelope = pow(max(
            sin(eventLocal * 3.14159265), 0.0), 8.0);
        float eventChance = smoothstep(0.34, 0.92,
            colossus_hash12(vec2(floor(eventClock), 83.1)));
        float stormDischarge = eventEnvelope * eventChance
            * energyVeins * (0.35 + energyFront * 0.65);

        vec3 deepColor;
        vec3 midColor;
        vec3 highColor;
        vec3 energyColor;
        vec3 atmosphereColor;
        colossus_time_palette(
            deepColor, midColor, highColor,
            energyColor, atmosphereColor);

        float lowToMid = smoothstep(
            0.07, 0.52, atmosphereLuma);
        float midToHigh = smoothstep(
            0.43, 0.88, atmosphereLuma + highClouds * 0.12);
        vec3 surfaceColor = mix(deepColor, midColor, lowToMid);
        surfaceColor = mix(surfaceColor, highColor, midToHigh);

        // Reconstruct a real spherical lighting response. The foreground sun
        // has a position inside the projected planet disc, so highlights and
        // the local atmospheric bloom point toward the visible sun instead of
        // being baked into the rotating texture.
        vec2 sunPlane = vec2(
            dot(sunDir, planetX), dot(sunDir, planetY))
            / max(angularRadiusSin, 0.001);
        vec3 localLightDirection = normalize(vec3(
            clamp(sunPlane.x, -1.0, 1.0),
            clamp(sunPlane.y, -1.0, 1.0), 0.88));
        float directLight = clamp(
            dot(surfaceNormal, localLightDirection), 0.0, 1.0);
        float limbDepth = pow(clamp(surfaceZ, 0.0, 1.0), 0.42);
        float sphericalLight = 0.30
            + directLight * 0.52 + limbDepth * 0.18;
        surfaceColor *= sphericalLight;
        surfaceColor *= 1.0 - pressureValleys * 0.14;
        surfaceColor += highColor * highClouds
            * (0.035 + energyFront * 0.035);

        float sunContactDistance = length(
            discCoord - sunPlane);
        float sunContactGlow = exp(-sunContactDistance
            * sunContactDistance / 0.075);
        surfaceColor += mix(highColor, energyColor, 0.55)
            * sunContactGlow * (0.10 + energyFront * 0.08);

        // A perfectly circular physical body remains stable. All violent
        // motion occurs in atmosphere and energy, while only the outer gaseous
        // radius breathes like a magnetosphere under solar pressure.
        float atmosphereBreath = 0.84
            + sin(clock * 0.143) * 0.10
            + sin(clock * 0.057 + atmosphereLuma * 4.0) * 0.06;
        float limbDistance = abs(radial - planetRadius);
        float limbWidth = max(
            0.030 * scale * atmosphereBreath, edgeAA * 2.0);
        float luminousLimb = exp(-limbDistance * limbDistance
            / max(limbWidth * limbWidth, 0.00002));
        float outsideDistance = max(
            radial - planetRadius, 0.0);
        float outerAura = exp(-outsideDistance * outsideDistance
            / max(pow(0.105 * scale * atmosphereBreath, 2.0), 0.0001));
        outerAura *= 1.0 - sphereMask;

        float limbAngle = atan(planetPlane.y, planetPlane.x);
        float arcWave = 0.5 + 0.5 * sin(
            limbAngle * 8.0 - clock * 0.39
                + sin(limbAngle * 3.0 + clock * 0.087) * 1.4);
        float magneticArcs = pow(max(arcWave, 0.0), 12.0)
            * outerAura;
        float knotWaveA = 0.5 + 0.5 * sin(
            limbAngle * 7.0 - clock * 0.54);
        float knotWaveB = 0.5 + 0.5 * cos(
            limbAngle * 11.0 + clock * 0.37 + 1.2);
        float limbKnots = (pow(max(knotWaveA, 0.0), 18.0)
            + pow(max(knotWaveB, 0.0), 22.0) * 0.70)
            * luminousLimb;

        float bodyAlpha = sphereMask * cinematicLowerReveal
            * (0.74 + directLight * 0.10 + highClouds * 0.045)
            * visibility;
        float atmosphereAlpha = clamp(
            (luminousLimb * 0.055 + outerAura * 0.026)
                * cinematicLowerReveal
                * visibility * ASTRAL_COLOSSUS_INTENSITY,
            0.0, 0.095);
        float alpha = clamp(
            1.0 - (1.0 - bodyAlpha) * (1.0 - atmosphereAlpha),
            0.0, 0.90);

        vec3 radiance = surfaceColor * bodyAlpha;
        radiance += energyColor
            * (energyVeins * (0.055 + energyFront * 0.10)
                + stormDischarge * 0.30)
            * sphereMask * cinematicLowerReveal
            * visibility * ASTRAL_COLOSSUS_INTENSITY;
        radiance += atmosphereColor
            * (luminousLimb * (0.11 + energyFront * 0.045)
                + outerAura * 0.045
                + magneticArcs * 0.15
                + limbKnots * 0.24)
            * cinematicLowerReveal
            * visibility * ASTRAL_COLOSSUS_INTENSITY;
        radiance += max(baseSkyColor, vec3(0.0))
            * atmosphereAlpha * 0.42;

        return vec4(max(radiance, vec3(0.0)), alpha);
    #endif
}

#endif // ASTRAL_COLOSSUS_PLANET_GLSL
