mountcontrol/cxx/tests/astrom_test.cpp
2025-07-04 12:24:34 +03:00

169 lines
6.2 KiB
C++

#include <fstream>
#include <iostream>
// #include "../mcc_coord.h"
#include "../mcc_traits.h"
#include "../mount_astrom.h"
template <typename T, typename U>
concept ncr_c =
std::is_lvalue_reference_v<T> && !std::is_const_v<T> && std::constructible_from<std::remove_reference_t<T>, U>;
template <typename T>
concept ae_c = requires(T t) {
requires mcc::traits::mcc_systime_c<mcc::traits::mcc_func_arg1_t<decltype(&T::greg2jul)>>;
requires mcc::traits::mcc_output_arg_c<mcc::traits::mcc_func_argN_t<decltype(&T::greg2jul), 2>>;
// requires ncr_c<mcc::traits::mcc_func_argN_t<decltype(&T::greg2jul), 2>, double>;
// requires std::constructible_from<mcc::traits::mcc_func_argN_t<decltype(&T::greg2jul), 2>, double>;
};
struct AE {
// int greg2jul(mcc::traits::mcc_systime_c auto t, double f) {}
int greg2jul(std::chrono::sys_time<std::chrono::nanoseconds> t, double& f)
{
return 0;
}
};
// BTA coords from maps.yandex.ru: 43.646711, 41.440732
int main(int argc, char* argv[])
{
int ecode = 0;
if (argc < 2) {
std::cerr << "Usage: " << argv[0] << " bulletinA-filename\n";
exit(1);
}
std::ifstream ist(argv[1]);
if (!ist) {
std::cerr << "Invalid filename!\n";
exit(1);
}
mcc::astrom::MccIersBulletinA bullA;
bullA.dump(std::cout);
std::cout << "\n\n";
mcc::astrom::MccLeapSeconds lps;
lps.dump(std::cout);
std::cout << "\n\n";
bool ok = bullA.load(ist);
bullA.dump(std::cerr);
ist.close();
std::cout << "\n\n\n------------------\n";
constexpr auto nnn = mcc::astrom::MccLeapSeconds::real_secs_t{std::numeric_limits<double>::quiet_NaN()};
auto now = std::chrono::system_clock::now();
// std::cout << "LEAP SECS for now: " << lps[now].value_or(std::numeric_limits<double>::quiet_NaN()) << "\n";
std::cout << "LEAP SECS for now: " << lps[now].value_or(nnn) << "\n";
// std::cout << "DUT1 for now: " << bullA.DUT1(now).value_or(std::numeric_limits<double>::quiet_NaN()) << " (" <<
// now
std::cout << "DUT1 for now: " << bullA.DUT1(now).value_or(nnn) << " (" << now << ")\n";
double mjd = 61077.4;
// std::cout << "DUT1 for MJD = " << mjd << ": " <<
// bullA.DUT1(mjd).value_or(std::numeric_limits<double>::quiet_NaN())
std::cout << "DUT1 for MJD = " << mjd << ": " << bullA.DUT1(mjd).value_or(nnn) << "\n";
auto st = mcc::astrom::mcc_julday(now, mjd);
std::cout << "MJD for now: " << std::setprecision(19) << mjd << " (" << now << ")\n";
// double pres = 786.6;
// double temp = 0.0;
// double hum = 0.8;
// double A, B;
// erfa::eraRefco(pres, temp, hum, 0.5, &A, &B);
// const double rr = 180.0 / std::numbers::pi * 60.0;
// std::cout << "A(arcmin) = " << A * rr << "; B(arcmin) = " << B * rr << "\n";
float alt_lim = 10.0;
std::string_view ra_str = "06:30:00.0", dec_str = "00:00:00.0";
std::cout << "\n\nTimes to object (RA = " << ra_str << ", DEC = " << dec_str << ") sets to given altitude ("
<< alt_lim << " degrees):\n";
using namespace std::chrono_literals;
// auto stm = mcc::astrom::mcc_time_to_alt_limit(alt_lim, ra_str, dec_str, 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";
// std::cout << "STM: " << std::chrono::duration_cast<std::chrono::minutes>(stm) << " minutes\n";
alt_lim = 85.0;
ra_str = "02:30:00.0", dec_str = "45:00:00.0";
std::cout << "\n\nTimes to object (RA = " << ra_str << ", DEC = " << dec_str << ") sets to given altitude ("
<< alt_lim << " degrees):\n";
auto stm =
mcc::astrom::mcc_time_to_alt({alt_lim, mcc::mcc_degrees}, {ra_str, mcc::mcc_hms}, dec_str, 43.646711_degs,
41.440732_degs, std::chrono::system_clock::now(), 0.041s, 32.184s, 37.0s);
// stm = mcc::astrom::mcc_time_to_alt_limit(alt_lim, ra_str, dec_str, 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";
std::cout << "STM: " << stm.first << ", " << stm.second << " seconds\n";
std::cout << "STM: " << std::chrono::duration_cast<std::chrono::minutes>(stm.first) << ", "
<< std::chrono::duration_cast<std::chrono::minutes>(stm.second) << " minutes\n";
std::cout << "\n\n\n";
const size_t N = 1000;
double mjds[N];
now = std::chrono::system_clock::now();
st = mcc::astrom::mcc_julday(now, mjds, std::chrono::milliseconds(100));
std::cout << std::format("comp time for {} MJDs = {}\n", N, std::chrono::system_clock::now() - now);
std::cout << std::format("MIN MJD: {:.16f}; \nMAX MJD: {:.16f}\n", mjds[0], mjds[N - 1]);
double elong = 43.646711 * std::numbers::pi / 180.0;
double phi = 41.440732 * std::numbers::pi / 180.0;
double aob, zob, hob, dob, rob, eo;
now = std::chrono::system_clock::now();
for (auto& el : mjds) {
st = erfa::eraAtco13(1.343523, 0.32352345, 0.0, 0.0, 0.0, 0.0, ERFA_DJM0, el, 0.041, elong, phi, 2100.0, 0.0,
0.0, 700.0, 10.0, 0.8, 0.5, &aob, &zob, &hob, &dob, &rob, &eo);
}
std::cout << std::format(
"comp time for {} coordinate transf = {}\n", N,
std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now() - now));
std::cout << "\n\n\nMccCoordinate class test:\n";
mcc::MccAngle ra("12:00:00", mcc::mcc_hms);
std::cout << "sin(12h) = " << std::sin(ra) << "\n";
std::cout << "cos(12h) = " << std::cos(ra) << "\n";
ra = "18:00:00"_hms;
std::cout << "sin(18h) = " << std::sin(ra) << "\n";
std::cout << "cos(18h) = " << std::cos(ra) << "\n";
ra = mcc::MccAngle{45.0, mcc::mcc_degrees};
std::cout << "sin(45) = " << std::sin(ra) << "\n";
std::cout << "cos(45) = " << std::cos(ra) << "\n";
std::cout << ra.sexagesimal(false, 4) << "\n";
std::cout << "\n\n\n";
std::cout << "ae_c = " << std::boolalpha << ae_c<AE> << "\n";
return ecode;
}