This commit is contained in:
Timur A. Fatkhullin 2025-03-22 00:48:52 +03:00
parent b9ee662850
commit b3fb445557
5 changed files with 71 additions and 40 deletions

View File

@ -79,7 +79,7 @@ endif()
# ERFA library
include(ExternalProject)
ExternalProject_Add(erfa_lib
ExternalProject_Add(erfalib
PREFIX ${CMAKE_BINARY_DIR}/erfa_lib
GIT_REPOSITORY "https://github.com/liberfa/erfa.git"
GIT_TAG "v2.0.1"
@ -89,13 +89,17 @@ ExternalProject_Add(erfa_lib
# SOURCE_DIR erfa
# INSTALL_DIR
LOG_CONFIGURE 1
CONFIGURE_COMMAND meson setup --reconfigure -Ddefault_library=static -Dbuildtype=release -Dc_args='-march=native' -Doptimization=3
CONFIGURE_COMMAND meson setup --reconfigure -Ddefault_library=static -Dbuildtype=release
-Dprefix=${CMAKE_BINARY_DIR}/erfa_lib -Dlibdir= -Dincludedir= -Ddatadir= <SOURCE_DIR>
# CONFIGURE_COMMAND meson setup --reconfigure -Ddefault_library=static -Dbuildtype=release -Dc_args='-march=native' -Doptimization=3
# -Dprefix=${CMAKE_BINARY_DIR}/erfa_lib -Dlibdir= -Dincludedir= -Ddatadir= <SOURCE_DIR>
BUILD_COMMAND ninja -C <BINARY_DIR>
INSTALL_COMMAND meson install -C <BINARY_DIR>
BUILD_BYPRODUCTS ${CMAKE_BINARY_DIR}/erfa_lib/liberfa.a
)
add_library(ERFA_LIB STATIC IMPORTED)
set_target_properties(ERFA_LIB PROPERTIES IMPORTED_LOCATION ${CMAKE_BINARY_DIR}/erfa_lib/liberfa.a)
add_dependencies(ERFA_LIB erfalib)
set(ERFA_INCLUDE_DIR ${CMAKE_BINARY_DIR}/erfa_lib)
include_directories(${ERFA_INCLUDE_DIR})

View File

@ -39,8 +39,8 @@ concept mcc_formattable =
requires(T v, std::format_context ctx) { std::formatter<std::remove_cvref_t<T>>().format(v, ctx); };
template <typename T>
// from https://stackoverflow.com/questions/74383254/concept-that-models-only-the-stdchrono-duration-types
template <typename T>
concept mcc_time_duration_c = requires {
[]<class Rep, class Period>(std::type_identity<std::chrono::duration<Rep, Period>>) {
@ -48,6 +48,12 @@ concept mcc_time_duration_c = requires {
};
template <typename T>
concept mcc_systime_c = requires {
[]<class DT>(std::type_identity<std::chrono::sys_time<DT>>) {}(std::type_identity<std::remove_cvref_t<T>>());
};
/* a callable concept and its signature traits */
template <typename T>

View File

@ -51,10 +51,10 @@ static constexpr double mcc_UT1_to_sideral_ratio = 1.002737909350795; // UT1/si
// modified Julian date (based on ERFA eraCal2jd)
template <traits::mcc_real_scalar_or_real_range_c ResT>
static int mcc_julday(const std::chrono::system_clock::time_point& start_time,
template <traits::mcc_real_scalar_or_real_range_c ResT, traits::mcc_time_duration_c DT = std::chrono::milliseconds>
static int mcc_julday(traits::mcc_systime_c auto const& start_time,
ResT& mjd,
const std::chrono::system_clock::duration& step = std::chrono::milliseconds(100))
const DT& step = std::chrono::milliseconds(100))
{
size_t mjd_size = 0;
if constexpr (std::ranges::range<ResT>) {
@ -124,7 +124,7 @@ static int mcc_julday(const std::chrono::system_clock::time_point& start_time,
if constexpr (std::ranges::random_access_range<ResT>) {
ptr += reg_size;
} else {
for (int k = 0; k < reg_size; ++k) {
for (size_t k = 0; k < reg_size; ++k) {
++ptr;
}
}
@ -156,7 +156,7 @@ double mcc_time_to_alt_limit(traits::mcc_real_or_char_range auto const& alt_limi
traits::mcc_real_or_char_range auto const& DEC,
traits::mcc_real_or_char_range auto const& LAT,
traits::mcc_real_or_char_range auto const& LON,
const std::chrono::system_clock::time_point& now,
traits::mcc_systime_c auto const& now,
traits::mcc_time_duration_c auto const& dut1, // UT1-UTC
traits::mcc_time_duration_c auto const& tt_tai, // TT-TAI
// TAI-UTC (leap seconds)
@ -166,51 +166,58 @@ double mcc_time_to_alt_limit(traits::mcc_real_or_char_range auto const& alt_limi
// HA = LST - RA
// cos(HA) = cos(LST)*cos(RA) + sin(LST)*sin(RA)
using AT = std::decay_t<decltype(alt_limit)>;
using RT = std::decay_t<decltype(RA)>;
using DT = std::decay_t<decltype(DEC)>;
using LT = std::decay_t<decltype(LAT)>;
using LGT = std::decay_t<decltype(LON)>;
// using AT = std::decay_t<decltype(alt_limit)>;
// using RT = std::decay_t<decltype(RA)>;
// using DT = std::decay_t<decltype(DEC)>;
// using LT = std::decay_t<decltype(LAT)>;
// using LGT = std::decay_t<decltype(LON)>;
double ra, dec, lat, lon, alt;
if constexpr (std::floating_point<AT>) {
alt = alt_limit * utils::deg2radCoeff;
} else {
alt = utils::parsAngleString(alt_limit);
alt *= utils::deg2radCoeff;
auto to_rads = [](const auto& v, bool hms = false) {
// using v_t = std::remove_cvref<decltype(v)>;
using v_t = std::remove_cvref_t<decltype(v)>;
double res;
if constexpr (!std::floating_point<v_t>) {
res = utils::parsAngleString(v, hms).value_or(std::numeric_limits<double>::quiet_NaN());
} else {
res = v;
}
if (!std::isfinite(res)) {
return res;
}
return res * utils::deg2radCoeff;
};
alt = to_rads(alt_limit);
if (!std::isfinite(alt)) {
return alt;
}
if (alt < 0.0) {
return std::numeric_limits<double>::quiet_NaN();
}
if constexpr (std::floating_point<RT>) {
ra = RA * utils::deg2radCoeff;
} else {
ra = utils::parsAngleString(RA, true);
ra *= utils::deg2radCoeff;
ra = to_rads(RA, true);
if (!std::isfinite(ra)) {
return ra;
}
if constexpr (std::floating_point<DT>) {
dec = DEC * utils::deg2radCoeff;
} else {
dec = utils::parsAngleString(DEC);
dec *= utils::deg2radCoeff;
dec = to_rads(DEC);
if (!std::isfinite(dec)) {
return dec;
}
if constexpr (std::floating_point<LT>) {
lat = LAT * utils::deg2radCoeff;
} else {
lat = utils::parsAngleString(LAT);
lat *= utils::deg2radCoeff;
lat = to_rads(LAT);
if (!std::isfinite(lat)) {
return lat;
}
if constexpr (std::floating_point<LGT>) {
lon = LON * utils::deg2radCoeff;
} else {
lon = utils::parsAngleString(LON);
lon *= utils::deg2radCoeff;
lon = to_rads(LON);
if (!std::isfinite(lon)) {
return lon;
}
if (lat >= 0.0) { // north hemisphere
@ -243,13 +250,18 @@ double mcc_time_to_alt_limit(traits::mcc_real_or_char_range auto const& alt_limi
return std::numeric_limits<double>::quiet_NaN();
}
double lst_now = erfa::eraGst06a(ERFA_DJM0, ut1_mjd, ERFA_DJM0, tt_mjd) + lon;
double lst_now = erfa::eraGst06a(ERFA_DJM0, ut1_mjd, ERFA_DJM0, tt_mjd);
lst_now += lon;
result = lst - lst_now;
if (result < 0.0) { // the next sideral day
result += 2.0 * std::numbers::pi;
}
if (result > std::numbers::pi) { // object is already below the limit
return 0.0;
}
result *= mcc_UT1_to_sideral_ratio; // to UT1 scale
return result;
@ -498,7 +510,6 @@ public:
return std::nullopt;
}
auto el_prev = _db.front();
for (auto const& el : _db) {
if (ymd <= el.ymd) {
return real_secs_t{el.dut1};

View File

@ -68,6 +68,12 @@ int main(int argc, char* argv[])
// const double rr = 180.0 / std::numbers::pi * 60.0;
// std::cout << "A(arcmin) = " << A * rr << "; B(arcmin) = " << B * rr << "\n";
using namespace std::chrono_literals;
double stm = mcc::astrom::mcc_time_to_alt_limit(10.0, "06:30:00.0", "00:00:00.0", 43.646711, 41.440732,
std::chrono::system_clock::now(), 0.041s, 32.184s, 37.0s);
// auto stm_d = mcc::astrom::mcc_chrono_radians{stm};
std::cout << "STM: " << stm * 12.0 / std::numbers::pi * 60.0 << " minutes\n";
return ecode;
}

View File

@ -49,6 +49,10 @@ static std::string_view trimSpaces(R&& r, TrimType type = TrimType::TRIM_BOTH)
return std::string_view(f1, f2);
}
static std::string_view trimSpaces(const char* r, TrimType type = TrimType::TRIM_BOTH)
{
return trimSpaces(std::string_view(r), type);
}
template <typename T, std::ranges::contiguous_range R>
std::optional<T> numFromStr(R&& r)