libTheSky
Routines to compute sky positions of Sun, Moon, planets and more
Loading...
Searching...
No Matches
Functions/Subroutines
thesky_datetime Module Reference

Date and time procedures. More...

Functions/Subroutines

subroutine set_date_and_time (year, month, day, hour, minute, second)
 Set global date/time variables (year, month, ..., minute, second in TheSky_local) to specified values.
 
subroutine get_date_and_time (year, month, day, hour, minute, second)
 Retrieve the current global date/time variables (year, month, ..., minute, second, stored in TheSky_local)
 
subroutine system_clock_2_ymdhms (year, month, day, hour, minute, second, tz)
 Return system-clock date and time in (year, month, ..., minute, second and tz)
 
subroutine set_date_and_time_to_system_clock ()
 Set global date/time variables (year, month, ..., minute, second in TheSky_local) to system clock.
 
subroutine set_date_and_time_to_jd2000 ()
 Set global date/time variables (year, month, ..., minute, second in TheSky_local) to JD2000.0.
 
subroutine calctime (ut, jd, jde)
 Compute UT, JD, JDE, DeltaT and TZ, using the date and (local) time and TZ stored in the module TheSky_local.
 
real(double) function calc_gmst (jd, deltat)
 Calculate Greenwich Mean Sidereal Time for any instant, in radians.
 
real(double) function gmst_meeus (jd)
 Calculate Greenwich Mean Sidereal Time for any instant, in radians, using Meeus.
 
subroutine localtime2jd (jd)
 Computes JD, DeltaT and TZ, from date/time variables in module TheSky_local.
 
real(double) function calc_deltat (jd, force_recompute)
 Compute DeltaT for a given JD.
 
real(double) function calc_deltat_ymd (y, m, d, force_recompute)
 Compute DeltaT for given y,m,d.
 
real(double) function find_deltat_in_range (y, y0)
 Find a precise value for DeltaT through linear interpolation of the two adjacent tabulated values.
 
real(double) function calc_deltat_approx (jd)
 Compute DeltaT for given JD, using a simple parabolic approximation.
 
subroutine jd2dtm (jd, yy, mm, d, h, m, s)
 Convert a Julian day (UT) to LOCAL date and time (h,m,s)
 
subroutine jd2dthm (jd, yy, mm, d, h, m)
 Convert a Julian day (UT) to LOCAL date and time (h,m - no seconds)
 
real(double) function jd2ltime (jd0)
 Convert a Julian day (UT) to a local time (LT, h)
 
subroutine printdate (jd, jde)
 Prints date/time of a given Julian day (UT) to standard output.
 
subroutine printdate1 (jd, jde)
 Prints date/time of a given Julian day (UT) to standard output, but without a newline.
 
subroutine dls (yr, jdb, jde)
 Find the two Julian days of the beginning and the end of daylight-savings time in the EU for a given year.
 
real(double) function gettz (jd, ltz0, ldsttp)
 Returns time zone: tz0 or tz0+1.
 
subroutine print_date_time_and_location (op, nlbef, nlaf, ut, jd, jde, locname, tzname)
 Display a banner with date, time and location of calculation. Computes and returns UT, JD and JDE.
 
integer function dow (jd0)
 Calculates day of week (0 - Sunday, ... 6)
 
integer function woy (jd)
 Calculate the week-of-year number.
 
subroutine easter_gauss (year, month, day)
 Calculate the date of Easter using Gauss' method.
 
subroutine passover_gauss (year, month, day)
 Calculate the date of Passover or Pesach (Jewish Easter) using Gauss' method.
 

Detailed Description

Date and time procedures.

Function/Subroutine Documentation

◆ calc_deltat()

real(double) function thesky_datetime::calc_deltat ( real(double), intent(in)  jd,
logical, intent(in), optional  force_recompute 
)

Compute DeltaT for a given JD.

Parameters
jdJulian day
force_recomputeForce recomputation of DeltaT, even if the year is the same as in the last call (done in calc_deltat_ymd())
Note
VERY SLOW, use calc_deltat_ymd() if y,m,d are known

Definition at line 340 of file date_time.f90.

