Atmospheric refraction affects to various degrees exoplanet transit, lunar eclipse, as well as stellar occultation observations. Exoplanet retrieval algorithms often use analytical expressions for the column abundance along a ray traversing the atmosphere as well as for the deflection of that ray, which are first order approximations valid for low densities in a spherically symmetric homogeneous isothermal atmosphere. We derive new analytical formulae for both of these quantities, which are valid for higher densities, and use them to refine and validate a new ray tracing algorithm which can be used for arbitrary atmospheric temperature-pressure profiles. We illustrate with simple isothermal atmospheric profiles the consequences of our model for different planets: temperate Earth-like and Jovian-like planets, as well as HD189733b, and GJ1214b. We find that, for both hot exoplanets, our treatment of refraction does not make much of a difference to pressures as high as 10 atmosphere, but that it is important to consider the variation of gravity with altitude for GJ1214b. However, we find that the temperate atmospheres have an apparent scale height significantly smaller than their actual density scale height at densities larger than 1 amagat, thus increasing the difficulty of detecting spectral features originating in these regions. These denser atmospheric regions form a refractive boundary layer where column abundances and ray deflection increases dramatically with decreasing impact parameter. This refractive boundary layer mimics a surface, and none of the techniques mentioned above can probe atmospheric regions denser than about 4 amagat on these temperate planets.