Skip to content

Commit

Permalink
testing tpdates
Browse files Browse the repository at this point in the history
  • Loading branch information
xanthospap committed Oct 17, 2023
1 parent 565afed commit 76d6112
Show file tree
Hide file tree
Showing 3 changed files with 123 additions and 2 deletions.
16 changes: 16 additions & 0 deletions src/datetime_tops.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -376,6 +376,22 @@ constexpr typename S::underlying_type max_days_allowed() {
S::max_in_day;
}

/** Transform any second type S to fractional days
*
* @warning This function assumes that a day is made up of exactly 86400 sec
* and is thus not able to represent a fractional day when the day at hand is
* on a leap second insertion.
*/
#if __cplusplus >= 202002L
template <gconcepts::is_sec_dt S>
#else
template <typename S, typename = std::enable_if_t<S::is_of_sec_type>>
#endif
double fractional_days(S nsec) noexcept {
const double sec = static_cast<double>(nsec.__member_ref__());
return sec / static_cast<double>(S::max_in_day);
}

/** Explicit cast of any second type to another second type.
*
* Cast an instance of any second type (aka any instance for which
Expand Down
31 changes: 29 additions & 2 deletions src/tpdate.hpp
Original file line number Diff line number Diff line change
@@ -1,3 +1,10 @@
/** @file
* A utility class to hold datetime instances, in a continuous time-scale
* (e.g. TT, TAI, etc). In construst to datetime<S>, this is not a template
* class and uses a storage method of two floating point numbers (one for
* MJD and one for time of day) to represent datetime epochs.
*/

#ifndef __DSO_DATETIME_TWOPARTDATES_HPP__
#define __DSO_DATETIME_TWOPARTDATES_HPP__

Expand Down Expand Up @@ -27,6 +34,13 @@ struct TwoPartJulianDate {
/** Given an MJD, turn it to a JD and return it split in a convinient way
* The way the (returned) JD is split is determined by the template parameter
* \p S.
*
* @param[in] mjd The Modified Julian Day (integral part)
* @param[in] fday The time of day as fractional days
* @return An TwoPartJulianDate instance; the two parts of the instance if
* added, give the JD corresponding to the passed in MJD (i.e. \p mjd).
* The method in which the JD is plit, is driven by the template
* parameter S.
*/
template <JdSplitMethod S = JdSplitMethod::J2000>
TwoPartJulianDate jd_split(double mjd, double fday) noexcept {
Expand All @@ -41,6 +55,15 @@ TwoPartJulianDate jd_split(double mjd, double fday) noexcept {
}
} /* namespace core */

/** A datetime class to represent epochs in any continuous system.
*
* A TwoPartDate instance conviniently splits a datetime into two numeric
* values:
* - the Modified Julian Day (i.e. anumeric value which only has an integral
* part), and
* - the time of day, which is stored in fractional days
*
*/
class TwoPartDate {
private:
double _mjd; /** Mjd */
Expand All @@ -55,7 +78,8 @@ class TwoPartDate {
#endif
TwoPartDate(const datetime<T> &d) noexcept
: _mjd((double)d.imjd().as_underlying_type()),
_fday(d.sec().fractional_days()) {
//_fday(d.sec().fractional_days()) {
_fday(fractional_days<T>(d.sec())) {
this->normalize();
}

Expand Down Expand Up @@ -83,7 +107,7 @@ class TwoPartDate {
* day
*
* This is not an 'actual date' but rather a datetime interval, but can be
* represented a TwoPartDate instance. If the calling instance is prior to
* represented by a TwoPartDate instance. If the calling instance is prior to
* the operand (i.e. d1-d2 with d2>d1) the interval is signed as 'negative'.
* This means that the number of days can be negative, but the fractional
* day will always be positive
Expand Down Expand Up @@ -237,6 +261,9 @@ class TwoPartDate {
_fday = std::modf(_fday, &extra);
_mjd += extra;
}
#ifdef DEBUG
assert(_fday >= 0e0 && _fday < 1e0);
#endif
// all done
return;
}
Expand Down
78 changes: 78 additions & 0 deletions test/unit_tests/tpdates1.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,78 @@
#include "dtcalendar.hpp"
#include "tpdate.hpp"
#include <cassert>
#include <random>

using namespace dso;

constexpr const long num_tests = 1'000'000;
using nsec = dso::nanoseconds;
typedef nsec::underlying_type SecIntType;
constexpr const SecIntType TOSEC = nsec::sec_factor<SecIntType>();
constexpr const double SECINDAY = nsec::max_in_day;

int main() {
/* Generators for random numbers ... */
std::random_device rd;
std::mt19937 gen(rd());
std::uniform_int_distribution<> ydstr(1972, 2050); /* range for years */
std::uniform_int_distribution<> mdstr(1, 12); /* range for months */
std::uniform_int_distribution<> ddstr(1, 31); /* range for day of month */
std::uniform_int_distribution<SecIntType> nsdstr(
0, nsec::max_in_day); /* range for day of month */

int testnr = 0, ok;
datetime<nsec> d1, d2, d3;
while (testnr < num_tests) {
/* do we have a valid date ? */
try {
d1 = datetime<nsec>{year(ydstr(gen)), month(mdstr(gen)),
day_of_month(ddstr(gen)), nsec(nsdstr(gen))};
d2 = datetime<nsec>{year(ydstr(gen)), month(mdstr(gen)),
day_of_month(ddstr(gen)), nsec(nsdstr(gen))};
/* d3 on same day as d2 */
d3 = datetime<nsec>{d2.imjd(), nsec(nsdstr(gen))};
ok = 1;
} catch (std::exception &) {
ok = 0;
}
if (ok) {
/* construct a TwoPartDate from a datetime */
TwoPartDate tpd1(d1);
assert(tpd1.imjd() == d1.imjd().as_underlying_type());
assert(tpd1.fday() == d1.sec().to_fractional_seconds()/SECINDAY);
/* constructor from MJD and fraction of day */
tpd1 = TwoPartDate(d1.imjd().as_underlying_type(),
d1.sec().to_fractional_seconds() / SECINDAY);
assert(tpd1.imjd() == d1.imjd().as_underlying_type());
assert(tpd1.fday() == d1.sec().to_fractional_seconds() / SECINDAY);

TwoPartDate tpd2(d2);
TwoPartDate tpd3(d3);

/* how its done using datetime_intervals */
const auto interval = d2 - d1;
auto d = d1 + interval;
assert(d == d2);
/* same thing using TwoPartDate's */
const auto i1 = tpd2 - tpd1;
auto tpd = tpd1 + i1;
assert(tpd == tpd2);

/* using datetime_interval ... */
const auto d32 = d3 - d2;
d = d2 + d32;
assert(d == d3);
/* same thing using TwoPartDate's */
const auto tpd32 = tpd3 - tpd2;
tpd = tpd2 + tpd32;
assert(tpd == tpd3);

++testnr;
if (testnr % 10)
printf("%d/%ld\r", testnr, num_tests);
}
}

return 0;
}

0 comments on commit 76d6112

Please sign in to comment.