diff --git a/SPG4Tool/PassPredict.cpp b/SPG4Tool/PassPredict.cpp new file mode 100644 index 0000000..a1148a9 --- /dev/null +++ b/SPG4Tool/PassPredict.cpp @@ -0,0 +1,402 @@ +#include +#include +#include +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +#include "SPG4Function.h" + +struct PassDetails +{ + libsgp4::DateTime aos; + libsgp4::DateTime los; + double max_elevation; +}; + +double FindMaxElevation( + const libsgp4::CoordGeodetic& user_geo, + libsgp4::SGP4& sgp4, + const libsgp4::DateTime& aos, + const libsgp4::DateTime& los) +{ + libsgp4::Observer obs(user_geo); + + bool running; + + double time_step = (los - aos).TotalSeconds() / 9.0; + libsgp4::DateTime current_time(aos); //! current time + libsgp4::DateTime time1(aos); //! start time of search period + libsgp4::DateTime time2(los); //! end time of search period + double max_elevation; //! max elevation + + running = true; + + do + { + running = true; + max_elevation = -99999999999999.0; + while (running && current_time < time2) + { + /* + * find position + */ + libsgp4::Eci eci = sgp4.FindPosition(current_time); + libsgp4::CoordTopocentric topo = obs.GetLookAngle(eci); + + if (topo.elevation > max_elevation) + { + /* + * still going up + */ + max_elevation = topo.elevation; + /* + * move time along + */ + current_time = current_time.AddSeconds(time_step); + if (current_time > time2) + { + /* + * dont go past end time + */ + current_time = time2; + } + } + else + { + /* + * stop + */ + running = false; + } + } + + /* + * make start time to 2 time steps back + */ + time1 = current_time.AddSeconds(-2.0 * time_step); + /* + * make end time to current time + */ + time2 = current_time; + /* + * current time to start time + */ + current_time = time1; + /* + * recalculate time step + */ + time_step = (time2 - time1).TotalSeconds() / 9.0; + } while (time_step > 1.0); + + return max_elevation; +} + +libsgp4::DateTime FindCrossingPoint( + const libsgp4::CoordGeodetic& user_geo, + libsgp4::SGP4& sgp4, + const libsgp4::DateTime& initial_time1, + const libsgp4::DateTime& initial_time2, + bool finding_aos) +{ + libsgp4::Observer obs(user_geo); + + bool running; + int cnt; + + libsgp4::DateTime time1(initial_time1); + libsgp4::DateTime time2(initial_time2); + libsgp4::DateTime middle_time; + + running = true; + cnt = 0; + while (running && cnt++ < 16) + { + middle_time = time1.AddSeconds((time2 - time1).TotalSeconds() / 2.0); + /* + * calculate satellite position + */ + libsgp4::Eci eci = sgp4.FindPosition(middle_time); + libsgp4::CoordTopocentric topo = obs.GetLookAngle(eci); + + if (topo.elevation > 0.0) + { + /* + * satellite above horizon + */ + if (finding_aos) + { + time2 = middle_time; + } + else + { + time1 = middle_time; + } + } + else + { + if (finding_aos) + { + time1 = middle_time; + } + else + { + time2 = middle_time; + } + } + + if ((time2 - time1).TotalSeconds() < 1.0) + { + /* + * two times are within a second, stop + */ + running = false; + /* + * remove microseconds + */ + int us = middle_time.Microsecond(); + middle_time = middle_time.AddMicroseconds(-us); + /* + * step back into the pass by 1 second + */ + middle_time = middle_time.AddSeconds(finding_aos ? 1 : -1); + } + } + + /* + * go back/forward 1second until below the horizon + */ + running = true; + cnt = 0; + while (running && cnt++ < 6) + { + libsgp4::Eci eci = sgp4.FindPosition(middle_time); + libsgp4::CoordTopocentric topo = obs.GetLookAngle(eci); + if (topo.elevation > 0) + { + middle_time = middle_time.AddSeconds(finding_aos ? -1 : 1); + } + else + { + running = false; + } + } + + return middle_time; +} + +std::list GeneratePassList( + const libsgp4::CoordGeodetic& user_geo, + libsgp4::SGP4& sgp4, + const libsgp4::DateTime& start_time, + const libsgp4::DateTime& end_time, + const int time_step) +{ + std::list pass_list; + + libsgp4::Observer obs(user_geo); + + libsgp4::DateTime aos_time; + libsgp4::DateTime los_time; + + bool found_aos = false; + + libsgp4::DateTime previous_time(start_time); + libsgp4::DateTime current_time(start_time); + + while (current_time < end_time) + { + bool end_of_pass = false; + + /* + * calculate satellite position + */ + libsgp4::Eci eci = sgp4.FindPosition(current_time); + libsgp4::CoordTopocentric topo = obs.GetLookAngle(eci); + + if (!found_aos && topo.elevation > 0.0) + { + /* + * aos hasnt occured yet, but the satellite is now above horizon + * this must have occured within the last time_step + */ + if (start_time == current_time) + { + /* + * satellite was already above the horizon at the start, + * so use the start time + */ + aos_time = start_time; + } + else + { + /* + * find the point at which the satellite crossed the horizon + */ + aos_time = FindCrossingPoint( + user_geo, + sgp4, + previous_time, + current_time, + true); + } + found_aos = true; + } + else if (found_aos && topo.elevation < 0.0) + { + found_aos = false; + /* + * end of pass, so move along more than time_step + */ + end_of_pass = true; + /* + * already have the aos, but now the satellite is below the horizon, + * so find the los + */ + los_time = FindCrossingPoint( + user_geo, + sgp4, + previous_time, + current_time, + false); + + struct PassDetails pd; + pd.aos = aos_time; + pd.los = los_time; + pd.max_elevation = FindMaxElevation( + user_geo, + sgp4, + aos_time, + los_time); + + pass_list.push_back(pd); + } + + /* + * save current time + */ + previous_time = current_time; + + if (end_of_pass) + { + /* + * at the end of the pass move the time along by 30mins + */ + current_time = current_time + libsgp4::TimeSpan(0, 30, 0); + } + else + { + /* + * move the time along by the time step value + */ + current_time = current_time + libsgp4::TimeSpan(0, 0, time_step); + } + + if (current_time > end_time) + { + /* + * dont go past end time + */ + current_time = end_time; + } + }; + + if (found_aos) + { + /* + * satellite still above horizon at end of search period, so use end + * time as los + */ + struct PassDetails pd; + pd.aos = aos_time; + pd.los = end_time; + pd.max_elevation = FindMaxElevation(user_geo, sgp4, aos_time, end_time); + pass_list.push_back(pd); + } + + return pass_list; +} + + +bool PassPredict(std::vector& aos, std::vector& los, std::string line1, std::string line2, double startTime, double predictDayLen, double lon, double lat, double ati) { + + libsgp4::CoordGeodetic geo(lat, lon, ati); + libsgp4::Tle tle("GALILEO-PFM (GSAT0101) ", + line1, + line2); + libsgp4::SGP4 sgp4(tle); + + std::vector result; + + // 将 double 时间戳转为 QDateTime + + QDateTime dateTime = QDateTime::fromSecsSinceEpoch(static_cast(startTime*1000)); + + + libsgp4::DateTime start_date( + dateTime.date().year(), + dateTime.date().month(), + dateTime.date().day(), + dateTime.time().hour(), + dateTime.time().minute(), + dateTime.time().second(), + dateTime.time().msec() * 1000); + + libsgp4::DateTime end_date(start_date.AddDays(predictDayLen)); + + std::list pass_list; + pass_list = GeneratePassList(geo, sgp4, start_date, end_date, 180); + + if (pass_list.begin() == pass_list.end()) + { + return false; + } + else + { + std::stringstream ss; + ss << std::right << std::setprecision(1) << std::fixed; + std::list::const_iterator itr = pass_list.begin(); + do + { + ss << "AOS: " << itr->aos + << ", LOS: " << itr->los + << ", Max El: " << std::setw(4) << libsgp4::Util::RadiansToDegrees(itr->max_elevation) + << ", Duration: " << (itr->los - itr->aos) + << std::endl; + + + // 构造AOS时间 + QDate aosDate(itr->aos.Year(), itr->aos.Month(), itr->aos.Day()); + QTime aosTimePart(itr->aos.Hour(), itr->aos.Minute(), itr->aos.Second(), itr->aos.Microsecond()); + QDateTime aosDateTime(aosDate, aosTimePart, Qt::UTC); // + + // 构造LOS时间 + QDate losDate(itr->los.Year(), itr->los.Month(), itr->los.Day()); + QTime losTimePart(itr->los.Hour(), itr->los.Minute(), itr->los.Second(), itr->los.Microsecond()); + QDateTime losDateTime(losDate, losTimePart, Qt::UTC); // + + aos.push_back(aosDateTime); + los.push_back(losDateTime); + + + + } while (++itr != pass_list.end()); + std::cout << ss.str(); + + return true; + + } + + return false; + +} + + + \ No newline at end of file diff --git a/SPG4Tool/SPG4Function.h b/SPG4Tool/SPG4Function.h index c9c8b81..4e99725 100644 --- a/SPG4Tool/SPG4Function.h +++ b/SPG4Tool/SPG4Function.h @@ -4,10 +4,11 @@ #include "SPG4Tool_global.h" #include "BaseConstVariable.h" #include - +#include //void RunTle(libsgp4::Tle tle, double start, double end, double inc); + SPG4TOOL_EXPORT double parseTLETimeOffset(const std::string& tle); /// /// 根据两行根数生成坐标与速度,请注意返回的单位是 km @@ -22,4 +23,7 @@ SPG4TOOL_EXPORT double parseTLETimeOffset(const std::string& tle); /// /// SPG4TOOL_EXPORT std::vector getGPSPoints(std::string line1, std::string line2, double start, double end, double inc, bool printfinfoflag = false, bool running = true, bool first_run = true); +SPG4TOOL_EXPORT bool PassPredict(std::vector& aos, std::vector& los, std::string line1, std::string line2, double startTime, double predictDayLen, double lon, double lat, double ati); + + #endif diff --git a/SPG4Tool/SPG4Tool.h b/SPG4Tool/SPG4Tool.h index 40103f3..cf443e1 100644 --- a/SPG4Tool/SPG4Tool.h +++ b/SPG4Tool/SPG4Tool.h @@ -8,6 +8,11 @@ #include #include #include "BaseConstVariable.h" - +#include +#include +#include +#include +#include + std::vector RunTle(libsgp4::Tle tle, double start, double end, double inc, bool printfinfoflag, bool running = true, bool first_run = true); diff --git a/SPG4Tool/SPG4Tool.vcxproj b/SPG4Tool/SPG4Tool.vcxproj index 878cdae..1b5bd90 100644 --- a/SPG4Tool/SPG4Tool.vcxproj +++ b/SPG4Tool/SPG4Tool.vcxproj @@ -115,6 +115,7 @@ + diff --git a/SPG4Tool/SPG4Tool.vcxproj.filters b/SPG4Tool/SPG4Tool.vcxproj.filters index e658568..d2d7c1f 100644 --- a/SPG4Tool/SPG4Tool.vcxproj.filters +++ b/SPG4Tool/SPG4Tool.vcxproj.filters @@ -80,6 +80,9 @@ Source Files + + Source Files +