File: parametrizationTemporalPolynomial.cpp

package info (click to toggle)
groops 0%2Bgit20250907%2Bds-1
  • links: PTS, VCS
  • area: non-free
  • in suites: forky, sid
  • size: 11,140 kB
  • sloc: cpp: 135,607; fortran: 1,603; makefile: 20
file content (120 lines) | stat: -rwxr-xr-x 3,388 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
/***********************************************/
/**
* @file parametrizationTemporalPolynomial.cpp
*
* @brief Legendre Polynomials in time intervals.
* @see ParametrizationTemporal
*
* @author Torsten Mayer-Guerr
* @date 2015-05-31
*
*/
/***********************************************/

#include "base/import.h"
#include "base/legendrePolynomial.h"
#include "config/config.h"
#include "classes/timeSeries/timeSeries.h"
#include "classes/parametrizationTemporal/parametrizationTemporal.h"
#include "parametrizationTemporalPolynomial.h"

/***********************************************/

ParametrizationTemporalPolynomial::ParametrizationTemporalPolynomial(Config &config)
{
  try
  {
    TimeSeriesPtr timeSeries;
    readConfig(config, "polynomialDegree", degree,          Config::MUSTSET, "0", "polynomial degree");
    readConfig(config, "interval",         timeSeries,      Config::DEFAULT, "",  "intervals of polynomials");
    readConfig(config, "includeLastTime",  includeLastTime, Config::DEFAULT, "0", "");
    if(isCreateSchema(config)) return;

    times      = timeSeries->times();
    isInterval = (times.size() != 0);
    if(!isInterval)
      times = {Time(), date2time(2500,1,1)};
    idxStart   = 0;
    idxEnd = times.size()-1;
  }
  catch(std::exception &e)
  {
    GROOPS_RETHROW(e)
  }
}

/***********************************************/

Bool ParametrizationTemporalPolynomial::setInterval(const Time &timeStart, const Time &timeEnd, Bool estimatePerArc)
{
  try
  {
    const UInt idxStartOld = idxStart;
    const UInt idxEndOld   = idxEnd;

    if(estimatePerArc && !isInterval)
      times = {timeStart, timeEnd};

    idxStart = 0;
    while((idxStart+1<times.size()) && (timeStart>=times.at(idxStart+1)))
      idxStart++;
    idxEnd = idxStart;
    while((idxEnd<times.size()-1) && (timeEnd>times.at(idxEnd)))
      idxEnd++;

    return (idxStartOld != idxStart) || (idxEndOld != idxEnd);
  }
  catch(std::exception &e)
  {
    GROOPS_RETHROW(e)
  }
}

/***********************************************/

void ParametrizationTemporalPolynomial::parameterName(std::vector<ParameterName> &name) const
{
  try
  {
    for(UInt i=idxStart; i<idxEnd; i++)
    {
      name.push_back(ParameterName("", "", "", times.at(i), times.at(i+1))); // constant part
      for(UInt n=1; n<=degree; n++)
        name.push_back(ParameterName("", "", "legendrePolynomial.n"+n%"%i"s, times.at(i), times.at(i+1)));
    }
  }
  catch(std::exception &e)
  {
    GROOPS_RETHROW(e)
  }
}

/***********************************************/

void ParametrizationTemporalPolynomial::factors(const Time &time, UInt startIndex, std::vector<UInt> &index, std::vector<Double> &factor) const
{
  try
  {
    if((time < times.at(idxStart)) || (time > times.at(idxEnd)) || (!includeLastTime && (time == times.at(idxEnd))))
      return;

    // find index (interval)
    UInt idx = idxStart;
    while((idx+1 < idxEnd) && (time >= times.at(idx+1)))
      idx++;

    const Double t  = 2*(time-times.at(idx)).mjd()/(times.at(idx+1)-times.at(idx)).mjd()-1; // t is [-1,+1]
    const Vector Pn = LegendrePolynomial::compute(t, degree);
    for(UInt n=0; n<=degree; n++)
    {
      index.push_back( (degree+1)*(idx-idxStart)+n+startIndex );
      factor.push_back( Pn(n) );
    }
  }
  catch(std::exception &e)
  {
    GROOPS_RETHROW(e)
  }
}

/***********************************************/