Notebook

これは日々の作業を通して学んだことや毎日の生活で気づいたことをを記録しておく備忘録である。

HTML ファイル生成日時: 2024/11/21 17:40:55.112 (台灣標準時)

sbpy を使って観測できる小惑星を探す方法

sbpy という Python のパッケージがある ことを知ったので、試してみたでござる。日時と観測地点を指定し、 JPL Horizons を使って小惑星の位置を調べ、高度と等級の条件を満たす小惑星が 見つかったら、その小惑星についての情報を書き出す、というプログラムを書 いてみたでござる。

書いてみたプログラムは以下の通りでござる。


#!/usr/pkg/bin/python3.10

#
# Time-stamp: <2023/06/23 21:51:33 (CST) daisuke>
#

# importing argparse module
import argparse

# importing sbpy module
import sbpy.data

# importing astropy module
import astropy.time
import astropy.units

# initialising a parser of argparse
parser = argparse.ArgumentParser (description='Finding observable asteroids')

# adding arguments
parser.add_argument ('-t1', '--time1', default='2000-01-01T00:00', \
                     help='start time in UT (default: 2000-01-01T00:00)')
parser.add_argument ('-t2', '--time2', default='2000-01-02T00:00', \
                     help='end time in UT (default: 2000-01-02T00:00)')
parser.add_argument ('-d', '--dt', type=float, default=1.0, \
                     help='time step in hour for ephemeris (default: 1.0)')
parser.add_argument ('-n', '--number', type=int, default=10, \
                     help='number of asteroids for search (default: 10)')
parser.add_argument ('-a', '--alt', type=float, default=30.0, \
                     help='acceptable altitude (default: 30.0 deg)')
parser.add_argument ('-s', '--site', default='D35', \
                     help='observatory site (default: D35)')
parser.add_argument ('-m', '--mag', type=float, default=18.0, \
                     help='limiting magnitude in V-band (default: 20.0)')
parser.add_argument ('-o', '--output', default='asteroids.list', \
                     help='output file name (default: asteroids.list)')

# parsing arguments
args = parser.parse_args ()

# given arguments
time_start  = args.time1
time_end    = args.time2
step_hr     = args.dt
n_asteroids = args.number
req_alt     = args.alt
obs_site    = args.site
limmag_V    = args.mag
file_output = args.output

# units
u_hr  = astropy.units.hour
u_deg = astropy.units.deg
u_mag = astropy.units.mag

# bar
bar1 = '-' * 71
bar2 = '=' * 78

# date/time in UT in astropy.time format
datetime_start_str = time_start
datetime_stop_str  = time_end
step               = step_hr * u_hr
datetime_start     = astropy.time.Time (datetime_start_str, scale='utc')
datetime_stop      = astropy.time.Time (datetime_stop_str, scale='utc')
datetime           = {'start': datetime_start, \
                      'stop':  datetime_stop, \
                      'step':  step}

# opening output file for writing
with open (file_output, 'w') as fh_output:
    # input parameters
    fh_output.write (f'# time_start  = {time_start}\n')
    fh_output.write (f'# time_end    = {time_end}\n')
    fh_output.write (f'# time step   = {step_hr}\n')
    fh_output.write (f'# n_asteroids = {n_asteroids}\n')
    fh_output.write (f'# req_alt     = {req_alt}\n')
    fh_output.write (f'# obs_site    = {obs_site}\n')
    fh_output.write (f'# limmag_V    = {limmag_V}\n')
    fh_output.write (f'# file_output = {file_output}\n')
    fh_output.write (f'{bar2}\n')
    
    # ephemeris
    for i in range (1, n_asteroids + 1):
        eph = sbpy.data.Ephem.from_horizons (i, \
                                             location=obs_site, \
                                             epochs=datetime)

        # a flag for visibility
        observable = 1

        # check of visibility
        for j in range (len (eph)):
            if (eph["EL"][j].value < req_alt):
                observable = 0

        # check of apparent V-band magnitude
        for j in range (len (eph)):
            if (eph["V"][j].value > limmag_V):
                observable = 0

        # printing ephemeris
        if (observable):
            fh_output.write (f'{eph["targetname"][0]}\n')
            fh_output.write (f'  Date/Time                  RA   Dec')
            fh_output.write (f'    Az    El     R Delta alpha V-mag\n')
            fh_output.write (f'  {bar1}\n')
            for j in range (len (eph)):
                time_jd    = eph["epoch"][j]
                time_isot  = time_jd.to_value ('isot')
                fh_output.write (f'  {time_isot}')
                fh_output.write (f' {eph["RA"][j].value:5.1f}')
                fh_output.write (f' {eph["DEC"][j].value:+5.1f}')
                fh_output.write (f' {eph["AZ"][j].value:5.1f}')
                fh_output.write (f' {eph["EL"][j].value:+5.1f}')
                fh_output.write (f' {eph["r"][j].value:5.1f}')
                fh_output.write (f' {eph["delta"][j].value:5.1f}')
                fh_output.write (f' {eph["alpha"][j].value:5.1f}')
                fh_output.write (f' {eph["V"][j].value:5.2f}\n')
            fh_output.write (f'{bar2}\n')

