Example: Venus Azimutal Distributions

This notebook computes the distribution of Venus rising and setting azimuths at dawn and dusk over millennia, as observed from Thebes. It distinguishes morning-star and evening-star apparitions and summarizes their preferred horizon directions.

If you are running this script in Google Colab you need first to install the package and create a directory for saving plots:

try:
    from google.colab import drive
    %pip install -Uq montu
except ImportError:
    print("Not running in Colab, skipping installation")
    import plotly.io as pio
    pio.renderers.default = "notebook_connected"
    %load_ext autoreload
    %autoreload 2
# Create folders for figures and temporal files
!mkdir -p ./gallery/ ./montu_dem/

Not running in Colab, skipping installation

Loading tools

We need to load the packages and the required data for it (star database, planet position database, etc.):

# Montu packages and data
%matplotlib inline
import montu
from montu import D2S, PRINTDF, TABLEDF, DEG, RAD

# Other libraries required
import pandas as pd
import copy
import matplotlib.pyplot as plt
import numpy as np

# Load allstars
allstars = montu.Stars()

MontuPython version 0.50.0. 𓇍𓇋𓇋𓏏𓅓𓊵 𓎛𓎡𓄿𓀭𓎛𓈖𓂝𓎡 (ii-ti m Htp, HkAx Hn'-k)
Loading stellar catalogue montu_stellar_catalogue_v38.csv

Venus appearances

Lets set observing conditions and initial time of exploration:

Tebas = montu.Observer(lon=33,lat=24,height=0)

Our target body:

venus = montu.Planet('Venus')
sun = montu.Sun()

Now let’s compute the conditions of observation of Venus for a given date:

mtime = montu.Time('-300-12-04 12:00:00')
#mtime = montu.Time()
print(f"Conditions for Venus observations at {mtime.readable.datemix}")
print('='*80)

# Compute Venus condition
venus.conditions_in_sky(at=mtime,observer=Tebas)

# Compute Sun conditions
sun.conditions_in_sky(at=mtime,observer=Tebas)
print("Sun rise time:",Tebas.get_local_time(sun.condition.rise_time))
print("Sun set time:",Tebas.get_local_time(sun.condition.set_time))
print()

dusk_time,dawn_time = montu.Sun.when_is_twilight(day=mtime,observer=Tebas,sunbelow=-18)
print("Time when dawn twilight ends: ",Tebas.get_local_time(dawn_time))
print("Time when dusk twilight starts: ",Tebas.get_local_time(dusk_time))
print()

print("Venus rise time:",Tebas.get_local_time(venus.condition.rise_time))
print("Venus set time:",Tebas.get_local_time(venus.condition.set_time))
print()

# Determine if it will be visible at dawn or dusk
mtime_obs = None
print("Visibility conditions:")
if venus.condition.set_time > sun.condition.set_time:
    print("\tVisible at dawn")

    # Determine if it will be visible when twilight has ended
    if venus.condition.set_time > dawn_time:
        print("\t\tVisible after dawn twilight ends")
        mtime_obs = montu.Time(dawn_time,format='jd')
        state = 'dawn'

elif venus.condition.rise_time < sun.condition.rise_time:
    print("\tVisible at dusk")

    # Determine if it will be visible when twilight has ended
    if venus.condition.rise_time < dusk_time:
        print("\t\tVisible before dusk twilight start")
        mtime_obs = montu.Time(dusk_time,format='jd')
        state = 'dusk'
else:
    print("Venus is transiting the Sun")

print()
if mtime_obs:
    venus.conditions_in_sky(at=mtime_obs,observer=Tebas)
    print(f"Venus position at {state}:")
    print(f"\tAz,el: {D2S(venus.position.az)}, {D2S(venus.position.el)}")
    print(f"\tMagnitude: {venus.condition.Vmag}")

Conditions for Venus observations at -300-12-08 12:00:00
================================================================================
Sun rise time: 06:35:06.991
Sun set time: 17:17:29.175

Time when dawn twilight ends:  18:37:41.647
Time when dusk twilight starts:  05:14:50.821

Venus rise time: 03:08:17.013
Venus set time: 14:43:12.573

Visibility conditions:
        Visible at dusk
                Visible before dusk twilight start

Venus position at dusk:
        Az,el: 114:08:45.447, 27:08:28.325
        Magnitude: -4.13

Let’s do it for many dates during a whole Synodic period:

mtime_initial = montu.Time('-3000-01-01 12:00:00')

venus_conditions = []
for dt in montu.PROGRESS(
    montu.Util.arange(0,4*montu.ALLPLANETS.loc['Venus','SynodicOrbit']*montu.YEAR,1*montu.DAY)
    ):

    mtime = mtime_initial + dt

    # Compute Sun conditions
    sun.conditions_in_sky(at=mtime,observer=Tebas)
    dusk_time,dawn_time = montu.Sun.when_is_twilight(day=mtime,observer=Tebas,sunbelow=-18)

    # Compute Venus condition
    venus.conditions_in_sky(at=mtime,observer=Tebas)

    # Conditions
    visible_at_dawn = False
    visible_at_dusk = False
    visible_before_dusk = False
    visible_after_dawn = False

    # Determine if it will be visible at dawn or dusk
    mtime_obs = None
    az_horizon = 361
    if venus.condition.set_time > sun.condition.set_time:
        visible_at_dawn = True
        az_horizon = venus.condition.set_az
        # Determine if it will be visible when twilight has ended
        if venus.condition.set_time > dawn_time:
            visible_after_dawn = True
            mtime_obs = montu.Time(dawn_time,format='jd')
            state = 'dawn'

    elif venus.condition.rise_time < sun.condition.rise_time:
        visible_at_dusk = True
        az_horizon = venus.condition.rise_az
        # Determine if it will be visible when twilight has ended
        if venus.condition.rise_time < dusk_time:
            visible_before_dusk = True
            mtime_obs = montu.Time(dusk_time,format='jd')
            state = 'dusk'
    else:
        pass

    az_twilight = None
    el_twilight = None
    if mtime_obs:
        venus.conditions_in_sky(at=mtime_obs,observer=Tebas)
        az_twilight = venus.position.az
        el_twilight = venus.position.el

    venus_conditions += [dict(
        jed = mtime.jed,
        sun_rise_time = sun.condition.rise_time,
        sun_set_time = sun.condition.set_time,
        dusk_time = dusk_time,
        dawn_time = dawn_time,
        venus_rise_time = venus.condition.rise_time,
        venus_set_time = venus.condition.set_time,
        visible_at_down = visible_at_dawn,
        visible_after_dawn = visible_after_dawn,
        visible_at_dusk = visible_at_dusk,
        visible_before_dusk = visible_before_dusk,
        az_horizon = az_horizon,
        az_twilight = az_twilight,
        el_twilight = el_twilight,
    )]
mtime_final = mtime.get_readable()
venus_conditions = pd.DataFrame(venus_conditions)

# Save data
suffix = f'{abs(mtime_initial.readable.year)}'
venus_conditions.to_csv(f'gallery/venus-conditions-year_{suffix}_bce.csv')

100%|██████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████████| 2336/2336 [00:00<00:00, 2935.97it/s]

Exclude problematic conditions:

cond = (venus_conditions.az_twilight == 361)|(venus_conditions.az_horizon == 361)
venus_conditions = venus_conditions[~cond]

Plot conditions:

with plt.style.context('dark_background'):
    fig,ax_el = plt.subplots(figsize=(10,5))
    ax_az = ax_el.twinx()

    ax_el.plot(venus_conditions.jed,venus_conditions.el_twilight,'c')

    ax_az.plot(venus_conditions.jed,venus_conditions.az_twilight,'y',label='Azimuth at twilight')
    ax_az.plot(venus_conditions.jed[::10],venus_conditions.az_horizon[::10],'y+',label='Azimuth at set')

    # Decorate
    montu.Time.set_time_ticks(ax_el,format='jd',timefmt='%Y-%m',fontsize=10)
    ax_el.set_xlabel('Date')
    ax_el.set_ylabel('Elevation at twilight [deg]',color='c')
    ax_az.set_ylabel('Azimuth [deg]',color='y')
    ax_az.legend(loc='center')

    ax_el.set_title(f"Venus appearances conditions between {mtime_initial.strftime('%Y-%m-%d')} and {mtime_final.strftime('%Y-%m-%d')}")
    ax_el.grid(alpha=0.3)
    montu.Util.montu_mark(ax_el);
    fig.savefig(f'gallery/venus-appearances-year_{suffix}_bce.png')

Histogram of azimuths:

fig,ax = plt.subplots()

hs,xs,bs = ax.hist(
    venus_conditions.az_twilight,
    bins=50,density=True,
    color='r',alpha=0.3,
    label = 'Azimuth at twilight'
    )

hs,xs,bs = ax.hist(
    venus_conditions.az_horizon,
    bins=50,density=True,
    color='b',alpha=0.2,
    label = 'Azimuth at set'
    )

hmax = hs.max()
# vertical bars
ax.axvline(90,color='k')
ax.text(90,1.08*hmax,'E',ha='center')
ax.axvline(180,color='k',alpha=0.3)
ax.text(180,1.08*hmax,'S',ha='center')
ax.axvline(270,color='k')
ax.text(270,1.08*hmax,'W',ha='center')

ax.set_xticks(np.arange(0,360,15))
ax.tick_params(axis='x',which='major',labelsize=8,rotation=90)
ax.set_yticks([])

ax.legend()
ax.grid(alpha=0.3)
ax.set_xlabel('Azimuth [deg]')
ax.set_ylabel('Frequency')
ax.set_title(f'Distribution of azimuths for Venus for {suffix}',x=0.5,y=1.05)

montu.Util.montu_mark(ax);
fig.savefig(f'gallery/venus-appearances-distribution-year_{suffix}_bce.png')

The end!

---

Powered by MontuPython. For more examples see MontuPython GitHub repo.

Jorge I. Zuluaga © 2023-present