Commit 329533c1 authored by Zhang Xin's avatar Zhang Xin
Browse files

add time segment and timepoint observed for input pointinglist;add calculat... 

add time segment and timepoint observed for input pointinglist;add calculat area pyfile; optimize survey model
parent f6bbdbae

setup.py

+6 −7
Original line number Diff line number Diff line 
"""

'''

Author: Zhang Xin zhangx@bao.ac.cn

Date: 2024-07-02 13:55:19

Date: 2024-09-26 19:15:51

LastEditors: Zhang Xin zhangx@bao.ac.cn

LastEditTime: 2024-09-26 15:43:32

LastEditTime: 2025-01-12 17:28:10

FilePath: /CSST_Survey/setup.py

Description: 这是默认设置,请设置`customMade`, 打开koroFileHeader查看配置 进行设置: https://github.com/OBKoro1/koro1FileHeader/wiki/%E9%85%8D%E7%BD%AE

"""

'''



from setuptools import setup, find_packages





setup(

    name="csst_survey_sim",

    version="0.0.1",
@@ -19,7 +17,7 @@ setup( 
    license="NAOC",

    packages=find_packages(),

    # python_requires=">=3.9",

    install_requires=["numpy", "jax"],

    install_requires=["numpy", "jax",  "astropy==5.3.4"],

    # requires=['numpy', 'scipy', 'astropy', 'drizzlepac', 'stwcs'],

    # long_description=read('README.rst'),

    # classifiers=[
@@ -33,5 +31,6 @@ setup( 
        "survey_sim.strategies.data": ["*"],

        "survey_sim.ephemeris": ["data/**/*"],

        "survey_sim.satOrbit": ["data/**/*"],

        "survey_sim.surveyStat": ["data/*"]

    },

)

survey_sim/config/infoOutput.py

+9 −5
Original line number Diff line number Diff line
@@ -2,7 +2,7 @@ 
Author: Zhang Xin zhangx@bao.ac.cn

Date: 2024-11-08 15:12:55

LastEditors: Zhang Xin zhangx@bao.ac.cn

LastEditTime: 2024-11-11 10:01:35

LastEditTime: 2025-01-14 09:38:41

FilePath: /CSST_Survey/survey_sim/config/infooutput.py

Description: 这是默认设置,请设置`customMade`, 打开koroFileHeader查看配置 进行设置: https://github.com/OBKoro1/koro1FileHeader/wiki/%E9%85%8D%E7%BD%AE

'''
@@ -26,8 +26,8 @@ class InfoOutput(object): 
        fh.setFormatter(formatter)

        self.logger.addHandler(fh)



        hdr1 = "# JDTime lonitude(ecliptic) latitude(ecliptic) RA Dec sun_x sun_y sun_z moon_x moon_y moon_z sat_x sat_y sat_z sat_vel_x sat_vel_y sat_vel_z isInDeep "

        fmt1 = "%15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %4d"

        hdr1 = "# JDTime lonitude(ecliptic) latitude(ecliptic) RA Dec sun_x sun_y sun_z moon_x moon_y moon_z sat_x sat_y sat_z sat_vel_x sat_vel_y sat_vel_z isInDeep area_wide area_deep isInSunSide expTime transAngle"

        fmt1 = "%15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %15.8f %4d  %15.8f  %15.8f  %15.8f  %15.8f  %15.8f"

        self.hdr = hdr1

        self.fmt = fmt1
@@ -53,10 +53,14 @@ class InfoOutput(object): 
            self.hdr += "\n"

        self.outInfo.write(self.hdr)



    def outInfo_add_obj(self, jdTime=2459766., p_lon_ecl=0., p_lat_ecl=0., p_ra=0, p_dec=0., sun=[0, 0, 0], moon=[0, 0, 0], sat=[0, 0, 0], sat_vel=[0, 0, 0], isInDeep=0):

    def close_outFile(self):

        self.outInfo.close()

        # self.logger.close()



    def outInfo_add_obj(self, jdTime=2459766., p_lon_ecl=0., p_lat_ecl=0., p_ra=0, p_dec=0., sun=[0, 0, 0], moon=[0, 0, 0], sat=[0, 0, 0], sat_vel=[0, 0, 0], isInDeep=0, areaW=0.0, areaD=0.0, isInSunSide=1, exp_time=150., trans_Angle=1.):



        line = self.fmt % (

            jdTime, p_lon_ecl, p_lat_ecl, p_ra, p_dec, sun[0], sun[1], sun[2], moon[0], moon[1], moon[2], sat[0], sat[1], sat[2], sat_vel[0], sat_vel[1], sat_vel[2], isInDeep)

            jdTime, p_lon_ecl, p_lat_ecl, p_ra, p_dec, sun[0], sun[1], sun[2], moon[0], moon[1], moon[2], sat[0], sat[1], sat[2], sat_vel[0], sat_vel[1], sat_vel[2], isInDeep, areaW, areaD, isInSunSide, exp_time, trans_Angle)

        # if not line.endswith("\n"):

        line += "\n"

        self.outInfo.write(line)

survey_sim/constraints/_utils.py

+2 −2
Original line number Diff line number Diff line
@@ -167,7 +167,7 @@ def dot_product_matrix(A, B): 
        if cols_A != rows_B:

            raise ValueError("矩阵 A 的列数必须等于矩阵 B 的行数")

        # 初始化结果矩阵 C，所有元素为 0

        C = [[0 for _ in range(cols_B)] for _ in range(rows_A)]

        C = [[0. for _ in range(cols_B)] for _ in range(rows_A)]



        # 执行矩阵点乘

        for i in range(rows_A):
@@ -183,7 +183,7 @@ def dot_product_matrix(A, B): 
        if cols_A != rows_B:

            raise ValueError("矩阵 A 的列数必须等于矩阵 B 的行数")

        # 初始化结果矩阵 C，所有元素为 0

        C = [0 for _ in range(rows_A)]

        C = [0. for _ in range(rows_A)]



        # 执行矩阵点乘

        for i in range(rows_A):

survey_sim/strategies/SurveySim_MSC.py

+35 −17
Original line number Diff line number Diff line
@@ -2,7 +2,7 @@ 
Author: Zhang Xin zhangx@bao.ac.cn

Date: 2024-09-19 15:53:27

LastEditors: Zhang Xin zhangx@bao.ac.cn

LastEditTime: 2024-12-10 21:31:50

LastEditTime: 2025-01-17 10:27:57

FilePath: /CSST_Survey/survey_sim/strategies/SurveySim_MSC.py

Description: 这是默认设置,请设置`customMade`, 打开koroFileHeader查看配置 进行设置: https://github.com/OBKoro1/koro1FileHeader/wiki/%E9%85%8D%E7%BD%AE