使い方は、以下の通りでござる。


% ./sbpy_search_asteroids.py -h
usage: sbpy_search_asteroids.py [-h] [-t1 TIME1] [-t2 TIME2] [-d DT]
                                [-n NUMBER] [-a ALT] [-s SITE] [-m MAG]
                                [-o OUTPUT]

Finding observable asteroids

options:
  -h, --help            show this help message and exit
  -t1 TIME1, --time1 TIME1
                        start time in UT (default: 2000-01-01T00:00)
  -t2 TIME2, --time2 TIME2
                        end time in UT (default: 2000-01-02T00:00)
  -d DT, --dt DT        time step in hour for ephemeris (default: 1.0)
  -n NUMBER, --number NUMBER
                        number of asteroids for search (default: 10)
  -a ALT, --alt ALT     acceptable altitude (default: 30.0 deg)
  -s SITE, --site SITE  observatory site (default: D35)
  -m MAG, --mag MAG     limiting magnitude in V-band (default: 20.0)
  -o OUTPUT, --output OUTPUT
                        output file name (default: asteroids.list)

実行結果の例は、以下の通りでござる。


% ./sbpy_search_asteroids.py -t1 2023-07-01T14:00:00 -t2 2023-07-01T18:00:00 \
? -d 1 -n 500 -a 35 -s D35 -m 17 -o asteroids_20230701.list
% ls -lF asteroids_20230701.list 
-rw-r--r--  1 daisuke  taiwan  4005 Jun 23 21:58 asteroids_20230701.list
% cat asteroids_20230701.list
# time_start  = 2023-07-01T14:00:00
# time_end    = 2023-07-01T18:00:00
# time step   = 1.0
# n_asteroids = 500
# req_alt     = 35.0
# obs_site    = D35
# limmag_V    = 17.0
# file_output = asteroids_20230701.list
==============================================================================
34 Circe (A855 GA)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 271.7 -14.9 148.8 +46.2   2.8   1.8   4.1 12.57
  2023-07-01T15:00:00.000 271.7 -14.9 169.6 +51.1   2.8   1.8   4.1 12.57
  2023-07-01T16:00:00.000 271.7 -14.9 192.8 +50.8   2.8   1.8   4.1 12.57
  2023-07-01T17:00:00.000 271.6 -14.9 213.1 +45.4   2.8   1.8   4.1 12.57
  2023-07-01T18:00:00.000 271.6 -14.9 228.2 +36.3   2.8   1.8   4.1 12.57
==============================================================================
153 Hilda (A875 VC)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 270.9 -15.5 150.2 +45.9   3.4   2.4   3.4 12.58
  2023-07-01T15:00:00.000 270.9 -15.5 170.9 +50.6   3.4   2.4   3.4 12.58
  2023-07-01T16:00:00.000 270.9 -15.5 193.6 +50.0   3.4   2.4   3.4 12.58
  2023-07-01T17:00:00.000 270.9 -15.5 213.4 +44.4   3.4   2.4   3.4 12.58
  2023-07-01T18:00:00.000 270.9 -15.5 228.2 +35.4   3.4   2.4   3.4 12.58
