File:Insolation of day 180 at 15 45 65n in earth 140000 0 bp 1.svg

1. 2. Plot ice age orbital parameters
   3. with Python climlab
   5. 19.10.2021 0000.0000

import numpy as np import matplotlib.pyplot as plt

from scipy.interpolate import make_interp_spline, BSpline

from climlab.solar.orbital import OrbitalTable from climlab.solar.insolation import daily_insolation from climlab import constants as const

yearsmax1=-140 numxlabels=8 smoothsteps=200

yearsmax2=(-1*yearsmax1)+1

years = np.linspace(yearsmax1, 0,yearsmax2) yearlabels=np.linspace(yearsmax1,0)

years2 = np.linspace(years.min(), years.max(), smoothsteps) yearlabels2=np.linspace(yearsmax1,0,smoothsteps)

xlabels=np.linspace(yearsmax1,0,numxlabels) xlabels2=xlabels*-1 xlabels2[numxlabels-1]=0

orb = OrbitalTable.interp(kyear=years)

print(orb)

kyear = orb['kyear'] ecc = orb['ecc'] long_peri = orb['long_peri']%360 obliquity = orb['obliquity'] precession = orb['precession'] insolation_65n=orb['65NJul']

1. print (obliquity)

S15_180 = daily_insolation(lat=15, day=180, orb=orb) S25_180 = daily_insolation(lat=25, day=180, orb=orb) S35_180 = daily_insolation(lat=35, day=180, orb=orb) S45_180 = daily_insolation(lat=45, day=180, orb=orb) S55_180 = daily_insolation(lat=55, day=180, orb=orb) S65_180 = daily_insolation(lat=65, day=180, orb=orb) S75_180 = daily_insolation(lat=75, day=180, orb=orb)

1. define spline

spl = make_interp_spline(years, S15_180, k=3) S15_180_smoothed = spl(years2) spl = make_interp_spline(years, S25_180, k=3) S25_180_smoothed = spl(years2) spl = make_interp_spline(years, S35_180, k=3) S35_180_smoothed = spl(years2) spl = make_interp_spline(years, S45_180, k=3) S45_180_smoothed = spl(years2) spl = make_interp_spline(years, S55_180, k=3) S55_180_smoothed = spl(years2) spl = make_interp_spline(years, S65_180, k=3) S65_180_smoothed = spl(years2) spl = make_interp_spline(years, S75_180, k=3) S75_180_smoothed = spl(years2)

1. 1. now plotting ...

fig, ax1 = plt.subplots()

1. ax2 = ax1.twinx()

plt.suptitle("Insolation of day 180, W/m²",size=18 )

ax1.set_xlabel("Thousands of years ago ",size=20 ) ax1.set_ylabel("Insolation of day 180, W/m²",color="#7f0000", size=20)

ax1.plot(years2, S15_180_smoothed, color="#ff0000", linestyle=":",label="15° N") ax1.plot(years2, S45_180_smoothed, color="#ff0000", linestyle="-.",label="45° N") ax1.plot(years2, S65_180_smoothed, color="#ffaf00", linewidth=6,label="65° N")

ax1.xaxis.set_tick_params(labelsize=16) ax1.yaxis.set_tick_params(labelsize=16)

ax1.set_xticks(xlabels) ax1.set_xticklabels(xlabels2)

ax1.legend(fontsize=16, loc=2)

plt.show()