//!Warning does not model intersection with the Planet itself in the middle of the atmosphere
float calculateHazeFactor(in float fragDepth, in vec3 viewDir)
{
    //hazeHeight and camPos are heights above sea level 
    bool aboveHaze = camPos>hazeHeight;
    if (aboveHaze&&viewDir.y>=0.0)
        return 1.0;

    float planetRadius = 6357.0;
    float hazeRadius = planetRadius+hazeHeight;
    float distToEarthCentre = planetRadius+camPos;

    float cos_thetaD = -distToEarthCentre*viewDir.y/length(viewDir);
    //float discriminant = cos_theta*cos_theta-1.0+(hazeRadius/distToEarthCentre)*(hazeRadius/distToEarthCentre);
    float discriminant = cos_thetaD*cos_thetaD-distToEarthCentre*distToEarthCentre+hazeRadius*hazeRadius;
    if (discriminant<0.0)
        return 1.0;

    float discrSqrt = sqrt(discriminant);
    float hazeLen = min(cos_thetaD+discrSqrt,fragDepth); //ray end
    if (aboveHaze)
        hazeLen -= min(cos_thetaD-discrSqrt,fragDepth); //optional rayStart

    return hazeLen;
}