==============================================================================
242 Kriemhild (A884 SA)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 277.9  -7.4 135.1 +48.8   3.2   2.2   5.0 14.04
  2023-07-01T15:00:00.000 277.9  -7.4 156.1 +56.6   3.2   2.2   5.0 14.04
  2023-07-01T16:00:00.000 277.9  -7.4 184.0 +59.1   3.2   2.2   5.0 14.04
  2023-07-01T17:00:00.000 277.9  -7.4 210.5 +54.9   3.2   2.2   5.0 14.04
  2023-07-01T18:00:00.000 277.9  -7.4 229.4 +45.9   3.2   2.2   5.0 14.04
==============================================================================
246 Asporina (A885 EA)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 278.6  +0.6 125.4 +54.0   2.4   1.5   9.8 11.86
  2023-07-01T15:00:00.000 278.6  +0.6 148.1 +63.6   2.4   1.5   9.8 11.86
  2023-07-01T16:00:00.000 278.5  +0.6 183.8 +67.1   2.4   1.5   9.8 11.86
  2023-07-01T17:00:00.000 278.5  +0.6 217.4 +62.1   2.4   1.5   9.7 11.86
  2023-07-01T18:00:00.000 278.5  +0.6 237.8 +51.8   2.4   1.5   9.7 11.86
==============================================================================
393 Lampetia (A894 VC)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 285.3  +6.7 110.6 +52.1   1.9   0.9  15.7 10.53
  2023-07-01T15:00:00.000 285.3  +6.7 127.4 +64.3   1.9   0.9  15.7 10.53
  2023-07-01T16:00:00.000 285.3  +6.7 162.4 +72.5   1.9   0.9  15.7 10.53
  2023-07-01T17:00:00.000 285.3  +6.7 210.9 +70.8   1.9   0.9  15.7 10.53
  2023-07-01T18:00:00.000 285.3  +6.7 238.9 +60.8   1.9   0.9  15.7 10.53
==============================================================================
417 Suevia (A896 JB)
  Date/Time                  RA   Dec    Az    El     R Delta alpha V-mag
  -----------------------------------------------------------------------
  2023-07-01T14:00:00.000 280.5 -12.7 137.2 +43.0   2.8   1.8   3.8 13.39
  2023-07-01T15:00:00.000 280.5 -12.7 155.5 +50.8   2.8   1.8   3.8 13.39
  2023-07-01T16:00:00.000 280.5 -12.7 179.1 +53.8   2.8   1.8   3.8 13.39
  2023-07-01T17:00:00.000 280.5 -12.7 203.0 +51.1   2.8   1.8   3.8 13.39
  2023-07-01T18:00:00.000 280.5 -12.7 221.8 +43.7   2.8   1.8   3.8 13.39
==============================================================================



Frequently accessed files

  1. Computer___Python/20220518_0.html
  2. Computer___Network/20230726_00.html
  3. Misc___Taiwan/20240207_00.html
  4. Computer___Network/20230516_00.html
  5. Computer___FreeBSD/20220621_0.html
  6. Computer___Python/20220715_0.html
  7. Computer___Network/20230508_00.html
  8. Food___Taiwan/20220429_0.html
  9. Computer___NetBSD/20220817_3.html
  10. Computer___Python/20220410_0.html
  11. Computer___Network/20240416_00.html
  12. Computer___Network/20240130_00.html
  13. Computer___Debian/20210223_1.html
  14. Computer___NetBSD/20230119_00.html
  15. Computer___Python/20210124_0.html
  16. Computer___Python/20221013_0.html
  17. Computer___NetBSD/20220818_1.html
  18. Computer___NetBSD/20220428_0.html
  19. Science___Math/20220420_0.html
  20. Computer___NetBSD/20240101_02.html
  21. Computer___NetBSD/20220808_0.html
  22. Computer___TeX/20230503_00.html
  23. Computer___NetBSD/20230515_00.html
  24. Science___Astronomy/20220503_0.html
  25. Computer___NetBSD/20210127_0.html
  26. Computer___Python/20240101_00.html
  27. Computer___Network/20220413_1.html
  28. Computer___Python/20220816_1.html
  29. Computer___NetBSD/20210204_0.html
  30. Travel___Taiwan/20220809_2.html


HTML file generated by Kinoshita Daisuke.