Impor semua paket yang dibutuhkan
from skyfield import api
from skyfield import almanac
from skyfield.nutationlib import iau2000b
from skyfield.units import Angle
from datetime import timedelta
import astropy.units as u
from datetime import datetime
from pytz import timezone
Buat timescale dan unduh ephemeris. Kita akan menggunakan ephemeris dari JPL NASA versi de421.bsp.
timescale = api.load.timescale()
ephemeris = api.load('de421.bsp')
Tentukan koodinat di permukaan Bumi di mana prediksi akan dilakukan, sebagai contoh kita gunakan Observatorium Bosscha. Koodinat yang digunakan dalam derajat.
# Observatorium Bosscha 107:37:00 BT; -6:49:50 LS
longitude = 107 + 37/60
latitude = -(6 + 49/60 + 50/3600)
Perhitungan posisi Bulan dan Matahari bergantung pada posisi pengamat di Bumi. Oleh karena itu kita perlu menentukan lokasi topografi Skyfield untuk koordinat di atas.
bosschaTopo = api.Topos(longitude=longitude, latitude=latitude)
Waktu Matahari dan Bulan Terbenam
Misalkan ingin menghitung waktu terbenam antara tanggal 5 April dan 6 April 2019.
t1 = timescale.utc(2019, 4, 5, 0)
t2 = timescale.utc(2019, 4, 6, 0)
Paket Skyfield telah menyediakan fungsi untuk menghitung waktu terbenam Matahari. Definisi matahari terbenam yang digunakan paket Skyfield sesuai dengan definisi Matahari terbenam dari USNO, yaitu ketika pusat piringan Matahari berada pada jarak zenit $90^o 50′$.
sunriset, sunBol = almanac.find_discrete(t1, t2, almanac.sunrise_sunset(ephemeris, bosschaTopo))
Hasilnya adalah sebuah array waktu, dan array boolean yang akan bernilai True jika matahari terbit dan bernilai False jika matahari terbenam.
Untuk menghitung waktu terbenam Bulan dengan definisi di atas kita dapat membuat fungsi sebagai berikut.
def moonrise_moonset(ephemeris, topos):
moon = ephemeris['moon']
topos_at = (ephemeris['earth'] + topos).at
def is_moon_up_at(t):
t._nutation_angles = iau2000b(t.tt)
return topos_at(t).observe(moon).apparent().altaz()[0].degrees > -50/60
is_moon_up_at.rough_period = 0.5
return is_moon_up_at
moonriset, moonBol = almanac.find_discrete(t1, t2, moonrise_moonset(ephemeris, bosschaTopo))
Keluarannya akan mirip dengan format waktu terbenam Matahari.
Waktu terbenam dalam Julian Date
sunriset[~sunBol]
<Time tt=[2458578.95413447]>
moonriset[~moonBol]
<Time tt=[2458578.96095373]>
Waktu terbenam dalam standar UTC
sunriset[~sunBol].utc_iso()
['2019-04-05T10:52:48Z']
moonriset[~moonBol].utc_iso()
['2019-04-05T11:02:37Z']
Posisi Matahari dan Bulan saat Matahari Terbenam
Definisikan Bumi, Bulan, dan Matahari dari ephemeris di atas.
earth = ephemeris['earth']
sun = ephemeris['sun']
moon = ephemeris['moon']
Lokasi pengamat adalah posisi pusat Bumi ditambahkan posisi lokal topografi.
bosscha = earth + bosschaTopo
Koordinat RA dan Dec saat Matahari terbenam
sunradec = bosscha.at(sunriset[~sunBol]).observe(sun).apparent()
moonradec = bosscha.at(sunriset[~sunBol]).observe(moon).apparent()
RA dan Dec epoch saat itu.
print(sunradec.radec(epoch='date'))
(<Angle 00h 56m 36.75s>, <Angle +06deg 03' 05.6">, <Distance [1.00036989] au>)
print(moonradec.radec(epoch='date'))
(<Angle 01h 04m 12.04s>, <Angle +01deg 56' 37.8">, <Distance [0.00265776] au>)
Koordinat dalam Altitud dan Azimut
print(sunradec.altaz())
(<Angle -00deg 49' 59.9">, <Angle 275deg 59' 42.8">, <Distance [1.00036989] au>)
print(moonradec.altaz())
(<Angle 01deg 32' 08.7">, <Angle 272deg 08' 33.2">, <Distance [0.00265776] au>)
Elongasi
sunradec.separation_from(moonradec)
<Angle 04deg 31' 21.0">
Fraksi Iluminasi Bulan saat Matahari terbenam
almanac.fraction_illuminated(ephemeris, 'moon', sunriset[~sunBol])[0]*100
0.19602608000200417
Waktu konjungsi
Waktu konjungsi dapat ditentukan dengan mencari waktu saat fase bulan bernilai 0.
t0 = timescale.utc(2019, 4, 4, 0)
t1 = timescale.utc(2019, 4, 6, 0)
t, y = almanac.find_discrete(t0, t1, almanac.moon_phases(ephemeris))
t[y==0].utc_iso()
['2019-04-05T08:50:29Z']
Fungsi
Untuk memudahkan menjalankan seluruh prosedur di atas, kita dapat membuat fungsi yang mencakup seluruh prosedur sekaligus membuat plot posisi Bulan dan Matahari.
%matplotlib inline
from matplotlib import pyplot as plt
plt.style.use('ggplot')
longitude = 107+37/60
latitude = -(6+49/60+50/3600)
bosschaTopo = api.Topos(longitude=longitude, latitude=latitude)
def ucupkiyah(date, plot=False, day=0):
wib = timezone('Asia/Jakarta')
date = wib.localize(date)
t0 = timescale.utc(date + timedelta(days=day))
t1 = timescale.utc(t0.astimezone(wib) - timedelta(days=28))
t2 = timescale.utc(t0.astimezone(wib) + timedelta(days=28))
t, y = almanac.find_discrete(t1, t2, almanac.moon_phases(ephemeris))
# t0 = t[y==0]
t1 = timescale.utc(t0.astimezone(wib) - timedelta(days=0))
t2 = timescale.utc(t0.astimezone(wib) + timedelta(days=1))
sunriset, sunBol = almanac.find_discrete(t1, t2, almanac.sunrise_sunset(ephemeris, bosschaTopo))
moonriset, moonBol = almanac.find_discrete(t1, t2, moonrise_moonset(ephemeris, bosschaTopo))
sunset = sunriset[~sunBol].astimezone(wib)[0]
moonset = moonriset[~moonBol].astimezone(wib)[0]
moonAge = sunriset[~sunBol] - t[y==0]
moonAge = moonAge[0]*24*3600
lagTime = (moonriset[~moonBol][0]-sunriset[~sunBol][0])*24*3600
earth = ephemeris['earth']
sun = ephemeris['sun']
moon = ephemeris['moon']
bosscha = earth + bosschaTopo
sunradec = bosscha.at(sunriset[~sunBol]).observe(sun).apparent()
moonradec = bosscha.at(sunriset[~sunBol]).observe(moon).apparent()
sunra, sundec, _ = sunradec.radec(epoch='date')
moonra, moondec, _ = moonradec.radec(epoch='date')
sunal, sunaz, _ = sunradec.altaz()
moonal, moonaz, _ = moonradec.altaz()
relal = moonal.degrees[0] - sunal.degrees[0]
relaz = moonaz.degrees[0] - sunaz.degrees[0]
elongation = sunradec.separation_from(moonradec)
illumination = almanac.fraction_illuminated(ephemeris, 'moon', sunriset[~sunBol])[0]*100
print('By the Name of Allah')
print('Ucupkiyah Calculator by Muhammad Yusuf')
print('')
print('G. Conjunction Time\t : {}'.format(t[y==0].astimezone(wib)[0]))
print('Sunset Time \t \t : {}'.format(sunset))
print('Moonset Time \t \t : {}'.format(moonset))
print('')
print('G. Moon Age\t \t : {:02.0f}h {:02.0f}m {:02.2f}s'.format(moonAge//3600, moonAge%3600//60, moonAge%60))
print('Lag Time \t \t : {:02.0f}h {:02.0f}m {:02.2f}s'.format(lagTime//3600, lagTime%3600//60, lagTime%60))
print('')
print('Sun RA \t \t \t : {}'.format(sunra))
print('Moon RA \t \t : {}'.format(moonra))
print('')
print('Sun Dec \t \t : {}'.format(sundec))
print('Moon Dec \t \t : {}'.format(moondec))
print('')
print('Sun Azimuth \t \t : {}'.format(sunaz))
print('Moon Azimuth \t \t : {}'.format(moonaz))
print('')
print('Sun Altitude \t \t : {}'.format(sunal))
print('Moon Altitude \t \t : {}'.format(moonal))
print('')
print('Relative Altitude \t : {}'.format(Angle(degrees=relal).dstr()))
print('Relative Azimuth \t : {}'.format(Angle(degrees=relaz).dstr()))
print('')
print('Moon Elongation \t : {}'.format(elongation))
print('Moon Illumination \t : {:02.2f}%'.format(illumination))
if plot:
lingkaranMatahari = plt.Circle((sunaz.to(u.deg)[0].value, sunal.to(u.deg)[0].value),
(0.5*u.deg).value,
fc='yellow')
lingkaranBulan = plt.Circle((moonaz.to(u.deg)[0].value, moonal.to(u.deg)[0].value),
(0.5*u.deg).value,
fc='black', alpha=.5)
plt.figure(figsize=(15, 15))
ax = plt.subplot(aspect=1)
ax.add_patch(lingkaranMatahari)
ax.add_patch(lingkaranBulan)
if sunaz.to(u.deg)[0] <= moonaz.to(u.deg)[0]:
xlimlow = sunaz.to(u.deg)[0].value
xlimhigh = moonaz.to(u.deg)[0].value
else:
xlimlow = moonaz.to(u.deg)[0].value
xlimhigh = sunaz.to(u.deg)[0].value
if sunal.to(u.deg)[0] <= moonal.to(u.deg)[0]:
ylimlow = sunal.to(u.deg)[0].value
ylimhigh = moonal.to(u.deg)[0].value
else:
ylimlow = moonal.to(u.deg)[0].value
ylimhigh = sunal.to(u.deg)[0].value
plt.xlim(xlimlow-2, xlimhigh+2)
plt.ylim(ylimlow-2, ylimhigh+2)
plt.text(sunaz.to(u.deg)[0].value-0.28, sunal.to(u.deg)[0].value-0.05, 'Matahari', fontsize=15)
plt.text(moonaz.to(u.deg)[0].value-0.2, moonal.to(u.deg)[0].value-0.05, 'Bulan', fontsize=15, color='white')
plt.axhline(y=0)
plt.plot([moonaz.to(u.deg)[0].value, sunaz.to(u.deg)[0].value],
[moonal.to(u.deg)[0].value, moonal.to(u.deg)[0].value],
label='$\Delta$ Azimut = {}'.format(Angle(degrees=relaz).dstr()))
plt.plot([sunaz.to(u.deg)[0].value, sunaz.to(u.deg)[0].value],
[moonal.to(u.deg)[0].value, sunal.to(u.deg)[0].value],
label='$\Delta$ Altitud = {}'.format(Angle(degrees=relal).dstr()))
plt.plot([moonaz.to(u.deg)[0].value, sunaz.to(u.deg)[0].value],
[moonal.to(u.deg)[0].value, sunal.to(u.deg)[0].value],
label='Elongasi = {}'.format(elongation))
plt.xlabel('Azimut (derajat)')
plt.ylabel('Altitud (derajat)')
plt.legend()
plt.title(sunset)
Prediksi Hilal
Syaban 1440H
ucupkiyah(datetime(2019, 4, 5), True)
By the Name of Allah
Ucupkiyah Calculator by Muhammad Yusuf
G. Conjunction Time : 2019-04-05 15:50:28.624000+07:00
Sunset Time : 2019-04-05 17:52:48.034000+07:00
Moonset Time : 2019-04-05 18:02:37.218000+07:00
G. Moon Age : 02h 02m 19.41s
Lag Time : 00h 09m 49.18s
Sun RA : 00h 56m 36.75s
Moon RA : 01h 04m 12.04s
Sun Dec : +06deg 03' 05.6"
Moon Dec : +01deg 56' 37.8"
Sun Azimuth : 275deg 59' 42.8"
Moon Azimuth : 272deg 08' 33.2"
Sun Altitude : -00deg 49' 59.9"
Moon Altitude : 01deg 32' 08.7"
Relative Altitude : 02deg 22' 08.6"
Relative Azimuth : -03deg 51' 09.6"
Moon Elongation : 04deg 31' 21.0"
Moon Illumination : 0.20%
Ramadhan 1440H
ucupkiyah(datetime(2019, 5, 5), True, 0)
By the Name of Allah
Ucupkiyah Calculator by Muhammad Yusuf
G. Conjunction Time : 2019-05-05 05:45:29.812000+07:00
Sunset Time : 2019-05-05 17:41:46.344000+07:00
Moonset Time : 2019-05-05 18:08:49.400000+07:00
G. Moon Age : 11h 56m 16.53s
Lag Time : 00h 27m 3.06s
Sun RA : 02h 48m 47.59s
Moon RA : 03h 13m 50.25s
Sun Dec : +16deg 14' 05.4"
Moon Dec : +13deg 52' 39.4"
Sun Azimuth : 286deg 15' 06.9"
Moon Azimuth : 284deg 43' 23.0"
Sun Altitude : -00deg 49' 59.9"
Moon Altitude : 05deg 28' 23.6"
Relative Altitude : 06deg 18' 23.5"
Relative Azimuth : -01deg 31' 43.9"
Moon Elongation : 06deg 29' 19.4"
Moon Illumination : 0.42%
Syawal 1440H
ucupkiyah(datetime(2019, 6, 3), True, 0)
By the Name of Allah
Ucupkiyah Calculator by Muhammad Yusuf
G. Conjunction Time : 2019-06-03 17:01:56.654000+07:00
Sunset Time : 2019-06-03 17:40:01.059000+07:00
Moonset Time : 2019-06-03 17:40:28.707000+07:00
G. Moon Age : 00h 38m 4.40s
Lag Time : 00h 00m 27.65s
Sun RA : 04h 44m 30.26s
Moon RA : 04h 43m 28.73s
Sun Dec : +22deg 18' 09.4"
Moon Dec : +19deg 30' 21.5"
Sun Azimuth : 292deg 21' 54.8"
Moon Azimuth : 289deg 33' 36.1"
Sun Altitude : -00deg 49' 59.9"
Moon Altitude : -00deg 43' 46.0"
Relative Altitude : 00deg 06' 13.9"
Relative Azimuth : -02deg 48' 18.8"
Moon Elongation : 02deg 48' 24.7"
Moon Illumination : 0.07%