What percentage of the sky could the sun be in?

Question

For every location on earth, there are parts of the sky that the sun can never be in (e.g., directly overhead anywhere outside the tropics).

I think the sun can be in about 40% of the sky (a hemisphere).

I'll withhold my reasoning/calculations for now, because I'd like the opinion(s) of others more-expert (most are) than me in this subject.

Please weigh in. I believe this "simple" geometric question will have been answered centuries ago, but I have not found it as an answer, nor even as a subject.

Answer 1

I'd like to expand more on Arvind H's answer. The Earth inclination is 23°, this means that the sun can have a maximum elevation above and below the celestial equator of $\alpha=23°$, as shown in the first figure.

Let's start by calculating the solid angle of the band on the celestial sphere between $\alpha$ and $-\alpha$ on the equator. In spherical coordinates where $\phi$ is the longitude and $\theta$ is the latitude:

$$\Omega = \int_0^{2\pi} \int_{-\alpha}^{+\alpha} \cos(\theta)d\theta d\phi = 4\pi \sin(\alpha)$$

Let's rapidly check that this is correct. If $\alpha=0$ this gives 0, and if $\alpha = \pi/2$ this gives $4\pi$, which is the whole sphere, as expected.

Now let's tackle the OP question. An observer in a particular location can only see half of the celestial sphere, that has non-trivial intersections with the area we have just calculated, as shown in the second figure.

Yet, the intersected area highlighted in the second figure is always exactly half of the band, i.e. $2 \pi \sin(\alpha)$. This is clearly true for the three cases I have drawn, but the symmetry of the problem implies that it is true in every case. A nice way to see it is by looking at the top and bottom pole cups where the Sun will never pass. The observer plane cuts them in a symmetrical way, so that if you combine the visible part of the top cup with the visible part of the bottom cup, you will always get a whole cup, as shown in the third figure.

This means that the percentage of the sky covered by the Sun at any given location is ${2\pi \sin\alpha \over 2\pi} = \sin\alpha \approx 40\%$

Answer 2

