genSurveyDFSByPointing.py

'''
Author: Zhang Xin zhangx@bao.ac.cn
Date: 2025-04-21 09:53:28
LastEditors: Zhang Xin zhangx@bao.ac.cn
LastEditTime: 2025-09-03 16:56:04
FilePath: /CSST_Survey/Users/zhangxin/Work/SurveyPlan/point/pointing_c9/genSurveyDFSByPointing.py
Description: 这是默认设置,请设置`customMade`, 打开koroFileHeader查看配置 进行设置: https://github.com/OBKoro1/koro1FileHeader/wiki/%E9%85%8D%E7%BD%AE
'''
import json
from datetime import datetime
from astropy.table import Table, join  # type: ignore
import numpy as np  # type: ignore
from astropy.time import Time
from datetime import datetime, timedelta
import pytz
import sys


def calculateTransTime(transAngle=1.0, shuterTime=1.5):
    angleVStime = np.array([[1, 20, 45, 180], [80, 127, 196, 581]])

    if transAngle < angleVStime[0, 0]:
        tTime = 70
    elif transAngle == angleVStime[0, 0]:
        tTime = angleVStime[1, 0]
    else:
        for i in np.arange(1, 4, 1):
            if (transAngle > angleVStime[0, i-1] and transAngle <= angleVStime[0, i]):
                tTime = angleVStime[1, i-1] * ((transAngle - angleVStime[0, i])) / (((angleVStime[0, i-1] - angleVStime[0, i]))) + \
                    angleVStime[1, i] * ((transAngle - angleVStime[0, i-1])) / \
                    (((angleVStime[0, i] - angleVStime[0, i-1])))
                break

    return tTime + shuterTime * 2.0


def JDTimeToBeijingTime(jd_time=2459766.):
    # 使用 astropy 转换到UTC时间
    t = Time(jd_time, format='jd', scale='utc')
    utc_time = t.datetime  # 转换为 datetime 对象
    # 转换为北京时间（UTC+8）
    beijing_time = utc_time + timedelta(hours=8)
    # 使用 pytz 添加时区信息（可选）
    beijing_tz = pytz.timezone('Asia/Shanghai')
    beijing_time = beijing_tz.localize(beijing_time)
    # 格式化输出
    formatted_time = beijing_time.strftime("%Y-%m-%d %H:%M:%S")
    return formatted_time, beijing_time


def dataInit():
    # 定义数据字典
    data = {
        "id": 1593,
        "opid": "101000001593",
        "obs_id": "10100001593",
        "instrument": "MSC",
        "obs_group": "None",
        "proposal_id": "None",
        # "project_id": "None",
        # "file_type": "SCI",

        "obs_type": "WIDE",
        "object": "m31",
        "mode": "OM1",
        "dataset": "None",
        "starttime": "2031-04-23 16:40:19",
        "endtime": "2031-04-23 16:43:37",
        # "duration": 150.0,
        "exptime": 150.0,
        "params": {
            "t0": "2031-04-23 16:41:04",
            "T1": 150.0,
            "bias_voltage": 0,
            "r_rate": 0,
            "h_rate": 0,
            "expstart": "2031-04-23 16:41:04",
            "nframes": 1,
            "cover": 0,
            "gainlvl": "",
            "shutter": "F",
            "submode": "OM1-1",
            "ledflag": "",
            "ledstat": "",
            "ledt01": 0.0,
            "ledt02": 0.0,
            "ledt03": 0.0,
            "ledt04": 0.0,
            "ledt05": 0.0,
            "ledt06": 0.0,
            "ledt07": 0.0,
            "ledt08": 0.0,
            "ledt09": 0.0,
            "ledt10": 0.0,
            "ledt11": 0.0,
            "ledt12": 0.0,
            "ledt13": 0.0,
            "ledt14": 0.0,
            "ra_obj": 146.4567,
            "dec_obj": 54.0354,
            "pa_obj": 25.0,
            "guide_window": 1,
            "guide_duration": 10,
            "guide_pos": [1.0, 2.0],
        },
        "simresult": {

            "jd": 2459766.56586422,
            "dec": 18.9566929133858,
            "ra": 55.716621253406,
            "sat_x": 2071.4244979379,
            "sat_y": 5303.28628499587,
            "sat_z": 3639.38736139007,
            "sun_x": -36018853.1282829,
            "sun_y": 135579783.188155,
            "sun_z": 58773491.2520049,
            "moon_x": -388299.686481602,
            "moon_y": -18964.2980441153,
            "moon_z": 20511.371245116,
            "deep_area": 0.0,
            "very_deep_area": 0.0,
            "ca": -1.0,
            "cb": 1.0,
            "exptime": 150.0,
            "turn_angle": 1.9229313453104,
            "is_in_sun": 1,
            "cmg": 0.0,
            "charge": 97200.0,
            "between_angle": 16.2831689171475,
            "saa_sec": 0.0,
            "skyid": 309072,
            "board_angle": 25.0,
            "quad_0": 0.999859206659046,
            "quad_i": 0.00671061788580071,
            "quad_j": -0.00768959137244524,
            "quad_k": -0.0133193208638005,
            "sky_type": 1.0,
            "id": 0,
            "sat_vx": -7.104108,
            "sat_vy": 2.912855,
            "sat_vz": 0.332436,
            "fov_lt_dec": 20.4566929133858,
            "fov_lt_ra": 47.6784104513669,
            "fov_lb_dec": 17.4566929133858,
            "fov_lb_ra": 47.6784104513669,
            "fov_rb_dec": 17.4566929133858,
            "fov_rb_ra": 47.9427837789041,
            "fov_rt_dec": 20.4566929133858,
            "fov_rt_ra": 47.9427837789041
        }
    }
    return data


