Commit 73eac811 authored by Alexander Lapshin's avatar Alexander Lapshin

golden section search invert interval fix; + bisection search

parent f104ee11
...@@ -127,7 +127,7 @@ void GFCSRotator::GetGcsToJ20006dConverter(const TimeJD& date, FrameConverter6D& ...@@ -127,7 +127,7 @@ void GFCSRotator::GetGcsToJ20006dConverter(const TimeJD& date, FrameConverter6D&
gcsToJ2000.Converter3D.Rotation.mul(w, gcsToJ2000.VRotation); gcsToJ2000.Converter3D.Rotation.mul(w, gcsToJ2000.VRotation);
} }
void GFCSRotator::GetJ2000ToNze6dConverter(const TimeJD& date, Vect3 siteGcsPos, FrameConverter6D& J2000toNze) const void GFCSRotator::GetJ2000ToNze6dConverter(const TimeJD& date, const Vect3& siteGcsPos, FrameConverter6D& J2000toNze) const
{ {
NZEFrame nzeFrame; NZEFrame nzeFrame;
nzeFrame.Init(siteGcsPos); nzeFrame.Init(siteGcsPos);
...@@ -143,7 +143,7 @@ void GFCSRotator::GetJ2000ToNze6dConverter(const TimeJD& date, Vect3 siteGcsPos, ...@@ -143,7 +143,7 @@ void GFCSRotator::GetJ2000ToNze6dConverter(const TimeJD& date, Vect3 siteGcsPos,
J2000toNze.Invert(); J2000toNze.Invert();
} }
void GFCSRotator::GetJ2000ToLocalJ20006dConverter(const TimeJD& date, Vect3 siteGcsPos, FrameConverter6D& J2000toLocalJ2000) const void GFCSRotator::GetJ2000ToLocalJ20006dConverter(const TimeJD& date, const Vect3& siteGcsPos, FrameConverter6D& J2000toLocalJ2000) const
{ {
TransMtx CJ2000toGFCS, CGFCStoJ2000; TransMtx CJ2000toGFCS, CGFCStoJ2000;
Rotate(date, CJ2000toGFCS, CGFCStoJ2000); Rotate(date, CJ2000toGFCS, CGFCStoJ2000);
...@@ -163,14 +163,14 @@ void GFCSRotator::GetJ2000ToLocalJ20006dConverter(const TimeJD& date, Vect3 site ...@@ -163,14 +163,14 @@ void GFCSRotator::GetJ2000ToLocalJ20006dConverter(const TimeJD& date, Vect3 site
J2000toLocalJ2000.Invert(); J2000toLocalJ2000.Invert();
} }
void GFCSRotator::J2000ToNze(const TimeJD& date, const Vect6& j2000, Vect3 siteGcsPos, Vect6& nze) const void GFCSRotator::J2000ToNze(const TimeJD& date, const Vect6& j2000, const Vect3& siteGcsPos, Vect6& nze) const
{ {
FrameConverter6D J2000toNze; FrameConverter6D J2000toNze;
GetJ2000ToNze6dConverter(date, siteGcsPos, J2000toNze); GetJ2000ToNze6dConverter(date, siteGcsPos, J2000toNze);
J2000toNze.Convert(j2000, nze); J2000toNze.Convert(j2000, nze);
} }
void GFCSRotator::J2000ToLocalJ2000(const TimeJD& date, const Vect6& j2000, Vect3 siteGcsPos, Vect6& localJ2000) const void GFCSRotator::J2000ToLocalJ2000(const TimeJD& date, const Vect6& j2000, const Vect3& siteGcsPos, Vect6& localJ2000) const
{ {
FrameConverter6D J2000toLocalJ2000; FrameConverter6D J2000toLocalJ2000;
GetJ2000ToLocalJ20006dConverter(date, siteGcsPos, J2000toLocalJ2000); GetJ2000ToLocalJ20006dConverter(date, siteGcsPos, J2000toLocalJ2000);
......
...@@ -15,15 +15,15 @@ public: ...@@ -15,15 +15,15 @@ public:
void Rotate(const TimeJD& date, TransMtx& CJ2000toGFCS, TransMtx& CGFCStoJ2000) const; void Rotate(const TimeJD& date, TransMtx& CJ2000toGFCS, TransMtx& CGFCStoJ2000) const;
void J2000ToGCS(const PhasePoint6D& j2000, PhasePoint6D& gcs) const; void J2000ToGCS(const PhasePoint6D& j2000, PhasePoint6D& gcs) const;
void J2000ToGCS(const TimeJD& date, const Vect6& j2000, Vect6& gcs) const; void J2000ToGCS(const TimeJD& date, const Vect6& j2000, Vect6& gcs) const;
void J2000ToLocalJ2000(const TimeJD& date, const Vect6& j2000, Vect3 siteGcsPos, Vect6& localJ2000) const; void J2000ToLocalJ2000(const TimeJD& date, const Vect6& j2000, const Vect3& siteGcsPos, Vect6& localJ2000) const;
void J2000ToNze(const TimeJD& date, const Vect6& j2000, Vect3 siteGcsPos, Vect6& nze) const; void J2000ToNze(const TimeJD& date, const Vect6& j2000, const Vect3& siteGcsPos, Vect6& nze) const;
void GCSToJ2000(const PhasePoint6D& gcs, PhasePoint6D& j2000) const; void GCSToJ2000(const PhasePoint6D& gcs, PhasePoint6D& j2000) const;
void GCSToJ2000(const TimeJD& date, const Vect6& gcs, Vect6& j2000) const; void GCSToJ2000(const TimeJD& date, const Vect6& gcs, Vect6& j2000) const;
private: private:
void GetGcsToJ20006dConverter(const TimeJD& date, FrameConverter6D& gcsToJ2000) const; void GetGcsToJ20006dConverter(const TimeJD& date, FrameConverter6D& gcsToJ2000) const;
void GetJ2000ToGcs6dConverter(const TimeJD& date, FrameConverter6D& j2000ToGcs) const; void GetJ2000ToGcs6dConverter(const TimeJD& date, FrameConverter6D& j2000ToGcs) const;
void GetJ2000ToNze6dConverter(const TimeJD& date, Vect3 siteGcsPos, FrameConverter6D& J2000toNze) const; void GetJ2000ToNze6dConverter(const TimeJD& date, const Vect3& siteGcsPos, FrameConverter6D& J2000toNze) const;
void GetJ2000ToLocalJ20006dConverter(const TimeJD& date, Vect3 siteGcsPos, FrameConverter6D& J2000toLocalJ2000) const; void GetJ2000ToLocalJ20006dConverter(const TimeJD& date, const Vect3& siteGcsPos, FrameConverter6D& J2000toLocalJ2000) const;
TransMtx BJ2000toGFCS; TransMtx BJ2000toGFCS;
TimeJD m_baseDate; TimeJD m_baseDate;
}; };
...@@ -14,8 +14,9 @@ double golden_section_search( ...@@ -14,8 +14,9 @@ double golden_section_search(
const double invphi = (std::sqrt(5.0) - 1) / 2.0; const double invphi = (std::sqrt(5.0) - 1) / 2.0;
const double invphi2 = (3 - std::sqrt(5.0)) / 2.0; const double invphi2 = (3 - std::sqrt(5.0)) / 2.0;
a = std::min(a, b); if (a > b) {
b = std::max(a, b); std::swap(a, b);
}
double h = b - a; double h = b - a;
...@@ -62,13 +63,14 @@ inline double golden_section_search_stdf( ...@@ -62,13 +63,14 @@ inline double golden_section_search_stdf(
double b, double b,
const double time_accuracy, const double time_accuracy,
const std::function<double(const double&)>& f const std::function<double(const double&)>& f
) )
{ {
const double invphi = (std::sqrt(5.0) - 1) / 2.0; const double invphi = (std::sqrt(5.0) - 1) / 2.0;
const double invphi2 = (3 - std::sqrt(5.0)) / 2.0; const double invphi2 = (3 - std::sqrt(5.0)) / 2.0;
a = std::min(a, b); if(a > b) {
b = std::max(a, b); std::swap(a, b);
}
double h = b - a; double h = b - a;
...@@ -110,3 +112,29 @@ inline double golden_section_search_stdf( ...@@ -110,3 +112,29 @@ inline double golden_section_search_stdf(
} }
} }
inline double bisection_search_stdf(
double a,
double b,
const double accuracy,
const std::function<bool(const double&)>& f
)
{
//f(a) should not be equal f(b) at the first step
if (a > b) {
std::swap(a, b);
}
while (fabs(b - a) > accuracy) {
const double c = (a + b) / 2.0;
if (f(a) != f(c)) {
b = c;
}
else {
a = c;
}
}
return (a + b) / 2.0;
}
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