References calc_deltat(), calc_deltat_ymd(), thesky_local::day, thesky_constants::deltat_forced, thesky_local::month, and thesky_local::year.

Referenced by calc_deltat(), thesky_planets::jupiterphys(), thesky_moon::moonpos_la(), thesky_planets::planet_position(), thesky_planets::planetelements(), and thesky_sun::sunpos_la().

Here is the call graph for this function:

◆ calc_deltat_approx()

real(double) function thesky_datetime::calc_deltat_approx ( real(double), intent(in)  jd)

Compute DeltaT for given JD, using a simple parabolic approximation.

Parameters
jdJulian day for computation
Return values
deltatThe value for DeltaT (= TDT-UT) in seconds.

Definition at line 505 of file date_time.f90.

References calc_deltat_approx(), and thesky_constants::jd1820.

Referenced by calc_deltat_approx().

Here is the call graph for this function:

◆ calc_deltat_ymd()

real(double) function thesky_datetime::calc_deltat_ymd ( integer, intent(in)  y,
integer, intent(in)  m,
real(double), intent(in)  d,
logical, intent(in), optional  force_recompute 
)

Compute DeltaT for given y,m,d.

Parameters
yYear
mMonth
dDay
force_recomputeForce recomputation of DeltaT, even if the year is the same as in the last call
Note
  • Faster than calc_deltat. Use this routine rather than calc_deltat() if y,m,d are known

Definition at line 394 of file date_time.f90.

References calc_deltat_ymd(), thesky_constants::deltat_0, thesky_constants::deltat_accel, thesky_constants::deltat_change, thesky_constants::deltat_forced, thesky_constants::deltat_maxyr, thesky_constants::deltat_minyr, thesky_constants::deltat_values, thesky_constants::deltat_years, and find_deltat_in_range().

Referenced by calc_deltat(), calc_deltat_ymd(), calctime(), and localtime2jd().

Here is the call graph for this function:

◆ calc_gmst()

real(double) function thesky_datetime::calc_gmst ( real(double), intent(in)  jd,
real(double), intent(in), optional  deltat 
)

Calculate Greenwich Mean Sidereal Time for any instant, in radians.

Parameters
jdJulian day of computation
DeltaTΔT (s) (optional - default: use DelTaT from TheSky_local)
Return values
calc_gmstGreenwich Mean Sidereal Time (radians)
See also
Explanatory Supplement to the Astronomical Almanac, 3rd edition, Eq. 6.66 (2012)

Definition at line 243 of file date_time.f90.

References calc_gmst(), and thesky_local::deltat.

Referenced by calc_gmst(), thesky_visibility::limmag_zenith_jd(), thesky_moon::moonpos_la(), thesky_planets::planet_position(), print_date_time_and_location(), and thesky_sun::sunpos_la().

Here is the call graph for this function:

◆ calctime()

subroutine thesky_datetime::calctime ( real(double), intent(out)  ut,
real(double), intent(out)  jd,
real(double), intent(out)  jde 
)

Compute UT, JD, JDE, DeltaT and TZ, using the date and (local) time and TZ stored in the module TheSky_local.

Parameters
utUniversal time (output)
jdJulian date (output)
jdeJulian date (output)
Note
  • Date and time are obtained from year, month, day, hour, minute, second, through the module TheSky_local, assuming LOCAL time
  • DeltaT and TZ are returned through the module TheSky_local
  • Computed JD is in UNIVERSAL time
Todo:
Remove recomputation of 'UT JD'!!!

Definition at line 203 of file date_time.f90.

References calc_deltat_ymd(), thesky_local::day, thesky_local::deltat, gettz(), thesky_local::hour, thesky_local::minute, thesky_local::month, thesky_local::second, thesky_local::tz, and thesky_local::year.

Referenced by print_date_time_and_location().

Here is the call graph for this function:

◆ dls()

subroutine thesky_datetime::dls ( integer, intent(in)  yr,
real(double), intent(out)  jdb,
real(double), intent(out)  jde 
)

Find the two Julian days of the beginning and the end of daylight-savings time in the EU for a given year.

Parameters
yrYear (CE)
jdbJulian day of beginning of DST (output)
jdeJulian day of end of DST (output)