I was surprised to see that the 40% confirmed in other answers does not depend on latitude (but certainly does depend on Earth's axial inclination!) and curious what the patterns looked like and needed to satisfy my "one script per day" minimum to keep my coding skills from getting rusty, so I wrote the script below.

I considered 0, 15, 30, 45, 60, 75, and 90 degrees latitudes and sampled one million equally spaced time points in the year 2022.

Sure enough it's about 40% for the Sun and in 2022¹ about 46% for the Moon.

---

¹The Moon's orbit slowly varies over timescales longer than one year.

---

below: boolean registry of "hits" for each degree elevation (vertical from -90 to +90 wrt horizon) and azimuth (horizontal) for one million equally spaced time samples in 2022 at 90, 75, 60, 45, 30, 15 and 0 degrees latitude on the prime meridian. There's still a few missing spots due to finite sampling but I won't rerun with ten million timepoints unless I figure out over how many years I'd need to sample to get a good understanding of the Moon's behavior.

---

Script is set to 100,000 time points so that it does't take to long to run, feel free to increase

from skyfield.api import Loader, Topos
import numpy as np
import matplotlib.pyplot as plt
from skyfield import searchlib

load = Loader('~/Documents/fishing/SkyData') # avoid multiple copies of large files
ts = load.timescale()
eph = load('de421.bsp')

earth, sun, moon = [eph[x] for x in ('earth', 'sun', 'moon')]

nlats = 7
latitudes = np.linspace(0, 90, nlats) # degrees

np.save('latitudes.npy', latitudes)

t0, t1 = ts.utc([2022, 2023], 1, 1)
n_points = 100000
times =  ts.linspace(t0, t1, n_points)

sun_data, moon_data = [], []

for lat in latitudes:
    location = Topos(latitude_degrees = lat,
                     longitude_degrees = 0.0)
    moon_astrometric = (earth + location).at(times).observe(moon)
    moon_alts, moon_azs, d = moon_astrometric.apparent().altaz()

    sun_astrometric = (earth + location).at(times).observe(sun)
    sun_alts, sun_azs, d = sun_astrometric.apparent().altaz()

    sun_data.append(np.vstack([sun_alts.degrees, sun_azs.degrees]))
    moon_data.append(np.vstack([moon_alts.degrees, moon_azs.degrees]))

nalt, naz = 180, 361  # nalt is even so there's no bin cetered on horizon
sun_data = np.rint(sun_data).astype(int)
moon_data = np.rint(moon_data).astype(int)

np.save('sun_data.npy', sun_data)
np.save('moon_data.npy', moon_data)

sun_plots, moon_plots = np.zeros((2, nlats, nalt+1, naz), dtype=bool)

for i, (sd, md) in enumerate(zip(sun_data, moon_data)):
    sun_plots[i, sd[0]+(nalt>>1), sd[1]] = True
    moon_plots[i, md[0]+(nalt>>1), md[1]] = True

halfpi, pi, twopi, fourpi = [f * np.pi for f in (0.5, 1, 2, 4)] 
altitude, azimuth = np.mgrid[-halfpi:halfpi:(nalt+1)*1j, 0:twopi:naz*1j] 

dalt, daz = np.radians(1), np.radians(1) ### implicit from .degrees() and .int()

d_Omega = np.cos(altitude) * dalt * daz
in_the_sky = altitude > 0

solid_angle_of_sky = d_Omega[in_the_sky].sum()
print('solid_angle_of_sky / fourpi: ', solid_angle_of_sky / fourpi)

sun_percents = 100 * (sun_plots * d_Omega * in_the_sky).sum(axis=(1, 2)) / twopi
moon_percents = 100 * (moon_plots * d_Omega * in_the_sky).sum(axis=(1, 2)) / twopi
np.save('sun_percents.npy', sun_percents)
np.save('moon_percents.npy', moon_percents)

big_sun = sun_plots.reshape(-1, sun_plots.shape[-1])
big_moon = moon_plots.reshape(-1, moon_plots.shape[-1])

np.save('big_sun.npy', big_sun)
np.save('big_moon.npy', big_moon)

if True:
    if False:
        import numpy as np
        import matplotlib.pyplot as plt
        big_sun = np.load('big_sun.npy')        
        big_moon = np.load('big_moon.npy')        
    fig, (ax_sun, ax_moon) = plt.subplots(1, 2, figsize=[12, 7.5])
    ax_sun.imshow(big_sun)
    ax_sun.set_title('Sun', fontsize=16)
    ax_moon.imshow(big_moon)
    ax_moon.set_title('Moon', fontsize=16)
    plt.show()

if True:
    if False:
        import numpy as np
        import matplotlib.pyplot as plt
        sun_percents = np.load('sun_percents.npy')
        moon_percents = np.load('moon_percents.npy')
        latitudes = np.load('latitudes.npy')
        nlats = 7
        latitudes = np.linspace(0, 90, nlats) # degrees
    fig, ax = plt.subplots(1, 1)
    ax.plot(latitudes, sun_percents, '-r', label='Sun')
    ax.plot(latitudes, sun_percents, 'ok')
    ax.plot(latitudes, moon_percents, '--g', label='Moon')
    ax.plot(latitudes, moon_percents, 'ok')
    ax.set_ylim(0, 100)
    ax.set_xlabel('latitude')
    ax.set_ylabel('approx percent of sky')
    plt.legend()
    plt.show()

Answer 3

Yes, you're right. The current inclination of the earth's axis of rotation to the plane of orbit around sun is 23.43708° degrees. The sun can effectively cover the area from latitudes +23.43708° to -23.43708° on a sphere. Ignoring the night side of the celestial sphere gives the required area as Sin(23.43708°) = 39.774% ~40%.