...
This commit is contained in:
parent
66c3c7e6c7
commit
09cf5f9c19
@ -251,35 +251,145 @@ public:
|
||||
|
||||
|
||||
// ICRS to observed
|
||||
// azimuth must be counted from the South through the West
|
||||
error_t icrs2obs(mcc_angle_c auto ra_icrs,
|
||||
mcc_angle_c auto dec_icrs,
|
||||
mcc_angle_c auto* ra_obs,
|
||||
mcc_angle_c auto* dec_obs,
|
||||
mcc_angle_c auto* ha_obs,
|
||||
mcc_angle_c auto* az,
|
||||
mcc_angle_c auto* zd,
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
// returned azimuth is counted from the South through the West
|
||||
error_t icrsToObs(mcc_angle_c auto const& ra_icrs,
|
||||
mcc_angle_c auto const& dec_icrs,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto* ra_obs,
|
||||
mcc_angle_c auto* dec_obs,
|
||||
mcc_angle_c auto* ha_obs,
|
||||
mcc_angle_c auto* az,
|
||||
mcc_angle_c auto* zd,
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
{
|
||||
return icrsTo(true, ra_icrs, dec_icrs, epoch, ra_obs, dec_obs, ha_obs, az, zd, obj_params);
|
||||
}
|
||||
|
||||
error_t icrs2obs(MccSkyRADEC_ICRS const& radec_icrs,
|
||||
MccSkyRADEC_OBS* radec_obs,
|
||||
MccSkyAZZD* azzd,
|
||||
mcc_angle_c auto* ha_obs) {
|
||||
error_t icrsToObs(MccSkyRADEC_ICRS const& radec_icrs,
|
||||
MccSkyRADEC_OBS* radec_obs,
|
||||
MccSkyAZZD* azzd,
|
||||
mcc_angle_c auto* ha_obs,
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
{
|
||||
double ra_obs, dec_obs, az, zd, ha;
|
||||
|
||||
MccCelestialCoordEpoch ep;
|
||||
ep.fromMJD((radec_obs->MJD()));
|
||||
|
||||
auto err = icrsToObs(radec_icrs.x(), radec_icrs.y(), ep, &ra_obs, &dec_obs, &ha, &az, &zd, obj_params);
|
||||
|
||||
if (!err) {
|
||||
if (radec_obs) {
|
||||
radec_obs->setX(ra_obs);
|
||||
radec_obs->setY(dec_obs);
|
||||
}
|
||||
|
||||
if (azzd) {
|
||||
azzd->setMJD(radec_obs->MJD());
|
||||
azzd->setX(az);
|
||||
azzd->setY(zd);
|
||||
}
|
||||
|
||||
if (ha_obs) {
|
||||
*ha_obs = ha;
|
||||
}
|
||||
}
|
||||
|
||||
return err;
|
||||
};
|
||||
|
||||
error_t obs2icrs(MccCoordPairKind obs_type,
|
||||
mcc_angle_c auto const& co_lon,
|
||||
mcc_angle_c auto const& co_lat,
|
||||
mcc_angle_c auto* ra_icrs,
|
||||
mcc_angle_c auto* dec_icrs)
|
||||
|
||||
// ICRS to apparent (in vacuo)
|
||||
// returned azimuth is counted from the South through the West
|
||||
error_t icrsToApp(mcc_angle_c auto const& ra_icrs,
|
||||
mcc_angle_c auto const& dec_icrs,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto* ra_app,
|
||||
mcc_angle_c auto* dec_app,
|
||||
mcc_angle_c auto* ha_app,
|
||||
mcc_angle_c auto* az,
|
||||
mcc_angle_c auto* zd, // should be interpretated as zenithal distance corrected for refraction
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
{
|
||||
return icrsTo(false, ra_icrs, dec_icrs, epoch, ra_app, dec_app, ha_app, az, zd, obj_params);
|
||||
}
|
||||
|
||||
|
||||
error_t obs2icrs(mcc_coord_pair_c auto const& xy_obs, MccSkyRADEC_ICRS* radec_icrs) {}
|
||||
error_t icrsToApp(MccSkyRADEC_ICRS const& radec_icrs,
|
||||
MccSkyRADEC_OBS* radec_app,
|
||||
MccSkyAZZD* azzd,
|
||||
mcc_angle_c auto* ha_app,
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
{
|
||||
double ra_app, dec_app, az, zd, ha;
|
||||
|
||||
MccCelestialCoordEpoch ep;
|
||||
ep.fromMJD((radec_app->MJD()));
|
||||
|
||||
auto err = icrsToApp(radec_icrs.x(), radec_icrs.y(), ep, &ra_app, &dec_app, &ha, &az, &zd, obj_params);
|
||||
|
||||
if (!err) {
|
||||
if (radec_app) {
|
||||
radec_app->setX(ra_app);
|
||||
radec_app->setY(dec_app);
|
||||
}
|
||||
|
||||
if (azzd) {
|
||||
azzd->setMJD(radec_app->MJD());
|
||||
azzd->setX(az);
|
||||
azzd->setY(zd);
|
||||
}
|
||||
|
||||
if (ha_app) {
|
||||
*ha_app = ha;
|
||||
}
|
||||
}
|
||||
|
||||
return err;
|
||||
}
|
||||
|
||||
|
||||
error_t obsToICRS(MccCoordPairKind obs_type,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto const& co_lon,
|
||||
mcc_angle_c auto const& co_lat,
|
||||
mcc_angle_c auto* ra_icrs,
|
||||
mcc_angle_c auto* dec_icrs)
|
||||
{
|
||||
return toICRS(true, obs_type, epoch, co_lon, co_lat, ra_icrs, dec_icrs);
|
||||
}
|
||||
|
||||
|
||||
error_t obsToICRS(mcc_coord_pair_c auto const& xy_obs, MccSkyRADEC_ICRS* radec_icrs)
|
||||
{
|
||||
double ra, dec;
|
||||
|
||||
MccCelestialCoordEpoch ep;
|
||||
ep.fromMJD((xy_obs->MJD()));
|
||||
|
||||
auto err = obsToICRS(xy_obs.pair_kind, ep, xy_obs.x(), xy_obs.y(), &ra, &dec);
|
||||
if (err) {
|
||||
return err;
|
||||
}
|
||||
|
||||
if (radec_icrs) {
|
||||
radec_icrs->setX(ra);
|
||||
radec_icrs->setY(dec);
|
||||
}
|
||||
|
||||
return err;
|
||||
}
|
||||
|
||||
|
||||
error_t appToICRS(MccCoordPairKind app_type,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto const& co_lon,
|
||||
mcc_angle_c auto const& co_lat,
|
||||
mcc_angle_c auto* ra_icrs,
|
||||
mcc_angle_c auto* dec_icrs)
|
||||
{
|
||||
return toICRS(false, app_type, epoch, co_lon, co_lat, ra_icrs, dec_icrs);
|
||||
}
|
||||
|
||||
|
||||
error_t equationOrigins(const mcc_julday_c auto& mjd, mcc_angle_c auto* eo)
|
||||
@ -415,6 +525,197 @@ protected:
|
||||
engine_state_t _currentState{};
|
||||
|
||||
std::unique_ptr<std::mutex> _stateMutex;
|
||||
|
||||
|
||||
error_t icrsTo(bool observed, // true - observed, false - apparent
|
||||
mcc_angle_c auto const& ra_icrs,
|
||||
mcc_angle_c auto const& dec_icrs,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto* ra,
|
||||
mcc_angle_c auto* dec,
|
||||
mcc_angle_c auto* ha,
|
||||
mcc_angle_c auto* az,
|
||||
mcc_angle_c auto* zd,
|
||||
obj_pars_t* obj_params = nullptr)
|
||||
{
|
||||
int err;
|
||||
double r, d, h, a, z, eo;
|
||||
double pressure = 0.0; // 0 for apparent coordinates type (see ERFA's refco.c: if pressure is zero then
|
||||
// refraction is also zero)
|
||||
|
||||
|
||||
error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;
|
||||
|
||||
std::lock_guard lock{*_stateMutex};
|
||||
|
||||
if (observed) {
|
||||
pressure = _currentState.meteo.pressure;
|
||||
}
|
||||
|
||||
auto dut1 = _currentState._bulletinA.DUT1(epoch.MJD());
|
||||
|
||||
if (!dut1.has_value()) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
|
||||
}
|
||||
|
||||
auto pol_pos = _currentState._bulletinA.polarCoords(epoch.MJD());
|
||||
if (!pol_pos.has_value()) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
|
||||
}
|
||||
|
||||
// const auto arcsec2rad = std::numbers::pi / 180 / 3600;
|
||||
const auto arcsec2rad = 1.0_arcsecs;
|
||||
pol_pos->x *= arcsec2rad;
|
||||
pol_pos->y *= arcsec2rad;
|
||||
|
||||
if (obj_params) {
|
||||
err = eraAtco13(ra_icrs, dec_icrs, obj_params->pm_RA, obj_params->pm_DEC, obj_params->parallax,
|
||||
obj_params->radvel, ERFA_DJM0, epoch.MJD(), dut1->count(), _currentState.lon,
|
||||
_currentState.lat, _currentState.elev, pol_pos->x, pol_pos->y, pressure,
|
||||
_currentState.meteo.temperature, _currentState.meteo.humidity, _currentState.wavelength, &a,
|
||||
&z, &h, &d, &r, &eo);
|
||||
} else {
|
||||
err = eraAtco13(ra_icrs, dec_icrs, 0.0, 0.0, 0.0, 0.0, ERFA_DJM0, epoch.MJD(), dut1->count(),
|
||||
_currentState.lon, _currentState.lat, _currentState.elev, pol_pos->x, pol_pos->y, pressure,
|
||||
_currentState.meteo.temperature, _currentState.meteo.humidity, _currentState.wavelength, &a,
|
||||
&z, &h, &d, &r, &eo);
|
||||
}
|
||||
|
||||
if (err == 1) {
|
||||
ret = MccCCTE_ERFAErrorCode::ERROR_DUBIOUS_YEAR;
|
||||
} else if (err == -1) {
|
||||
ret = MccCCTE_ERFAErrorCode::ERROR_UNACCEPTABLE_DATE;
|
||||
}
|
||||
|
||||
|
||||
if (ra) {
|
||||
*ra = r;
|
||||
}
|
||||
|
||||
if (dec) {
|
||||
*dec = d;
|
||||
}
|
||||
|
||||
if (az) {
|
||||
// NOTE: according to definition of astronomical azimuth it is counted from the South through the West, but
|
||||
// in the ERFA the azimuth is counted from the North through the East!!!
|
||||
//
|
||||
*az = MccAngle(a + std::numbers::pi).normalize<MccAngle::NORM_KIND_0_360>();
|
||||
}
|
||||
|
||||
if (zd) {
|
||||
*zd = z;
|
||||
}
|
||||
|
||||
return ret;
|
||||
}
|
||||
|
||||
error_t toICRS(bool observed, // true - observed, false - apparent
|
||||
MccCoordPairKind pair_type,
|
||||
mcc_coord_epoch_c auto const& epoch,
|
||||
mcc_angle_c auto const& co_lon,
|
||||
mcc_angle_c auto const& co_lat,
|
||||
mcc_angle_c auto* ra_icrs,
|
||||
mcc_angle_c auto* dec_icrs)
|
||||
{
|
||||
error_t ret = MccCCTE_ERFAErrorCode::ERROR_OK;
|
||||
|
||||
// check coordinate pair consistency
|
||||
if (mcc_is_app_coordpair(pair_type) && observed) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
|
||||
}
|
||||
if (mcc_is_obs_coordpair(pair_type) && !observed) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
|
||||
}
|
||||
|
||||
std::lock_guard lock{*_stateMutex};
|
||||
|
||||
auto dut1 = _currentState._bulletinA.DUT1(epoch.MJD());
|
||||
|
||||
if (!dut1.has_value()) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
|
||||
}
|
||||
|
||||
auto pol_pos = _currentState._bulletinA.polarCoords(epoch.MJD());
|
||||
if (!pol_pos.has_value()) {
|
||||
return MccCCTE_ERFAErrorCode::ERROR_BULLETINA_OUT_OF_RANGE;
|
||||
}
|
||||
|
||||
// const auto arcsec2rad = std::numbers::pi / 180 / 3600;
|
||||
const auto arcsec2rad = 1.0_arcsecs;
|
||||
pol_pos->x *= arcsec2rad;
|
||||
pol_pos->y *= arcsec2rad;
|
||||
|
||||
std::string type;
|
||||
double x, y, ra, dec;
|
||||
double pressure = 0.0;
|
||||
|
||||
if (observed) {
|
||||
pressure = _currentState.meteo.pressure;
|
||||
}
|
||||
|
||||
switch (pair_type) {
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_AZZD:
|
||||
// NOTE: according to definition of astronomical azimuth it is counted from the South through the West,
|
||||
// but in the ERFA the azimuth is counted from the North through the East!!!
|
||||
//
|
||||
x = co_lon + std::numbers::pi;
|
||||
y = co_lat;
|
||||
type = "A";
|
||||
break;
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_AZALT:
|
||||
// NOTE: according to definition of astronomical azimuth it is counted from the South through the West,
|
||||
// but in the ERFA the azimuth is counted from the North through the East!!!
|
||||
//
|
||||
x = co_lon + std::numbers::pi;
|
||||
y = MccCCTE_ERFA::PI_2 - co_lat; // altitude to zenithal distance
|
||||
type = "A";
|
||||
break;
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_HADEC_OBS:
|
||||
type = "H";
|
||||
x = co_lon;
|
||||
y = co_lat;
|
||||
break;
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_RADEC_OBS:
|
||||
type = "R";
|
||||
x = co_lon;
|
||||
y = co_lat;
|
||||
break;
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_HADEC_APP:
|
||||
type = "H";
|
||||
x = co_lon;
|
||||
y = co_lat;
|
||||
break;
|
||||
case mcc::MccCoordPairKind::COORDS_KIND_RADEC_APP:
|
||||
type = "R";
|
||||
x = co_lon;
|
||||
y = co_lat;
|
||||
break;
|
||||
default:
|
||||
return MccCCTE_ERFAErrorCode::ERROR_UNSUPPORTED_COORD_PAIR;
|
||||
};
|
||||
|
||||
int err =
|
||||
eraAtoc13(type.c_str(), x, y, ERFA_DJM0, epoch.MJD(), dut1->count(), _currentState.lon, _currentState.lat,
|
||||
_currentState.elev, pol_pos->x, pol_pos->y, pressure, _currentState.meteo.temperature,
|
||||
_currentState.meteo.humidity, _currentState.wavelength, &ra, &dec);
|
||||
|
||||
if (err == 1) {
|
||||
ret = MccCCTE_ERFAErrorCode::ERROR_DUBIOUS_YEAR;
|
||||
} else if (err == -1) {
|
||||
ret = MccCCTE_ERFAErrorCode::ERROR_UNACCEPTABLE_DATE;
|
||||
}
|
||||
|
||||
if (ra) {
|
||||
*ra_icrs = ra;
|
||||
}
|
||||
|
||||
if (dec) {
|
||||
*dec_icrs = dec;
|
||||
}
|
||||
|
||||
return ret;
|
||||
}
|
||||
};
|
||||
|
||||
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user