Definition at line 761 of file date_time.f90.

References dow().

Here is the call graph for this function:

◆ dow()

integer function thesky_datetime::dow ( real(double), intent(in)  jd0)

Calculates day of week (0 - Sunday, ... 6)

Parameters
jd0Julian day (UT)
Note
Input in UT, output in local time
Todo:
Switch using dow(jd) to dow_ut(jd+tz/24.d0) to make it general

Definition at line 984 of file date_time.f90.

References dow(), and thesky_local::tz.

Referenced by dls(), dow(), gettz(), print_date_time_and_location(), and woy().

Here is the call graph for this function:

◆ easter_gauss()

subroutine thesky_datetime::easter_gauss ( integer, intent(in)  year,
integer, intent(out)  month,
integer, intent(out)  day 
)

Calculate the date of Easter using Gauss' method.

Parameters
yearJulian/Gregorian year to compute date of Easter for.
monthJulian/Gregorian month of year of Easter Sunday.
dayJulian/Gregorian day of month of Easter Sunday.
Note

Definition at line 1039 of file date_time.f90.

◆ find_deltat_in_range()

real(double) function thesky_datetime::find_deltat_in_range ( integer, intent(in)  y,
real(double), intent(in)  y0 
)

Find a precise value for DeltaT through linear interpolation of the two adjacent tabulated values.

Parameters
yCurrent year
y0Current date, as fractional year

Definition at line 470 of file date_time.f90.

References thesky_constants::deltat_n, thesky_constants::deltat_values, thesky_constants::deltat_years, and find_deltat_in_range().

Referenced by calc_deltat_ymd(), and find_deltat_in_range().

Here is the call graph for this function:

◆ get_date_and_time()

subroutine thesky_datetime::get_date_and_time ( integer, intent(out)  year,
integer, intent(out)  month,
integer, intent(out)  day,
integer, intent(out)  hour,
integer, intent(out)  minute,
real(double), intent(out)  second 
)

Retrieve the current global date/time variables (year, month, ..., minute, second, stored in TheSky_local)

Parameters
yearYear (output)
monthMonth (output)
dayDay (output)
hourHour (output)
minuteMinute (output)
secondSecond (output)

Definition at line 96 of file date_time.f90.

References thesky_local::day, thesky_local::hour, thesky_local::minute, thesky_local::month, thesky_local::second, and thesky_local::year.

◆ gettz()

real(double) function thesky_datetime::gettz ( real(double), intent(in)  jd,
real(double), intent(in), optional  ltz0,
integer, intent(in), optional  ldsttp 
)

Returns time zone: tz0 or tz0+1.

Parameters
jdJulian day (UT?)
ltz0Default time zone for the current location ('winter time')
ldsttpDaylight-savings time rules to use: 0: no DST (e.g. UT), 1: EU, 2: USA/Canada
Note
  • currently implemented for no DST (dsttp=0), EU (dsttp=1) and USA/Canada >2007 (dsttp=2) only
  • tz0 and dsttp can be provided through the module TheSky_local, or using the optional arguments. Note that using the latter will update the former!

Definition at line 802 of file date_time.f90.

References dow(), thesky_local::dsttp, gettz(), thesky_local::tz, and thesky_local::tz0.

Referenced by thesky_riset::best_obs_date(), thesky_visibility::best_planet_xsmag(), calctime(), gettz(), localtime2jd(), and print_date_time_and_location().

Here is the call graph for this function:

◆ gmst_meeus()

real(double) function thesky_datetime::gmst_meeus ( real(double), intent(in)  jd)

Calculate Greenwich Mean Sidereal Time for any instant, in radians, using Meeus.

Parameters
jdJulian day of computation
Return values
gmst_meeusGreenwich Mean Sidereal Time in radians
See also
Meeus, Sect. 12: Siderial time at Greenwich; Eq.12.4

Definition at line 278 of file date_time.f90.

References gmst_meeus().

Referenced by gmst_meeus().

Here is the call graph for this function:

◆ jd2dthm()

subroutine thesky_datetime::jd2dthm ( real(double), intent(in)  jd,
integer, intent(out)  yy,
integer, intent(out)  mm,
integer, intent(out)  d,
integer, intent(out)  h,
integer, intent(out)  m 
)

