16 #include "NumericalApproximation.h"
30 init(targetCode, observerCode,
false);
63 bool swapObserverTarget) {
64 init(targetCode, observerCode, swapObserverTarget);
80 const bool &swapObserverTarget) {
106 m_swapObserverTarget = swapObserverTarget;
110 if ( m_swapObserverTarget ) {
180 QString abcorr(correction);
181 abcorr.remove(QChar(
' '));
182 abcorr = abcorr.toUpper();
184 validAbcorr <<
"NONE" <<
"LT" <<
"LT+S" <<
"CN" <<
"CN+S"
185 <<
"XLT" <<
"XLT+S" <<
"XCN" <<
"XCN+S";
187 if (validAbcorr.indexOf(abcorr) >= 0) {
191 QString msg =
"Invalid abberation correction [" + correction +
"]";
290 QString msg =
"A SpicePosition cache has already been created";
294 if(startTime > endTime) {
295 QString msg =
"Argument startTime must be less than or equal to endTime";
299 if((startTime != endTime) && (size == 1)) {
300 QString msg =
"Cache size must be more than 1 if startTime endTime differ";
311 for(
int i = 0; i < size; i++) {
366 QString msg =
"A SpicePosition cache has already been created";
399 "Invalid value for CacheType keyword in the table "
406 for (
int r = 0; r < table.
Records(); r++) {
411 else if (rec.
Fields() == 4) {
415 QString msg =
"Expecting four or seven fields in the SpicePosition table";
419 std::vector<double> j2000Coord;
420 j2000Coord.push_back((
double)rec[0]);
421 j2000Coord.push_back((
double)rec[1]);
422 j2000Coord.push_back((
double)rec[2]);
427 std::vector<double> j2000Velocity;
428 j2000Velocity.push_back((
double)rec[3]);
429 j2000Velocity.push_back((
double)rec[4]);
430 j2000Velocity.push_back((
double)rec[5]);
440 std::vector<double> coeffX, coeffY, coeffZ;
442 for (
int r = 0; r < table.
Records() - 1; r++) {
448 coeffX.push_back((
double)rec[0]);
449 coeffY.push_back((
double)rec[1]);
450 coeffZ.push_back((
double)rec[2]);
454 double baseTime = (double)rec[0];
455 double timeScale = (double)rec[1];
456 double degree = (double)rec[2];
515 Table table(tableName, record);
519 for (
int i = 0; i < (int)
p_cache.size(); i++) {
520 record[inext++] =
p_cache[i][0];
521 record[inext++] =
p_cache[i][1];
522 record[inext++] =
p_cache[i][2];
549 record += spacecraftX;
550 record += spacecraftY;
551 record += spacecraftZ;
553 Table table(tableName, record);
555 for(
int cindex = 0; cindex <
p_degree + 1; cindex++) {
566 record[2] = (double) p_degree;
575 "Cannot create Table, no Cache is loaded.",
591 QString tabletype =
"";
594 tabletype =
"Linear";
597 tabletype =
"HermiteSpline";
600 tabletype =
"PolyFunction";
641 QString msg =
"Only cached positions can be returned as a line cache of positions and time";
645 return Cache(tableName);
668 QString msg =
"The SpicePosition has not yet been fit to a function";
687 for (std::vector<double>::size_type pos = 0; pos <
p_cacheTime.size(); pos++) {
727 if(cacheTimeSize == 1)
return Cache(tableName);
739 QString msg =
"A SpicePosition polynomial function has not been created yet";
767 double extremumXtime = -b1 / (2.*c1) +
p_baseTime;
769 if(extremumXtime < firstTime) {
770 extremumXtime = firstTime;
772 if(extremumXtime > lastTime) {
773 extremumXtime = lastTime;
778 double extremumYtime = -b2 / (2.*c2) +
p_baseTime;
780 if(extremumYtime < firstTime) {
781 extremumYtime = firstTime;
783 if(extremumYtime > lastTime) {
784 extremumYtime = lastTime;
789 double extremumZtime = -b3 / (2.*c3) +
p_baseTime;
791 if(extremumZtime < firstTime) {
792 extremumZtime = firstTime;
794 if(extremumZtime > lastTime) {
795 extremumZtime = lastTime;
819 for(
int i = 0; i < (int)
p_cacheTime.size(); i++) {
821 std::vector <double> time;
840 return Cache(tableName);
851 std::vector<double> XC, YC, ZC;
861 else if (
p_cache.size() == 2) {
880 std::vector<double> time;
903 for(
int cIndex = 0; cIndex < 3; cIndex++) {
905 slope[cIndex] = posline.
Slope();
908 XC.push_back(intercept[0]);
909 XC.push_back(slope[0]);
910 YC.push_back(intercept[1]);
911 YC.push_back(slope[1]);
912 ZC.push_back(intercept[2]);
913 ZC.push_back(slope[2]);
921 for(std::vector<double>::size_type pos = 0; pos <
p_cacheTime.size(); pos++) {
974 const std::vector<double>& YC,
975 const std::vector<double>& ZC,
1023 std::vector<double>& YC,
1024 std::vector<double>& ZC) {
1087 SpicePosition::PartialType partialVar,
int coeffIndex) {
1090 std::vector<double> coordinate(3, 0);
1093 int coordIndex = partialVar;
1128 SpicePosition::PartialType partialVar,
int coeffIndex) {
1131 std::vector<double> dvelocity(3, 0);
1134 int coordIndex = partialVar;
1137 double derivative = 0.;
1140 double Epsilon = 1.0e-15;
1141 if (fabs(time) <= Epsilon) time = 0.0;
1145 derivative = coeffIndex * pow(time, coeffIndex-1) /
p_timeScale;
1148 dvelocity[coordIndex] = derivative;
1164 double derivative = 1.;
1167 if(coeffIndex > 0 && coeffIndex <=
p_degree) {
1168 derivative = pow(time, coeffIndex);
1170 else if(coeffIndex == 0) {
1175 " exceeds degree of polynomial";
1191 QString msg =
"No velocity vector available";
1227 std::vector<double>::iterator pos;
1239 if(cacheIndex < 0) cacheIndex = 0;
1244 std::vector<double> p2 =
p_cache[cacheIndex+1];
1245 std::vector<double> p1 =
p_cache[cacheIndex];
1254 p_velocity[0] = (p2[0] - p1[0]) * mult + p1[0];
1255 p_velocity[1] = (p2[1] - p1[1]) * mult + p1[1];
1256 p_velocity[2] = (p2[2] - p1[2]) * mult + p1[2];
1292 vector<double> xvalues(
p_cache.size());
1293 vector<double> xhermiteYValues(
p_cache.size());
1294 vector<double> yhermiteYValues(
p_cache.size());
1295 vector<double> zhermiteYValues(
p_cache.size());
1297 vector<double> xhermiteVelocities(
p_cache.size());
1298 vector<double> yhermiteVelocities(
p_cache.size());
1299 vector<double> zhermiteVelocities(
p_cache.size());
1301 for(
unsigned int i = 0; i <
p_cache.size(); i++) {
1302 vector<double> &cache =
p_cache[i];
1306 xhermiteYValues[i] = cache[0];
1307 yhermiteYValues[i] = cache[1];
1308 zhermiteYValues[i] = cache[2];
1312 xhermiteVelocities[i] = cacheVelocity[0];
1313 yhermiteVelocities[i] = cacheVelocity[1];
1314 zhermiteVelocities[i] = cacheVelocity[2];
1374 std::vector<double> rtime;
1415 std::vector<double> hermiteVelocity =
p_velocity;
1418 for (
int index = 0; index < 3; index++) {
1477 "Source type is not Memcache, cannot convert.",
1490 vector <int> inputIndices;
1491 inputIndices.push_back(0);
1492 inputIndices.push_back(n / 2);
1493 inputIndices.push_back(n);
1496 vector <int> indexList =
HermiteIndices(tolerance, inputIndices);
1499 for(
int i = n; i >= 0; i--) {
1500 if(!binary_search(indexList.begin(), indexList.end(), i)) {
1533 unsigned int n = indexList.size();
1539 for(
unsigned int i = 0; i < indexList.size(); i++) {
1559 for(
unsigned int i = indexList.size() - 1; i > 0; i--) {
1565 for(
int line = indexList[i-1] + 1; line < indexList[i]; line++) {
1573 if(xerror > tolerance || yerror > tolerance || zerror > tolerance) {
1579 if(xerror < tolerance && yerror < tolerance && zerror < tolerance) {
1585 indexList.push_back((indexList[i] + indexList[i-1]) / 2);
1589 if(indexList.size() > n) {
1590 sort(indexList.begin(), indexList.end());
1650 for(
int icoef =
p_degree + 1; icoef <= degree; icoef++) {
1651 coefX.push_back(0.);
1652 coefY.push_back(0.);
1653 coefZ.push_back(0.);
1660 std::vector<double> coefX(degree + 1),
1664 for(
int icoef = 0; icoef <= degree; icoef++) {
1699 double cacheSlope = 0.0;
1732 double velocity = 0.;
1734 for (
int icoef=1; icoef <=
p_degree; icoef++) {
1757 double diffTime = timeEt -
p_et;
1758 std::vector<double> extrapPos(3, 0.);
1760 extrapPos[1] =
p_coordinate[1] + diffTime*p_velocity[1];
1761 extrapPos[2] =
p_coordinate[2] + diffTime*p_velocity[2];
1777 "Hermite coordinates only available for PolyFunctionOverHermiteConstant",
1786 return hermiteCoordinate;
1871 const QString &refFrame,
1872 const QString &abcorr,
1873 double state[6],
bool &hasVelocity,
1874 double &lightTime)
const {
1879 spkez_c((SpiceInt) target, (SpiceDouble) et, refFrame.toAscii().data(),
1880 abcorr.toAscii().data(), (SpiceInt) observer, state, &lightTime);
1885 SpiceBoolean spfailure = failed_c();
1888 hasVelocity =
false;
1890 spkezp_c((SpiceInt) target, (SpiceDouble) et, refFrame.toAscii().data(),
1891 abcorr.toAscii().data(), (SpiceInt) observer, state, &lightTime);
1922 const bool &hasVelocity) {
1928 if ( hasVelocity ) {
1942 if (m_swapObserverTarget) {