Skip to content
GitLab
Switch branch/tag
Find file Normal viewHistoryPermalink
shao_test.py 3.3 KB
Newer Older
Fang Yuedong's avatar
add orbit position inputs, add astrometric modeling  
Fang Yuedong committed
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75 
import os
import sys

from ObservationSim.Astrometry.Astrometry_util import on_orbit_obs_position

def readFits2List(fn):
    from astropy.io import fits
    if not os.path.exists(fn):
        print("Can not find "+fn); return False
    
    hdul=fits.open(fn)
    data = hdul[1].data
    ra_list = data['RA'].tolist()
    dec_list = data['Dec'].tolist()
    pmra_list = data['pmra'].tolist()
    pmdec_list = data['pmdec'].tolist()
    parallax_list = data['parallax'].tolist()
    rv_list = [0.0 for i in range(len(ra_list))]
    hdul.close()
    return ra_list, dec_list, pmra_list, pmdec_list, rv_list, parallax_list
    
def usageExample(targets_fn, results_fn, x, y, z, vx, vy, vz, date_str, time_str):
    if not os.path.exists(targets_fn):
        print("Can not find " + targets_fn)
        sys.exit(1)
    ra_list, dec_list, pmra_list, pmdec_list, rv_list, parallax_list = readFits2List(targets_fn)
    output_ra_list, output_dec_list = on_orbit_obs_position(ra_list, dec_list, pmra_list, \
                                      pmdec_list, rv_list, parallax_list, len(ra_list), \
                                      x, y, z, vx, vy, vz, "J2015.5", date_str, time_str, lib_path='./libshao.so')
    f = open(results_fn, "w")
    ll = list()
    ll.append("n,date,time,ra,dec,pm_ra,pm_dec,rv,parallax,obs_ra,obs_dec\n")
    for i in range(len(output_ra_list)):
        ra_str = str(ra_list[i]) 
        dec_str = str(dec_list[i]) 
        pmra_str = str(pmra_list[i])
        pmdec_str = str(pmdec_list[i])
        rv_str = str(rv_list[i])
        parallax_str = str(parallax_list[i])
        ora_str = str(output_ra_list[i])
        odec_str = str(output_dec_list[i])      
        l = str(i)+date_str+","+time_str+","+ra_str+","+dec_str+","+pmra_str+","+pmdec_str+","+rv_str+","+parallax_str+","+ora_str+","+odec_str+"\n"
        ll.append(l)
    f.writelines(ll)
    f.close()
    print("Process finished. Results save to "+results_fn)
    
def usageTips():
    print("Usage 1: python3 ./shao_test.py")
    print("Usage 2: python3 ./shao_test.py yyyy-MM-dd HH:mm:ss.ms Positon1_KM Positon2_KM Positon3_KM Velocity1_MK/S Velocity2_MK/S Velocity3_MK/S $PATH1/MMW_Gaia_Cluster_D20_SS_astrometry.fits $PATH2/results.csv")
    print("Caution: Do no include space in path; Unit is KM or KM/S")    
    print("Example: python3 ./shao_test.py 2025-03-05 22:20:15.12 2347.766100 5132.421392 3726.591334 5.282357 4.644825 -3.074722 ./MMW_Gaia_Cluster_D20_SS_astrometry.fits ./results.csv")
    
if __name__ == "__main__":
    args_value = sys.argv
    if len(args_value) == 1:
        usageExample("./MMW_Gaia_Cluster_D20_SS_astrometry.fits", "./results.csv", 2347.766100, 5132.421392, 3726.591334, 5.282357, 4.644825, -3.074722, "2025-03-05", "22:20:15.12"); sys.exit(1)
    elif len(args_value) == 11:
        date_str = args_value[1]
        time_str = args_value[2]
        x = float(args_value[3])
        y = float(args_value[4])
        z = float(args_value[5])
        vx = float(args_value[6])
        vy = float(args_value[7])
        vz = float(args_value[8])
        targets_fn = args_value[9]
        results_fn = args_value[10]
        if not os.path.exists(targets_fn):
            print("Can not find " + targets_fn)
            sys.exit(1)
        usageExample(targets_fn, results_fn, x, y, z, vx, vy, vz, date_str, time_str)
    else:
        usageTips()