Commit 15d98c3b authored by Alexander Lapshin's avatar Alexander Lapshin

.

parent 22ec10db
......@@ -54,7 +54,7 @@ inline bool extract_value(const jsonval& parent, const std::string& key, bool ex
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<std::string>& vec, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<std::string>& vec, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -95,7 +95,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, std::vector<
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<int>& vec, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<int>& vec, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -120,7 +120,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, std::vector<
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<double>& vec, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, std::vector<double>& vec, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -146,7 +146,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, std::vector<
}
template <typename T>
bool loadjson(const jsonval& parent, std::string key, std::vector<T>& vec, bool exist = true, bool* loaded = 0)
bool loadjson(const jsonval& parent, std::string key, std::vector<T>& vec, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -169,7 +169,7 @@ bool loadjson(const jsonval& parent, std::string key, std::vector<T>& vec, bool
}
template <typename T, typename Tid>
bool loadjson(const jsonval& parent, std::string key, std::map<Tid, T>& map, bool exist = true, bool* loaded = 0)
bool loadjson(const jsonval& parent, std::string key, std::map<Tid, T>& map, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -192,7 +192,7 @@ bool loadjson(const jsonval& parent, std::string key, std::map<Tid, T>& map, boo
}
template <typename T>
bool loadjson(const jsonval& parent, std::string key, T& item, void (*f)(const jsonval&, T&), bool exist = true, bool* loaded = 0)
bool loadjson(const jsonval& parent, std::string key, T& item, void (*f)(const jsonval&, T&), const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -211,7 +211,7 @@ bool loadjson(const jsonval& parent, std::string key, T& item, void (*f)(const j
}
template <typename T>
bool loadjson(const jsonval& parent, const std::string& key, T& item, bool exist = true, bool* loaded = 0)
bool loadjson(const jsonval& parent, const std::string& key, T& item, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -273,7 +273,7 @@ void loadjson(const jsonval& parent, std::map<T2, T>& map)
}
}
static bool loadjson(const jsonval& parent, const std::string& key, double& value, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, double& value, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -306,7 +306,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, double& valu
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, int& value, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, int& value, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -339,7 +339,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, int& value,
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, bool& value, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, bool& value, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -375,7 +375,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, bool& value,
return true;
}
static bool loadjson(const jsonval& parent, const std::string& key, std::string& value, bool exist = true, bool* loaded = 0)
static bool loadjson(const jsonval& parent, const std::string& key, std::string& value, const bool exist = true, bool* loaded = nullptr)
{
jsonval::ConstMemberIterator itr;
if (!extract_value(parent, key, exist, loaded, itr)) {
......@@ -400,7 +400,7 @@ static bool loadjson(const jsonval& parent, const std::string& key, std::string&
}
template <typename T>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const T& object, void (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
void PutValue(jsonalloc& alc, jsonval& ss, const std::string& conName, const T& object, void (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
{
jsonval obj(rapidjson::kObjectType);
f(alc, obj, object);
......@@ -410,13 +410,14 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const T& object,
}
template <typename T>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vector<T>& vec, void (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
void PutValue(jsonalloc& alc, jsonval& ss, const std::string& conName, const std::vector<T>& vec, bool (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
{
jsonval jsvec(rapidjson::kArrayType);
for (typename std::vector<T>::const_iterator it = vec.begin(); it != vec.end(); ++it) {
for (auto it = vec.begin(); it != vec.end(); ++it) {
jsonval child(rapidjson::kObjectType);
f(alc, child, *it);
jsvec.PushBack(child, alc);
if (f(alc, child, *it)) {
jsvec.PushBack(child, alc);
}
}
jsonval jname(conName, alc);
......@@ -424,13 +425,14 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vecto
}
template <typename T>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vector<const T*>& vec, void (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T*))
void PutValue(jsonalloc& alc, jsonval& ss, const std::string conName, const std::vector<T*>& vec, bool (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
{
jsonval jsvec(rapidjson::kArrayType);
for (typename std::vector<const T*>::const_iterator it = vec.begin(); it != vec.end(); ++it) {
for (typename std::vector<T*>::const_iterator it = vec.begin(); it != vec.end(); ++it) {
jsonval child(rapidjson::kObjectType);
f(alc, child, *it);
jsvec.PushBack(child, alc);
if (f(alc, child, *(*it))) {
jsvec.PushBack(child, alc);
}
}
jsonval jname(conName, alc);
......@@ -438,19 +440,19 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vecto
}
template <typename T, typename TK>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::map<TK, T>& vec, void (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
void PutValue(jsonalloc& alc, jsonval& ss, const std::string conName, const std::map<TK, T>& vec, bool (*f)(rapidjson::Document::AllocatorType&, jsonval&, const T&))
{
jsonval jsvec(rapidjson::kObjectType);
for (typename std::map<TK, T>::const_iterator it = vec.begin(); it != vec.end(); ++it) {
jsonval child(rapidjson::kObjectType);
f(alc, child, it->second);
std::stringstream stream;
stream << it->first;
std::string key = stream.str();
if (f(alc, child, it->second)) {
std::stringstream stream;
stream << it->first;
std::string key = stream.str();
jsonval jkey(key, alc);
jsvec.AddMember(jkey, child, alc);
jsonval jkey(key, alc);
jsvec.AddMember(jkey, child, alc);
}
}
jsonval jname(conName, alc);
......@@ -458,7 +460,7 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::map<T
}
template <typename T>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vector<T>& vec)
void PutValue(jsonalloc& alc, jsonval& ss, const std::string conName, const std::vector<T>& vec)
{
jsonval jsvec(rapidjson::kArrayType);
for (typename std::vector<T>::const_iterator it = vec.begin(); it != vec.end(); ++it) {
......@@ -472,7 +474,7 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const std::vecto
}
template <typename T>
void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const T& object)
void PutValue(jsonalloc& alc, jsonval& ss, const std::string conName, const T& object)
{
jsonval obj(rapidjson::kObjectType);
PutValue(alc, obj, object);
......@@ -481,7 +483,7 @@ void PutValue(jsonalloc& alc, jsonval& ss, std::string conName, const T& object)
ss.AddMember(jname, obj, alc);
}
static void PutValue(jsonalloc& alc, jsonval& ss, const std::string& conName, size_t value)
static void PutValue(jsonalloc& alc, jsonval& ss, const std::string& conName, const size_t value)
{
jsonval jname(conName, alc);
ss.AddMember(jname, value, alc);
......@@ -519,7 +521,7 @@ static void PutValue(jsonalloc& alc, jsonval& ss, const std::string& conName, co
ss.AddMember(jname, jvalue, alc);
}
static void PutValue(jsonalloc& alc, jsonval& ss, double value)
static void PutValue(jsonalloc& alc, jsonval& ss, const double value)
{
ss.PushBack(value, alc);
}
......@@ -530,7 +532,7 @@ static void PutValue(jsonalloc& alc, jsonval& ss, const std::string& value)
ss.PushBack(jvalue, alc);
}
static void flushjson(const jsondoc& outroot, const std::string& outputfile, bool pretty = true)
static void flushjson(const jsondoc& outroot, const std::string& outputfile, const bool pretty = true)
{
if (!outputfile.length()) {
if (pretty) {
......@@ -603,3 +605,12 @@ static void loadFileJson(const std::string& filename, rapidjson::Document& docum
throw Exp() << "json parse error " << filename << ", code: " << ok.Code() << ", offset: " << ok.Offset();
}
}
static void loadStrJson(const std::string& data, rapidjson::Document& document)
{
rapidjson::ParseResult ok = document.Parse(data.c_str());
if (!ok) {
throw Exp() << "json parse error " << ", code: " << ok.Code() << ", offset: " << ok.Offset();
}
}
......@@ -3,36 +3,72 @@
#include <string>
#include <sstream>
class Converter
{
// enum STR2INT_ERROR { XSUCCESS, XOVERFLOW, XUNDERFLOW, XINCONVERTIBLE };
//
// STR2INT_ERROR str2int(int& i, char const* s, int base = 0)
// {
// char* end;
// long l;
// errno = 0;
// l = strtol(s, &end, base);
// if ((errno == ERANGE && l == LONG_MAX) || l > INT_MAX) {
// return XOVERFLOW;
// }
// if ((errno == ERANGE && l == LONG_MIN) || l < INT_MIN) {
// return XUNDERFLOW;
// }
// if (*s == '\0' || *end != '\0') {
// return XINCONVERTIBLE;
// }
// i = l;
// return XSUCCESS;
// }
class Converter {
public:
static bool ToNumber(const std::string& str, unsigned short int& val)
{
std::istringstream ss(str);
ss >> val;
return !ss.fail();
// std::istringstream ss(str);
// ss >> val;
// return !ss.fail();
// int ival;
// bool res = str2int(ival, str.c_str()) == XSUCCESS;
// val = ival;
val = std::stoi(str, nullptr);
return true;
}
static bool ToNumber(const std::string& str, int& val)
{
std::istringstream ss(str);
ss >> val;
return !ss.fail();
//std::istringstream ss(str);
//ss >> val;
//return !ss.fail();
val = std::stoi(str, nullptr);
return true;
}
static bool ToNumber(const std::string& str, double& val)
{
std::istringstream ss(str);
ss >> val;
return !ss.fail();
// std::istringstream ss(str);
// ss >> val;
// return !ss.fail();
val = std::stod(str, nullptr);
return true;
}
static bool ToNumber(const std::string& str, float& val)
{
std::istringstream ss(str);
ss >> val;
return !ss.fail();
// std::istringstream ss(str);
// ss >> val;
// return !ss.fail();
val = std::stof(str, nullptr);
return true;
}
static bool ToNumber(const std::string& str, bool& val)
......@@ -41,14 +77,17 @@ public:
val = false;
return true;
}
else if (str == "true" || str == "TRUE") {
if (str == "true" || str == "TRUE") {
val = true;
return true;
}
else {
std::istringstream ss(str);
ss >> val;
return !ss.fail();
}
//std::istringstream ss(str);
//ss >> val;
//return !ss.fail();
val = std::stoi(str, nullptr);
return true;
}
};
......@@ -356,3 +356,11 @@ std::string implode(std::string delim, const std::vector<T>& vec)
return ss.str();
}
template <typename T>
bool in_array(const std::vector<T>& src, T subject)
{
return std::find(std::begin(src), std::end(src), subject) != std::end(src);
}
......@@ -2,7 +2,6 @@
#include "Converter.h"
#include "TM.h"
#include <vector>
#include <string>
#include <iostream>
#include <iomanip>
......@@ -91,6 +90,13 @@ static bool StrToDate(const std::string& date, TimeGD& gd)
}
}
static TimeJD StrToDate(const std::string& sdate)
{
TimeGD date;
StrToDate(sdate, date);
return date;
}
static int isleap(int year)
{
return (year % 4 == 0 && year % 100 != 0) || (year % 400 == 0);
......
#ifndef XMLFUNCTIONS_H
#define XMLFUNCTIONS_H
#pragma once
#include <myXML_Document.hpp>
#include "MException.h"
#include <iostream>
#include <fstream>
#include <vector>
#include <string>
#include "TM.h"
......@@ -18,7 +15,6 @@ bool FindFirstContainer(myXML::Contain* parent, const std::string& path, myXML::
bool FindXMLContainer(const std::map<std::string, myXML::Contain*>& items, const std::string& name, myXML::Contain*& result);
bool FindXMLContainer(myXML::Contain* parent, const std::string& path, myXML::Contain*& result);
template <typename T>
void ReadXML(T storage, TimeGD& result);
template <typename T>
......@@ -36,7 +32,6 @@ void ReadXML(T storage, std::string subname, std::vector<T2>& result);
template <typename T, typename T2, typename T3>
void ReadXML(T storage, std::string subname, std::map<T2, T3>& result);
template <typename T>
void ReadXMLVector(myXML::Contain* con, std::string subname, std::vector<T>& result)
{
......@@ -364,6 +359,14 @@ void PutValue(std::ostream& ss, std::string conName, const T& object)
ss << "</" << conName << ">";
}
template <typename T>
void PutValue(std::ostream& ss, std::string conName, const T* object)
{
ss << "<" << conName << ">";
object->ToXml(ss);
ss << "</" << conName << ">";
}
template <typename T>
void PutValueAttr(std::ostream& ss, std::string conName, const T& object)
{
......@@ -432,8 +435,7 @@ static void PutValue(std::ostream& ss, std::string conName, std::string subName,
static void PutValue(std::ostream& ss, const std::string& conName, const double value, const std::string& attributes = "")
{
int precision = 15;
ss << std::setprecision(precision);
ss << std::setprecision(std::numeric_limits<double>::digits10 + 2);
if (attributes.empty()) {
ss << "<" << conName << ">" << value << "</" << conName << ">";
}
......@@ -466,4 +468,3 @@ static void LoadXmlString(const std::string& str, myXML::Document& doc)
}
}
#endif
#pragma once
struct eh
{
double prf;
double xpond_xmit_delay;
double xpond_utc_off;
double xpond_osc_drift;
double cofm_corr;
int cospar_id;
int sic;
int norad;
int eph_seq;
int sttpyear;
int sttpmonth;
int sttpday;
int sttphour;
int sttpmin;
int sttpsec;
int endpyear;
int endpmonth;
int endpday;
int endphour;
int endpmin;
int endpsec;
int ephsep;
int compat;
int version;
int epyear;
int epmon;
int epday;
int ephour;
int tar_type;
int ref_frame;
int rot_type;
int cofm_app;
int ltcorr;
int atrk_ro_0;
int xtrk_ro_0;
int rtrk_ro_0;
int atrk_ro_1;
int xtrk_ro_1;
int rtrk_ro_1;
int atrk_ro_2;
int xtrk_ro_2;
int rtrk_ro_2;
char eph_source[4];
char tar_name[10];
char eph_notes[10];
};
// class export_ei {
// std::string jds;
// std::string x;
// std::string y;
// std::string z;
// std::string vx;
// std::string vy;
// std::string vz;
// };
struct ei
class ei
{
public:
ei();
~ei();
double* jdi;
double* jdf;
double* outv[6];
double* c12;
double* c23;
int* leapflag;
double* jds;
int nv;
};
int cephinit(char* str);
int cinterp_pvac(ei& ephi, double jdint, double seceph, double outvec[6], int fvel, int* ircode);
int cinterp_pvac(const ei& ephi, double jdint, double seceph, double outvec[6], int fvel, int* ircode);
#include <stdio.h>
#include <math.h>
#include <cmath>
#include "ephi.h"
int debug = 0;
int bs_lower_bound(double a[], int n, double x) {
int bs_lower_bound(const double a[], const int n, const double x)
{
int l = 0;
int h = n; // Not n - 1
while (l < h) {
int mid = (l + h) / 2;
const int mid = (l + h) / 2;
if (x <= a[mid]) {
h = mid;
}
......@@ -21,6 +21,30 @@ int bs_lower_bound(double a[], int n, double x) {
void hermite(int index, int ityp, const double* x, const double* y, const double* z, int nmax, int nval, double xp, double* yp, double* zp, int* ircode);
ei::ei()
:
jdi(nullptr),
jdf(nullptr),
outv(),
jds(nullptr),
nv(-1)
{
for (auto& el : outv) {
el = nullptr;
}
}
ei::~ei()
{
delete[] jdi;
delete[] jdf;
delete[] jds;
for (auto& el : outv) {
delete[] el;
}
}
/*-----------------------------------------------------------------------
void cinterp_pvac (jdint, seceph, outvec, backvec,
aberout,aberback, corr12, corr23, ircode)
......@@ -43,22 +67,21 @@ Date: 1 May, 2003
History:
-----------------------------------------------------------------------*/
int cinterp_pvac(
ei& ephi,
double jdint,
const ei& ephi,
double jdint,
double seceph,
double outvec[6],
int fvel,
double outvec[6],
const int fvel,
int* ircode)
{
double itime, z, vz; /* z is dummy for derivative field */
int j;
double z, vz; /* z is dummy for derivative field */
const double itime = ((seceph / 86400.e0 - ephi.jdf[0]) + (jdint - ephi.jdi[0])) * 86400.e0;
itime = ((seceph / 86400.e0 - ephi.jdf[0]) + (jdint - ephi.jdi[0])) * 86400.e0;
/* 3 or 6 elements? Handle velocity if present */
int pindex = bs_lower_bound((double*)ephi.jds, ephi.nv, itime);
const int pindex = bs_lower_bound((double*)ephi.jds, ephi.nv, itime);
for (j = 0; j < 6; j++) {
for (int j = 0; j < 6; j++) {
/* Out */
/* Position */
if (j < 3) {
......@@ -71,8 +94,8 @@ int cinterp_pvac(
hermite(pindex, -fvel, ephi.jds, &ephi.outv[j][0], &z, ephi.nv, 10, itime, &outvec[j], &outvec[j + 3], ircode);
}
}
/* Velocity */
/* Interpolate if they exist */
/* Velocity */
/* Interpolate if they exist */
else if (fabs(ephi.outv[j][0]) > 1.e-10) {
hermite(pindex, 1, ephi.jds, &ephi.outv[j][0], &z, ephi.nv, 10, itime, &outvec[j], &vz, ircode);
}
......@@ -84,4 +107,3 @@ int cinterp_pvac(
return 0;
}
......@@ -4,8 +4,7 @@
#include <algorithm>
template <bool dir, typename vartype>
struct cmppos : public std::binary_function<vartype, vartype, bool>
{
struct cmppos : public std::binary_function<vartype, vartype, bool> {
bool operator()(vartype s1, vartype s2) const
{
const TimeJD& v1 = s1.get_date();
......@@ -15,23 +14,25 @@ struct cmppos : public std::binary_function<vartype, vartype, bool>
}
};
static double get_step(double period, double ecc)
static double get_step(const double period, const double ecc)
{
if (ecc < 0.1) {
return std::max((double)(period / 100.0 * 60), 1.0);
}
else if (ecc < 0.2) {
if (ecc < 0.2) {
return std::max((double)(period / 200.0 * 60), 1.0);
}
else if (ecc < 0.4) {
if (ecc < 0.4) {
return std::max((double)(period / 400.0 * 60), 1.0);
}
else if (ecc < 0.6) {
if (ecc < 0.6) {
return std::max((double)(period / 600.0 * 60), 1.0);
}
else {
return std::max((double)(period / 1000.0 * 60), 1.0);
}
return std::max((double)(period / 1000.0 * 60), 1.0);
}
static double get_step(const SInitOrbit& orbit)
......@@ -41,31 +42,22 @@ static double get_step(const SInitOrbit& orbit)
return get_step(period, ecc);
}
XInterpolator::XInterpolator()
: m_ephi(0)
{
}
XInterpolator::XInterpolator() = default;
XInterpolator::XInterpolator(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end)
:
m_ephi(0),
m_orbit(orbit)
{
create(orbit, beg, end);
}
XInterpolator::~XInterpolator()
XInterpolator::XInterpolator(const std::vector<IRes>& vec)
{
if (m_ephi) {
delete m_ephi;
m_ephi = 0;
}
create(vec);
}
void XInterpolator::propagate_dir(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, int dir, std::vector<IRes>& result)
void XInterpolator::propagate_dir(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, const int dir, std::vector<IRes>& result)
{
Propagator prop(orbit, false);
double step = get_step(orbit);
const double step = get_step(orbit);
bool exit = false;
TimeJD date = beg;
......@@ -80,12 +72,13 @@ void XInterpolator::propagate_dir(const SInitOrbit& orbit, const TimeJD& beg, co
PhasePoint6D pos;
prop.Propagate(date, pos);
result.push_back(IRes(pos.CoordsVel, date, step, age, orbit.GetId()));
result.emplace_back(pos.CoordsVel, date, age, orbit.GetId());
date = ShiftDate(date, dir * step);
}
}
void XInterpolator::propagate(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, int dir, std::vector<IRes>& result)
void XInterpolator::propagate(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, int dir,
std::vector<IRes>& result)
{
if (dir != 0) {
propagate_dir(orbit, beg, end, dir, result);
......@@ -95,19 +88,19 @@ void XInterpolator::propagate(const SInitOrbit& orbit, const TimeJD& beg, const
propagate_dir(orbit, ShiftDate(orbit.GetDate(), -get_step(orbit) / 2.0), beg, -1, result);
}
std::sort(result.begin(), result.end(), [](const IRes& a, const IRes& b)
{
std::sort(result.begin(), result.end(), [](const IRes& a, const IRes& b) {
return a.get_date() < b.get_date();
});
}
void XInterpolator::calc_positions(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, std::vector<IRes>& result)
void XInterpolator::calc_positions(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end,
std::vector<IRes>& result)
{
int n = 100;
int step = get_step(orbit);
const int n = 100;
const int step = get_step(orbit);
TimeJD beg_ext = ShiftDate(beg, -n * step);
TimeJD end_ext = ShiftDate(end, +n * step);
const TimeJD beg_ext = ShiftDate(beg, -n * step);
const TimeJD end_ext = ShiftDate(end, +n * step);
if (orbit.GetDate() <= beg_ext) {
propagate(orbit, beg_ext, end_ext, 1, result);
......@@ -120,17 +113,8 @@ void XInterpolator::calc_positions(const SInitOrbit& orbit, const TimeJD& beg, c
}
}
static void init(ei& ephi, int size)
static void init(ei& ephi, const int size)
{
int j;
// int idir;
// int mjd;
// double sod;
// double dgast = 0.e0;
// char ephstr[256];
/* initialize all the arrays, in case there is nothing to read */
ephi.jdi = new double[size]();
memset(ephi.jdi, 0, size * sizeof(double));
......@@ -140,71 +124,17 @@ static void init(ei& ephi, int size)
ephi.jds = new double[size]();
memset(ephi.jds, 0, size * sizeof(double));
// ephi.offsetjdi = new double[size]();
// memset(ephi.offsetjdi, 0, size*sizeof(double));
// ephi.offsetjdf = new double[size]();
// memset(ephi.offsetjdf, 0, size*sizeof(double));
for (j = 0; j < 6; j++) {
ephi.outv[j] = new double[size]();
memset(ephi.outv[j], 0, size * sizeof(double));
//ephi.backv[j] = new double[size]();
//memset(ephi.backv[j], 0, size*sizeof(double));
}
for (j = 0; j < 3; j++) {
//ephi.aberoutv[j] = new double[size]();
//memset(ephi.aberoutv[j], 0, size*sizeof(double));
//ephi.aberbackv[j] = new double[size]();
//memset(ephi.aberbackv[j], 0, size*sizeof(double));
//ephi.offsetoutv[j] = new double[size]();
//memset(ephi.offsetoutv[j], 0, size*sizeof(double));
//ephi.offsetbackv[j] = new double[size]();
//memset(ephi.offsetbackv[j], 0, size*sizeof(double));
//ephi.rotangv[j] = new double[size]();
//memset(ephi.rotangv[j], 0, size*sizeof(double));
for (auto& j : ephi.outv) {
j = new double[size]();
memset(j, 0, size * sizeof(double));
}
ephi.c12 = new double[size]();
memset(ephi.c12, 0, size * sizeof(double));
ephi.c23 = new double[size]();
memset(ephi.c23, 0, size * sizeof(double));
ephi.leapflag = new int[size]();
memset(ephi.leapflag, 0, size * sizeof(int));
// ephi.oscrelcorr = new double[size]();
// memset(ephi.oscrelcorr, 0, size*sizeof(double));
// ephi.gast = new double[size]();
// memset(ephi.gast, 0, size*sizeof(double));
// ephi.eopxy[0] = new double[size]();
// memset(ephi.eopxy[0], 0, size*sizeof(double));
// ephi.eopxy[1] = new double[size]();
// memset(ephi.eopxy[1], 0, size*sizeof(double));
// ephi.eopdutc = new double[size]();
// memset(ephi.eopdutc, 0, size*sizeof(double));
ephi.nv = -1;
// ephi.nvoff = -1;
// ephi.nvrot = -1;
// ephi.nveop = -1;
}
void fill(const Vect6& vec, double jdi, double jdf, ei& ephi)
void fill(const Vect6& vec, const double jdi, const double jdf, ei& ephi)
{
ephi.nv++;
ephi.leapflag[ephi.nv] = 0;
ephi.outv[0][ephi.nv] = vec.x * 1E6;
ephi.outv[1][ephi.nv] = vec.y * 1E6;
......@@ -220,47 +150,15 @@ void fill(const Vect6& vec, double jdi, double jdf, ei& ephi)
ephi.jdi[ephi.nv] += 1.e0;
ephi.jdf[ephi.nv] -= 1.e0;
}
/* In case there is not back vector (e.g., satellites)... */
// ephi.backv[0][ephi.nv] = -ephi.outv[0][ephi.nv];
// ephi.backv[1][ephi.nv] = -ephi.outv[1][ephi.nv];
// ephi.backv[2][ephi.nv] = -ephi.outv[2][ephi.nv];
// ephi.backv[3][ephi.nv] = -ephi.outv[3][ephi.nv];
// ephi.backv[4][ephi.nv] = -ephi.outv[4][ephi.nv];
// ephi.backv[5][ephi.nv] = -ephi.outv[5][ephi.nv];
}
void finalize(ei& ephi)
{
/* create a small usable time tag */
(ephi.nv)++;
for (int i = 0; i < ephi.nv; i++) {
ephi.jds[i] = ((ephi.jdf[i] - ephi.jdf[0]) +
(ephi.jdi[i] - ephi.jdi[0])) * 86400.e0;
ephi.jds[i] = ((ephi.jdf[i] - ephi.jdf[0]) + (ephi.jdi[i] - ephi.jdi[0])) * 86400.e0;
}
// (ephi.nvoff)++;
/*for (int i = 0; i < ephi.nvoff; i++){
ephi.offsetjds[i] = ((ephi.offsetjdf[i] - ephi.offsetjdf[0]) +
(ephi.offsetjdi[i] -
ephi.offsetjdi[0])) * 86400.e0;
}*/
// (ephi.nvrot)++;
/*for (int i = 0; i < ephi.nvrot; i++){
ephi.rotjds[i] = ((ephi.rotjdf[i] - ephi.rotjdf[0]) +
(ephi.rotjdi[i] - ephi.rotjdi[0])) * 86400.e0;
}*/
// (ephi.nveop)++;
/*for (int i = 0; i < ephi.nveop; i++){
ephi.eopjds[i] = ((ephi.eopjdf[i] - ephi.eopjdf[0]) +
(ephi.eopjdi[i] - ephi.eopjdi[0])) * 86400.e0;
}*/
}
/*void Interpolator::filter_addons()
......@@ -298,32 +196,36 @@ void finalize(ei& ephi)
void XInterpolator::create(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end)
{
assert(m_ephi == 0);
std::vector<IRes> positions;
calc_positions(orbit, beg, end, positions);
calc_positions(orbit, beg, end, m_positions);
create(positions);
}
m_ephi = new ei;
m_beg = m_positions.front().get_date();
m_end = m_positions.back().get_date();
void XInterpolator::create(const std::vector<IRes>& positions)
{
m_beg = positions.front().get_date();
m_end = positions.back().get_date();
init(*m_ephi, m_positions.size());
init(m_ephi, positions.size());
for (std::vector<IRes>::const_iterator it = m_positions.begin(); it != m_positions.end(); ++it) {
Vect6 pos = it->get_pos();
TimeJD date = it->get_date();
for (const auto& p : positions) {
Vect6 pos = p.get_pos();
TimeJD date = p.get_date();
PhasePoint6D itrf(static_cast<const IFrame&>(J2000Frame()), pos, date);
//itrf.Convert(GreenwichFrame());
//double jdi = date.GetDays();
//double jdf = date.GetFraction();
fill(itrf.CoordsVel, date.GetDays(), date.GetFraction(), *m_ephi);
fill(itrf.CoordsVel, date.GetDays(), date.GetFraction(), m_ephi);
}
finalize(*m_ephi);
finalize(m_ephi);
m_positions = positions;
}
Vect6 XInterpolator::get_pos(const TimeJD& date) const
Vect6 XInterpolator::get_pos(const TimeJD& date, const bool tks) const
{
if (!m_ephi) {
if (!m_ephi.jdi) {
throw Exp() << "interpolation not initialised\n";
}
......@@ -332,22 +234,33 @@ Vect6 XInterpolator::get_pos(const TimeJD& date) const
}
double outvec[6];
int fvel = 1;
const int fvel = 1;
int ircode;
double jdi = date.GetDays();
double sec = date.GetFraction() * 86400;
cinterp_pvac(*m_ephi, jdi, sec, outvec, fvel, &ircode);
const double jdi = date.GetDays();
const double sec = date.GetFraction() * 86400;
cinterp_pvac(m_ephi, jdi, sec, outvec, fvel, &ircode);
if (ircode != 0) {
throw Exp() << "interpolation fail\n";
}
if (!tks) {
return {
outvec[0],
outvec[1],
outvec[2],
outvec[3],
outvec[4],
outvec[5]
};
}
return Vect6(
outvec[0],
outvec[1],
outvec[2],
outvec[3],
outvec[4],
outvec[5]
);
return {
outvec[0] / 1E6,
outvec[1] / 1E6,
outvec[2] / 1E6,
outvec[3] / 1E3,
outvec[4] / 1E3,
outvec[5] / 1E3
};
}
......@@ -12,21 +12,19 @@ class XInterpolator
public:
XInterpolator();
XInterpolator(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end);
~XInterpolator();
Vect6 get_pos(const TimeJD& date) const;
int get_size() const { return m_positions.size(); }
XInterpolator(const std::vector<IRes>& vec);
Vect6 get_pos(const TimeJD& date, bool tks = false) const;
TimeJD get_beg() const { return m_beg; }
TimeJD get_end() const { return m_end; }
const SInitOrbit& get_orbit() const { return m_orbit; }
static void calc_positions(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, std::vector<IRes>& result);
void export_data() const;
std::vector<IRes> m_positions;
private:
void create(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end);
void create(const std::vector<IRes>& positions);
static void propagate(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, int dir, std::vector<IRes>& result);
static void propagate_dir(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, int dir, std::vector<IRes>& result);
TimeJD m_beg;
TimeJD m_end;
ei* m_ephi;
std::vector<IRes> m_positions;
SInitOrbit m_orbit;
ei m_ephi;
};
#include "ires.h"
#include <utility>
IRes::IRes(Vect6 vec, const TimeJD& date, int step, double age, TOrbitId orbitid, bool addon)
IRes::IRes(const Vect6& vec, const TimeJD& date, const double age, std::string orbitid, const bool addon)
:
m_orbitid(std::move(orbitid)),
m_vec(vec),
m_date(date),
m_step(step),
m_age(age),
m_orbitid(orbitid),
m_addon(addon)
{
}
......@@ -2,7 +2,7 @@
#include "DVectors.h"
#include "TM.h"
#include "OrbitId.h"
#include <string>
class IRes
{
......@@ -10,17 +10,15 @@ public:
IRes()
{
};
IRes(Vect6 vec, const TimeJD& date, int step, double age, TOrbitId orbitid, bool addon = false);
IRes(const Vect6& vec, const TimeJD& date, double age, std::string orbitid, bool addon = false);
TimeJD get_date() const { return m_date; }
Vect6 get_pos() const { return m_vec; }
int get_step() const { return m_step; }
bool is_addon() const { return m_addon; }
double get_age() const { return m_age; }
private:
TOrbitId m_orbitid;
std::string m_orbitid;
Vect6 m_vec;
TimeJD m_date;
int m_step;
double m_age;
bool m_addon;
};
#include "orbBlock.h"
#include <utility>
#include "TextFunctions.h"
OrbBlock::OrbBlock()
:
m_ip(nullptr)
{
}
OrbBlock::OrbBlock(
const SInitOrbit& orbit,
const TimeJD& beg,
......@@ -23,16 +30,16 @@ OrbBlock::OrbBlock(
const SInitOrbit& orbit,
const TimeJD& beg,
const TimeJD& end,
double mind,
const std::string& next_orbitid
const double mind,
std::string next_orbitid
)
:
m_orbit(orbit),
m_beg(beg),
m_end(end),
m_mind(mind),
m_next_orbitid(next_orbitid),
m_ip(0)
m_next_orbitid(std::move(next_orbitid)),
m_ip(nullptr)
{
#ifdef _DEBUG
m_beg_s = DateToStr(beg);
......@@ -49,3 +56,36 @@ void OrbBlock::interpolate(const TimeJD& beg, const TimeJD& end)
{
m_ip = new XInterpolator(m_orbit, beg, end);
}
void OrbBlock::load(const jsonval& parent)
{
std::string orbitid;
loadjson(parent, "id", orbitid);
loadjson(parent, "beg", m_beg);
loadjson(parent, "end", m_end);
loadjson(parent, "mind_km", m_mind, false);
std::string vecs;
loadjson(parent, "vecs", vecs, false);
static std::vector<std::string> arr = splitstr(vecs, ",");
std::vector<IRes> vres;
for (size_t i = 0; i < arr.size(); i += 8) {
int days;
Converter::ToNumber(arr[i], days);
double fraction;
Converter::ToNumber(arr[i + 1], fraction);
double x, y, z, vx, vy, vz;
Converter::ToNumber(arr[i + 2], x);
Converter::ToNumber(arr[i + 3], y);
Converter::ToNumber(arr[i + 4], z);
Converter::ToNumber(arr[i + 5], vx);
Converter::ToNumber(arr[i + 6], vy);
Converter::ToNumber(arr[i + 7], vz);
vres.emplace_back(Vect6(x, y, z, vx, vy, vz), TimeJD(days, fraction), 0, orbitid);
}
}
......@@ -3,15 +3,13 @@
#include "TM.h"
#include "StateVector.h"
#include "interpolate.h"
#include "json_functions.h"
class OrbBlock
{
class OrbBlock {
public:
OrbBlock()
{
};
OrbBlock();;
OrbBlock(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end);
OrbBlock(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, double mind, const std::string& next_orbitid);
OrbBlock(const SInitOrbit& orbit, const TimeJD& beg, const TimeJD& end, double mind, std::string next_orbitid);
~OrbBlock();
TimeJD get_beg() const { return m_beg; }
......@@ -21,6 +19,7 @@ public:
std::string get_next_orbitid() const { return m_next_orbitid; }
void interpolate(const TimeJD& beg, const TimeJD& end);
const XInterpolator* get_interpolator() const { return m_ip; }
void load(const jsonval& parent);
private:
SInitOrbit m_orbit;
TimeJD m_beg;
......
#include "orbBlocks.h"
#include "TimeFunctions.h"
#include "Propagator.h"
#include "xalgorithm.h"
#include "orb_date.h"
......@@ -18,13 +17,12 @@ void OrbBlocksStore::init(
{
SInitOrbits sordet_orbits = orbits;
std::sort(sordet_orbits.begin(), sordet_orbits.end(), [](const SInitOrbit& a, const SInitOrbit& b)
{
std::sort(sordet_orbits.begin(), sordet_orbits.end(), [](const SInitOrbit& a, const SInitOrbit& b) {
return a.GetDate() < b.GetDate();
});
m_beg_before_filter = beg ? *beg : sordet_orbits.front().GetDate();
m_end_before_filter = end ? *end : sordet_orbits.back().GetDate();
// m_beg_before_filter = beg ? *beg : sordet_orbits.front().GetDate();
// m_end_before_filter = end ? *end : sordet_orbits.back().GetDate();
if (filter_orbits) {
sordet_orbits = filter(sordet_orbits);
......@@ -103,6 +101,16 @@ SInitOrbits OrbBlocksStore::filter(const SInitOrbits& orbits) const
return result;
}
void loadjson(const jsonval& parent, OrbBlock& item)
{
item.load(parent);
}
void OrbBlocksStore::load(const jsonval& parent)
{
loadjson(parent, "records", m_blocks);
}
void OrbBlocksStore::init_min_dist(const SInitOrbits& orbits, const TimeJD& beg, const TimeJD& end, bool forced)
{
TimeJD left = beg;
......@@ -113,7 +121,7 @@ void OrbBlocksStore::init_min_dist(const SInitOrbits& orbits, const TimeJD& beg,
const SInitOrbit& orb1 = orbits[i];
if (i == orbits.size() - 1) {
m_blocks.push_back(OrbBlock(orb1, left, end));
m_blocks.emplace_back(orb1, left, end);
break;
}
......@@ -132,7 +140,7 @@ void OrbBlocksStore::init_min_dist(const SInitOrbits& orbits, const TimeJD& beg,
//no intersection
if (orb2beg < left) {
//orbit1 if out of left and orbit2 has left border inside, so orbit1 no longer needed
throw (Exp() << "assert point 0");
throw Exp() << "assert point 0";
}
TimeJD right;
......@@ -172,13 +180,12 @@ void OrbBlocksStore::init_min_dist(const SInitOrbits& orbits, const TimeJD& beg,
}
}
class MinDistEvaluator
{
class MinDistEvaluator {
public:
MinDistEvaluator(const XInterpolator& ip1, const XInterpolator& ip2, const TimeJD& base) :
m_base(base),
m_ip1(ip1),
m_ip2(ip2),
m_base(base)
m_ip2(ip2)
{
}
......@@ -186,8 +193,8 @@ public:
{
TimeJD date = ShiftDate(m_base, time * 86400);
Vect6 p1 = m_ip1.get_pos(date);
Vect6 p2 = m_ip2.get_pos(date);
const Vect6 p1 = m_ip1.get_pos(date);
const Vect6 p2 = m_ip2.get_pos(date);
return sqrt(
sqr(p1.x - p2.x) +
......@@ -201,8 +208,7 @@ private:
const XInterpolator& m_ip2;
};
void OrbBlocksStore::calc_min_dist(const SInitOrbit& orb1, const SInitOrbit& orb2, const TimeJD& beg, const TimeJD& end,
TimeJD& result_jd, double& result_dist) const
void OrbBlocksStore::calc_min_dist(const SInitOrbit& orb1, const SInitOrbit& orb2, const TimeJD& beg, const TimeJD& end, TimeJD& result_jd, double& result_dist) const
{
XInterpolator ip1(orb1, beg, end);
XInterpolator ip2(orb2, beg, end);
......@@ -242,30 +248,29 @@ void OrbBlocksStore::calc_min_dist(const SInitOrbit& orb1, const SInitOrbit& orb
Vect6 OrbBlocksStore::get_pos(const TimeJD& date) const
{
for (OrbBlocks::const_iterator it = m_blocks.begin(); it != m_blocks.end(); ++it) {
if (date >= it->get_beg() && date <= it->get_end()) {
return it->get_interpolator()->get_pos(date);
for (const auto& block : m_blocks) {
if (date >= block.get_beg() && date <= block.get_end()) {
return block.get_interpolator()->get_pos(date);
}
}
throw (Exp() << "no orbit found");
throw Exp() << "no orbit found";
}
const SInitOrbit& OrbBlocksStore::get_orbit(const TimeJD& date) const
{
for (OrbBlocks::const_iterator it = m_blocks.begin(); it != m_blocks.end(); ++it) {
if (date >= it->get_beg() && date <= it->get_end()) {
return it->get_interpolator()->get_orbit();
for (const auto& block : m_blocks) {
if (date >= block.get_beg() && date <= block.get_end()) {
return block.get_orbit();
}
}
throw (Exp() << "no orbit found");
throw Exp() << "no orbit found";
}
SInitOrbits OrbBlocksStore::calc_orbits(std::vector<TimeJD> dates, bool variates) const
{
std::sort(dates.begin(), dates.end(), [](const TimeJD& a, const TimeJD& b)
{
std::sort(dates.begin(), dates.end(), [](const TimeJD& a, const TimeJD& b) {
return a < b;
});
......
......@@ -4,6 +4,7 @@
#include <vector>
#include "orbBlock.h"
#include "StateVector.h"
#include "json_functions.h"
typedef std::vector<OrbBlock> OrbBlocks;
......@@ -18,13 +19,14 @@ public:
const SInitOrbit& get_orbit(const TimeJD& date) const;
const OrbBlocks& get_blocks() const { return m_blocks; }
OrbBlocks& get_blocks() { return m_blocks; }
const TimeJD& get_beg_before_filter() const { return m_beg_before_filter; }
const TimeJD& get_end_before_filter() const { return m_end_before_filter; }
//const TimeJD& get_beg_before_filter() const { return m_beg_before_filter; }
//const TimeJD& get_end_before_filter() const { return m_end_before_filter; }
void load(const jsonval& parent);
private:
void init_min_dist(const SInitOrbits& orbits, const TimeJD& beg, const TimeJD& end, bool forced);
void calc_min_dist(const SInitOrbit& orb1, const SInitOrbit& orb2, const TimeJD& beg, const TimeJD& end, TimeJD& result_jd, double& result_dist) const;
SInitOrbits filter(const SInitOrbits& orbits) const;
OrbBlocks m_blocks;
TimeJD m_beg_before_filter;
TimeJD m_end_before_filter;
//TimeJD m_beg_before_filter;
//TimeJD m_end_before_filter;
};
......@@ -62,7 +62,7 @@ void propagate_date(Propagator& prop, const TimeJD& date, IRes& result)
{
PhasePoint6D pos;
prop.Propagate(date, pos);
result = IRes(pos.CoordsVel, date, 0, DaysInterval(prop.GetInitOrbit().GetDate(), date), prop.GetInitOrbit().GetId());
result = IRes(pos.CoordsVel, date, DaysInterval(prop.GetInitOrbit().GetDate(), date), prop.GetInitOrbit().GetId());
}
void propagate_date(Propagator& prop, const TimeJD& date, SInitOrbit& result)
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment