DailyInsolation | R Documentation |
Computes the daily average insolation, converted and adapted from the Eisenmann and Huybers,2006 Matlab Code.
This script contains orbital parameter data for the past 50000 years from Berger and Loutre (1991).
Detailed description of calculation: Values for eccentricity, obliquity, and longitude of perihelion for the past 5 Myr are taken from Berger and Loutre 1991 (data from ncdc.noaa.gov). If using calendar days, solar longitude is found using an approximate solution to the differential equation representing conservation of angular momentum (Kepler's Second Law) (day_type=1), or an iterative solution of Keplers Equation (day_type=3).
Given the orbital parameters and solar longitude, daily average insolation is calculated exactly following Berger 1978.
Orbital parameters from Berger and Loutre 1991.
Calculation of daily average insolation following Berger A. (1978).
Authors of original MATLAB version : Ian Eisenman and Peter Huybers, Harvard University, August 2006 eisenman@fas.harvard.edu This file is available online at http://deas.harvard.edu/~eisenman/downloads
Translation into R and extension to arbitray calendars; Thomas Laepple, 2008
DailyInsolation( kyear, lat, day, day_type = 1, fast = TRUE, T.alpha = 80, alpha = 0 )
kyear |
Thousands of years before present (0 to 5000). |
lat |
Latitude in degrees (-90 to 90). |
day |
Indicator of time of year; calendar day by default. |
day_type |
Convention for specifying time of year (1,2,3) [optional], see details |
fast |
TRUE = use precalculated values of orbital parameters in 100yr steps |
T.alpha |
for custom calendar (day_type=3) |
alpha |
for custom calendar (day_type=3) |
list(Fsw=Fsw,ecc=ecc,obliquity=obliquity,long_perh=long_perh,lambda=lambda/2/pi*360) Fsw = Daily average solar radiation in W/m^2. ... and all orbital parameters
day_type=1 (default): day input is calendar day (1-365.24), where day 1 is January first. The calendar is referenced to the vernal equinox which always occurs at day 80.
day_type=2: day input is solar longitude (0-360 degrees). Solar longitude is the angle of the Earth's orbit measured from spring equinox (21 March). Note that calendar days and solar longitude are not linearly related because, by Kepler's Second Law, Earth's angular velocity varies according to its distance from the sun.
day_type=3: Custom calendar day input is calendar day (1-365.24), where day 1 is January first. The calendar is referenced to the solar longitude alpha (0..2pi, angle of the earth orbit from spring equinox) which always occurs at day T.alpha. This uses an iterative solution of Keplers equation and is therefore much slower than the day_type=1 which uses a trigeometric approximation of Kepler's equation.
Thomas Laepple
Berger A. and Loutre M.F. (1991). Insolation values for the climate of the last 10 million years. Quaternary Science Reviews, 10(4), 297-317.
Berger A. (1978). Long-term variations of daily insolation and Quaternary climatic changes. Journal of Atmospheric Science, 35(12), 2362-2367.
Grossman, N.: Orbits under the Inverse Square Law, in The Sheer Joy of Celestial Mechanics, edited by N. Grossman, pp. 41–61, Birkhäuser Boston, Boston, MA., 1996.
par(mfcol = c(2, 1)) plot(0:1000, DailyInsolation(0:1000, 40, 180)$Fsw, main = "day 180 (= peak summer) insolation at 40degN", xlab = "kyr BP", ylab = "W/m2", type = "l") plot(1:365, DailyInsolation(0, 40, 1:365)$Fsw, main = "modern daily insolation at 40degN", xlab = "day of the year", ylab = "W/m2", type = "l")
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.