Convert a Julian day (UT) to LOCAL date and time (h,m - no seconds)

Parameters
jdJulian day (UT)
yyYear (CE, LT) (output)
mmMonth (LT) (output)
dDay (LT) (output)
hHour (LT) (output)
mMinute (LT) (output)

Definition at line 601 of file date_time.f90.

References thesky_local::tz.

◆ jd2dtm()

subroutine thesky_datetime::jd2dtm ( real(double), intent(in)  jd,
integer, intent(out)  yy,
integer, intent(out)  mm,
integer, intent(out)  d,
integer, intent(out)  h,
integer, intent(out)  m,
real(double), intent(out)  s 
)

Convert a Julian day (UT) to LOCAL date and time (h,m,s)

Parameters
jdJulian day (UT)
yyYear (CE, LT) (output)
mmMonth (LT) (output)
dDay (LT) (output)
hHour (LT) (output)
mMinute (LT) (output)
sSecond (+ fraction, LT) (output)

Definition at line 532 of file date_time.f90.

References thesky_local::tz.

◆ jd2ltime()

real(double) function thesky_datetime::jd2ltime ( real(double), intent(in)  jd0)

Convert a Julian day (UT) to a local time (LT, h)

Parameters
jd0Julian day (UT)

Definition at line 658 of file date_time.f90.

References jd2ltime(), and thesky_local::tz.

Referenced by jd2ltime().

Here is the call graph for this function:

◆ localtime2jd()

subroutine thesky_datetime::localtime2jd ( real(double), intent(out)  jd)

Computes JD, DeltaT and TZ, from date/time variables in module TheSky_local.

Parameters
JDJulian date (output)
Note
  • Date and time are obtained from year, month, day, hour, minute, second, through the module TheSky_local
  • DeltaT and TZ are returned through the module TheSky_local

Definition at line 308 of file date_time.f90.

References calc_deltat_ymd(), thesky_local::day, thesky_local::deltat, gettz(), thesky_local::hour, thesky_local::minute, thesky_local::month, thesky_local::second, thesky_local::tz, and thesky_local::year.

Here is the call graph for this function:

◆ passover_gauss()

subroutine thesky_datetime::passover_gauss ( integer, intent(in)  year,
integer, intent(out)  month,
integer, intent(out)  day 
)

Calculate the date of Passover or Pesach (Jewish Easter) using Gauss' method.

Parameters
yearJulian/Gregorian year to compute date of Passover for.
monthJulian/Gregorian month of year of the DAY of Passover.
dayJulian/Gregorian day of month of the DAY of Passover.
Note
  • Based on Gauss, C.F., 1802. Berechnung des Judischen Osterfestes. Werke, pp.80-81, as described in Meeus, Ch.9.
  • Passover is always on 15 Nisan. Note that this is the DAY of Passover, the day AFTER the evening of the (first) Seder.
  • Valid for the Julian and Gregorian calendar.
  • The Jewish year is the Julian/Gregorian year + 3760.
  • For the years -4000 - +4000, this corresponds EXACTLY to Har’El, Z., 2005. Gauss Formula for the Julian Date of Passover. A∆(A), 1, p.0. Har’El also computes the day of week, which for the same years corresponds exactly to dow(cal2jd(year,month,day)) with cal2jd from libSUFR/date_and_time, dow() from libTheSky/datetimeSky, and year,month,day as in the current subroutine. Note that Passover day can only fall on Sun, Tue, Thu or Sat.

Definition at line 1101 of file date_time.f90.

◆ print_date_time_and_location()

subroutine thesky_datetime::print_date_time_and_location ( integer, intent(in), optional  op,
integer, intent(in), optional  nlbef,
integer, intent(in), optional  nlaf,
real(double), intent(out), optional  ut,
real(double), intent(out), optional  jd,
real(double), intent(out), optional  jde,
character, dimension(*), intent(in), optional  locname,
character, dimension(*), intent(in), optional  tzname 
)

Display a banner with date, time and location of calculation. Computes and returns UT, JD and JDE.