'''
@@ -27,6 +27,8 @@ from survey_sim.skyField import skyMap 
from survey_sim.skyField.csstSkyCoord import csstSkyCoordGeocentricEquator



from survey_sim.config.infoOutput import InfoOutput



from survey_sim.surveyStat.obsAreaTools import initAreaStatLib, calculateAccArea, freeMemForAreaStat

import time

import os

import pickle
@@ -70,6 +72,8 @@ class SurveySim(object): 
        self.outInfo = InfoOutput(

            dir=dir, logger_filename=logName, info_out_filename=outFn)



        self.areaCalLib = initAreaStatLib(hplayer=11)



    def initSkyMap(self, surveyCons=None):



        skydat_path = resource_filename("survey_sim", "strategies/data")
@@ -105,11 +109,10 @@ class SurveySim(object): 
            smF = open(skyMap_init_data, "wb")

            pickle.dump(self.skymap, smF)



    def do_observe_plan(self, jump_time=60):

    def do_observe_plan(self, jump_time_norm=60, jump_time_unObs=120):

        self.betaCons.get_survey_time_segment_MSC()

        curPSky = None

        orbitDataId = 0

        jump_time = jump_time

        pNum = 0

        for start_end_t in self.betaCons.MSC_time:

            sec_times = np.arange(
@@ -138,16 +141,19 @@ class SurveySim(object): 
            curTime = start_end_t[0]

            while curTime <= start_end_t[1]:

                print("Time:", curTime)

                curTime, obsTime, curPSky, sun, moon, sat, satV, orbitDataId, _ = self.searchOnePointing(

                    curTime=curTime, curPSky=curPSky, orbitDataId=orbitDataId, jump_time=jump_time, endTime=start_end_t[1], satOrbit_interp=[sat_xi, sat_yi, sat_zi], subPoint_interp=[subPoint_xi, subPoint_yi])

                curTime, obsTime, curPSky, sun, moon, sat, satV, orbitDataId, isInSunSide, exp_time, trans_Angle = self.searchOnePointing(

                    curTime=curTime, curPSky=curPSky, orbitDataId=orbitDataId, jump_time_norm=jump_time_norm, jump_time_unObs=jump_time_unObs, endTime=start_end_t[1], satOrbit_interp=[sat_xi, sat_yi, sat_zi], subPoint_interp=[subPoint_xi, subPoint_yi])

                if sun is not None:

                    print("obsTime-------------------:", obsTime)

                    areaWide, areaDeep = calculateAccArea(

                        curPSky.ecl_lat, curPSky.ecl_lon, self.areaCalLib)

                    self.outInfo.outInfo_add_obj(jdTime=obsTime, p_lon_ecl=curPSky.ecl_lon, p_lat_ecl=curPSky.ecl_lat, p_ra=curPSky.ra,

                                                 p_dec=curPSky.dec, sun=sun, moon=moon, sat=sat, sat_vel=satV, isInDeep=curPSky.inDeepFlag)

                                                 p_dec=curPSky.dec, sun=sun, moon=moon, sat=sat, sat_vel=satV, isInDeep=curPSky.inDeepFlag, areaW=areaWide, areaD=areaDeep, isInSunSide=isInSunSide, exp_time=exp_time, trans_Angle=trans_Angle)

                    pNum = pNum + 1

                if pNum > 20:

                if pNum > 3:

                    break



            self.outInfo.close_outFile()

        self.areaCalLib.freeObject()

    # 搜索一个是否满足的节点

    # 返回 returnTime,otime, self.skymap.skymap[final_index], sun_r1, moon_r1, sat_r1, satVel_r1, orbitDataId_r1

    #    下一次开始时间，开始曝光时间，返回指向的天区（或None），太阳位置（或None），月亮位置（或None），卫星位置（或None），卫星速度（或None），轨道开始搜索的位置
@@ -158,7 +164,8 @@ class SurveySim(object): 
        curTime=2459766.0,

        curPSky=None,

        orbitDataId=0,

        jump_time=60,

        jump_time_norm=60,

        jump_time_unObs=120,

        endTime=2460766.0,

        satOrbit_interp=None,

        subPoint_interp=None
@@ -205,6 +212,7 @@ class SurveySim(object): 
            curTime), subPoint_interp[1](curTime)])

        # print(underStarPoint - underStarPoint1)

        battery_q = self.bettery_level

        jump_time = jump_time_norm

        while saa_constraint.IsInSAA(uderStarPoint=underStarPoint):

            lastSAATime = curTime

            curTime = curTime + jump_time / 86400.0
@@ -238,7 +246,8 @@ class SurveySim(object): 
            battery_q = (

                battery_q

                + sun_front_time * 86400.0 *

                energy_constraint.getPlaneEnergy(1.0, yr)

                energy_constraint.getPlaneEnergy(

                    1.0, yr, sConstraint=self.surveyCons)

                - jump_time * self.surveyCons.POWER_CONSUMPTION

            )
@@ -298,6 +307,8 @@ class SurveySim(object): 
        curTime_yr = (curTime - self.surveyCons.CSST_START_TIME) / 365.25

        # 将当前的儒略日转换为【0~10】中的年份



        unObs_continue_flag = 0



        if isInSunSide == 1:

            # skynum = len(self.skymap.skymap_id)

            for s_id in self.skymap.skymap_id:
@@ -747,6 +758,12 @@ class SurveySim(object): 
        returnTime = 0



        if final_min_value >= MAX_VALUE:

            if unObs_continue_flag == 0:

                jump_time = jump_time_norm

            else:

                jump_time = jump_time_unObs



            unObs_continue_flag = unObs_continue_flag + 1

            returnTime = curTime + jump_time / 86400.0



            e_t1 = 0  # 1 orbit, sun time left
@@ -769,8 +786,9 @@ class SurveySim(object): 
                battery_q = self.surveyCons.BATTERY_LOW

                print("ERROR:-----Energy not balance")

            self.bettery_level = battery_q

            return returnTime, returnTime, curPSky, None, None, None, None, orbitDataId, isInSunSide

            return returnTime, returnTime, curPSky, None, None, None, None, orbitDataId, isInSunSide, None, None

        else:

            unObs_continue_flag = 0

            returnTime = curTime + (final_use_time) / 86400.0

            self.skymap.skymap[final_index].flag = self.skymap.skymap[final_index].flag + 1

            otime = curTime+(final_use_time-final_exp_time) / \
@@ -804,7 +822,7 @@ class SurveySim(object): 


            # gc.collect()



            return returnTime, otime, self.skymap.skymap[final_index], sun_r1, moon_r1, sat_r1, satVel_r1, orbitDataId_r1, isInSunSide

            return returnTime, otime, self.skymap.skymap[final_index], sun_r1, moon_r1, sat_r1, satVel_r1, orbitDataId_r1, isInSunSide, final_exp_time, final_Angle



        # locateSun(infp, otime, sun_r1);

        # locateMoon(infp, otime, moon_r1);
@@ -1096,7 +1114,7 @@ class SurveySim(object): 
    def calculateExposureTime(self, obseveSkyUnit=None, p_to_zenith_cosV=0.5):

        exTime = self.surveyCons.EXPTIME  # 曝光时间

        extraTime = 5  # s 宽场额外增加5s，深场额外增加12s再天顶30-40°直接（两边70-80°之间）

        if obseveSkyUnit.inDeepFlag >= 0:  # 深场

        if obseveSkyUnit.inDeepFlag > 0:  # 深场

            exTime = self.surveyCons.EXPTIME_DEEP

            extraTime = 12

        if p_to_zenith_cosV <= 1:  # 指向地球亮边
@@ -1110,10 +1128,10 @@ class SurveySim(object): 
        return exTime





ss = SurveySim(time_start=2459768, time_end=2459775.0,

               bettery_level=97200., cmgCons=None, dir="/Users/zhangxin/Work/SurveyPlan/CSST_Survey", logName="MSC_survey_20241108_log", outFn="MSC_survey_20241108_log_outInfo")

ss = SurveySim(time_start=2461958., time_end=2461959.,

               bettery_level=97200., cmgCons=None, dir="/Users/zhangxin/Work/SurveyPlan/CSST_Survey", logName="MSC_survey_20250113_log.dat", outFn="MSC_survey_20250113_log_outInfo.dat")



t1 = time.time()

tt = ss.do_observe_plan(jump_time=300)

tt = ss.do_observe_plan(jump_time_norm=60, jump_time_unObs=150)

t2 = time.time()

print("ust time %f minute" % ((t1-t2)/60))
 
print("ust time %f minute" % ((t2-t1)/60))

survey_sim/surveyStat/.DS_Store

0 → 100644
+10 KiB

File added.

No diff preview for this file type.

View file