From 012ac08252f1d939a9879029475f646658bec43b Mon Sep 17 00:00:00 2001 From: Fan Dongwei Date: Fri, 19 Dec 2025 21:36:16 +0800 Subject: [PATCH 1/3] =?UTF-8?q?out=5Fpath=5Fparam.wokr=5Fdir=E6=9B=B4?= =?UTF-8?q?=E6=AD=A3=E4=B8=BAwork=5Fdir=EF=BC=9BSurveySim=E7=B1=BB?= =?UTF-8?q?=E5=A2=9E=E5=8A=A0=E4=B8=80=E4=B8=AAorbit=5Fdata=E5=88=9D?= =?UTF-8?q?=E5=A7=8B=E5=8C=96=E5=8F=82=E6=95=B0=EF=BC=8C=E4=BB=A5=E5=87=8F?= =?UTF-8?q?=E5=B0=91=E6=98=9F=E5=8E=86=E5=88=9D=E5=A7=8B=E5=8C=96?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- config/config_surveysim.yaml | 2 +- config/config_surveysim_fdw.yaml | 2 +- run_survey_sim.py | 2 +- survey_sim/strategies/SurveySim_MSC_mpi.py | 13 +++++++------ 4 files changed, 10 insertions(+), 9 deletions(-) diff --git a/config/config_surveysim.yaml b/config/config_surveysim.yaml index bfceb9d..b803be5 100644 --- a/config/config_surveysim.yaml +++ b/config/config_surveysim.yaml @@ -1,5 +1,5 @@ out_path_param: - wokr_dir: /Users/zhangxin/Work/SurveyPlan/CSST_Survey # define the folder of output file + work_dir: /Users/zhangxin/Work/SurveyPlan/CSST_Survey # define the folder of output file out_filename_prefix: test_2461876_1 # output file prefix, 3 files: prefix_plan.json, prefix_info.dat, prefix_log.log # input_survey_data ---- whether to use input survey data; if yes, set the file name (including path), otherwise set null diff --git a/config/config_surveysim_fdw.yaml b/config/config_surveysim_fdw.yaml index ca4fc34..bdb16fe 100644 --- a/config/config_surveysim_fdw.yaml +++ b/config/config_surveysim_fdw.yaml @@ -1,5 +1,5 @@ out_path_param: - wokr_dir: D:/Projects/CSST/dev/csst-planning/upload/sim/output # define the folder of output file + work_dir: D:/Projects/CSST/dev/csst-planning/upload/sim/output # define the folder of output file out_filename_prefix: test_2461876_1 # output file prefix, 3 files: prefix_plan.json, prefix_info.dat, prefix_log.log # input_survey_data ---- whether to use input survey data; if yes, set the file name (including path), otherwise set null diff --git a/run_survey_sim.py b/run_survey_sim.py index cb11b07..6727fe1 100644 --- a/run_survey_sim.py +++ b/run_survey_sim.py @@ -158,6 +158,6 @@ if __name__ == '__main__': if config_data.get("input_survey_data") is not None: input_sim_data = config_data.get("input_survey_data") sim_from_input(jsonfn=input_sim_data, - outDir=config_data["out_path_param"].get("wokr_dir"), info_out_filename_pre=config_data["out_path_param"].get("out_filename_prefix")) + outDir=config_data["out_path_param"].get("work_dir"), info_out_filename_pre=config_data["out_path_param"].get("out_filename_prefix")) else: run_survey_sim(config=config_data) diff --git a/survey_sim/strategies/SurveySim_MSC_mpi.py b/survey_sim/strategies/SurveySim_MSC_mpi.py index 11f889f..c5e885c 100644 --- a/survey_sim/strategies/SurveySim_MSC_mpi.py +++ b/survey_sim/strategies/SurveySim_MSC_mpi.py @@ -267,13 +267,13 @@ def calulateSearchPointWeight(skymap, sky_attr_tab, idx, sun_grcs_car, sun_ecl_ class SurveySim(object): - def __init__(self, time_start=2459767.0, time_end=2463416.0, bettery_level=97200., cmgCons=None, dir=None, outFn=None, obsStartNum=1, config=None, beta_flag=True): + def __init__(self, time_start=2459767.0, time_end=2463416.0, bettery_level=97200., cmgCons=None, dir=None, outFn=None, obsStartNum=1, config=None, beta_flag=True, orbit_data=None): if config is not None: self.plan_startTime = config["sim_param"].get("start_time") self.plan_endTime = config["sim_param"].get("end_time") self.bettery_level = config["sim_param"].get("bettery_level") - self.outDir = config["out_path_param"].get("wokr_dir") + self.outDir = config["out_path_param"].get("work_dir") self.info_out_filename_pre = config["out_path_param"].get( "out_filename_prefix") self.config = config @@ -298,10 +298,11 @@ class SurveySim(object): self.initSkyMap(surveyCons=self.surveyCons) # print("memory analysis 2 : %f" % ( # (psutil.Process(os.getpid()).memory_info().rss / 1024 / 1024 / 1024))) - tmpData = loadSatOrbitDat() - ids1 = tmpData[:, 0] > self.plan_startTime - 0.2 - # ids2 = tmpData[:, 0] < self.plan_endTime + 2 - self.orbitDat = tmpData[ids1] + if orbit_data is None: + orbit_data = loadSatOrbitDat() + ids1 = orbit_data[:, 0] > self.plan_startTime - 0.2 + # ids2 = orbit_data[:, 0] < self.plan_endTime + 2 + self.orbitDat = orbit_data[ids1] # self.orbitDat = loadSatOrbitDat() # print("memory analysis 3 : %f" % ( # (psutil.Process(os.getpid()).memory_info().rss / 1024 / 1024 / 1024))) -- GitLab From c47ea711e92747dcc17cc431749e5666e3d888b6 Mon Sep 17 00:00:00 2001 From: Fan Dongwei Date: Sat, 20 Dec 2025 14:50:29 +0800 Subject: [PATCH 2/3] =?UTF-8?q?moonAgle=E5=BA=94=E4=B8=BAmoonAngle?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- survey_sim/constraints/surveyConstraint.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/survey_sim/constraints/surveyConstraint.py b/survey_sim/constraints/surveyConstraint.py index 29c5452..9867b2f 100644 --- a/survey_sim/constraints/surveyConstraint.py +++ b/survey_sim/constraints/surveyConstraint.py @@ -5,13 +5,13 @@ class surveyConstraint(object): def __init__( self, sunAngle=50, - moonAgle=40, + moonAngle=40, earthLLimbAngle_limit=80, earthLLimbAngle_plus=70, earthDLimgAngle=30, ): self.sun_los_angle_cos = math.cos(sunAngle * math.pi / 180.0) - self.moon_los_angle_cos = math.cos(moonAgle * math.pi / 180.0) + self.moon_los_angle_cos = math.cos(moonAngle * math.pi / 180.0) # 约为70°是根据400km轨道，地球平均半径计算的，切线和卫星矢量夹角 self.satZenith_angle_light_min_cos = math.cos( (180 - 70 - earthLLimbAngle_limit) * math.pi / 180.0 -- GitLab From 9dd22166383df944591be21c3dfe16a69e66c43c Mon Sep 17 00:00:00 2001 From: Fan Dongwei Date: Sat, 20 Dec 2025 22:07:47 +0800 Subject: [PATCH 3/3] =?UTF-8?q?=E4=BD=BF=E7=94=A8=E4=BA=8C=E5=88=86?= =?UTF-8?q?=E6=9F=A5=E6=89=BE=E6=89=BE=E6=97=B6=E9=97=B4?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- survey_sim/config/infoOutput.py | 2 +- survey_sim/satOrbit/locate_sat.py | 198 +++++++++++++++++++----------- 2 files changed, 129 insertions(+), 71 deletions(-) diff --git a/survey_sim/config/infoOutput.py b/survey_sim/config/infoOutput.py index 07b5268..fbcd418 100644 --- a/survey_sim/config/infoOutput.py +++ b/survey_sim/config/infoOutput.py @@ -183,7 +183,7 @@ class InfoOutput(object): data['simresult']['is_in_sun'] = isInSunSide data['simresult']['cmg'] = cmg data['simresult']['charge'] = energy - data['simresult']['skyid'] = sky_ID + data['simresult']['skyid'] = int(sky_ID) # data['simresult']['board_angle'] = d1['boardAngle'] data['simresult']['quad_0'] = sat_q[0] data['simresult']['quad_i'] = sat_q[1] diff --git a/survey_sim/satOrbit/locate_sat.py b/survey_sim/satOrbit/locate_sat.py index aa765c7..088eb52 100644 --- a/survey_sim/satOrbit/locate_sat.py +++ b/survey_sim/satOrbit/locate_sat.py @@ -68,6 +68,59 @@ def loadSatOrbitDat(datDir=None): return oData +def findTimeIndex(time, time_col, startId=0): + """ + 查找时间在轨道数据中的索引位置（使用二分查找优化，并利用 startId） + + 参数: + time: 目标时间 + time_col: 时间列数组 + startId: 起始搜索索引（用于优化连续调用） + + 返回: + int: 找到的索引，如果未找到返回 -1 + """ + orbDataLen = len(time_col) + + if orbDataLen < 2: + return -1 + + # 如果 startId 有效，先检查 startId 附近是否包含目标时间 + if startId > 0 and startId < orbDataLen - 1: + # 检查 startId 附近的时间范围（前后各检查一定范围） + search_range = 100 # 检查范围 + start_check = max(0, startId - search_range) + end_check = min(orbDataLen - 1, startId + search_range) + + # 检查 startId 附近的时间是否包含目标时间 + t_start = time_col[start_check].item() + t_end = time_col[end_check].item() + + if t_start <= time <= t_end: + # 在 startId 附近使用线性搜索（因为范围小，线性搜索更快） + for i in range(start_check, end_check): + if i >= orbDataLen - 1: + break + t1 = time_col[i].item() + t2 = time_col[i + 1].item() + if time >= t1 and time < t2: + return i + + # 如果 startId 附近没找到，使用全局二分查找 + # 使用 np.searchsorted 进行二分查找 + idx = np.searchsorted(time_col, time, side='right') + + # 检查是否找到有效的时间区间 + if idx > 0 and idx < orbDataLen: + i = idx - 1 + t1 = time_col[i].item() + t2 = time_col[i + 1].item() + if time >= t1 and time < t2: + return i + + return -1 + + def locateSat(time=2459766.0, OrbitData=None, startId=0, orbDataLen=-1): if orbDataLen <= 0 and OrbitData is not None: @@ -75,78 +128,83 @@ def locateSat(time=2459766.0, OrbitData=None, startId=0, orbDataLen=-1): satPos = np.zeros(3, dtype=np.float64) satVel = np.zeros(3, dtype=np.float64) nextSid = startId - for i in np.arange(startId, orbDataLen - 1, 1): - t1 = OrbitData[i, 0] - t2 = OrbitData[i + 1, 0] - if time >= t1 and time < t2: - if (t2 - t1) > 130.0 / 86400.0: # should be 120s, for error, set 120+10s - print("Warning: time interval is not 120s!!!!!!!!!") - # return satPos-1, 0 - # break - x1 = OrbitData[i, 1] - y1 = OrbitData[i, 2] - z1 = OrbitData[i, 3] - x2 = OrbitData[i + 1, 1] - y2 = OrbitData[i + 1, 2] - z2 = OrbitData[i + 1, 3] - - line1 = np.sqrt(x1 * x1 + y1 * y1 + z1 * z1) - line2 = np.sqrt(x2 * x2 + y2 * y2 + z2 * z2) - theta = np.arccos((x1 * x2 + y1 * y2 + z1 * z2) / (line1 * line2)) - theta1 = (time - t1) / (t2 - t1) * theta - theta2 = theta - theta1 - line = (t2 - time) / (t2 - t1) * line1 + \ - (time - t1) / (t2 - t1) * line2 - # // double ef = sin(theta2)/sin(theta1); - - x0 = np.sin(theta2) * x1 / line1 + np.sin(theta1) * x2 / line2 - y0 = np.sin(theta2) * y1 / line1 + np.sin(theta1) * y2 / line2 - z0 = np.sin(theta2) * z1 / line1 + np.sin(theta1) * z2 / line2 - line_ = np.sqrt(x0 * x0 + y0 * y0 + z0 * z0) - - satPos[0] = x0 * line / line_ - satPos[1] = y0 * line / line_ - satPos[2] = z0 * line / line_ - - xdot1 = OrbitData[i, 4] - ydot1 = OrbitData[i, 5] - zdot1 = OrbitData[i, 6] - - xdot2 = OrbitData[i + 1, 4] - ydot2 = OrbitData[i + 1, 5] - zdot2 = OrbitData[i + 1, 6] - - ldot1 = np.sqrt(xdot1 * xdot1 + ydot1 * ydot1 + zdot1 * zdot1) - ldot2 = np.sqrt(xdot2 * xdot2 + ydot2 * ydot2 + zdot2 * zdot2) - thetadot = math.acos( - (xdot1 * xdot2 + ydot1 * ydot2 + zdot1 * zdot2) / (ldot1 * ldot2) - ) - - ldot = (t2 - time) / (t2 - t1) * ldot1 + \ - (time - t1) / (t2 - t1) * ldot2 - - if time != t1: - thetadot1 = (time - t1) / (t2 - t1) * thetadot - thetadot2 = thetadot - thetadot1 - efdot = math.sin(thetadot2) / math.sin(thetadot1) - - xdot0 = efdot * xdot1 / ldot1 + xdot2 / ldot2 - ydot0 = efdot * ydot1 / ldot1 + ydot2 / ldot2 - zdot0 = efdot * zdot1 / ldot1 + zdot2 / ldot2 - else: - xdot0 = xdot1 / ldot1 - ydot0 = ydot1 / ldot1 - zdot0 = zdot1 / ldot1 - - ldot_ = np.sqrt(xdot0 * xdot0 + ydot0 * ydot0 + zdot0 * zdot0) - - satVel[0] = xdot0 * ldot / ldot_ - satVel[1] = ydot0 * ldot / ldot_ - satVel[2] = zdot0 * ldot / ldot_ + # 使用独立的查找函数 + time_col = OrbitData[:, 0] + i = findTimeIndex(time, time_col, startId) + + # 如果找到匹配的时间区间，进行插值计算 + if i >= 0: + t1 = OrbitData[i, 0].item() + t2 = OrbitData[i + 1, 0].item() + if (t2 - t1) > 130.0 / 86400.0: # should be 120s, for error, set 120+10s + print("Warning: time interval is not 120s!!!!!!!!!") + # return satPos-1, 0 + # break + + x1 = OrbitData[i, 1] + y1 = OrbitData[i, 2] + z1 = OrbitData[i, 3] + x2 = OrbitData[i + 1, 1] + y2 = OrbitData[i + 1, 2] + z2 = OrbitData[i + 1, 3] + + line1 = np.sqrt(x1 * x1 + y1 * y1 + z1 * z1) + line2 = np.sqrt(x2 * x2 + y2 * y2 + z2 * z2) + theta = np.arccos((x1 * x2 + y1 * y2 + z1 * z2) / (line1 * line2)) + theta1 = (time - t1) / (t2 - t1) * theta + theta2 = theta - theta1 + line = (t2 - time) / (t2 - t1) * line1 + \ + (time - t1) / (t2 - t1) * line2 + # // double ef = sin(theta2)/sin(theta1); + + x0 = np.sin(theta2) * x1 / line1 + np.sin(theta1) * x2 / line2 + y0 = np.sin(theta2) * y1 / line1 + np.sin(theta1) * y2 / line2 + z0 = np.sin(theta2) * z1 / line1 + np.sin(theta1) * z2 / line2 + line_ = np.sqrt(x0 * x0 + y0 * y0 + z0 * z0) + + satPos[0] = x0 * line / line_ + satPos[1] = y0 * line / line_ + satPos[2] = z0 * line / line_ + + xdot1 = OrbitData[i, 4] + ydot1 = OrbitData[i, 5] + zdot1 = OrbitData[i, 6] + + xdot2 = OrbitData[i + 1, 4] + ydot2 = OrbitData[i + 1, 5] + zdot2 = OrbitData[i + 1, 6] + + ldot1 = np.sqrt(xdot1 * xdot1 + ydot1 * ydot1 + zdot1 * zdot1) + ldot2 = np.sqrt(xdot2 * xdot2 + ydot2 * ydot2 + zdot2 * zdot2) + thetadot = math.acos( + (xdot1 * xdot2 + ydot1 * ydot2 + zdot1 * zdot2) / (ldot1 * ldot2) + ) + + ldot = (t2 - time) / (t2 - t1) * ldot1 + \ + (time - t1) / (t2 - t1) * ldot2 + + if time != t1: + thetadot1 = (time - t1) / (t2 - t1) * thetadot + thetadot2 = thetadot - thetadot1 + efdot = math.sin(thetadot2) / math.sin(thetadot1) + + xdot0 = efdot * xdot1 / ldot1 + xdot2 / ldot2 + ydot0 = efdot * ydot1 / ldot1 + ydot2 / ldot2 + zdot0 = efdot * zdot1 / ldot1 + zdot2 / ldot2 + else: + xdot0 = xdot1 / ldot1 + ydot0 = ydot1 / ldot1 + zdot0 = zdot1 / ldot1 + + ldot_ = np.sqrt(xdot0 * xdot0 + ydot0 * ydot0 + zdot0 * zdot0) + + satVel[0] = xdot0 * ldot / ldot_ + satVel[1] = ydot0 * ldot / ldot_ + satVel[2] = zdot0 * ldot / ldot_ + + nextSid = i - nextSid = i - break return satPos, satVel, nextSid -- GitLab