if __name__ == "__main__":
    if len(sys.argv) < 4:
        print('input_file_name  output_file_name datasec')
    obs_group = "None"
    if len(sys.argv) == 5:
        obs_group = sys.argv[4]
    subMod = "OM1-1"
    if len(sys.argv) == 6:
        obs_group = sys.argv[4]
        subMod = sys.argv[5]
    fn1 = 'E17.5_b17.5_beta_11.6_opt_transtime_1_CMG_1_dp_2_0.25_da_10_Texp_1.5_DEC60_500_0.1_800_1000_+5deg.obsid.csv'
    # fn2 = 'pointing_50_5_n.dat'
    fn2 = sys.argv[1]
    outFn = sys.argv[2]
    t1 = Table.read(fn1, format='ascii')
    t2 = Table.read(fn2, format='ascii')
    t1['obs_id'] = (t1['obs_id']//1e11*1e8+t1['obs_id'] % 1e8).astype(np.int64)

    t2.rename_column('id', 'obs_id')
    inputData = join(t1, t2, keys='obs_id')

    # dfn = 'pointing_50_5_n_obsid.dat'
    # inputData = Table.read(dfn, format='ascii')

    all_data = []
    for i in np.arange(len(inputData)):
        d1 = inputData[i]
        data = dataInit()
        data['id'] = int(d1['id'])
        data['obs_id'] = d1['obs_id']
        data['opid'] = (str(d1['obs_id'])[0:3]+'0'+str(d1['obs_id'])[3:])
        data['instrument'] = 'MSC'
        data['obs_group'] = obs_group
        # data['file_type'] = 'SCI'
        if d1['sky_type'] == 2:
            data['obs_type'] = 'DEEP'
        else:
            data['obs_type'] = 'WIDE'
        data['object'] = d1['skyid']
        data['dataset'] = sys.argv[3]
        data['starttime'] = JDTimeToBeijingTime(
            d1['jd'])[1].isoformat()
        data['endtime'] = JDTimeToBeijingTime(
            d1['jd']+d1['exposure']/86400.)[1].isoformat()
        data['exptime'] = d1['exposure']

        # param
        data['params']['ra_obj'] = d1['ra_2']
        data['params']['dec_obj'] = d1['dec_2']
        data['params']['pa_obj'] = d1['pos_angle']
        data['params']['T1'] = d1['exposure']
        data['params']['t0'] = JDTimeToBeijingTime(
            d1['jd'])[1].isoformat()
        data['params']['expstart'] = data['params']['t0']
        data['params']['cover'] = 1
        data['params']['shutter'] = "T"
        data['params']['submode'] = subMod

        # simresult
        data['simresult']['jd'] = d1['jd']
        data['simresult']['dec'] = d1['decH']
        data['simresult']['ra'] = d1['raH']
        data['simresult']['sat_x'] = d1['sat_x_1']
        data['simresult']['sat_y'] = d1['sat_y_1']
        data['simresult']['sat_z'] = d1['sat_z_1']
        data['simresult']['sun_x'] = d1['sun_x_1']
        data['simresult']['sun_y'] = d1['sun_y_1']
        data['simresult']['sun_z'] = d1['sun_z_1']
        data['simresult']['moon_x'] = d1['moon_x_1']
        data['simresult']['moon_y'] = d1['moon_y_1']
        data['simresult']['moon_z'] = d1['moon_z_1']
        data['simresult']['deep_area'] = d1['deepArea']
        data['simresult']['very_deep_area'] = d1['veryDeepArea']
        data['simresult']['ca'] = float(d1['ca'])
        data['simresult']['cb'] = float(d1['cb'])
        data['simresult']['exptime'] = d1['exposure']
        data['simresult']['turn_angle'] = d1['turnAngle']
        data['simresult']['is_in_sun'] = d1['isInSun']
        data['simresult']['cmg'] = d1['cmg']
        data['simresult']['charge'] = d1['charge']
        data['simresult']['skyid'] = int(d1['skyid'])
        data['simresult']['board_angle'] = d1['boardAngle']
        data['simresult']['quad_0'] = d1['quad_0']
        data['simresult']['quad_i'] = d1['quad_i']
        data['simresult']['quad_j'] = d1['quad_j']
        data['simresult']['quad_k'] = d1['quad_k']
        data['simresult']['sky_type'] = float(d1['sky_type'])
        data['simresult']['id'] = int(d1['id'])
        data['simresult']['sat_vx'] = d1['sat_vx']
        data['simresult']['sat_vy'] = d1['sat_vy']
        data['simresult']['sat_vz'] = d1['sat_vz']

        data['simresult']['fov_lt_dec'] = 0.
        data['simresult']['fov_lt_ra'] = 0.
        data['simresult']['fov_lb_dec'] = 0.
        data['simresult']['fov_lb_ra'] = 0.
        data['simresult']['fov_rb_dec'] = 0.
        data['simresult']['fov_rb_ra'] = 0.
        data['simresult']['fov_rt_dec'] = 0.
        data['simresult']['fov_rt_ra'] = 0.

        all_data.append(data)
        # break

    json_data = json.dumps(all_data, indent=4, default=str)

    # 打印JSON格式的字符串
    # print(json_data)

    with open(outFn, "w") as file:
        file.write(json_data)