Parameters
opOutput unit (optional, default 6)
nlbefNumber of newlines before output (optional, default 0)
nlafNumber of newlines after output (optional, default 0)
utUT: Universal Time (optional)
jdJD: Julian day (optional)
jdeJDE: Apparent Julian day (optional)
locnameLocation name (optional)
tznameTime-zone name (optional)
Note
  • uses the module local to obtain the variables year, month, etc.
  • calls calctime(), gettz()

Definition at line 910 of file date_time.f90.

References calc_gmst(), calctime(), thesky_local::day, thesky_local::deltat, dow(), gettz(), thesky_local::height, thesky_local::hour, thesky_local::lat0, thesky_local::lon0, thesky_local::minute, thesky_local::month, thesky_local::second, thesky_local::tz, and thesky_local::year.

Here is the call graph for this function:

◆ printdate()

subroutine thesky_datetime::printdate ( real(double), intent(in)  jd,
real(double), intent(in), optional  jde 
)

Prints date/time of a given Julian day (UT) to standard output.

Parameters
jdJulian day (UT)
jdeJulian day ephemeris time (dynamical time)
Note
  • calls printdate1()

Definition at line 690 of file date_time.f90.

References printdate1().

Here is the call graph for this function:

◆ printdate1()

subroutine thesky_datetime::printdate1 ( real(double), intent(in)  jd,
real(double), intent(in), optional  jde 
)

Prints date/time of a given Julian day (UT) to standard output, but without a newline.

Parameters
jdJulian day (UT)
jdeJulian day ephemeris time (dynamical time)

Definition at line 712 of file date_time.f90.

Referenced by thesky_moonroutines::moon_phase_next(), and printdate().

◆ set_date_and_time()

subroutine thesky_datetime::set_date_and_time ( integer, intent(in)  year,
integer, intent(in)  month,
integer, intent(in)  day,
integer, intent(in)  hour,
integer, intent(in)  minute,
real(double), intent(in)  second 
)

Set global date/time variables (year, month, ..., minute, second in TheSky_local) to specified values.

Parameters
yearYear
monthMonth
dayDay
hourHour
minuteMinute
secondSecond

Definition at line 65 of file date_time.f90.

References thesky_local::day, thesky_local::hour, thesky_local::minute, thesky_local::month, thesky_local::second, and thesky_local::year.

Referenced by set_date_and_time_to_jd2000(), and set_date_and_time_to_system_clock().

◆ set_date_and_time_to_jd2000()

subroutine thesky_datetime::set_date_and_time_to_jd2000

Set global date/time variables (year, month, ..., minute, second in TheSky_local) to JD2000.0.

Definition at line 179 of file date_time.f90.

References set_date_and_time().

Here is the call graph for this function:

◆ set_date_and_time_to_system_clock()

subroutine thesky_datetime::set_date_and_time_to_system_clock

Set global date/time variables (year, month, ..., minute, second in TheSky_local) to system clock.

Definition at line 159 of file date_time.f90.

References set_date_and_time(), and thesky_local::tz.

Here is the call graph for this function:

◆ system_clock_2_ymdhms()

subroutine thesky_datetime::system_clock_2_ymdhms ( integer, intent(out), optional  year,
integer, intent(out), optional  month,
integer, intent(out), optional  day,
integer, intent(out), optional  hour,
integer, intent(out), optional  minute,
real(double), intent(out), optional  second,
real(double), intent(out), optional  tz 
)

Return system-clock date and time in (year, month, ..., minute, second and tz)

Parameters
yearYear (output)
monthMonth of year (output)
dayDay of month (output)
hourHour of day (output)
minuteMinute (output)
secondSecond (output)
tzTime zone w.r.t. Greenwich - >0 = east (output)

Definition at line 130 of file date_time.f90.

◆ woy()

integer function thesky_datetime::woy ( real(double), intent(in)  jd)

Calculate the week-of-year number.

Parameters
jdJulian day
Todo:
  • CHECK: int or floor?
  • Depends on dow(), which depends on module local -> switch to dow_ut()

Definition at line 1010 of file date_time.f90.

References dow(), and woy().

Referenced by woy().

Here is the call graph for this function: