Fix docs

docs/src/find_equinox_perihelion_dates.jl

@@ -8,13 +8,28 @@ FT = Float64
 include(joinpath(pkgdir(Insolation), "parameters", "create_parameters.jl"))
 param_set = create_insolation_parameters(FT)
 
+epoch_string = "2000-01-01T11:58:56.816"
+date0 = DateTime(epoch_string, dateformat"y-m-dTHH:MM:SS.s")
+
 # Difference in NH and SH zenith angles at time x in given year
 function zdiff(x, year, od)
     date = xtomarchdate(x, year)
-    theta_s, dist =
-        daily_zenith_angle(date, od, -45.0, param_set, milankovitch = true)
-    theta_n, dist =
-        daily_zenith_angle(date, od, 45.0, param_set, milankovitch = true)
+    theta_s, dist = daily_zenith_angle(
+        date,
+        date0,
+        od,
+        -45.0,
+        param_set,
+        milankovitch = true,
+    )
+    theta_n, dist = daily_zenith_angle(
+        date,
+        date0,
+        od,
+        45.0,
+        param_set,
+        milankovitch = true,
+    )
     return theta_n - theta_s
 end
 
@@ -28,7 +43,8 @@ end
 # Earth-Sun distance
 function edist(x, year, od)
     date = xtojandate(x, year)
-    _, dist = daily_zenith_angle(date, od, 0.0, param_set, milankovitch = true)
+    _, dist =
+        daily_zenith_angle(date, date0, od, 0.0, param_set, milankovitch = true)
     return dist / IP.orbit_semimaj(param_set)
 end
 

docs/src/plot_diurnal_cycle.jl

@@ -15,11 +15,14 @@ function diurnal_cycle(lat, lon, date, od, timezone, filename)
     hours = collect(range(0, stop = 24, length = nhours))
     insol = zeros(nhours)
     sza = zeros(nhours)
+    epoch_string = "2000-01-01T11:58:56.816"
+    date0 = DateTime(epoch_string, dateformat"y-m-dTHH:MM:SS.s")
+
     for (i, hr) in enumerate(hours)
         h = Int(round(hr + timezone))
         m = Int(round((hr + timezone - h) * 60))
         datetime = date + Dates.Hour(h) + Dates.Minute(m)
-        S, mu = solar_flux_and_cos_sza(datetime, od, lon, lat, param_set)
+        S, mu = solar_flux_and_cos_sza(datetime, date0, od, lon, lat, param_set)
         insol[i] = S * mu
         sza[i] = rad2deg(acos(mu))
     end

docs/src/plot_insolation.jl

@@ -22,12 +22,15 @@ function calc_day_lat_insolation(
     d_arr = Array{I}(round.(collect(range(0, stop = 365, length = n_days))))
     l_arr = collect(range(-90, stop = 90, length = n_lats))
     F_arr = zeros(n_days, n_lats)
+    epoch_string = "2000-01-01T11:58:56.816"
+    date0 = DateTime(epoch_string, dateformat"y-m-dTHH:MM:SS.s")
     # loop over days
     for (i, d) in enumerate(d_arr)
         for (j, lat) in enumerate(l_arr)
             date = Dates.DateTime(2000, 1, 1) + Dates.Day(d)
             θ, dist = daily_zenith_angle(
                 date,
+                date0,
                 od,
                 lat,
                 param_set,